# Electronic properties of twisted multilayer graphene

V. Hung Nguyen1, Trinh X. Hoang2 and J.-C. Charlier1 1Institute of Condensed
Matter and Nanosciences, Université catholique de Louvain (UCLouvain), B-1348
Louvain-la-Neuve, Belgium
2Institute of Physics, Vietnam Academy of Science and Technology, Ba Dinh,
Hanoi 11108, Vietnam

###### Abstract

Twisted bilayer graphene displays many fascinating properties that can be
tuned by varying the relative angle (also called twist angle) between its
monolayers. As a remarkable feature, both the electronic flat bands and the
corresponding strong electron localization have been obtained at a specific
“magic” angle ($\sim 1.1^{\circ}$), leading to the observation of several
strongly correlated electronic phenomena. Such a discovery has hence inspired
the creation of a novel research field called twistronics, i.e., aiming to
explore novel physical properties in vertically stacked 2D structures when
tuning the twist angle between the related layers. In this paper, a
comprehensive and systematic study related to the electronic properties of
twisted multilayer graphene (TMG) is presented based on atomistic
calculations. The dependence of both the global and the local electronic
quantities on the twist angle and on the stacking configuration are analyzed,
fully taking into account atomic reconstruction effects. Consequently, the
correlation between structural and electronic properties are clarified,
thereby highlighting the shared characteristics and differences between
various TMG systems as well as providing a comprehensive and essential
overview. On the basis of these investigations, possibilities to tune the
electronic properties are discussed, allowing for further developments in the
field of twistronics.

##  I Introduction

Twisted multilayer graphene (TMG) refers to structures of vertically stacked
graphene layers when one or a few layers are rotated relative to the others by
arbitrary angles wu2014 ; chen2016 ; mogera2017 ; shen2020 ; liu2020 ;
park2021 ; zhang2021 ; chen2021e ; zeyu2021 ; maria2021 ; jeong2021 , thus
creating the so-called moiré superlattices. The first member of TMG systems,
twisted bilayer graphene (TBLG), has recently attracted great interest from
the scientific community since it exhibits fascinating properties, leading to
the exploration of several novel features as surprising as superconductivity
that can not be observed in conventional graphene systems andrei2020 . In
particular and most remarkably, electronic properties of TBLG can be tuned by
varying the twist angle $\theta$ wang2019 ; carr2020 , especially leading to
the observation of localized flat bands near the Fermi level at a critical
angle ($\simeq 1.1^{\circ}$) called the magic angle rafi2011 ; utama2021 . At
this critical angle, many-body effects are therefore significantly enlarged
choi2018 ; kerelsky2019 ; choi2021 ; gadelha2021 ; gadelha2021b ; barbosa2022
and charge carriers do not possess enough kinetic energy to escape from their
strong mutual interactions, thus forming novel strongly correlated electronic
states. Such a strong electronic correlation in magic-angle TBLG is the key
ingredient for observing several novel phenomena such as superconductivity,
correlated insulating states, magnetism, and even quantized anomalous Hall
effects cao2018a ; cao2018b ; yankowitz2019 ; xie2019 ; lu2019 ; choi2019 ;
jiang2019 ; aaron2019 ; serlin2020 ; leon2020 .

The observations mentioned above have been considered as real scientific
breakthroughs, creating a nascent field of twistronics carr2017 . Indeed,
after the discovery of magic-angle TBLG, the research of other twisted systems
has become an emerging topic. In addition to TMG, the explored systems include
also twisted superlattices of other 2D layered materials such as graphene on
hexagonal boron nitride, graphene on transition-metal dichalcogenide layers,
van der Waals moiré superlattices of transition-metal dichalcogenide layers,
and so on (e.g., see Refs. liu2014 ; kang2017 ; naik2018 ; xian2019 ;
yaping2020 ; kennes2021 ; xu2021 ; tao2022  as well as the recent review
liu2022  and references therein). These explorations have also inspired
researchers to extendedly apply the twist approach to other physical systems,
e.g., magnetic 2D materials tong2018 ; wang2020 ; tiancheng2021 ; luo2021  and
photonic layered structures wang2020p ; mao2021 ; zhang2021p ; chen2021 . In
overall, the properties of those systems are definitely dependent of their
constituting layers and importantly, can be tuned by varying their twist
angle, similar to those observed in TBLG. In particular, in the TMG cases,
electronic flat bands and accordingly, the strongly correlated electronic
phenomena have been similarly observed around their magic angle (e.g., see in
Refs. khalaf2019 ; haddadi2020 ; jeong2021 ; shen2020 ; youngju2020 ;
chen2021e ). Moreover, more tunable possibilities can be observed in TMG
bandwidth ; yankowitz2019 ; liu2020 ; shen2020 ; chen2021e ; zeyu2021 ;
park2021 ; rickhaus2021 , due to a larger number of layers and the plenty of
stacking configurations georgios2020 , compared to TBLG. In this context, a
comprehensive understanding of the electronic properties of TMG is highly
desirable. Some new TMGs zeyu2021 ; park2021 ; jeong2021  have been very
recently produced experimentally, and hence their corresponding electronic
properties are needed to be more exhaustively clarified. In addition, a good
understanding of those twisted systems could be also helpful for research of
twisted structures based on other 2D materials.

On the theoretical side, though a large number of works devoted to the
electronic properties of some TMG has been reported (e.g., in carr2017 ;
carr2019 ; khalaf2019 ; carr2020 ; tritsaris2020 ; haddadi2020 ; shin2021 ), a
systematical study on this subject is however still missing. Remarkably, it
has been shown that atomic reconstruction (i.e., the structure relaxation)
significantly occurs in TBLG at low angles ($\lesssim 1.1^{\circ}$), thus
strongly influencing their electronic properties, i.e., especially the
formation of electronic flat bands as well as the electronic localization
pictures gargiulo2017 ; gadelha2021 ; nguyen2021 . The structural and
electronic properties of TBLG also exhibit a strong correlation nguyen2021 ,
that explains essentially the corresponding electronic localization features.
Fully taking into account the atomic reconstruction effects is therefore
mandatory to compute accurately the electronic properties of TMG, which has
not been properly considered yet in all works previously published in the
literature.

This paper aims to investigate systematically the most common and essential
electronic properties of TMG by atomistic calculations, which fully take the
atomic reconstruction into account. First, twisted bilayer graphene is
investigated to illustrate the overall structural and electronic properties
induced by the twisting effects. Some of these properties are commonly
observed in other multilayer structures. Afterward, the study focuses on
several typical twisted multilayer (more than 2 graphene layers) systems
around their magic angles, at which electronic flat bands and strong electron
localization are also observed. The correlation between the obtained
electronic and structural properties is also clarified. Electronic (both
global and local) quantities obtained depending on the twist angle and the
stacking configuration are analyzed in details, thus highlighting the shared
characteristics and differences between TMGs. Finally, on the basis of these
investigations, possibilities to tune these electronic properties are
discussed.

##  II Modeling Methodology

Figure 1: Twisted multilayer graphene systems investigated in this work. SM-
structures are twisted systems created by putting and rotating a graphene
monolayer on top of a multilayer (n layers) graphene that is either in AB
stacked or in ABC stacked. These structures are hereafter abbreviated as
SM-(AB)n and SM-(ABC)n, respectively. AS-structures are alternatively twisted
systems of n graphene monolayers, abbreviated as ASn. DB is the abbreviation
of twisted double bilayer graphene superlattices. At last, AMn structures (n =
1,2,3,4,5) are alternatively twisted systems containing, at least, one Bernal
stacking graphene bilayer.

As mentioned earlier, the atomic reconstruction is significant and hence
taking into account these structural effects is mandatory to model accurately
the electronic properties of twisted graphene systems, especially, at small
twist angles around and below the magic one. In this work, we used classical
force-fields to compute the atomic structure relaxation. More specifically,
intralayer forces are computed using the optimized Tersoff and Brenner
potentials lindsay2010 , while interlayer forces are modeled using the
Kolmogorov-Crespi potentials kolmogorov2005 ; leven2016 . The atomic structure
was optimized until all force components are smaller than 0.5 meV/atom and the
electronic calculations were then performed. The atomic structure relaxation
obtained by this method is shown to be in good agreement with experiments
brihuega2012 ; kerelsky2019 ; yoo2019 ; zhang2020  as well as the results
obtained by other (both quantum and semi-classical) methods uchida2014 ;
gargiulo2017 ; cantele2020 ; lamparski2020 , except that the intralayer
equilibrium C-C bond distance is slightly overestimated by $\sim 1\%$
leven2016 .

Similarly to previous works gadelha2021 ; nguyen2021 , the electronic
properties of TMG (schematized in Fig.1) were then computed using the $p_{z}$
tight-binding (TB) Hamiltonians after the structural optimization has been
fully achieved. Hopping energies $t_{nm}$ between carbon C-sites are
determined by the standard Slater-Koster formula

$t_{nm}=\cos^{2}\phi_{nm}V_{pp\sigma}(r_{nm})+\sin^{2}\phi_{nm}V_{pp\pi}(r_{nm})$
(1)

where the direction cosine of $\vec{r}_{nm}=\vec{r}_{m}-\vec{r}_{n}$ along Oz
axis is $\cos\phi_{nm}=z_{nm}/r_{nm}$. The distance-dependent Slater-Koster
parameters are determined as trambly12

$\displaystyle
V_{pp\pi}(r_{nm})=V_{pp\pi}^{0}\exp\left[q_{\pi}\left(1-\frac{r_{nm}}{a_{0}}\right)\right]F_{c}(r_{nm})$
$\displaystyle
V_{pp\sigma}(r_{nm})=V_{pp\sigma}^{0}\exp\left[q_{\sigma}\left(1-\frac{r_{nm}}{d_{0}}\right)\right]F_{c}(r_{nm})$

with a smooth cutoff function
$F_{c}(r_{nm})=\left[1+\exp\left(\frac{r_{nm}-r_{c}}{\lambda_{c}}\right)\right]^{-1}$.

In this work, the TB Hamiltonian was further simplified from the model
considered in nguyen2021 , i.e., only couplings limited by $r_{inplane}\leq
2.16$ Å were taken into account. This approximation presents an advantage that
the matrix of TB Hamiltonian is very sparse so as to make calculations
feasible for small twist angles and/or multilayer systems in which the number
of C atoms could be huge. Within this approximation, the TB parameters
nguyen2021  are slightly adjusted as follows:

$\displaystyle
V_{pp\pi}^{0}=-2.7\,~{}\mathrm{eV},\,\,\,V_{pp\sigma}^{0}=367.5\,~{}\mathrm{meV},$
$\displaystyle\frac{q_{\pi}}{a_{0}}=\frac{q_{\sigma}}{d_{0}}=22.18\,~{}\mathrm{nm}^{-1},$
$\displaystyle
a_{0}=0.1439\,~{}\mathrm{nm},\,\,\,d_{0}=0.33\,~{}\mathrm{nm},\,\,\,r_{c}=0.614\,~{}\mathrm{nm},\,\,\,\lambda_{c}=0.0265\,~{}\mathrm{nm}$

Indeed, such an adjusted TB Hamiltonian fairly accurately reproduces the
electronic properties of TBLG as obtained in nguyen2021  as well as in
experiments, e.g., in cao2018a ; gadelha2021 ; barbosa2022 . As a clear
illustration, the almost flat bands near the Fermi level is perfectly
reproduced in the relaxed TBLG in the vicinity of the magic angle (i.e., $\sim
1.1^{\circ}$, see Fig. 3 below) as experimentally observed (e.g., see in Ref.
cao2018a ).

##  III Twisted bilayer graphene

Figure 2: Stacking structure of TBLG imaged by the spatial variation of
stacking vector t (as defined in the inset of bottom-left panel) when varying
the twist angle $\theta$. Four top images are 2D real-space representation of
$|$t$|$ while the bottom-left one is a plot of $|$t$|$ extracted along the
line connecting BA - SP (soliton line) - AB - AA stacking regions. The bottom-
right diagram is a simple description of the stacking properties of TBLG in
the different $\theta$-regimes (from large to small twist-angles).

As the simplest TMG structure, TBLG exhibits most clearly essential electronic
properties induced by the twisting effects. These properties as well as their
correlation with the structural properties, that generally also appear in
other TMGs (as examples, see in Refs. nguyen2021 ; haddadi2020  as well as
illustration in Figs.14-15 in section V), will be therefore revised and
systematized in this section.

Typical images illustrating the evolution of the stacking structure in TBLG
when varying the twist angle $\theta$ are presented in Fig.2. The stacking
vector t, as illustrated in the bottom-left panel, indicates the smallest
displacement applied locally to one layer to recover the AA stacking at the
considered position. In general, TBLG superlattice includes three typical
stacking regions: AA stacking, AB/BA stacking regions and soliton lines
between AB/BA ones (called SP for $saddle$ $point$). The presented images
illustrate that the magic angle $\sim 1.1^{\circ}$ separates TBLG into two
regimes: $\theta>1.1^{\circ}$ (large angles) and $\theta<1.1^{\circ}$ (small
ones) nguyen2021 . At large angles, the moiré superlattice evolves smoothly
and both AA and AB/BA stacking regions are continuously enlarged when
decreasing the twist angle. However, AA stacking regions get saturation when
$\theta$ decreasingly approaches $1.1^{\circ}$. At small angles, whereas the
size of AA stacking regions is unchanged, AB/BA regions are continuously
enlarged when decreasing $\theta$, according to the enlargement of moiré
cells. Thus, in addition to the critical angle (i.e., magic angle) $\sim
1.1^{\circ}$, there are two TBLG classes, separated by such angle and
exhibiting distinct stacking structures, as summarized in the bottom-right
diagram of Fig.2.

Figure 3: Electronic bandstructure of TBLG for the twist angle
$\theta<1.1^{\circ}$ (left in blue), $\theta\simeq 1.1^{\circ}$ (center in
red) and $\theta>1.1^{\circ}$ (right in black). The corresponding density of
states (DOS) are presented on the extreme right with the same color scale for
the three moiré structures.

The structural properties discussed above crucially govern the electronic
features of TBLG, as illustrated by the typical electronic structures in Fig.3
as well as the local densities of states (LDOS) images in Fig.4. In large-
angle TBLG, Dirac fermions are still preserved at low energies but a
remarkable renormalization of their Fermi velocity is obtained yin2015 , as
seen in Fig.3 for $\theta\simeq 2.28^{\circ}$. In addition, saddle points
emerge at the crossing of Dirac cones (i.e., at the M-point in the
superlattice Brillouin region), yielding van Hove singularities in the density
of states (DOS) at energies that are lowered down when decreasing $\theta$.
When $\theta$ approaches $1.1^{\circ}$, two (electron and hole) van Hove
singularities merge at zero energy and accordingly the low energy bands are
flattened, resulting in the observation of an extremely high DOS peak at that
energy (see Fig.3 for $\theta\simeq 1.08^{\circ}$). The angle $1.1^{\circ}$
has been thus called the “magic angle” rafi2011 . When further decreasing
$\theta$, the electronic bands become more dispersive again (see Fig.3 for
$\theta\simeq 0.76^{\circ}$). Remarkably, it has been shown that fully taking
the structure relaxation into account, the electronic flat band (accordingly,
the localized DOS peak at zero energy) as at $1.1^{\circ}$ is no longer
observed for small angles $<1.1^{\circ}$ as discussed in details in Ref.
nguyen2021 . Moreover, even though it is not present in the global electronic
quantities, the strong electronic localization is still observed locally in AA
stacking regions, which has been experimentally confirmed in Ref. gadelha2021
.

Figure 4: Real-space LDOS images in TBLG at different twist angles. Energy is
fixed where the strongest electronic localization is obtained, i.e., at the
van Hove singularity point in DOS (see Fig.3) for $\theta>1.1^{\circ}$ and at
zero energy for $\theta\lesssim 1.1^{\circ}$.

These distinct electronic properties of TBLG at large and small angles are
visually illustrated by the LDOS pictures in Fig.4 as well as their
correlation with the structural properties discussed in Fig.2. Actually, the
LDOS in Fig.4 are represented in each case at the strongest localization
point. First, the strong electronic localization is observed in the AA
stacking region for all twist angles. Remarkably, the real-space electronic
localization exhibits the similar $\theta$-dependence, compared to that of
stacking structure presented in Fig.2. Indeed, while it is continuously
enlarged when reducing $\theta$ above $1.1^{\circ}$, the real-space LDOS peak
obtained in AA stacking regions is maximum at $\theta\simeq 1.1^{\circ}$ and
saturates in the regime $<1.1^{\circ}$. These results demonstrate clearly a
strong correlation between the structural properties and the electronic
localization features in TBLG.

More remarkably, the observation of electronic flat bands is found to be
essentially related to the maximum electronic localization in AA stacking
regions and their maximum contribution to the global electronic properties of
the system (see also discussions in section V). These features are
concurrently obtained at $\theta\simeq 1.1^{\circ}$. To clarify such
observation, we notice additionally here that the strongest electronic
localization in energy is favorable (unfavorable) in AA stacking regions
(other stacking regions) nguyen2021 . At large angles $>1.1^{\circ}$, the
electronic localization in AA stacking regions (see Fig.4 for $\theta\simeq
4.41^{\circ}$) and consequently the global localization in energy (see the
extreme right panel of Fig.3) are relatively weak. When decreasing $\theta$ in
this regime, the electronic localization in AA stacking regions and their
contribution are continuously enhanced. Accordingly, global localization
(i.e., DOS) peaks increase and their positions concurrently move towards zero
energy. At $\theta\simeq 1.1^{\circ}$, all these effects are maximized and
consequently the electronic flat bands are obtained. When further decreasing
$\theta$ in the regime below $1.1^{\circ}$, even though a strong electronic
localization is still locally observed in the saturated AA stacking region,
the increasing contribution of AB/BA stacking ones gradually diminishes the
global electronic localization and hence flat bands can not be observed again
nguyen2021 , thus precluding the emergence of other (i.e., small) magic angles
predicted by other theoretical calculations rafi2011 .

Finally, it is worth noting that even though they could be quantitatively
affected by the presence of a larger number of layers and/or by a more
complicatedly stacked configuration, the structural and electronic features
similar to those presented above are still observed in other TMGs, simply
because all these features are essentially governed by twisting effects.

##  IV Twisted multilayer graphene

In this section, the electronic properties of TMGs schematized in Fig.1 are
investigated. The creation of these systems could be simply considered as
adding one or a few graphene layers to a TBLG one. Since these TMG share a
similar $\theta$-dependence (e.g., see Refs. haddadi2020 ; nguyen2021 ) with
TBLG, the present study is mainly focused on the magic-angle cases where low-
energy bands are most flat and thus strongly correlated electronic phenomena
can be observed.

Figure 5: Electronic bandstructures of different SM-systems at their magic
angle. On top, the results in SM-(AB)2 (i.e., twisted monolayer-bilayer)
structure are presented. Other systems including SM-(AB)n and SM-(ABC)n (see
the definition in Fig.1) at their magic angle $\sim 1.2^{\circ}$ are computed
and presented in the second and third rows, respectively, with n = 3,4,5,6
(from left to right). On top-right, DOS versus energy in SM-(AB)2 structure is
also superimposed.

We first consider TMG obtained when one rotated graphene monolayer is placed
on top of a graphene multilayer, which is either in an AB-stacked form or an
ABC-stacked one (i.e., SM-systems in Fig.1). Fig.5 presents the electronic
bandstructures of different SM-systems in both two cases. First of all, the
SM-(AB)2 (i.e., twisted monolayer-bilayer) system at its magic angle exhibits
a bandstructure (see on top of Fig.5) with isolated flat bands, similar to
that obtained in the magic-angle TBLG. The main difference between two cases
is that these bands (accordingly, zero-energy DOS peak) are observed to be
less flat (lower) than those of TBLG. Note that continuum Hamiltonian
calculations youngju2020  have shown that the bandwidth of flat bands in
SM-(AB)2 system is generally narrower than that obtained in magic-angle TBLG.
This is indeed also confirmed by our calculations for systems at the same
twist angle in the regime $\theta>>1.1^{\circ}$. The results in the magic
angle case however exhibit an opposite behavior, i.e., the bandwidth of the
magic-angle

Figure 6: Total LDOS and layer decomposed LDOS in SM-(AB)2 system in Fig.5 at
zero energy.

SM-(AB)2 system is larger than that obtained in magic-angle TBLG. This is
actually consistent with the results presented in Refs. choi2021 ; bandwidth
for atomic reconstructed superlattices. These results could be understood as a
consequence of the fact, as discussed in the previous section, that the band
flatness has a strong correlation with the real-space electronic localization
that is basically stronger in AA-stacking regions nguyen2021 . In SM-(AB)2
systems, the corresponding regions are obtained in an AAB-stacked form, i.e.,
a mixture of AA- and AB-stacking configurations while the latter is
unfavorable for the electronic localization. Such mixture diminishes the
electronic localization in the SM-(AB)2 system. This is visibly illustrated by
comparing the total LDOS in Fig.6 with that obtained in magic-angle TBLG (see
Fig.4 for $\theta\simeq 1.08^{\circ}$) and is further demonstrated in Fig.14
in section V, thus explaining essentially the less flat bands obtained in the
SM-(AB)2 case.

In addition, another interesting feature is found in Fig.6 that flat bands in
the magic-angle SM-(AB)2 system exhibit a coexistence of layer-resolved
localized-delocalized behaviors, i.e., they are spatially localized on the
twisted side whereas delocalized on the bilayer side. This unique feature has
been also experimentally confirmed in a recent work tong2021  and the
delocalized correlated electronic states have been similarly explored in
twisted double bilayer graphene zhang2021b . Note that electronic states
delocalized either in the graphene plane or in different layers are generally
obtained in twisted systems consisting of more than 2 layers as the SM-(AB)2
system here and others considered below, where certain layers are strongly
influenced by the moiré superlattice interactions whereas the weak effects
take place in other layers. Moreover, the observed coexistence of localized-
delocalized electronic states essentially has a tight relationship with the
weak real-space localization in the SM-(AB)2 system discussed above.

Figure 7: Layer decomposed DOS (LD-DOS) in some SM-systems considered in
Fig.5.

When the number of graphene layers of the multilayer counterpart increases,
both systems SM-(AB)n and SM-(ABC)n exhibit clearly distinct electronic
properties (see Fig.5). In particular, while the electronic structure of
SM-(AB)n systems presents both flat bands and dispersive ones at low energies,
flat bands in SM-(AB)n systems are still clearly isolated from high-energy
bands for $n\leq 5$. However, the energy gaps between the mentioned bands in
the latter case is reduced when increasing n and is closed for $n>5$. Such
distinct properties could be explained by the difference in electronic
properties of AB- and ABC-stacked graphene multilayer systems, i.e., while
ABC-stacked systems host a single band that is partially flat around zero
energy, there is a mixture of several (increasing with n) parabolic and linear
bands in the AB-stacked ones koshino2010 .

To further clarify the difference between these two systems, we analyze their
layer decomposed DOS (LD-DOS) in function of energy in Fig.7. For the SM-(AB)2
system, besides the coexistence of layer-resolved localized-delocalized states
presented in Fig.6, LD-DOS still exhibits the usual behaviours similar to that
observed in magic-angle TBLG. However, other distinct features are found for
SM-systems with a larger number of graphene layers. In particular, the
electron localization is mostly observed in three top layers of SM-(AB)n
structures where the moiré superlattice interactions take place. However, the
localization is obtained in both top and bottom outermost layers of SM-(ABC)3
and SM-(ABC)4 systems (see illustration in Fig.7 for SM-(ABC)4 one). Moreover,
when the number of layers n further increases, such electron localization in
SM-(ABC)n surprisingly changes, i.e., the effect in the top outermost layers
is diminished while it is enhanced in the bottom ones. When $n$ is large
enough, the electron localization is mostly observed in the bottom outermost
layers as illustrated in Fig.7 for the SM-(ABC)6 system.

Figure 8: Electronic bandstructures of alternatively twisted systems of
graphene monolayers schematized on top-left (i.e., the ASn systems specified
in Fig.1) with $\theta$ near their magic angle $\theta_{n}^{MA}$ (see text).
On top-right, DOS versus energy in AS2-structure is also superimposed.

Next, TMGs created by alternatively rotating ($0,\theta,0,\theta,...$)
graphene monolayers (i.e., ASn structures in Fig.1) are considered in Fig.8.
The magic angle of these ASn structures can be estimated by the simple
expression $\theta_{n}^{MA}=2\theta_{n}^{MA}\cos(\pi/(n+1)))$ with
$\theta_{2}^{MA}\simeq 1.1^{\circ}$ as proposed in Refs. khalaf2019 ;
jeong2021 . This is indeed demonstrated in Fig.8 where (almost) flat bands are
obtained in ASn structures near $\theta_{n}^{MA}$. Different from the results
suggested by the effective continuum models khalaf2019 ; park2021 ; zeyu2021 ;
jeong2021 , the TB bandstructure here presents four degenerate flat bands
(more visibly seen in section V when a vertical electric field is applied) and
n-2 additional linear bands at low energies. Moreover, saddle points
(accordingly, van Hove singularities in DOS) occur in dispersive bands for
n$\geqslant$4\. The low-energy bands in these ASn systems are also shown to be
generally flatter than those obtained for the SM-structures presented in
Fig.5, thus confirming again our discussion above about the property of
electron localization in the presence of Bernal stacking bilayer graphene.
Note that another difference between systems without and with Bernal stacking
bilayer graphene has been also discussed in Ref. choi2021 . Indeed, in this
study, significant sublattice polarization has been observed in the latter
case, leading to weaker electron-phonon couplings than those obtained in
twisted systems of graphene monolayers only.

Figure 9: Layer decomposed DOS (LD-DOS) in some AS-systems considered in
Fig.8. On bottom-left, a zoom-combined picture of the results in AS5 structure
is presented.

To understand more deeply their electronic properties, the layer decomposed
DOS obtained in these ASn structures is presented in Fig.9. Remarkably, the
strong electronic localization takes place in their centered layers, in
contrast to those observed in SM-ones presented in Fig.7. Moreover, LD-DOS in
the centered layers presents behaviors that are very similar to those obtained
in magic-angle TBLG. In particular, the strongly localized DOS peak is likely
isolated from high-energy electronic states by finite energy gaps.
Contrastingly, LD-DOS in other layers (especially, the outermost ones) exhibit
some signatures of dispersive bands, i.e., the characteristics of linear bands
and van Hove singularities induced by the twist as visibly shown in the zoom-
image on the bottom of Fig.9. Thus, the flat bands and dispersive ones of ASn
structures exhibit oppositely spatial properties, i.e., whereas the former is
localized in the system center, the latter is observed in outer layers. This
implies that there is aslo a coexistence of localized-delocalized electronic
states in the ASn systems of large n.

Figure 10: Electronic bandstructures of ASn systems when one or a few layers
are additionally rotated by $60^{\circ}$ (see the illustration on top). In
(a,b), the top layer of AS3 at $\theta\simeq 1.30^{\circ}$ and AS4 at
$\theta\simeq 1.54^{\circ}$, respectively, is additionally rotated by
$60^{\circ}$. In (c), such additional $60^{\circ}$-rotation is applied to two
top layers of AS4 at $\theta\simeq 1.35^{\circ}$. Figure 11: The structural
and electronic properties of two AS3 structures: AS3 at $\theta\simeq
1.54^{\circ}$ in Fig.8 (see in the left) and AS3 at $\theta\simeq
1.30^{\circ}$ considered in Fig.10a (see in the right). Images in the first
and second rows present the spatial variation of out-of-plane displacements in
the top and middle layers, respectively. LDOS of two systems at zero energy is
presented in the bottom images.

When ASn systems are created, another remarkable possibility could in
principle be obtained, in particular, one or a few graphene layers are rotated
by $\theta_{60}=\theta$+60∘ (instead of $\theta$) as illustrated on top of
Fig.10, similar to that reported in Ref. finney2019  for graphene/hBN
structures. Actually, such $\theta_{60}$-rotations do not change the lattice
orientation of the graphene layer, compared to the simple $\theta$-rotated
systems considered in Figs.8-9. However, three AAA-, ABA- and BAB-stacking
regions obtained in the latter case are changed to AAB-, ABB- and ABC-stacking
ones in the ASn systems with $\theta$+60∘ rotations, leading to significant
changes in their electronic structure. Indeed, the flat electronic bands are
observed at another twist angle as shown in Fig.10, compared to similar
structures (i.e., with the same number of layers) presented in Figs.8-9. In
addition, dispersive bands at low energies are generally no longer linear
ones. These changes can be clarified by analyzing the structural properties of
two considered types of ASn systems (e.g., see in Fig.11 for the AS3 cases).
In particular, when compared to the first type (i.e., with $\theta$ rotations
only), the top layer of the AS3 system in Fig.10a presents a totally different
buckling pattern (see the top panels in Fig.11). Similarly, the middle layer
of the latter system is no longer flat as that observed in the former one (see
the panels in the second row of Fig.11). Accordingly to these different
structural properties, these two systems exhibit totally different LDOS images
as shown in the bottom of Fig.11. These results thus clarify the difference in
their electronic structures discussed above.

Figure 12: Electronic bandstructures of TMGs including DB- (i.e., twisted
double bilayer graphene) and AM-structures (alternatively twisted structures
containing, at least, one Bernal stacking bilayer, as schematized in Fig.1).
Calculations were performed at their magic angles where low-energy bands are
most flat.

At last, the electronic bandstructures of other TMGs schematized in Fig.1 at
their magic angles are presented in Fig.12. First, the isolated flat bands are
only obtained in magic-angle DB-structure (i.e., twisted double bilayer
graphene) haddadi2020 , similar to those obtained in TBLG and in SM-(AB)2 one.
AMn-systems (n=1,2,3,4,5) are actually alternatively twisted structures but
different from the AS ones, containing at least one Bernal stacking bilayer.
Most remarkably, the results obtained in those systems show that the flat
bands can be easily obtained if the middle layer is a graphene monolayer
(i.e., AM1 and AM3 structures), whereas the low-energy bands are generally
less flat if it is a Bernal stacking bilayer. Moreover, similar to the
SM-(AB)n and ASn structures, dispersive bands are always additionally obtained
at low energies, which is likely an inherent property of alternatively twisted
systems. Once more, the low-energy bands of those systems in the presence of
Bernal stacking bilayers are shown to be less flat than those presented above
for systems of graphene monolayers only.

Figure 13: Electronic bandstructures of moiré of moiré graphene superlattices
created by putting one rotated (by $\theta_{23}$) graphene monolayer on top of
magic-angle TBLG (i.e., $\theta_{12}\simeq 1.08^{\circ}$). Two typical cases,
particularly, large angle ($\theta_{23}\simeq 23.96^{\circ}$) and small angle
($\theta_{23}\simeq 3.25^{\circ}$), are presented.

For completeness, we finally present in Fig.13 the electronic bandstructures
of moiré-of-moiré graphene superlattices, which are here trilayer graphene
systems created by vertically stacking magic-angle TBLG (i.e.,
$\theta_{12}\simeq 1.08^{\circ}$) and one rotated (by the angle $\theta_{23}$)
graphene monolayer. We particularly consider two cases: a large angle
$\theta_{23}\simeq 23.96^{\circ}$ and a small one $\theta_{23}\simeq
3.25^{\circ}$. Actually, when $\theta_{23}$ is large, the interaction between
the graphene monolayer and the magic angle TBLG has no significant effect at
low energies, i.e., the obtained low-energy bands are approximately a simple
combination of the bands of those two systems. Actually, this feature is very
similar to those obtained in TBLGs in the large angle regime when their bands
can be considered as a combination of the bands of two graphene monolayers.
However, when $\theta_{23}$ is small (i.e., $3.25^{\circ}$ considered here),
the mentioned interaction has significant effects on the electronic structure,
in particular, the bandwidth of flat bands is significantly enlarged, compared
to that obtained in magic-angle TBLG. These results suggest that besides
twisted graphene systems investigated above, other more complicated twisted
multilayer ones, forming moiré-of-moiré superlattices, shall generally exhibit
more complicated electronic structures with less opportunity for observing the
flat bands and accordingly, strong electron localization.

##  V Flattening electronic bands and tunability

As emphasized in section III, the electronic flat bands in reconstructed TBLG
is essentially obtained in the situation when the real-space localization in
AA stacking regions and their contribution to the global electronic properties
are concurrently maximized at the magic angle. In fact, as it will be
demonstrated below, such feature is commonly observed in TMGs and hence could
be considered as the practical origin of their observed flat bands.

Figure 14: Spatial variation of LDOS in TBLG, AS3 and SM-(AB)2 systems when
varying the twist angle around their magic one. The presented LDOS is
extracted along the line connecting regions where AA, SP, AB/BA stacking
configurations occur. Calculations are performed at energy where the strongest
electronic localization is obtained, i.e., at the electron van Hove
singularity point for large $\theta$ and at zero energy for $\theta$ below the
magic angle. Figure 15: Electronic structures of TMGs modulated by strain and
vertical pressure. TBLG at $\theta\simeq 1.30^{\circ}$ is considered in (a)
for the pristine case, in (b) with a biaxial strain of 9$\%$, and in (c) with
a vertical pressure $P\simeq 1.45$ GPa. The corresponding LDOS obtained at the
strongest electronic localization point is presented in (d). Electronic
bandstructures of AS3 at $\theta\simeq 1.89^{\circ}$ with $P=0$ and $\simeq$
2.27 GPa is displayed in (e) and (f), respectively. The results of SM-(AB)2 at
$\theta\simeq 1.47^{\circ}$ with $P=0$ and $\simeq$ 2.27 GPa are in (g) and
(h), respectively.

First, it is worth reminding that based on continuum models rafi2011 , the
first magic angle of TBLG can be estimated using the equation

$\theta_{MA}\approx 3w/(2\pi t)$ (2)

where w and t stand for interlayer and intralayer tunneling strengths,
respectively. Indeed, $\theta_{MA}\approx 1.1^{\circ}$ is obtained, in very
good agreement with experiments cao2018a ; yankowitz2019 ; gadelha2021 ;
barbosa2022 , with $w=0.109$ eV and $t=2.7$ eV. Calculations using continuum
models rafi2011 ; khalaf2019  have additionally predicted a discrete series of
other magic angles $<1.1^{\circ}$. However, in agreement with the existing
experiments (most clearly, see recent articles gadelha2021 ; barbosa2022 ),
atomistic calculations have demonstrated in Ref. nguyen2021  and as discussed
in section III that these small magic angles can not be practically obtained
in reconstructed TBLGs. This is essentially due to the fact that while it gets
maximum at $\theta\simeq 1.1^{\circ}$, the global electronic localization is
progressively reduced, according to the increasing contribution of AB/BA
stacking regions when reducing $\theta$ below $<1.1^{\circ}$ as seen in Fig.4
and more visibly in the left-panel of Fig.14. Interestingly, all these
features are similarly observed in other TMGs as seen in two examples, i.e.,
AS3 and SM-(AB)2 systems, presented in Fig.14. Indeed, their electronic flat
bands (i.e., magic angle) are exactly obtained when the real-space
localization in AAA and AAB stacking regions, respectively, and their
contribution to the global electronic properties of the system are
concurrently maximized. Moreover, the commonness of such observation is
further emphasized as it is also valid in TMGs when in-plane strains and/or
vertical pressure are applied to tune the magic angle value (see Fig.15 below,
especially, the full illustration in Figs.15(a-d)). On the basis of those
results, we thus conclude that the maximization of real-space electronic
localization in regions where the AA stacking configuration is presented could
be considered as the essential origin of the flat electronic bands observed in
TMGs.

In addition, the results in Fig.14 are another illustration of feature
discussed above that at their magic angles, the real-space electronic
localization in the presence of Bernal bilayer graphene (i.e., the SM-(AB)2
structure here) is significantly weaker than that in the systems of graphene
monolayers only. This feature explained why less flat bands (accordingly,
smaller localized DOS peak in energy) are obtained in magic-angle systems
containing Bernal bilayer graphene, compared to other structures.

Figure 16: Effects of vertical electric field $E_{\bot}$ on the electronic
bandstructure of TBLG at $\theta\simeq 1.08^{\circ}$ and DM-structure at
$\theta\simeq 1.25^{\circ}$.

As mentioned, the dependence of magic-angle $\theta_{MA}$ on interlayer and
intralayer tunneling strengths described in Eq.(2) suggests that $\theta_{MA}$
can be tuned by modifying $w$ and $t$, for instance, by applying in-plane
strains and vertical pressure nguyen2015 ; huder2018 ; yankowitz2019 ; lin2020
; liu2020s ; wu2021 ; balint2021  as illustrated in Fig.15. In particular, the
intralayer hopping energies are reduced by a tensile in-plane strain while the
interlayer tunneling strength increases when a vertical pressure is applied
pressure , thus increasing $\theta_{MA}$ as shown. Even though they could be
quantitatively different, these tunable possibilities are commonly obtained in
all TMGs as illustrated in Fig.15 for AS3 and SM-(AB)2 structures.

Another tunability has been actually explored bandwidth ; liu2020 ; shen2020 ;
chen2021e ; rickhaus2021 , i.e., by a vertical electric field (E-field) as
illustrated in Figs.16-18. First, Fig.16 presents typical pictures showing the
E-field effects in magic-angle twisted double bilayer graphene (i.e., DB-
structure), compared to those in the TBLG structure. It shows that the E-field
generally presents stronger effects in DB-structure than in TBLG one. This
could be simply understood as a direct consequence of the fact that E-field
can modulate the bandstructure and open a bandgap in Bernal stacking bilayer
graphene, which is not observed for monolayer graphene. These results thus
emphasize that vertical E-field is a good degree of freedom to control the
flat band and strongly correlated electronic phenomena in twisted systems
containing Bernal stacking bilayer graphene, indeed as reported in Refs.
bandwidth ; liu2020 ; shen2020 ; chen2021e ; rickhaus2021 .

Figure 17: Electric field ($E_{\bot}$) effects on the electronic bandstructure
of AS3 at $\theta\simeq 1.54^{\circ}$ (a,e), AS4 at $\theta\simeq
1.70^{\circ}$ (b,f), AS5 at $\theta\simeq 1.79^{\circ}$ (c,g), and AS3
considered in Fig.10a at $\theta\simeq 1.30^{\circ}$ (d,h). Figure 18: Layer
decomposed DOS (LD-DOS) in systems considered in Figs.17(e-f) when an electric
field $E_{\perp}=10$ meV/Å is applied.

Alternatively twisted multilayer systems of graphene monolayers (i.e., AS
structures) is another example that has shown a good tunability by vertical
E-field zeyu2021 ; park2021 ; jeong2021 . In particular, when an E-field is
applied, a distinct feature observed in AS structures is that the flat bands
can strongly hybridize with the Dirac bands, leading to the controllability of
the bandwidth and interaction strength in the flat bands. In Fig.17, the
electronic bandstructure of ASn (n=3,4,5) is presented for two different
E-fields. In addition to the mentioned hybridization, some other interesting
effects are found. Indeed, while two flat bands hybridize with two Dirac
bands, two other flat bands are still observed in the AS3 case. For the AS4
systems, four Dirac bands are observed at low energies (see Fig.8) and hence
the hybridization takes place for all four flat bands with these four Dirac
bands when a large E-field is applied. Moreover, the hybridized flat bands can
be separated in energy by the E-field. Similar features occur for ASn systems
with $n\geqslant 5$, except that besides four Dirac bands hybridizing with
flat bands, the remaining Dirac bands are still observed near the Fermi level
and a larger separation of hybridized flat bands is obtained.

For the ASn structures investigated in Fig.10 when a $\theta_{60}$ rotation is
applied, the E-field even presents different effects, compared to the
conventional cases in Fig.8. In particular, instead of the band hybridization,
flat bands in such systems are increasingly isolated from other dispersive
bands when an E-field is applied, as seen in Fig.17(d,h) for the AS3 system of
Fig.10(a). These results could be simply understood by its different
structural properties, compared to the conventional AS3 one of Fig.8 as
analyzed in Fig.11.

Finally, to clarify more completely those observed electronic features of ASn
systems, we analyze their layer decomposed DOS, when an E-field is applied, in
Fig.18. In the AS3 case, even though the band hybridization occurs, changes in
LD-DOS are relatively weak at the considered E-fields. This could be because
besides two hybridized flat bands, other two flat bands are still preserved
near the Fermi level as mentioned above. Stronger E-field effects on LD-DOS
are found in the cases of a larger number of layers. In particular, the single
localized DOS peak at zero E-field is generally separated into two peaks when
an E-field is applied. Interestingly, those two peaks are shown to be likely
separated in different layers. Indeed, the electron localized peak is mostly
observed in the bottom layers whereas the hole peak in the top ones in the
case of AS4 and AS5 systems in Fig.18. These results suggest possibilities to
control the spatial dependence of electron localization and therefore of
strongly correlated electronic phenomena in those twisted systems.

##  VI Summary

Using atomistic calculations, the electronic properties of twisted multilayer
graphene systems are systematically studied. First, it is shown that both
their structural and electronic properties exhibit a common feature that there
is a unique magic angle at which the electronic localization in regions
containing AA stacking configuration and its contribution to the global
properties of the system are concurrently maximized. Flat electronic bands and
strongly correlated electronic phenomena are therefore obtained at this
critical case. In addition, the structural and electronic (particularly,
electronic localization) properties are shown to exhibit a strong correlation.
Consequently, the magic angle separates TMG systems into two classes, i.e., in
small and large angle regimes, exhibiting distinct properties as clearly
illustrated by the presented results obtained in TBLG nguyen2021 . Moreover,
the strong correlation between the maximization of electronic localization in
AA stacking region and the observation of flat bands presents an unknown
physical view about the origin of magic angle. The electronic properties of
magic-angle TMGs containing more than 2 graphene layers are shown to exhibit
some other distinct properties, due to the presence of a large number of
layers as well as various stacking configurations. In particular, the flat
bands and dispersive ones at low energies can be concurrently obtained in some
structures such as SM-(AB)n and alternatively twisted ones. The coexistence of
localized-delocalized electronic states spatially separated are observed. At
last, these TMGs present good possibilities for tuning electronic properties
by external fields. Finally, the magic angle is also shown to be nicely
tunable by the effects of strain and/or vertical pressure. Our study thereby
provides a comprehensive overview as well as highlights the essential and
outstanding electronic properties of twisted graphene systems that could be
helpful for further developments in twistronics research.

## Acknowledgments

V.-H.N. and J.-C.C. acknowledge financial support from the Fédération
Wallonie-Bruxelles through the ARC Grants (N∘ 16/21-077 and N∘ 21/26-116),
from the European Union’s Horizon 2020 Research Project and Innovation Program
— Graphene Flagship Core3 (N∘ 881603), from the Flag-Era JTC project “TATTOOS”
(N∘ R.8010.19), from the EOS project “CONNECT” (N∘ 40007563) and from the
Belgium F.R.S.-FNRS through the research projects (N∘ T.0051.18 and N∘
T.029.22F). Computational resources have been provided by the CISM
supercomputing facilities of UCLouvain and the CÉCI consortium funded by
F.R.S.-FNRS of Belgium (N∘ 2.5020.11). T.X.H. acknowledges support from
International Centre of Physics at the Institute of Physics, VAST, under grant
number ICP.2022.05.

## References

  * [1] Jiang-Bin Wu, Xin Zhang, Mari Ijäs, Wen-Peng Han, Xiao-Fen Qiao, Xiao-Li Li, De-Sheng Jiang, Andrea C. Ferrari, and Ping-Heng Tan.  Resonant raman spectroscopy of twisted multilayer graphene.  Nat. Commun., 5:5309, 2014.
  * [2] Xu-Dong Chen, Wei Xin, Wen-Shuai Jiang, Zhi-Bo Liu, Yongsheng Chen, and Jian-Guo Tian.  High-precision twist-controlled bilayer and trilayer graphene.  Adv. Mater., 28:2563–2570, 2016.
  * [3] Umesha Mogera, Sunil Walia, Bharath Bannur, Murali Gedda, and Giridhar U. Kulkarni.  Intrinsic nature of graphene revealed in temperature-dependent transport of twisted multilayer graphene.  J. Phys. Chem. C, 121:13938–13943, 2017.
  * [4] Cheng Shen, Yanbang Chu, QuanSheng Wu, Na Li, Shuopei Wang, Yanchong Zhao, Jian Tang, Jieying Liu, Jinpeng Tian, Kenji Watanabe, Takashi Taniguchi, Rong Yang, Zi Yang Meng, Dongxia Shi, Oleg V. Yazyev, and Guangyu Zhang.  Correlated states in twisted double bilayer graphene.  Nat. Phys., 16:520–525, 2020.
  * [5] Xiaomeng Liu, Zeyu Hao, Eslam Khalaf, Jong Yeon Lee, Yuval Ronen, Hyobin Yoo, Danial Haei Najafabadi, Kenji Watanabe, Takashi Taniguchi, Ashvin Vishwanath, and Philip Kim.  Tunable spin-polarized correlated states in twisted double bilayer graphene.  Nature, 583:221–225, 2020.
  * [6] Jeong Min Park, Yuan Cao, Kenji Watanabe, Takashi Taniguchi, and Pablo Jarillo-Herrero.  Tunable strongly coupled superconductivity in magic-angle twisted trilayer graphene.  Nature, 590:249–255, 2021.
  * [7] Xi Zhang, Kan-Ting Tsai, Ziyan Zhu, Wei Ren, Yujie Luo, Stephen Carr, Mitchell Luskin, Efthimios Kaxiras, and Ke Wang.  Correlated insulating states and transport signature of superconductivity in twisted trilayer graphene superlattices.  Phys. Rev. Lett., 127:166802, 2021.
  * [8] Shaowen Chen, Minhao He, Ya-Hui Zhang, Valerie Hsieh, Zaiyao Fei, K. Watanabe, T. Taniguchi, David H. Cobden, Xiaodong Xu, Cory R. Dean, and Matthew Yankowitz.  Electrically tunable correlated and topological states in twisted monolayer–bilayer graphene.  Nat. Phys., 17:374–380, 2021.
  * [9] Zeyu Hao, A. M. Zimmerman, Patrick Ledwith, Eslam Khalaf, Danial Haie Najafabadi, Kenji Watanabe, Takashi Taniguchi, Ashvin Vishwanath, and Philip Kim.  Electric field–tunable superconductivity in alternating-twist magic-angle trilayer graphene.  Science, 371:1133–1138, 2021.
  * [10] Maria Brzhezinskaya, Oleg Kononenko, Victor Matveev, Aleksandr Zotov, Igor I. Khodos, Vladimir Levashov, Vladimir Volkov, Sergey I. Bozhko, Sergey V. Chekmazov, and Dmitry Roshchupkin.  Engineering of numerous moire superlattices in twisted multilayer graphene for twistronics and straintronics applications.  ACS Nano, 15:12358–12366, 2021.
  * [11] Jeong Min Park, Yuan Cao, Liqiao Xia, Shuwen Sun, Kenji Watanabe, Takashi Taniguchi, and Pablo Jarillo-Herrero.  Magic-angle multilayer graphene: A robust family of moire superconductors.  arXiv:2112.10760, 2021.
  * [12] Eva Y. Andrei and Allan H. MacDonald.  Graphene bilayers with a twist.  Nat. Mater., 19:1265–1275, 2020.
  * [13] J. Wang, X. Mu, L. Wang, and M. Sun.  Properties and applications of new superlattice: twisted bilayer graphene.  Mater. Today Phys., 9:100099, 2019.
  * [14] Stephen Carr, Shiang Fang, and Efthimios Kaxiras.  Electronic-structure methods for twisted moire layers.  Nat. Rev. Mater., 5:748–763, 2020.
  * [15] Rafi Bistritzer and Allan H. MacDonald.  Moire bands in twisted double-layer graphene.  PNAS, 108:12233–12237, 2011.
  * [16] M. Iqbal Bakti Utama, Roland J. Koch, Kyunghoon Lee, Nicolas Leconte, Hongyuan Li, Sihan Zhao, Lili Jiang, Jiayi Zhu, Kenji Watanabe, Takashi Taniguchi, Paul D. Ashby, Alexander Weber-Bargioni, Alex Zettl, Chris Jozwiak, Jeil Jung, Eli Rotenberg, Aaron Bostwick, and Feng Wang.  Visualization of the flat electronic band in twisted bilayer graphene near the magic angle twist.  Nat. Phys., 17:184–188, 2021.
  * [17] Young Woo Choi and Hyoung Joon Choi.  Strong electron-phonon coupling, electron-hole asymmetry, and nonadiabaticity in magic-angle twisted bilayer graphene.  Phys. Rev. B, 98:241412, 2018.
  * [18] Alexander Kerelsky, Leo J. McGilly, Dante M. Kennes, Lede Xian, Matthew Yankowitz, Shaowen Chen, K. Watanabe, T. Taniguchi, James Hone, Cory Dean, Angel Rubio, and Abhay N. Pasupathy.  Maximized electron interactions at the magic angle in twisted bilayer graphene.  Nature, 572:95–100, 2019.
  * [19] Young Woo Choi and Hyoung Joon Choi.  Dichotomy of electron-phonon coupling in graphene moire flat bands.  Phys. Rev. Lett., 127:167001, 2021.
  * [20] Andreij C. Gadelha, Douglas A. A. Ohlberg, Cassiano Rabelo, Eliel G. S. Neto, Thiago L. Vasconcelos, João L. Campos, Jessica S. Lemos, Vinícius Ornelas, Daniel Miranda, Rafael Nadas, Fabiano C. Santana, Kenji Watanabe, Takashi Taniguchi, Benoit van Troeye, Michael Lamparski, Vincent Meunier, Viet-Hung Nguyen, Dawid Paszko, Jean-Christophe Charlier, Leonardo C. Campos, Luiz G. Cançado, Gilberto Medeiros-Ribeiro, and Ado Jorio.  Localization of lattice dynamics in low-angle twisted bilayer graphene.  Nature, 590:405–409, 2021.
  * [21] Andreij C Gadelha, Viet-Hung Nguyen, Eliel G Neto, Fabiano Santana, Markus B Raschke, Michael Lamparski, Vincent Meunier, Jean-Christophe Charlier, and Ado Jorio.  Electron-phonon coupling in a magic-angle twisted-bilayer graphene device.  arXiv:2110.14916, 2021.
  * [22] Tiago C Barbosa, Andreij C Gadelha, Douglas A A Ohlberg, Kenji Watanabe, Takashi Taniguchi, Gilberto Medeiros-Ribeiro, Ado Jorio, and Leonardo C Campos.  Raman spectra of twisted bilayer graphene close to the magic angle.  2D Mater., 9:025007, 2022.
  * [23] Yuan Cao, Valla Fatemi, Shiang Fang, Kenji Watanabe, Takashi Taniguchi, Efthimios Kaxiras, and Pablo Jarillo-Herrero.  Unconventional superconductivity in magic-angle graphene superlattices.  Nature, 556:43–50, 2018.
  * [24] Yuan Cao, Valla Fatemi, Ahmet Demir, Shiang Fang, Spencer L. Tomarken, Jason Y. Luo, Javier D. Sanchez-Yamagishi, Kenji Watanabe, Takashi Taniguchi, Efthimios Kaxiras, Ray C. Ashoori, and Pablo Jarillo-Herrero.  Correlated insulator behaviour at half-filling in magic-angle graphene superlattices.  Nature, 556:80–84, 2018.
  * [25] Matthew Yankowitz, Shaowen Chen, Hryhoriy Polshyn, Yuxuan Zhang, K. Watanabe, T. Taniguchi, David Graf, Andrea F. Young, and Cory R. Dean.  Tuning superconductivity in twisted bilayer graphene.  Science, 363:1059–1064, 2019.
  * [26] Yonglong Xie, Biao Lian, Berthold Jãck, Xiaomeng Liu, Cheng-Li Chiu, Kenji Watanabe, Takashi Taniguchi, B. Andrei Bernevig, and Ali Yazdani.  Spectroscopic signatures of many-body correlations in magic-angle twisted bilayer graphene.  Nature, 572:101–105, 2019.
  * [27] Xiaobo Lu, Petr Stepanov, Wei Yang, Ming Xie, Mohammed Ali Aamir, Ipsita Das, Carles Urgell, Kenji Watanabe, Takashi Taniguchi, Guangyu Zhang, Adrian Bachtold, Allan H. MacDonald, and Dmitri K. Efetov.  Superconductors, orbital magnets and correlated states in magic-angle bilayer graphene.  Nature, 574:653–657, 2019.
  * [28] Youngjoon Choi, Jeannette Kemmer, Yang Peng, Alex Thomson, Harpreet Arora, Robert Polski, Yiran Zhang, Hechen Ren, Jason Alicea, Gil Refael, Felix von Oppen, Kenji Watanabe, Takashi Taniguchi, and Stevan Nadj-Perge.  Electronic correlations in twisted bilayer graphene near the magic angle.  Nat. Phys., 15:1174–1180, 2019.
  * [29] Yuhang Jiang, Xinyuan Lai, Kenji Watanabe, Takashi Taniguchi, Kristjan Haule, Jinhai Mao, and Eva Y. Andrei.  Charge order and broken rotational symmetry in magic-angle twisted bilayer graphene.  Nature, 573:91–95, 2019.
  * [30] Aaron L. Sharpe, Eli J. Fox, Arthur W. Barnard, Joe Finney, Kenji Watanabe, Takashi Taniguchi, M. A. Kastner, and David Goldhaber-Gordon.  Emergent ferromagnetism near three-quarters filling in twisted bilayer graphene.  Science, 365:605–608, 2019.
  * [31] M. Serlin, C. L. Tschirhart, H. Polshyn, Y. Zhang, J. Zhu, K. Watanabe, T. Taniguchi, L. Balents, and A. F. Young.  Intrinsic quantized anomalous hall effect in a moire heterostructure.  Science, 367:900–903, 2020.
  * [32] Leon Balents, Cory R. Dean, Dmitri K. Efetov, and Andrea F. Young.  Superconductivity and strong correlations in moiré flat bands.  Nat. Phys., 16:725–733, 2020.
  * [33] Stephen Carr, Daniel Massatt, Shiang Fang, Paul Cazeaux, Mitchell Luskin, and Efthimios Kaxiras.  Twistronics: Manipulating the electronic properties of two-dimensional layered structures through their twist angle.  Phys. Rev. B, 95:075420, 2017.
  * [34] Kaihui Liu, Liming Zhang, Ting Cao, Chenhao Jin, Diana Qiu, Qin Zhou, Alex Zettl, Peidong Yang, Steve G. Louie, and Feng Wang.  Evolution of interlayer coupling in twisted molybdenum disulfide bilayers.  Nat. Commun., 5:4966, 2014.
  * [35] Peng Kang, Wan-Ting Zhang, Vincent Michaud-Rioux, Xiang-Hua Kong, Chen Hu, Guang-Hua Yu, and Hong Guo.  Moiré impurities in twisted bilayer black phosphorus: Effects on the carrier mobility.  Phys. Rev. B, 96:195406, 2017.
  * [36] Mit H. Naik and Manish Jain.  Ultraflatbands and shear solitons in moire patterns of twisted bilayer transition metal dichalcogenides.  Phys. Rev. Lett., 121:266401, 2018.
  * [37] Lede Xian, Dante M. Kennes, Nicolas Tancogne-Dejean, Massimo Altarelli, and Angel Rubio.  Multiflat bands and strong correlations in twisted bilayer boron nitride: Doping-induced correlated insulator and superconductor.  Nano Lett., 19:4934–4940, 2019.
  * [38] Yaping Yang, Jidong Li, Jun Yin, Shuigang Xu, Ciaran Mullan, Takashi Taniguchi, Kenji Watanabe, Andre K. Geim, Konstantin S. Novoselov, and Artem Mishchenko.  In situ manipulation of van der waals heterostructures for twistronics.  Sci. Adv., 6:eabd3655, 2020.
  * [39] Dante M. Kennes, Martin Claassen, Lede Xian, Antoine Georges, Andrew J. Millis, James Hone, Cory R. Dean, D. N. Basov, Abhay N. Pasupathy, and Angel Rubio.  Moire heterostructures as a condensed-matter quantum simulator.  Nat. Phys., 17:155–163, 2021.
  * [40] Qiaoling Xu, Yuzheng Guo, and Lede Xian.  Moiré flat bands in twisted 2d hexagonal vdW materials.  2D Mater., 9(1):014005, 2021.
  * [41] Shengdan Tao, Xuanlin Zhang, Jiaojiao Zhu, Pimo He, Shengyuan A. Yang, Yunhao Lu, and Su-Huai Wei.  Designing ultra-flat bands in twisted bilayer materials at large twist angles: Theory and application to two-dimensional indium selenide.  J. Am. Chem. Soc., 144:3949–3956, 2022.
  * [42] Mengya Liu, Liping Wang, and Gui Yu.  Developing graphene-based moire heterostructures for twistronics.  Adv. Sci., 9:2103170, 2022.
  * [43] Qingjun Tong, Fei Liu, Jiang Xiao, and Wang Yao.  Skyrmions in the moire of van der waals 2d magnets.  Nano Lett., 18:7194–7199, 2018.
  * [44] Chong Wang, Yuan Gao, Hongyan Lv, Xiaodong Xu, and Di Xiao.  Stacking domain wall magnons in twisted van der waals magnets.  Phys. Rev. Lett., 125:247201, 2020.
  * [45] Tiancheng Song, Qi-Chao Sun, Eric Anderson, Chong Wang, Jimin Qian, Takashi Taniguchi, Kenji Watanabe, Michael A. McGuire, Rainer Stöhr, Di Xiao, Ting Cao, Jörg Wrachtrup, and Xiaodong Xu.  Direct visualization of magnetic domains and moire magnetism in twisted 2d magnets.  Science, 374:1140–1144, 2021.
  * [46] Xi-Wang Luo and Chuanwei Zhang.  Spin-twisted optical lattices: Tunable flat bands and larkin-ovchinnikov superfluids.  Phys. Rev. Lett., 126:103201, 2021.
  * [47] Wenhui Wang, Wenlong Gao, Xiaodong Chen, Fulong Shi, Guixin Li, Jianwen Dong, Yuanjiang Xiang, and Shuang Zhang.  Moire fringe induced gauge field in photonics.  Phys. Rev. Lett., 125:203901, 2020.
  * [48] Xin-Rui Mao, Zeng-Kai Shao, Hong-Yi Luan, Shao-Lei Wang, and Ren-Min Ma.  Magic-angle lasers in nanostructured moire superlattice.  Nat. Nanotechnol., 16:1099–1105, 2021.
  * [49] Qing Zhang, Guangwei Hu, Weiliang Ma, Peining Li, Alex Krasnok, Rainer Hillenbrand, Andrea Alu, and Cheng-Wei Qiu.  Interface nano-optics with van der waals polaritons.  Nature, 597:187–195, 2021.
  * [50] Jialin Chen, Xiao Lin, Mingyuan Chen, Tony Low, Hongsheng Chen, and Siyuan Dai.  A perspective of twisted photonic structures.  Appl. Phys. Lett., 119:240501, 2021.
  * [51] Eslam Khalaf, Alex J. Kruchkov, Grigory Tarnopolsky, and Ashvin Vishwanath.  Magic angle hierarchy in twisted graphene multilayers.  Phys. Rev. B, 100:085109, 2019.
  * [52] Fatemeh Haddadi, QuanSheng Wu, Alex J. Kruchkov, and Oleg V. Yazyev.  Moiré flat bands in twisted double bilayer graphene.  Nano Lett., 20:2410–2415, 2020.
  * [53] Youngju Park, Bheema Lingam Chittari, and Jeil Jung.  Gate-tunable topological flat bands in twisted monolayer-bilayer graphene.  Phys. Rev. B, 102:035411, 2020.
  * [54] Young Woo Choi and Hyoung Joon Choi.  Intrinsic band gap and electrically tunable flat bands in twisted double bilayer graphene.  Phys. Rev. B, 100:201402, 2019.
  * [55] Peter Rickhaus, Folkert K. de Vries, Jihang Zhu, Elías Portoles, Giulia Zheng, Michele Masseroni, Annika Kurzmann, Takashi Taniguchi, Kenji Watanabe, Allan H. MacDonald, Thomas Ihn, and Klaus Ensslin.  Correlated electron-hole state in twisted double-bilayer graphene.  Science, 373:1257–1260, 2021.
  * [56] Georgios A Tritsaris, Stephen Carr, Ziyan Zhu, Yiqi Xie, Steven B Torrisi, Jing Tang, Marios Mattheakis, Daniel T Larson, and Efthimios Kaxiras.  Electronic structure calculations of twisted multi-layer graphene superlattices.  2D Mater., 7:035028, 2020.
  * [57] Stephen Carr, Shiang Fang, Ziyan Zhu, and Efthimios Kaxiras.  Exact continuum model for low-energy electronic states of twisted bilayer graphene.  Phys. Rev. Research, 1:013001, 2019.
  * [58] Georgios A Tritsaris, Stephen Carr, Ziyan Zhu, Yiqi Xie, Steven B Torrisi, Jing Tang, Marios Mattheakis, Daniel T Larson, and Efthimios Kaxiras.  Electronic structure calculations of twisted multi-layer graphene superlattices.  2D Mater., 7(3):035028, 2020.
  * [59] Jiseon Shin, Bheema Lingam Chittari, and Jeil Jung.  Stacking and gate-tunable topological flat bands, gaps, and anisotropic strip patterns in twisted trilayer graphene.  Phys. Rev. B, 104:045413, 2021.
  * [60] Fernando Gargiulo and Oleg V Yazyev.  Structural and electronic transformation in low-angle twisted bilayer graphene.  2D Mater., 5:015019, 2017.
  * [61] V. Hung Nguyen, D. Paszko, M. Lamparski, B. van Troeye, V. Meunier, and J.-C. Charlier.  Electronic localization in small-angle twisted bilayer graphene.  2D Mater., 8:035046, 2021.
  * [62] L. Lindsay and D. A. Broido.  Optimized tersoff and brenner empirical potential parameters for lattice dynamics and phonon thermal transport in carbon nanotubes and graphene.  Phys. Rev. B, 81:205441, 2010.
  * [63] Aleksey N. Kolmogorov and Vincent H. Crespi.  Registry-dependent interlayer potential for graphitic systems.  Phys. Rev. B, 71:235415, 2005.
  * [64] Itai Leven, Tal Maaravi, Ido Azuri, Leeor Kronik, and Oded Hod.  Interlayer potential for graphene/h-bn heterostructures.  J. Chem. Theory Comput., 12:2896–2905, 2016.
  * [65] I. Brihuega, P. Mallet, H. González-Herrero, G. Trambly de Laissardière, M. M. Ugeda, L. Magaud, J. M. Gómez-Rodríguez, F. Ynduráin, and J.-Y. Veuillen.  Unraveling the intrinsic and robust nature of van hove singularities in twisted bilayer graphene by scanning tunneling microscopy and theoretical analysis.  Phys. Rev. Lett., 109:196802, 2012.
  * [66] Hyobin Yoo, Rebecca Engelke, Stephen Carr, Shiang Fang, Kuan Zhang, Paul Cazeaux, Suk Hyun Sung, Robert Hovden, Adam W. Tsen, Takashi Taniguchi, Kenji Watanabe, Gyu-Chul Yi, Miyoung Kim, Mitchell Luskin, Ellad B. Tadmor, Efthimios Kaxiras, and Philip Kim.  Atomic and electronic reconstruction at the van der waals interface in twisted bilayer graphene.  Nat. Mater., 18(5):448–453, 2019.
  * [67] Shuai Zhang, Aisheng Song, Lingxiu Chen, Chengxin Jiang, Chen Chen, Lei Gao, Yuan Hou, Luqi Liu, Tianbao Ma, Haomin Wang, Xi-Qiao Feng, and Qunyang Li.  Abnormal conductivity in low-angle twisted bilayer graphene.  Sci. Adv., 6(47):eabc5555, 2020.
  * [68] Kazuyuki Uchida, Shinnosuke Furuya, Jun-Ichi Iwata, and Atsushi Oshiyama.  Atomic corrugation and electron localization due to moire patterns in twisted bilayer graphenes.  Phys. Rev. B, 90:155451, 2014.
  * [69] Giovanni Cantele, Dario Alfè, Felice Conte, Vittorio Cataudella, Domenico Ninno, and Procolo Lucignano.  Structural relaxation and low-energy properties of twisted bilayer graphene.  Phys. Rev. Research, 2:043127, 2020.
  * [70] Michael Lamparski, Benoit Van Troeye, and Vincent Meunier.  Soliton signature in the phonon spectrum of twisted bilayer graphene.  2D Mater., 7:025050, 2020.
  * [71] G. Trambly de Laissardière, D. Mayou, and L. Magaud.  Numerical studies of confined states in rotated bilayers of graphene.  Phys. Rev. B, 86:125413, 2012.
  * [72] Long-Jing Yin, Jia-Bin Qiao, Wen-Xiao Wang, Wei-Jie Zuo, Wei Yan, Rui Xu, Rui-Fen Dou, Jia-Cai Nie, and Lin He.  Landau quantization and fermi velocity renormalization in twisted graphene bilayers.  Phys. Rev. B, 92:201408, 2015.
  * [73] Ling-Hui Tong, Qingjun Tong, Li-Zhen Yang, Yue-Ying Zhou, Qilong Wu, Yuan Tian, Li Zhang, Lijie Zhang, and Long-Jing Yin Zhihui Qin.  Spectroscopic visualization of flat bands in magic-angle twisted monolayer-bilayer graphene: localization-delocalization coexisting electronic states.  arXiv:2109.13528, 2021.
  * [74] Canxun Zhang, Tiancong Zhu, Salman Kahn, Shaowei Li, Birui Yang, Charlotte Herbig, Xuehao Wu, Hongyuan Li, Kenji Watanabe, Takashi Taniguchi, Stefano Cabrini, Alex Zettl, Michael P. Zaletel, Feng Wang, and Michael F. Crommie.  Visualizing delocalized correlated electronic states in twisted double bilayer graphene.  Nat. Commun., 12:2516, 2021.
  * [75] Mikito Koshino.  Interlayer screening effect in graphene multilayers with $aba$ and $abc$ stacking.  Phys. Rev. B, 81:125304, 2010.
  * [76] Nathan R. Finney, Matthew Yankowitz, Lithurshanaa Muraleetharan, K. Watanabe, T. Taniguchi, Cory R. Dean, and James Hone.  Tunable crystal symmetry in graphene–boron nitride heterostructures with coexisting moiré superlattices.  Nat. Nanotechnol., 14:1029–1034, 2019.
  * [77] V Hung Nguyen and P Dollfus.  Strain-induced modulation of dirac cones and van hove singularities in a twisted graphene bilayer.  2D Mater., 2:035005, 2015.
  * [78] Loïc Huder, Alexandre Artaud, Toai Le Quang, Guy Trambly de Laissardière, Aloysius G. M. Jansen, Gérard Lapertot, Claude Chapelier, and Vincent T. Renard.  Electronic spectrum of twisted graphene layers under heterostrain.  Phys. Rev. Lett., 120:156405, 2018.
  * [79] Xianqing Lin, Haotian Zhu, and Jun Ni.  Pressure-induced gap modulation and topological transitions in twisted bilayer and twisted double bilayer graphene.  Phys. Rev. B, 101:155405, 2020.
  * [80] Yi-Wen Liu, Ying Su, Xiao-Feng Zhou, Long-Jing Yin, Chao Yan, Si-Yu Li, Wei Yan, Sheng Han, Zhong-Qiu Fu, Yu Zhang, Qian Yang, Ya-Ning Ren, and Lin He.  Tunable lattice reconstruction, triangular network of chiral one-dimensional states, and bandwidth of flat bands in magic angle twisted bilayer graphene.  Phys. Rev. Lett., 125:236102, 2020.
  * [81] Zewen Wu, Xueheng Kuang, Zhen Zhan, and Shengjun Yuan.  Magic angle and plasmon mode engineering in twisted trilayer graphene with pressure.  Phys. Rev. B, 104:205104, 2021.
  * [82] Bálint Szentpéteri, Peter Rickhaus, Folkert K. de Vries, Albin Márffy, Bálint Fülöp, Endre Tóvári, Kenji Watanabe, Takashi Taniguchi, Andor Kormányos, Szabolcs Csonka, and Péter Makk.  Tailoring the band structure of twisted double bilayer graphene with pressure.  Nano Lett., 21:8777–8784, 2021.
  * [83] To model the effects of vertical pressure P, the average interlayer distance d is fixed and modified (in relaxation simulations) with respect to its equilibrium value $d_{0}$. The same approach is applied to superlattice vectors to model in-plane strains. P is estimated as in Phys. Rev. B 98, 085144 (2018), i.e., $P=A(e^{-B\epsilon}-1)$ with A = 5.73 GPa, B = 9.54, and $\epsilon=1-(d/d_{0})$.


position and in velocity. Galaxies in a normal font match the catalogued
optical counterpart in position, but with no prior catalogued redshift. For
these galaxies our data give the first redshift measurements. For the galaxies
in an italic font, there is a spatial match with a previously uncatalogued
optical counterpart. The “Galaxies with uncertain parameters” are only
partially included in the data presented here. The central velocities listed
are approximate.

##  9 Summary

We have presented a description of the MHONGOOSE survey, a large, deep, high-
resolution H i survey on MeerKAT of 30 nearby gas-rich spiral and dwarf
galaxies. The total survey comprises 1650h of observing time, with 55h
allocated per galaxy. The main science goal of MHONGOOSE is to characterise
the low-column density H i in these galaxies, linking it with various cold gas
accretion scenarios. Other science goals include probing in detail the link
between H i and star formation, the kinematics of the galaxies, constraining
the presence and distribution of dark matter, as well as the environments of
the targets. The sample covers $\sim 3$ orders of magnitude in H i mass, from
low-mass $M_{\rm H{\hskip 0.56917pt\tt I}}\sim 10^{7}$ $M_{\odot}$ as
typically found in dwarf galaxies to $M_{\rm H{\hskip 0.56917pt\tt I}}\sim
10^{10}$ $M_{\odot}$ as found in large spiral galaxies.

In this paper, we present the survey, the sample, and a summary of the
observations. In addition we describe the data reduction procedure and data
products.

Each galaxy will eventually be observed for 55h in 10 tracks of 5.5h, and we
produce full-depth data cubes at a velocity resolution of 1.4 km s-1 and a
range in angular resolutions from $\sim 8^{\prime\prime}$ as the highest
resolution to $\sim 90^{\prime\prime}$ as our lowest resolution. The resulting
column density sensitivity ($3\sigma$ over 16 km s-1) ranges from $4.3\cdot
10^{19}$ cm-2 at the highest resolution to $4.2\cdot 10^{17}$ cm-2 at the
lowest resolution. The $3\sigma$, 50 km s-1 mass-detection limit is $5.4\cdot
10^{5}$ $M_{\odot}$ at 10 Mpc, the median distance of the sample.

We present first results based on the full-depth data of a limited number of
galaxies and single-track observations of the complete sample (the latter
observations comprising 10% of the survey data). Comparing full-depth, single-
track and single-dish data, we show that MeerKAT is excellent at recovering
the total fluxes of our sample galaxies. We also show that the full-depth data
detect on average $\sim 2\%$ more H i than the single-track data. While this
increase may seem modest, the additional H i mass it represents is sufficient
to allow for the existence of a $\sim 10^{18}$ cm-2 low-column density H i
component in the disks of our galaxies. This shows the difficulty of
investigating the presence of low-column density H i using total H i mass
measurements.

We also investigate the link between the local H i surface density and the
velocity “dispersion” in all of our sample galaxies. We recover the well-known
upper limit on H i surface density (measured on kpc scales), and find a lower
limit to the measured second-moment value of $\sim 8.5$ km s-1. We find that
the H i surface density has two distinct locations in a surface density-
velocity dispersion diagram. One at $\sim 5\cdot 10^{20}$ cm-2, representing
the gas in the star-forming disks, and one at $\sim 5\cdot 10^{19}$ cm-2
representing gas not participating in star formation. The two regions are
distinct and separated at an H i surface density of $2\cdot 10^{20}$ cm-2.
Previous studies (e.g. Kanekar et al. 2011) have shown that at this surface
density a rapid change in spin temperature occurs, and that this is the column
density where the formation of a cool and dense component of the ISM starts.

We give an overview of all additional galaxies detected in a the single-track
data. These are defined as the sources detected in H i that were not targeted
by the survey, but were nevertheless found within the $1.5^{\circ}\times
1.5^{\circ}$ area and spanning 1000 km s-1 in velocity range around the target
galaxies. We find that half of these 49 galaxies have optical and catalogued
counterparts that coincide in position and velocity. Just under half have
catalogued optical counterparts that coincide in position but that had until
now no previously measured redshifts. Five of our additional galaxies have
optical counterparts that have not been previously catalogued. We have not
found any H i detections that do not have optical counterparts.

Studies of the full-depth data will be presented in future papers. However,
based on the results shown here, it is already clear that the MHONGOOSE data
is of excellent quality, showing great promise for future work on the low-
column density H i in nearby galaxies.

Data: The single-track moment maps presented in this paper are available for
download at the MHONGOOSE website (https://mhongoose.astron.nl) or through a
DOI (https://zenodo.org/doi/10.5281/zenodo.10907079).

###### Acknowledgements.

It is a pleasure to thank SARAO staff members past and present for their help
in making MHONGOOSE a reality. The MeerKAT telescope is operated by the South
African Radio Astronomy Observatory, which is a facility of the National
Research Foundation, an agency of the Department of Science and Innovation.
This work has received funding from the European Research Council (ERC) under
the European Union’s Horizon 2020 research and innovation programme (grant
agreement No. 882793 “MeerGas”). We thank SURF, and particularly Dr. Raymond
Oonk for the support in the early stages of the data reduction. KAO is
supported by a Royal Society Dorothy Hodgkin Fellowship, by STFC through grant
ST/T000244/1, and by the European Research Council (ERC) through Advanced
Investigator grant to C.S. Frenk, DMIDAS (GA 786910). DJP, NZ and SK are
supported through the South African Research Chairs Initiative of the
Department of Science and Technology and National Research Foundation (Grant
number 77825). EA and AB thank the CNES for financial support. LC acknowledges
funding from the Chilean Agencia Nacional de Investigación y Desarrollo (ANID)
through Fondo Nacional de Desarrollo Científico y Tecnológico (FONDECYT)
Regular Project 1210992. The work at AIRUB is supported by the German Federal
Ministry of Education and Research (BMBF) Verbundforschung grant 05A20PC4
(Verbundprojekt D-MeerKAT-II). KS acknowledges support from the Natural
Sciences and Engineering Research Council of Canada (NSERC). LVM and AS
acknowledge financial support from the grant CEX2021-001131-S funded by
MCIN/AEI/10.13039/501100011033 and from the grant PID2021-123930OB-C21 funded
by MCIN/AEI/10.13039/501100011033, by “ERDF A way of making Europe” and by the
European Union. BN acknowledges financial support from the Severo Ochoa grant
CEX2021-001131-S funded by MCIN/AEI/10.13039/501100011033. This research has
made use of the NASA/IPAC Extragalactic Database (NED), which is funded by the
National Aeronautics and Space Administration and operated by the California
Institute of Technology. This research has made use of the SIMBAD database,
operated at CDS, Strasbourg, France.

## References

  * Anand et al. (2021) Anand, G. S., Lee, J. C., Van Dyk, S. D., et al. 2021, MNRAS, 501, 3621
  * Auld et al. (2006) Auld, R., Minchin, R. F., Davies, J. I., et al. 2006, MNRAS, 371, 1617
  * Barnes et al. (2001) Barnes, D. G., Staveley‐Smith, L., de Blok, W. J. G., et al. 2001, MNRAS, 322, 486
  * Bigiel et al. (2010) Bigiel, F., Leroy, A., Walter, F., et al. 2010, AJ, 140, 1194
  * Bigiel et al. (2008) Bigiel, F., Leroy, A., Walter, F., et al. 2008, AJ, 136, 2846
  * Bigiel et al. (2011) Bigiel, F., Leroy, A. K., Walter, F., et al. 2011, ApJ, 730, L13
  * Braun (1997) Braun, R. 1997, ApJ, 484, 637
  * Braun & Thilker (2004) Braun, R. & Thilker, D. A. 2004, A&A, 417, 421
  * Camilo (2018) Camilo, F. 2018, Nature Astronomy, 2, 594
  * Crain et al. (2016) Crain, R. A., Bahé, Y. M., Lagos, C. d. P., et al. 2016, MNRAS, 464, 4204
  * Davé et al. (1999) Davé, R., Hernquist, L., Katz, N., & Weinberg, D. H. 1999, ApJ, 511, 521
  * de Blok et al. (2016) de Blok, W. J. G., Adams, E. A. K., Amram, P., et al. 2016, Proceedings of MeerKAT Science: On the Pathway to the SKA — PoS(MeerKAT2016), 007
  * de Blok et al. (2020) de Blok, W. J. G., Athanassoula, E., Bosma, A., et al. 2020, A&A, 643, A147
  * de Blok et al. (2014) de Blok, W. J. G., Keating, K. M., Pisano, D. J., et al. 2014, A&A, 569, A68
  * de Blok et al. (2018) de Blok, W. J. G., Walter, F., Ferguson, A. M. N., et al. 2018, ApJ, 865, 26
  * Dey et al. (2019) Dey, A., Schlegel, D. J., Lang, D., et al. 2019, AJ, 157, 168
  * Di Teodoro & Fraternali (2014) Di Teodoro, E. M. & Fraternali, F. 2014, A&A, 567, A68
  * Di Teodoro & Fraternali (2015) Di Teodoro, E. M. & Fraternali, F. 2015, MNRAS, 451, 3021
  * Elagali et al. (2019) Elagali, A., Staveley-Smith, L., Rhee, J., et al. 2019, MNRAS, 487, 2797
  * Espada et al. (2011a) Espada, D., Muñoz-Mateos, J. C., Gil de Paz, A., et al. 2011a, ApJ, 736, 20
  * Espada et al. (2011b) Espada, D., Verdes-Montenegro, L., Huchtmeier, W. K., et al. 2011b, A&A, 532, A117
  * Fraternali (2017) Fraternali, F. 2017, arXiv [1612.00477]
  * Fraternali et al. (2001) Fraternali, F., Oosterloo, T., Sancisi, R., & Moorsel, G. v. 2001, ApJ, 562, L47
  * Fraternali & Tomassetti (2012) Fraternali, F. & Tomassetti, M. 2012, MNRAS, 426, 2166
  * Haynes et al. (2018) Haynes, M. P., Giovanelli, R., Kent, B. R., et al. 2018, ApJ, 861, 49
  * Heald et al. (2016) Heald, G., de Blok, W. J. G., Lucero, D., et al. 2016, MNRAS, 462, 1238
  * Heald et al. (2011) Heald, G., Józsa, G., Serra, P., et al. 2011, A&A, 526, A118
  * Hess et al. (2015) Hess, K. M., Jarrett, T. H., Carignan, C., Passmoor, S. S., & Goedhart, S. 2015, MNRAS, 452, 1617
  * Hunter et al. (2012) Hunter, D. A., Ficut-Vicas, D., Ashley, T., et al. 2012, AJ, 144, 134
  * Ianjamasimanana et al. (2017) Ianjamasimanana, R., de Blok, W. J. G., & Heald, G. H. 2017, AJ, 153, 213
  * Irwin et al. (2009) Irwin, J. A., Hoffman, G. L., Spekkens, K., et al. 2009, ApJ, 692, 1447
  * Jarrett et al. (2023) Jarrett, T. H., Cluver, M. E., Taylor, E. N., et al. 2023, ApJ, 946, 95
  * Jonas (2018) Jonas, J. 2018, Proceedings of MeerKAT Science: On the Pathway to the SKA — PoS(MeerKAT2016), 001
  * Jones et al. (2017) Jones, M. G., Espada, D., Verdes-Montenegro, L., et al. 2017, A&A, 609, A17
  * Józsa et al. (2020) Józsa, G. I., Andati, L. A., de Blok, W., et al. 2020, Astronomical Society of the Pacific Conference Series, 532, 447
  * Kamphuis et al. (2022) Kamphuis, P., Jütte, E., Heald, G. H., et al. 2022, A&A, 668, A182
  * Kanekar et al. (2011) Kanekar, N., Braun, R., & Roy, N. 2011, ApJ, 737, L33
  * Karachentsev et al. (2002) Karachentsev, I. D., Sharina, M. E., Dolphin, A. E., et al. 2002, A&A, 385, 21
  * Kennicutt (1998) Kennicutt, R. C. 1998, ApJ, 498, 541
  * Kereš et al. (2005) Kereš, D., Katz, N., Weinberg, D. H., & Davé, R. 2005, MNRAS, 363, 2
  * Koribalski et al. (2004) Koribalski, B. S., Staveley-Smith, L., Kilborn, V. A., et al. 2004, AJ, 128, 16
  * Koribalski et al. (2020) Koribalski, B. S., Staveley-Smith, L., Westmeier, T., et al. 2020, Ap&SS, 365, 118
  * Koribalski et al. (2018) Koribalski, B. S., Wang, J., Kamphuis, P., et al. 2018, MNRAS, 478, 1611
  * Kourkchi & Tully (2017) Kourkchi, E. & Tully, R. B. 2017, ApJ, 843, 16
  * Kregel et al. (2004) Kregel, M., Kruit, P. C. V. D., & de Blok, W. J. G. 2004, MNRAS, 352, 768
  * Krumholz et al. (2009) Krumholz, M. R., McKee, C. F., & Tumlinson, J. 2009, ApJ, 693, 216
  * Leon et al. (2008) Leon, S., Verdes-Montenegro, L., Sabater, J., et al. 2008, A&A, 485, 475
  * Leroy et al. (2019) Leroy, A. K., Sandstrom, K. M., Lang, D., et al. 2019, ApJS, 244, 24
  * Leroy et al. (2008) Leroy, A. K., Walter, F., Brinks, E., et al. 2008, AJ, 136, 2782
  * Leroy et al. (2013) Leroy, A. K., Walter, F., Sandstrom, K., et al. 2013, AJ, 146, 19
  * Lisenfeld et al. (2007) Lisenfeld, U., Verdes-Montenegro, L., Sulentic, J., et al. 2007, A&A, 462, 507
  * Liu et al. (2023) Liu, Y., Zhu, M., Yu, H., et al. 2023, MNRAS, 523, 3905
  * Lucero et al. (2015) Lucero, D. M., Carignan, C., Elson, E. C., et al. 2015, MNRAS, 450, 3935
  * Maccagni et al. (2022) Maccagni, F. M., Healy, J., de Blok, W. J. G., & Serra, P. 2022, MeerKAT Technical Report, DOI: 10.48479/bhpj-nz95
  * Maccagni et al. (2017) Maccagni, F. M., Morganti, R., Oosterloo, T. A., Geréb, K., & Maddox, N. 2017, A&A, 604, A43
  * Maccagni et al. (2021) Maccagni, F. M., Serra, P., Gaspari, M., et al. 2021, A&A, 656, A45
  * Makhathini (2018) Makhathini, S. 2018, PhD thesis, Rhodes University, South Africa
  * Mao et al. (2021) Mao, Y.-Y., Geha, M., Wechsler, R. H., et al. 2021, ApJ, 907, 85
  * Marasco et al. (2012) Marasco, A., Fraternali, F., & Binney, J. J. 2012, MNRAS, 419, 1107
  * Marasco et al. (2019) Marasco, A., Fraternali, F., Heald, G., et al. 2019, A&A, 631, A50
  * Mauch et al. (2020) Mauch, T., Cotton, W. D., Condon, J. J., et al. 2020, ApJ, 888, 61
  * McKee & Krumholz (2010) McKee, C. F. & Krumholz, M. R. 2010, ApJ, 709, 308
  * Meurer et al. (2006) Meurer, G. R., Hanish, D. J., Ferguson, H. C., et al. 2006, ApJS, 165, 307
  * Meurer et al. (1998) Meurer, G. R., Staveley‐Smith, L., & Killeen, N. E. B. 1998, MNRAS, 300, 705
  * Meurer et al. (2013) Meurer, G. R., Zheng, Z., & de Blok, W. J. G. 2013, MNRAS, 429, 2537
  * Meyer et al. (2004) Meyer, M. J., Zwaan, M. A., Webster, R. L., et al. 2004, MNRAS, 350, 1195
  * Moster et al. (2013) Moster, B. P., Naab, T., & White, S. D. M. 2013, MNRAS, 428, 3121
  * Offringa et al. (2014) Offringa, A. R., McKinley, B., Hurley-Walker, N., et al. 2014, MNRAS, 444, 606
  * Oosterloo et al. (2007) Oosterloo, T., Fraternali, F., & Sancisi, R. 2007, AJ, 134, 1019
  * Pingel et al. (2018) Pingel, N. M., Pisano, D. J., Heald, G., et al. 2018, ApJ, 865, 36
  * Pisano (2014) Pisano, D. J. 2014, AJ, 147, 48
  * Popping et al. (2009) Popping, A., Davé, R., Braun, R., & Oppenheimer, B. D. 2009, A&A, 504, 15
  * Pritzl et al. (2003) Pritzl, B. J., Knezek, P. M., III, J. S. G., et al. 2003, ApJ, 596, L47
  * Putman et al. (2012) Putman, M., Peek, J., & Joung, M. 2012, ARA&A, 50, 491
  * Ramesh et al. (2023) Ramesh, R., Nelson, D., & Pillepich, A. 2023, MNRAS, 518, 5754
  * Sabater et al. (2012) Sabater, J., Verdes-Montenegro, L., Leon, S., Best, P., & Sulentic, J. 2012, A&A, 545, A15
  * Saintonge et al. (2013) Saintonge, A., Lutz, D., Genzel, R., et al. 2013, ApJ, 778, 2
  * Sancisi et al. (2008) Sancisi, R., Fraternali, F., Oosterloo, T., & Hulst, T. v. d. 2008, A&A Rev., 15, 189
  * Sardone et al. (2021) Sardone, A., Pisano, D. J., Pingel, N. M., et al. 2021, ApJ, 910, 69
  * Schruba et al. (2017) Schruba, A., Leroy, A. K., Kruijssen, J. M. D., et al. 2017, ApJ, 835, 278
  * Serra et al. (2023) Serra, P., Maccagni, F. M., Kleiner, D., et al. 2023, A&A, 673, A146
  * Serra et al. (2012) Serra, P., Oosterloo, T., Morganti, R., et al. 2012, MNRAS, 422, 1835
  * Serra et al. (2015) Serra, P., Westmeier, T., Giese, N., et al. 2015, MNRAS, 448, 1922
  * Shapiro & Field (1976) Shapiro, P. R. & Field, G. B. 1976, ApJ, 205, 762
  * Sorgho et al. (2018) Sorgho, A., Carignan, C., Pisano, D. J., et al. 2018, MNRAS, 482, 1248
  * Spekkens et al. (2014) Spekkens, K., Urbancic, N., Mason, B. S., Willman, B., & Aguirre, J. E. 2014, ApJ, 795, L5
  * Sternberg et al. (2014) Sternberg, A., Petit, F. L., Roueff, E., & Bourlot, J. L. 2014, ApJ, 790, 10
  * Sun et al. (2020) Sun, J., Leroy, A. K., Schinnerer, E., et al. 2020, ApJ, 901, L8
  * Sun et al. (2018) Sun, J., Leroy, A. K., Schruba, A., et al. 2018, ApJ, 860, 172
  * Tacconi et al. (2018) Tacconi, L. J., Genzel, R., Saintonge, A., et al. 2018, ApJ, 853, 179
  * Tully & Fisher (1977) Tully, R. & Fisher, J. 1977, A&A, 54, 661
  * Tully et al. (2013) Tully, R. B., Courtois, H. M., Dolphin, A. E., et al. 2013, AJ, 146, 86
  * Tully et al. (2016) Tully, R. B., Courtois, H. M., & Sorce, J. G. 2016, AJ, 152, 50
  * Tumlinson et al. (2017) Tumlinson, J., Peeples, M. S., & Werk, J. K. 2017, ARA&A, 55, 1
  * Verheijen & Sancisi (2001) Verheijen, M. A. W. & Sancisi, R. 2001, A&A, 370, 765
  * Veronese et al. (2023) Veronese, S., de Blok, W. J. G., & Walter, F. 2023, A&A, 672, A55
  * Wakker & Woerden (1997) Wakker, B. P. & Woerden, H. v. 1997, ARA&A, 35, 217
  * Walter et al. (2008) Walter, F., Brinks, E., de Blok, W. J. G., et al. 2008, AJ, 136, 2563
  * Walter et al. (2020) Walter, F., Carilli, C., Neeleman, M., et al. 2020, ApJ, 902, 111
  * Westmeier et al. (2021) Westmeier, T., Kitaeff, S., Pallot, D., et al. 2021, MNRAS, 506, 3962
  * Wolfe et al. (2016) Wolfe, S. A., Lockman, F. J., & Pisano, D. J. 2016, ApJ, 816, 81
  * Wolfe et al. (2013) Wolfe, S. A., Pisano, D. J., Lockman, F. J., McGaugh, S. S., & Shaya, E. J. 2013, Nature, 497, 224
  * Wolfire et al. (2003) Wolfire, M. G., McKee, C. F., Hollenbach, D., & Tielens, A. G. G. M. 2003, ApJ, 587, 278
  * Wong et al. (2016) Wong, O. I., Meurer, G. R., Zheng, Z., et al. 2016, MNRAS, 460, 1106
  * Wong & Blitz (2002) Wong, T. & Blitz, L. 2002, ApJ, 569, 157
  * Xu et al. (2022) Xu, C. K., Cheng, C., Appleton, P. N., et al. 2022, Nature, 610, 461
  * Zhu & Putman (2023) Zhu, J. & Putman, M. E. 2023, arXiv [2303.00763]

##  Appendix A Atlas of single-track moment maps

In Fig. A.18 we show the r05_t00 single-track moment maps of the entire
sample. A short description of the panels that are shown (from left to right)
for each galaxy follows.

The left-most panel shows a combined optical _grz_ color image from DECaLS DR9
(except for J1303–17b where we use a _giz_ image from DR10). The next panel
shows the zeroth-moment or integrated H i map. The lowest column density
contour shown corresponds to a S/N of 3 (see Sect. 5.9 for a description).
This value is listed at the bottom of the zeroth-moment panel. Subsequent
contours then increase by a factor $2^{n}$, where $n=0,2,4,...$.

In the third panel we show the first-moment map or intensity-weighted velocity
field. Here we indicate the central velocity (as listed in Table 1) with a
thick contour. The central velocity is defined from the integrated H i profile
as the velocity halfway between the velocities where the flux density equals
20% of the peak value (see Sect. 6.2). Other contours are spaced by 10 or 20
km s-1 with respect to this velocity, as listed in the panel. The color map
indicates whether velocities are receding (red) or approaching (blue).

The second-moment maps are shown using the same intensity scale for all
galaxies. For all galaxies, the colour map shows the same velocity range from
0 (light-blue) to 30 (red) km s-1. Contours start at 12 km s-1, and increase
in steps of 12 km s-1. For reference, the 24 km s-1 contour is shown in black.
Both the first- and second-moment maps have been masked such that only values
with ${\rm S/N}>3$ in the zeroth-moment map are shown. While we do not present
a full analysis of the first- and second-moment maps in this paper, we do
caution that for the edge-on galaxies, these have to be interpreted with care
because of the strong projection effects.

##  Appendix B The effect of beam smearing and inclination on second-moment
values

The column density-second moment distributions shown in Figs. 12 and 13 were
derived using the r15_t00 standard resolution with a beam size of $\sim
30^{\prime\prime}$. This beam size, in combination with the angular sizes of
the galaxies, could potentially lead to inflated second-moment values due to
beam smearing (especially in the centers of the galaxies). In addition, long
line-of-sight effects in the high-inclination galaxies could also contribute
to higher second-moment values. Here we address these potential issues and
show that they do not significantly affect the distribution shown in Fig. 13.

We test the effect of beam smearing by creating a model galaxy similar to one
of our sample galaxies, convolving it to the observed resolution and inserting
it in a data cube with a noise level equal to that found for the observations.
We then treat it exactly as the sample galaxies for the derivation of the
moment maps and creation of the moment distributions.

We use as the basis of our model the rotation curve and H i radial profile of
NGC1566 as presented in Elagali et al. (2019). We approximate the rotation
curve with a linear, steep rise from 0 to 180 km s-1 in the inner 5 kpc. At
larger radii we assume a constant rotation velocity of 180 km s-1 (cf. Fig. 10
in Elagali et al. 2019). For the H i distribution we assume a central value of
$\log(N_{\rm H{\hskip 0.56917pt\tt I}}/M_{\odot}\,{\rm pc}^{-2})=1$, with
$\log(N_{\rm H{\hskip 0.56917pt\tt I}})$ linearly decreasing with radius until
$\log(N_{\rm H{\hskip 0.56917pt\tt I}}/M_{\odot}\,{\rm pc}^{-2})=-0.5$ at 50
kpc (cf. Fig. 12 in Elagali et al. 2019). For the purposes of this exercise we
extend the radial range of the model to 100 kpc by assuming that the rotation
curve stays flat, and that the $\log(N_{\rm H{\hskip 0.56917pt\tt I}})$ values
keep decreasing at the same rate as inside 50 kpc. We adopt a value of 30
degrees for the inclination (cf. Fig 10 in Elagali et al. 2019). For the
velocity dispersion of the H i we use a constant value of 10 km s-1. The H i
disk is assumed to be thin, with a scale height of 100 pc. We assume the same
distance and systemic velocity as for the real NGC 1566 (17.6 Mpc and 1496 km
s-1, respectively; cf. Tab. 1). With these parameters we create a model cube
using the GALMOD section of the 3DBAROLO software package (Di Teodoro &
Fraternali 2015) with a beam size equal to that of the r15_t00 standard
resolution. We insert this model in a cube filled with an rms noise of $0.48$
mJy (again at the r15_t00 resolution), consistent with the noise values
typically found in the single-track cubes.

Zeroth- and second-moment maps of the model are then created with SoFiA-2 in
an identical manner as for the real galaxies. We plot the resulting column
density-second moment distribution in Fig. B.18 (right panel), comparing it
with the observed NGC 1566 distribution (left panel). Several noteworthy
features can be seen. Firstly, in the model we see a plume towards high
second-moment values at the highest column densities. This is caused by beam
smearing due to the slope of the inner rotation curve. Secondly, on the low-
column density end of the distribution we see the familiar tail towards low
second-moment value due to masking (as discussed in Sect. 7.3). Finally, in
the intermediate column density range, we see the constant second-moment
values. The limited resolution has slightly increased the average value from
its original 10 km s-1, but only by $\sim 1$ km s-1 or less. These therefore
are a good reflection of the velocity dispersion in the model disk.

In summary, only the high second-moment values at the highest column densities
(corresponding to the innermost part of the disk where the rotation curve
rises) can be attributed to beam smearing. Higher second-moment values at
lower column densities (i.e., further out in the disk) reflect higher velocity
dispersions, or the presence of multiple components along the line of sight.

We can compare the distribution shown in Fig. B.18 (right panel) with the much
more complex distribution for the real NGC 1566 in Fig. B.18 (left panel). The
latter shows higher second-moment values over a much larger range in column
densities, indicating that a significant fraction of these higher values are
not due to beam smearing. Inspection of Fig. 12 shows that these are rather
caused by high-second moment components in galaxies with overall low second-
moment values (e.g., J0419–54) as well as galaxies that have overall higher
second-moment values (e.g., J1318–21).

We can test this further by taking the inclinations of the galaxies into
account. In general, the broadening of profiles due to beam smearing and
length of the line of sight through the disk increases dramatically towards
high inclinations, with the effect strongest in edge-on galaxies. There is
therefore a possibility that high-inclination galaxies contribute to the
higher second-moment values we observe at intermediate column densities (i.e.,
around $\sim 10^{20}$ cm-2). In Fig. B.19 (left panel) we show again the
distribution from Fig. 13, which includes all galaxies in the sample, at all
inclinations, including edge-ons. To remove the effect that high-inclination
galaxies may have on the distribution, we show in the centre panel the
distribution limited to galaxies with an inclination $i<60^{\circ}$ (16
galaxies) and in the right panel the distribution limited to inclination
$i<70^{\circ}$ (22 galaxies). Comparison of these three distributions shows
that omitting the high-inclination galaxies mainly affects the width of the
second-moment distribution in the low-column density tail. The effects in the
intermediate and high-column density ranges are much less pronounced, and the
distribution remains essentially the same in these column density regimes.

Figure A.18: Single-track moment maps for the MHONGOOSE galaxies using the
r05_t00 resolution. From left to right: _(i):_ Combined _grz_ -colour image
from DECaLS. _(ii):_ primary-beam corrected zeroth-moment or integrated H i
intensity map. Contours as indicated in the Figure. The lowest contour
represents S/N = 3, with subsequent contour levels increasing by factors of
two. _(iii):_ First-moment map or intensity-weighted velocity field. Red
colours indicate the receding side, blue colours the approaching side. The
central velocity (listed in Tab. 1) is indicated by the thick contour. Other
contours are spaced by 10 or 20 km s-1, as indicated in the Figure. _(iv):_
Second-moment map: the colours show the range from 0 (light-blue) to 30 (red)
km s-1. The lowest contour shows the 12 km s-1 level, and subsequent contours
are spaced by 12 km s-1. The 24 km s-1 contour is shown in black. See text for
a more extensive description.

Figure A.18 (Cont.): Single-track moment maps of MHONGOOSE galaxies.

Figure A.18 (Cont.): Single-track moment maps of MHONGOOSE galaxies.

Figure A.18 (Cont.): Single-track moment maps of MHONGOOSE galaxies.

Figure A.18 (Cont.): Single-track moment maps of MHONGOOSE galaxies.

Figure B.18: Comparison of the column density-second moment distribution of a
model galaxy built to resemble NGC 1566 with that of NGC 1566 itself. The
observed distribution is shown in the left panel; the model in the right
panel. For the model we use a single, thin H i disk with a velocity dispersion
of 10 km s-1. The model galaxy is “observed” and treated just like the single-
track MHONGOOSE galaxies. Note the plume due to beam smearing at the high
column density end, and the masking tail at the low column density end. The
model does not show the high second-moment values at intermediate column
densities that are present in the observations.

Figure B.19: Comparison of the column density-second moment distribution for a
number of inclination limits. The left panel shows the distribution as shown
in Fig. 13 without any inclination limits. The centre panel shows the
distribution when limited to galaxies with inclination $i<60^{\circ}$. The
right panel shows the same for $i<70^{\circ}$.