SciELO - Scientific Electronic Library Online

 
vol.67 número5Colloidal soft matter physicsLondon superconductor and time-varying mesoscopic LC circuits índice de autoresíndice de materiabúsqueda de artículos
Home Pagelista alfabética de revistas  

Servicios Personalizados

Revista

Articulo

Indicadores

Links relacionados

  • No hay artículos similaresSimilares en SciELO

Compartir


Revista mexicana de física

versión impresa ISSN 0035-001X

Rev. mex. fis. vol.67 no.5 México sep./oct. 2021  Epub 28-Mar-2022

https://doi.org/10.31349/revmexfis.67.050102 

Research

Reviews

Electronic properties of 2D materials and its heterostructures: a minimal review

G. G. Naumisa 

aDepartamento De Sistemas Complejos, Instituto de Física, Universidad Nacional Autónoma de México, Apartado Postal 20-364, 01000, Ciudad de Mexico01000, México. e-mail: naumis@fisica.unam.mx


Abstract

An overview of two-dimensional (2D) materials electronic properties is presented, including research in multilayered heterostructures. An emphasis is made on simple models that contain the representative physical features seen among 2D materials while presenting different and important perspectives that have been ignored or overlooked in other reviews. Starting with a short section on the crystallographic and diffraction properties, the review continues with a discussion of the theoretical models needed to describe the electronic properties. A special emphasis is made on the rise of the Dirac equation in terms of the electronic wavefunctions’ frustration due to the underlying triangular symmetry of graphene. Then a new method to deal with such problems in other systems is presented. Also, a section concerning the less known graphene’s free-electron bands is presented, which is important to describe interactions with metals and liquids as water. These bands are explained in terms of the electron interaction with its charge image, resulting in an effective Hydrogen model leading to a Rydberg series. We also discuss the effects of the disorder, flexural modes, strain, and electromagnetic waves, using novel techniques developed in collaborations with other groups in Mexico. Using all of the previous techniques, other exotic matter phases are studied like Kekule and Moiré patterns, flat bands, topological insulators, and time-dependent topological states. Finally, heterostructures made by stacking layers of 2D materials are studied. A special section is devoted to the latest discovered superconductivity of graphene over graphene at magic angles, including our latest reduction of the problem onto a simple 2×2 Hamiltonian, which describes the phenomena. Moreover, any other stacking of graphene layers like trilayer graphene, can be reduced using such method.

Keywords: 2D Materials; graphene; silicene; borophene; phoshorene; transition metals dichalcogenides; electronic properties; optical properties

PACS: 73.22.Pr; 72.80.Vp; 78.67.Wj; 63.22.Rc; 61.48.Gh

1. Introduction and outline

It should not be a controversial undertaking to change the dimensionality of a physical system. On the one hand and when things go well, it leads to suggestive new viewpoints. On the other hand, new polemics and paradoxes arise. In the classic book “A Brief History of Time”, Stephen Hawkings wrote on the impossibility of two-dimensional (2D) life-forms by showing that digestion will be next to impossible. Yet, there are recent claims that 2D worlds are not necessarily “too simple” for life to exist [1]. If Phil Anderson’s manifesto “More Is Different” brings out the importance of having more elements to increase the complexity and have new emerging laws [2], the discovery of 2D materials showed the other side of the coin [3]: less is also different. And when looking at its physical properties, 2D materials can result in a quite different manifesto: sometimes less is more.

More amazing electronic, optical, and mechanical properties, like the material with the highest known electrical and thermal conductivity, highest known electronic moblility, and elasticity [4], yet many of these 2D materials can be strained, bent, and wrinkled as a soft material [5]. The 2D materials revolution started with graphene, a one-atom layer of carbon extracted from graphite [3,6]. A 2D crystal was thought to be impossible to make, as according to the Mermin-Wigner theorem there is no order with short-range interactions at a temperature different from zero. Eventually, it turns out that such theorem was literally too rigid to be applied in a realworld situation; graphene is a 2D membrane that lives in a 3D world, and it vibrates like a drum pad. Such vibrations, known as flexural modes, allow 2D materials to evade the fate predicted by many eminent theoretical physicists. In a record time of 6 years from graphene’s discovery, A. Geim and K. Novoselov were awarded the 2010 Nobel prize. As a lesson, this shows that still in the XXI century, there is a place for simple tools as adhesive tape, graphite, and an optical microscope if good ideas, creativity, and hard work is present. A second Nobel prize was given in 2016 to D.J. Thouless, F.D.M. Haldane and M. Kosterlitz for another 2D related discovery: topological phase transitions and topological phases of matter. This work was made nearly 40 years before the graphene’s discovery but pointed out in the same direction: order is not impossible in D< 3 at finite temperature, it is just a bit different.

Since the discovery of graphene, the main feature of any 2D material is now clear: what makes them special is their all surface nature. This leads to novel paths in physics as interactions are maximized. Afterward, new 2D materials were discovered, and now they are classified in families [7,8]. A short-list of these main 2D families is the following,

  1. Almost flat structures. Graphene is the most representative case. Other examples are Silicene, Germanene, hexagonal Boron Nitrate (h-BN), Borophene, Stanene. Most of the early research was made in this area.

  2. Puckered structures that are not completely flat. The most famous examples are transition metal dichalcogenides (TMDs) with a general formula MX2, with M a transition metal (Mo, W, Ti, Nb, etc.). Atom X can be a chalcogen element (Se, S, or Te), or group III-VI or IV-VI elements, producing materials as InSe, GaSe, GeSe, SnSe, SnS2, SnSe2.

  3. 2D transition metal carbides and nitrides, also known as MXenes. Their general formula is Mn+1Xn. Here M represents a transition metal, and X represents C or N (where n is 1, 2, or 3) with a surface terminated by O, OH, or F atoms.

The other three families are in the process of being investigated, 2D oxides and hydroxides, group-IV monochalcogenides (MX with M = Ge, Sn, and X = S, Se, Te), and 2D organic materials as pentacene or αT 3 graphene. From there, one can stack layers to form heterostructures as Moiré patterns [9], twisted bilayer and trilayer graphene [1’], etc. At this moment, there is a lot of research in this area due to the recent discovery of superconducting phases [11]. This results in a quantum phase diagram akin to those seen in highTc superconductors [1-12]. An important area is the study of how strain affects the electronic and optical [7,8,13-17]. The reason is twofold. On the one hand, the substrates induce lattice deformations, and on the other hand, the control of strain allows to modify the characteristics of a material. This results in a field known as straintonics. Similarly, there are other ways to control the optical and electronic properties: twistronics [11], valleytronics [18], origami and kirigami, and spintronics [19].

An important feature of 2D materials is the use of relativistic quantum equations as the Dirac and Weyl equations [20]. Moreover, strain can be included as pseudo magnetic fields with the notable property of producing fields much higher than real ones [20]. As expected, the marriage of solid-state physics with high energy approach has been very fruitful as allowed to produce in the lab effects that were impossible to be seen in otherwise much experimentally demanding systems [21]. We can cite the Klein effect [20], effects of the Dirac equation in curved spaces [22], and even analogies to black holes [23].

As expected, many reviews and books are covering all these topics [24,25]. Now the question is why to make another one. There are many answers to this question. The first is that the present review is focused on developing the tools to understand heterostructures and, from there, allow the reader to dig into complex quantum phase diagrams. It also covers topological properties and time-dependent topological properties. Also, we try to focus on questions that our group in Mexico pursued systematically, even a few months after the discovery of graphene in 2004; for example, why do Dirac cones appear, a question that the reader will probably not find in other reviews as somehow is bypassed as a legitimate question and eclipsed by the mathematics. Other questions posed during such period with Ph. D. students at the Instituto de Física, UNAM were the effects of the disorder, electromagnetic radiation and time-dependent strain without using perturbative methods, a subject that now is mainstream but was starting at that moment. After many years we developed tools and results that have been proved to be useful in the study of other 2D materials. And more importantly, in this review, we stick with the idea that most of the complex properties are contained in simple models that give a lot of insights into physics, and eventually give a physical sense to mathematics. In particular, here we include our new results on producing the simplest hamiltonian that describes twisted graphene at magic angles, bringing to the forefront its topological nature and the three main physical driving forces: quantum confinement, frustration balance and non-Abelian gauge fields.

The layout of this work is as follows. We start with a section devoted to defining the structure and diffraction properties. Then we study Dirac materials and their variants, multilayered 2D materials, and finally, we discuss superlattice properties.

2. Structure and diffraction of 2D materials

Compared with its 3D counterparts, the study of 2D materials structure and diffraction may seem to be simple, but in fact, it has its particularities and difficulties. If we stick with family one of our list, i.e., almost flat structures, we describe its structure using the point lattice,

rn1,n2=n1a1+n2a2, (1)

where a1 ,2 are the Bravais unitary vectors (see Fig. 1) and n 1 ,n 2 integers. The most simple case is graphene, in which the point lattice is given by a trigonal one that we will call A sublattice with vectors,

a1=a23,3,   a2=a2-3,3, (2)

where a = 1.42 Å is the distance between C atoms [26].˚ Also, we need to add a basis vector b j that points to the different atoms of the unit cell; in other words, it describes the so-called decoration of the unit cell. In graphene, if C atoms are at sites r n1 ,n 2, the others are made by a simple translation δ1 of these points, forming the B sublattice. As seen in Fig. 1, the lattice turns out to be a honeycomb one. For further reference, it is useful to define the vectors,

δ1=a23,1,δ2=a2-3,1,δ3=a0,-1, (3)

FIGURE 1 a) Graphene honeycomb lattice showing the unit cell (shaded), the Bravais vectors a 1 and a 2, which define a trigonal lattice. The first-neighbor vectors δ1, δ2, and δ3 are also indicated. The bipartite sublattices A and B, which define two trigonal lattices, are shown. b) First Brillouin zone (shaded) of graphene and the high-symmetry points K+ and K. The Fermi level and the Dirac points coincide withthe inequivalent high symmetry points K+ and K 

The vectors δ2 and δ3 can be written using the basis and Bravais lattice vectors,

δ2=δ1+a2-a1,δ3=δ1-a1. (4)

An important property is that A sublattice atoms only have as neighbors B lattice atoms and viceversa. Such lattices are called bipartite and are fundamental as they give an extra symmetry that proves relevant in all of its properties. Figure 1b) presents the reciprocal lattice vectors:

G1=2π3a3,1,G2=2π3a-3,1, (5)

and the first Brillouin zone (1BZ). Two high-symmetry points K±=±4π/33a,0 are observed, each related by timereversal symmetry. Other authors label them as K and K’; however, when considering strain this notation leads to confusion as we have explained in other works [14,16]. Points K+ and K are at intersections of the diffraction Bragg lines 2k·G1 = ±|G1|2,2k·G2 = ±|G2|2 and 2k·(G1+G2) = ±|G1 +G2|2. As usual, the physics is dominated by such an intersection as the Fermi surface falls exactly there, leading to a singularity in the density of electronic states (DOS) due to the vanishing electron group velocity [27]. The singularity is known as the Dirac cone.

The diffraction can be obtained from the form and atomic structure factors. A simple example is to consider Dirac delta-function scatterers centered at each lattice site:

Vr=V0lδr-rl, (6)

where r l are the atom positions. The diffraction pattern is given by the norm of the Fourier transform of V (r),

V~k=SVreikrdS. (7)

As an example, we consider a decorated version of the honeycomb lattice as in h-BN. Here the A sublattice can be occupied by Nitrogen atoms while the other sublattice is occupied by Boron atoms. These atoms each have a scattering potential V 0 and V 1, from where

V~k=GV0+V1eikδ1δk-G. (8)

As seen in Fig. 2, the diffraction is made from spots at the trigonal lattice reciprocal vectors. A spot in location lG1 + hG2 has intensity,

Fl,h=V0-V12+4V0V1cos2π32l-h, (9)

FIGURE 2 A typical diffraction pattern for a 2D materials, in this case graphene. The Miller indexes l and h label each diffraction spot position and the amplitude is given by F(l,h).  

where l,h are the Miller indexes. If we now make the atoms the same as in graphene V 0 = V 1, we recover graphene’s case,

Fgpl,h=4V02cos2π32l-h. (10)

with two possible intensities 4V02 or V02 for the diffraction peaks. For V 1 = 0, the spots have the same amplitude and therefore correspond to the diffraction of a pure trigonal lattice. A comparison with the experimental results is found in Ref. [28]. Observe that experiments are performed in 3D, so the spots are, in fact, rods. The flatness of graphene, as well as the number of layers, can be tested by looking at modulations of the rods’ intensities [28]. This simple picture is modified by several factors:

  • Strain and flexural deformations. The recipe is to add a strain vectorial field au(r).

  • Stacking of multiple layers. Diffraction is just the sum of each layer diffraction.

The strain and flexural deformation field can be periodic, random, or quasiperiodic and has been studied in a previous review [29]. For periodic strain, satellites peaks appear close to the Bragg lines, while for disordered cases, the spots broaden. Also, the Brillouin zone can be folded back as seen in Kekule bond patterns [18,30].´

In the multilayered case, with rotated layers and if we suppose that the lattices are not deformed, the diffraction pattern is a superposition of the rotated reciprocal lattices corresponding to each layer [9]. Rotations and displacements lead to different space groups [31] (known in 2D as wallpaper groups). A very powerful method, originally developed to treat quasicrystals, allows describing the lattice positions and diffraction using higher-dimensional lattices, even providing extra information as the coincidence lattice, which determines many physical properties of heterostructures [32].

3. Electronic properties of Dirac materials

The experimental values of the electronic properties of graphene and other 2D materials have been extensively discussed in several reviews [20,25]. Here we only made few useful remarks. In particular, what matters for the electronic and optical properties are the most energetic electrons, i.e., those with energies within a band of width k B T (being T the temperature and k B the Boltzmann constant) around the Fermi energy (E F ). Their velocity, concentration, and mean free path will control such properties.

Although in graphene’s we use relativistic quantum mechanics to describe electrons, in reality, the Fermi velocity v F is around cv F /300 with c the speed of light, which is by no means unusual and is comparable to other 3D material, as for example, Silicium. Moreover, as we will see below, the charge density is very low. This points out to the most two remarkable properties of graphene, the large mean free path even at room temperatures due to the absence of backscattering [24] and the possibility of changing the carrier concentration by electric doping. In other words, while semiconductors are doped chemically with p or n donors, in graphene, such doping is tuned dynamically by an external field to generate holes or electrons. This property has not been fully used in electronic devices, but it stands out as unique.

The electronic properties of graphene are surprisingly well described by a one-band tight-binding Hamiltonian [26,33]. As what matters is the Fermi level, a low-energy approach leading to an effective Dirac Hamiltonian turns out to be very useful. As we will see, this idea permeates to other 2D materials. Below we present a minimal approach to these models.

3.1. Dirac materials: graphene

Carbon has four electrons in valence orbitals, but three of them are used to make in-plane sp 2 covalent bonds with its three neighbors [34]. The remaining electron half-fills each Carbon π− hybrid orbital [34]. Neglecting the three valence electrons as their energetic contributions are far from the Fermi level, this leads to a one-band tight-binding (TB) Hamiltonian [26]:

H0=-t0rjn=13arjbrj+δn+H.c., (11)

where r j indicates atoms in the A sites positions. The hopping or transfer integral is t 0 ≈ 2.7 eV by fitting experimental data [26]. arj and brj+δn are creation and annihilation operators on the A sublattice at position r j , and on the B sublattice, respectively. To reproduce the energy dispersion of graphene in the whole Brillouin zone using a TB Hamiltonian, it is required the use of second and third nearest neighbors interactions [25,35]. These corrections are important to study the transport properties [36-38] and to describe strain and surface effects [39]. Equation (11) is further reduced to a 2 × 2 operator using the following Fourier transform,

arj=1Nkeikrjak, (12)

where k is a wavevector, and a similar transformation holds for br+δ n , resulting in,

H0=-t0kn=13e-ikδnakbk+H.c. (13)

The leads to an effective 2 × 2 Hamiltonian H(k). The Schrödinger equation is now,

0HABkHAB*k0akbk=Ekakbk, (14)

where H AB (k) = −t 0 f(k), and f(k) is a complex function:

fk=n=13e-ikδn. (15)

The energy dispersion is readily found from the eigenvalues of Eq. (14):

E±k=±t0fk2=±t0n,s=13e-ikδn-δs2 (16)

=±t03+2cos3kxa+4cos3kxa/2cos3kya/2. (17)

For a system with N sites, the eigenfunctions in real space are,

Ψk,±rj=eikrjN1±e-ikδ1 (18)

Figure 3 presents the surface E(k) obtained from Eq. (17). First, notice the symmetry around E = 0, known as the particle-hole symmetry. This is a consequence of the bipartite property and can be understood from the CyrotLackmann theorem as a result of the absence of odd-rings in the lattice [40]. Without charge doping by external fields, the π orbitals are half-filled, and therefore, the Fermi energy (E F ) lies at E = 0. As indicated in Fig. 3, Eq. (17) displays a conical dispersion near E = 0. The condition E = 0 produces two special points labeled by K D , known as Dirac points. For pristine graphene, K D coincides with K±. For strained graphene, these points do not coincide, and there have been many confusions in the literature about this point [14,16]. The existence of two inequivalent Dirac points leads to the concept of the valley to classify the regions around each cone tip [26]. Let us see how Equation (17) produces cones: for a crystal momentum k = K D + q with qa1, the energy dispersion Eq. (17) can be developed in powers of qa:

Ek=Eq=±vFq, (19)

where v F is the Fermi velocity:

vF=3t0a2. (20)

FIGURE 3 Energy dispersion obtained from Eq. (17). To the left, we show its corresponding DOS (ρ(E)). A zoom-in at the Fermi energy showing a cone is displayed as well. The vertices of the cones touch at the Dirac point, here at the high-symmetry points K± with linear dispersion. The Fermi energy is at E = 0. At the Dirac cones, the DOS behave as ρ(E) ∼ E, while the saddle points of Eq. (17) produce the Van Hove singularities, which are the two peaks in the DOS. They correspond to highly degenerated dimer states and signal the energy where frustration plays a role, therefore depleting the DOS and producing the Dirac cone.  

This leads to a linear DOS:

ρ0E=2Eπ2vF2, (21)

and to the following carrier density:

n0E=sgnE2E2π2vF2. (22)

Dirac cones are disorder protected by topology [41]. Graphene nanoribbons are studied using similar techniques leading to gaps that depend upon the width allowing to design electronic devices [42].

3.1.1. Rise of Dirac cones as a frustration effect in triangular lattices

An important and legitimate question is, after all, why do physically Dirac cones appear. Not many people ask this important question, so here we discuss how this is a consequence of the bipartite nature of the lattice and, more importantly, to the wavefunction frustration due to the triangular symmetry [38,43]. To see this, consider the squared Hamiltonian obtained from Eq. (13):

H02=HAB2k00HAB2kakbk=E2kakbk. (23)

Now the wave function components on the A and B sublattices are decoupled. Therefore H 2 0 describes a triangular lattice as it removes one bipartite sublattice [38,43]. As illustrated in Fig. 4, the spectrum of H 2 0 is obtained by folding around E = 0 the H 0 spectrum.

FIGURE 4 Sketch of the Hamiltonian eigenvalue renormalization from a graphene hexagonal lattice into a triangular one by the transformation H 2: the graphene’s density of states ρ(E) is transformed into ρ(E 2), resulting in a folding around E = 0 that is indicated by arrows. Band edges, central states, and phase differences among sites are represented by ± signs. Central states at E = 0 have a zero amplitude in one sublattice [37]. When one of the sublattices is renormalized, states near E = 0 result in edge band states with an antibonding nature in a triangular lattice [38,43], as indicated in the triangle that appears inside the hexagon. Due to frustration, states are pushed to higher energies, leading to van Hove’s degeneracies seen as peaks.  

The eigenstates of H 0 at E = 0 have zero amplitude in one bipartite lattice (see Fig. 4]. A state with E = 0 is a ground state of H 2 0. Nearby states have an antibonding nature in the sense phases differences are maximal for neighbors. We can see this by writing Eq. (17) without any reference to the first neighbor vectors δ n by using Eq. (4),

E2k=t023+2coska1+2coska2+2coska1-a2. (24)

Equation (24) gives the energy in terms of one of the triangular sublattice amplitudes and phases. The cosine terms in Eq. (24) give the bond contribution due to triangles while the first term is the self-energy of H 2 0, which is always the local coordination [40], in this case, Z = 3.

The ground state of Eq. (24) can be estimated from a variational procedure by considering a wave function with a maximal phase difference to decrease frustration, such as,

ka1=ka2=π, (25)

from where

ka1-a2=0. (26)

While in Eq. (24), the first two cosine terms will give each a contribution −2t 0, the third cosine increases the energy by 2t 0, and is thus frustrated, as shown in Fig. 4. This results in E(k) = ±t 0, which coincides with the Van Hove singularities of graphene. They are made from localized uncoupled dimers, disconnected from the lattice by having neighbors with zero amplitude. The degeneracy is given by the number of dimers, which is always N/2.

To further reduce the energy below E 2 = 1, the frustration needs also to be reduced. To see this, again, we use a variational wave function by proposing in Eq. (24) the following phases,

ka1=ϕ,ka2=-ϕka1-a2=2ϕ, (27)

where ϕ = 2π/3, and finally we obtain E(k)2 = 0. It turns out that this special wavevector is precisely the point k = K+. We get the other Dirac point k = K by choosing the other possible combination of phases. As for E = 0, the states are pseudo-spin polarized, the minimal frustrated state is four-fold degenerated. Then, as the energy goes from E 2 = t 2 0 to E 2 = 0, the DOS goes from a massive degenerated state to a four-fold degenerated state. For disordered systems, the probability of finding regions without frustration decreases exponentially, and a Lifshitz tail appears. In graphene, it leads to the Dirac cone where the DOS goes linearly to zero towards E → 0. Once graphene’s periodicity is broken by impurities or edges, confined zero-energy modes appear.

These E = 0 flat-band modes appear due to disorder or boundaries and are associated with topological properties [37,40]. Zero-energy modes are strictly localized and confined [40,44]. For doped graphene, the number of zeroenergy modes is obtained by summing disordered configurations [37]. Such flat bands are especially prominent at graphene on top of twisted graphene, inducing superconductivity [11]. The Lifshitz tail produces a pseudo mobility edge near the Dirac cones [43,45], as confirmed in ARPES experiments [46]. Resonant states are seen near the Fermi energy [45,47,48], and the wavefunctions are multifractal [38].

To further quantify the frustration, especially for the continuous models, notice how the electron dispersion depends only on the phase difference between neighboring bonds. Using Eq. (18), this suggests defining a function to measure frustration,

gkrj=±akrjbk*rj=±eikδ1N, (28)

where the label δ1 was dropped in b k(r j + δ1) as on to each site in A, we assign only one neighbor site B. The minus sign is now irrelevant in the squared Hamiltonian. Now we sum all bonds contributions and observe that,

fk2=12rjs=16gk*rjgkrj+as+3. (29)

In the previous expression, a s denotes second neighbors in H 0, which are first-neighbors in H02, i.e., a s = ±a1 ,±a2 ,±(a1 − a2). The last term in Eq. (29) is the self-energy obtained from the the interaction of the phases at the same site. Notice that it corresponds to the coordination Z = 3, and for quasiperiodic or disordered lattices is the local coordination, acting as a space-dependent repulsive potential [40,49]. Using the energy dispersion,

E2k-3t02=t02Fk, (30)

where the bonding (frustration) plus antibonding contributions are,

Fk=12rjs=16gk*rjgkrj+as. (31)

The ground state E = 0 is obtained when the frustration and antibonding contributions balance the self-energy in such a way that F(K±) = −3. The frustration increases as the momentum depart from the K± points.

3.1.2. Low-energy approximation: Dirac equation

As what matters most is electrons around the Fermi energy, it is very useful to find an effective hamiltonian by considering k = K± + q. The expansion of the hamiltonian Eq. (14) up to first order in q and around K+ is [50,51]:

HK+=vF0qx-iqyqx+iqy0=vFσp, (32)

where σ = (σ x y ) are Pauli matrices. The momentum is p = ℏq. A similar expansion can be made around K to give,

HK-=HK+t=vFσ*p, (33)

with σ = (σ x ,σ y ). Now replace q → −i ℏ ∇ in correspondence to the k · p approximation. The hamiltonians (32) and (33) are two-dimensional Dirac-like hamiltonian for massless fermions [52]. Notice that we used the word Diraclike as in a strict sense, there is no Dirac equation in; despite this, most people refer to this 2D version as the Dirac hamiltonian. The description given so far was in terms of Pauli matrices. They operate on spinors that here describe the wave function on each sublattice instead of a real spin. The term pseudospin is thus used to call this sublattice degree of freedom. In fact, before we kept both valleys separately, but the full structure of the hamiltonian taken into account both valleys is,

H=vF0Π00Π000000Π00Π0, (34)

where Π = ℏ (p x +ip y ), and H operates into a bispinor wavefunction. Their components now describe pseudospin and valley,

(|ΨK+,A,|ΨK+,B,|ΨK-,A,|ΨK-,B)T. (35)

The full hamiltonian structure is sketched out in Fig. 5. One interesting property is chirality. To see this, write Eq. (32) and Eq. (33) as,

HK+=vFph^2cmHK-=vFphT^, (36)

where h^ is the helicity operator,

h^=σ^p^p. (37)

FIGURE 5 Energy dispersion obtained from Eq. (17) considering all degrees of freedom: spin, pseudo spin, and valley. The arrows indicate the direction of the pseudospin along the with momentum on each Dirac cone and each band. For a given Dirac cone, the helicity is inverted on each band.  

It represents the projection of pseudospin on the momentum, and as h^ commutes with H, is a conserved quantity. This is understood by writing p = |p|exp[ p ], with θ p = tan−1(p y /p x ). Then Eq. (32) and Eq. (33) are written as

Hζp=vFph^ζ, (38)

and the valley helicity (chirality) operator is

h^ζ=0e-iζθpeiζθp0. (39)

The index ζ = 1 indicates valley K + and ζ = −1 valley K−. ζ is known as the valley degree of freedom. The helicity operator has eigenvalues h 1 = 1 and h −1 = −1 with eigenvectors,

|Ψζ,s=121seiζθp. (40)

Figure 5 presents a short summary of how the pseudospin, band, and valley are related to each other in the band structure.

3.2. Measuring frustration in continuous models

Let us show how frustration ideas translate into continuous low-energy models [53]. Consider first the Dirac low-energy model. To evaluate frustration, we look for the least frustrated state, which occurs for K±. Then we expand F(k+) using k = K+ + q where qa 0 ≪ 1,

Fk=FK++vF/t02q2. (41)

The Fermi velocity can be thus interpreted as the rate of frustration growth with q. For other 2D materials, it is useful to introduce a technique to measure frustration in such a way that although the energy dispersion may not be known, the ground state is reachable by using a variational procedure [53].

The method relies on noticing that still for continuous models, the phases between neighbors are given by the components of the spinor, and it follows that gkr=ψkrχk*r where the site index j was dropped as we deal with a continuous. Then we compare phases with the renormalized lattice first-neighbors a s , i.e., with g k(r+a s ). Using a Taylor series and by summing over neighbors, gkr+a¯s results in,

FkAd2rgk*rgkr+12s=1Z2a¯sgkr, (42)

where Z 2 is the number of neighbors and the integral is made in the unitary cell. The ground state is thus the one with ∇g k(r) = 0.

3.2.1. Graphene’s nearly free electron bands: charge images near a conducting plane

While many of graphene’s electronic properties are extremely well described by the Dirac approach or the TB approximation, it turns out that not all properties are well described under such approximations. Figure 6 presents a numerical calculation of graphene’s bands showing parabolic nearlyfree electron bands 3 eV above the Fermi energy [54-57]. The explanation of why this happens comes from a simple physical mechanism. If we think of graphene as a conducting plane, electronic clouds outside the plane will induce image charges [54,58-60]. When graphene interacts with substrates, adatoms, free molecules, etc., such states are important. As an example, graphene can become metallic with water, as the dipolar water layer lowers the nearly-free electron charge until it touches the Fermi energy [61]. This mechanism is akin to the flat band seen in twisted graphene bilayers, and we believe that such mechanism is in play for superconductivity in doped graphene laminates [61-64]. Further evidence is provided by recent reports of high-temperature superconductivity for water-treated graphite [65,66]. By using a simple model of an electron interacting with its charged image, a model akin to an effective Hydrogen atom is found [59]. The energy dispersion is,

Enk=k22m*+En, (43)

here k|| is the wave-vector parallel to the graphene plane, m is the effective mass.

FIGURE 6 Graphene bands evaluated numerically using the full potential linearized plane wave method (FP-LAPW) and the localdensity approximation (LDA). The lines are single graphene bands, while the gray background corresponds to the continuous spectrum. Each band is classified according to the symmetry group representations. The nearly parabolic free-electron bands are indicated with labels B 2u . Their shape is well approximated by Eq. (43) and Eq. (44) using n = 0 and n = 1. Reprinted figure with permission from Ref. [55]. Copyright (2012) by the American Physical Society.  

E n is found by solving the Schrödinger equation of an electron attracted by its image charge [67]. This leads to the Rydberg like series of states seen in Fig. 6. To better treat the regions close to the plane, LDA calculations are needed and matched with the analytical solutions, resulting in [54,58],

En=-0.85eV/n+a2, (44)

where a 2 is the matching distance of the analytical and numeric solution. This Rydberg series has been confirmed by photon photoemission spectroscopy [68]. A very important warning is in play here because the long tail of the image charge Coulomb potential is not well described in DFT calculations [60], a point ignored by many workers who model interactions of graphene with other materials.

3.2.2. Disorder effects

The coupling of the valley and pseudospin degrees of freedom is important in charge transport. As seen in Fig. 5, conduction electrons near a given valley, say the K valley, have momenta anti-parallel to the pseudospin. For back-scattering to occur, the momenta need to be reversed, i.e, to change from p to −p. This requires to invert the pseudospin and therefore valley, but valleys are far away in momentum, so backscattering is not possible, explaining the graphene’s long coherence lengths [24].

Disorder effects are understood through the study of the Dirac equation and the induced broken symmetries [69,70]. Despite this, the Dirac equation does not capture all effects of the disorder. Vacancies and impurities with enough selfenergy can produce backscattering [43,46]. Based on the frustration picture, a quasi-mobility edge was predicted in doped graphene [43] and was experimentally confirmed [46]. This topic is still controversial, as according to the orthodox theory, no mobility edges are seen in 2D systems. However, graphene has many symmetries, and disorder induces multifractal power-law localized states [37], for which the general theory used may not be valid.

3.2.3. Electromagnetic, strain and flexural fields

To understand the interaction of 2D materials with electromagnetic fields, there are two paths. One is to use perturbation techniques in which the field does not induce changes in the band structure. The first approach allows to calculate the transmittance, optical, and ac electrical conductivity for weak fields [19,71]. One can also include secondnearest-neighbors in the optical conductivity calculations using the non-equilibrium Green’s function formalism [72]. For graphene, the transmittance is almost frequency independent.

Here we do not follow these paths are they can be found elsewhere [71,72]. Instead, we focus on using a Floquet approach that fully includes the space-time periodicity of the external field and material [73,74].

Also, related to this topic, we add that strain and flexural deformations enter as a pseudo-electromagnetic field in the Dirac equation, and therefore, many results found for the electromagnetic theory can be applied to strain fields [20]. Despite its similarities, a real electromagnetic field breaks the time-reversal symmetry while a pseudo electromagnetic field does not. From a formal point of view, this only changes some signs on each valley; however, the resulting currents induced by the pseudo fields on each valley are canceled out. This not happens for electromagnetic fields, as the currents have the same sign on each valley [73,74].

We remark here that in 2013 we showed that the basic equations used in the continuous low-energy Dirac approximation had a problem as they were not able to solve even the simplest case, an isotropic uniform expansion [14]. The reason was found in Ref. [14]; when making the low-energy approximation in the tight-binding problem, one needs first find the Dirac points and then perform the approximation, as strain moves the Dirac points. Other groups done the approximation in the inverse order, and thus many results were wrong. This problem was fixed in Ref. [14] by producing the proper form of the hamiltonian for uniform strain. From there, a solution for general strain fields was constructed [16].

The analogy between strain and electromagnetic fields is very useful; for example, one can find solutions for Dirac fermions in magnetic fields [75] and translate them to strain. Strain leads to many effects, for example, a light-induced Faraday effect [76], tunable dichroism [17,77]. Moreover, strain and real magnetic fields can be combined to build new kinds of devices; for example, there is a theory that proposed quantum engine based on the gauge fields formalism [78] or valley-polarization of currents by nanobubbles [79].

We started in 2007 working in non-perturbative methods for electromagnetic waves [73,74]. Afterward, some results were extended for deformations. For example, ripples in graphene can be studied by using the minimal coupling in the Dirac hamiltonian as if it were an electromagnetic field, [73,74,80],

H=vFσηp^-ηA(r,t)+Vr,t, (45)

where η = ±1 labels the K + and K + Dirac points. A and V are pseudo-vector and pseudo-scalar potentials [20],

Vr,t=gεxx+εyy, (46)

Ar,t=Ax,Ay=β2accεxx-εyy,-2εxy, (47)

where the strain is characterized by the 2×2 tensor ε µν with components µ = x,y and ν = x,y, defined below. The dimensionless Grüneisen coefficient β3.0 gives the magnitude of the coupling between deformations and hopping parameter, and the parameter g ranges from 0 to 20 eV. Consider a deformation wave made of an out of plane displacement h and an in-plane displacement u. The strain tensor is,

εμν=12μhνh+12μuν+νuμ. (48)

The out-of plane deformation is,

hr,t=h0cosϕ, (49)

with a phase of the wave given by,

ϕ=Qr-Ωt. (50)

h(r,t) measures the deformation with respect to the normal of flat graphene and Q is the deformation wave vector with components (Q 1 ,Q 2). The time-frequency of the strain is given by Ω, which is related with Q = |Q| as Ω/Q = v s , where v s is the deformation propagation velocity, usually the velocity of sound. Then in-plane strain can be written is a similar way,

ur,t=ϵ0r+uccosϕ. (51)

The electron dynamics is now dictated by the equation,

itΨη(r,t)=&#91;vFση(p^-ηA(r,t)+V(r,t)&#93;Ψη(r,t). (52)

The method to find its solution is the same that we used to solve the analogous of Eq. (52) for electromagnetic fields [73,74]. The idea is to combine a Floquet theory with a Volkov approach, i.e., the solution is made from a plane wave and a function that depends only on the phase of the field. This is equivalent to solve the Dirac equation in the reference frame of the wave [73,74],

Ψρ=expikr-iEtΦρϕ, (53)

where Φ ρ (ϕ) is a function to be determined for each sublattice ρ = A,B, and E = v F ~|k|. This ansatz transforms the partial differential equations (52) into an ordinary secondorder differential equation. In all cases that we studied, the resulting effective equation turned out to be a Mathieu equation that describes a classical parametric pendulum [73,80]. Here we offer a simplified version for the case of an electron with momentum parallel to the wave propagation,

d2Γρϕdϕ2+a±-2qcos2ϕΓρϕ=0, (54)

where Γρϕ=eik~ϕΦρϕ and the signs ± labels each valley equation. The vector k~=k~x,k~y, defined as k~=k/Q, is projected into parallel and perpendicular directions of the propagating corrugation. The parameters a ± and q depend upon the momentum and field intensities [73] but can be better understood in terms of the parametric pendulum. a ± represents the ratio between the frequencies of the pendulum and forcing, while q is the amplitude of forcing. In Fig. 7, we present the stability regions of the Mathieu equation. Here represents the allowed spectrum of periodic solutions and is made from resonance gaps known as Arnold tongues. The solutions are the Mathieu special functions.

FIGURE 7 (Color online) Spectrum of the Mathieu equation showing the allowed and forbidden (shaded area) regions. The order of the gaps is indicated. This structure is known as the Arnold tongues. This diagram represents the stability of solutions for a parametric pendulum, where a ± is the ratio between the period of forcing and the natural pendulum frequency, and q is the magnitude of the external forcing.  

Whenever the frequency ratio is an integer, a gap is open in the spectrum, and resonance occurs. In the present approach, the role of the natural frequency is played by the energy of an electron with a given momentum, and the forcing frequency is given by the external field. q is proportional to the field intensity; whenever q = 0, we recover the case of graphene without a field. An algebraic approach to similar problems leads to the same conclusion [81].

The important result here is that light and strain open gaps. Light induces currents and even changes the sign of the current depending on the type of polarization [73]. Moreover, we found a strong non-linear response in the sense that an ac field induces a huge content of odd harmonics allowing graphene to work as a rectifier [73]. It is important to remark that for some angles of light incidence, q can be complex, and thus the diagram looks different [74]. An interesting effect is an electron focusing by sound as gaps are open in some preferred directions [21,80,82].

3.3. Other 2D materials effective low-energy hamiltonians

Graphene paved the way to look for other effective hamiltonians. The procedure is to solve a tight-binding or perform a DFT simulation and then use an approximation around the Fermi energy. Typical examples are Silicene and Borophene. For the borophene phase with a space group 8-Pmmn, the model hamiltonian is [83]

H^=ζvxkxσ^x+vykyσy^+vtkyσ^0. (55)

Notice that here there are three different velocities along each coordinate axis, where by {v x ,v y ,v t } = {0.86,0.69,0.32} in units of v F = 106 m/s, and the valley index is ζ = ±1 degree of freedom, (σ x y ) and σ 0 is the identity Pauli matrix. The resulting bands are,

Eλ,kζ=ζvtky+λvx2kx2+vy2ky2, (56)

with eigenfunctions,

ψλ,kζ= ζei k r2 ( 1 λeiΘ) (57)

λ = ±1 denotes the valence or conduction band as in graphene and Θ = tan−1(v y k y /v x k x ) is a measure of the anisotropy. In Fig. 8, we plot Eq. (57). Two features are seen here, i) the cone is titled and ii) there is a strong anisotropy. The tilting is produced by the v t σ 0 term, while v y /v x controls the anisotropy. Both effects imply many restrictions for optical transitions, producing a very transparent material and dynamical gap opening by electromagnetic radiation [84]. To study this system under strong electromagnetic fields as well as any other periodically time-driven quantum system, we developed a very efficient and simple monodromy approach [85], which numerically is an improvement over our original analytic algebraic method [81]. The method also reproduces the results of linearly polarized light at all field intensities, for which we were able to find an exact solution [86]. Finally, this kind of hamiltonian holds in general for other anisotropic honeycomb lattices, including mechanical, acoustic, microwave, etc., analogous systems [77].

FIGURE 8 Energy dispersion in borophene resulting from the hamiltonian (55) near one of the Dirac points. The cone is tilted and anisotropic.  

Interesting models are obtained by decorating graphene. An example is the α − T3 graphene, in which an additional molecule or atom enters at the center of each hexagon, coupling with just one bipartite sublattice [87-89]. Interestingly, this is the minimal model that presents flat bands coexisting with Dirac cone states. The α − T3 Hamiltonian low energy approximation is [90],

H^=0t0Cαfk0t0Cαf*k0t0Sαfk0t0Sαf*k0, (58)

with Cα=1/1+α2, Sα=α/1+α2, and the parameter α is in the interval 0 ≤ α ≤ 1. The model contains two inequivalent Dirac points and thus can be written as,

H^ξ0=vFkS^, (59)

where S^=ξS^x,S^y with ξ = ±1 is the valley index. The pseudo-spin operators S^x and S^yare,

S^x=0Cα0Cα0Sα0Sα0, (60)

S^y=0-iCα0iCα0-iSα0iSα0. (61)

When considering the optical and electronic properties, this model behaves as a three-level Rabi problem [91].

2D Materials can also be periodically modulated when interacting with substrates. An example is the Kekule pattern seen in graphene [92]. The model can be written as [93],

H=v0kσΔ0kx-ikyσ0Δ0kx+ikyσ0kσ, (62)

or H = v 0(k·σ)⊗τ 0+v F 0 σ 0⊗(k·τ), where τ = (τ x y ) is a second pair of Pauli matrices, τ 0 the unitary matrix, and ∆0 is the energetic parameter that measures the strength of the bond modulation. The corresponding spectrum is

ϵkαβ=αsβv0k, (63)

where α = ±1 labels the conduction and valence bands. β = ±1 is a label used to define two velocities, s β = (1 + β0).

The energy dispersion folds each graphene’s valleys into the Γ point of the Brillouin zone. This results in two cones that have the same tip, but one cone has a fast Fermi velocity v F (1 + ∆0), and the other a slow velocity v F (1 − ∆0). The eigenfunctions are [93,94],

Ψαα'k=Ψαk|Ψα'k, (64)

where |Ψ α i is the graphene’s single-valley eigenvector. This model has several important features. As both valleys have the same Dirac point in momentum space, it allows to perform valleytronics by strain [18], a fact recently confirmed experimentally [95]. Also, clear signatures are seen in the optical/electronic conductivities [94,96] and plasmons due to the Moiré interference between the two-electron gases with slightly different velocities [30]. The inclusion of secondneighbor interactions modifies the picture as a gap opens in one of the cones [97].

4. Multilayered 2D materials

Multilayered materials are made by the successive stacking of 2D layers [98,99]. Figure 9 presents the results of stacking two and three graphenes. In Fig. 9b)), the B atoms of layer 2 (B2) lie on top of layer 1 A atoms of layer 1 [100,101]. Figure 9c) presents the structure of an ABA trilayer. We can identify in Fig. 9 that the electronic properties of 2D stacked materials are ruled by:

  • Intralayer interactions and geometry on each monolayer.

  • The kind of interactions between monolayers, known as interlayer interactions.

  • The stacking geometry, defined by a translation and rotation.

FIGURE 9 Unstrained lattice structure (upper panels) and a sketch of their corresponding energy dispersion (lower panels) for a) monolayer graphene, b) Bernal stacked bilayer graphene, and c) Bernal stacked trilayer graphene. Blue atoms belong to the A bipartite lattice, while red atoms belong to the B lattice. The arrows indicate the different kinds of interactions that appear in a TB calculation, parametrized by t 0 for intra-layer interaction, and γ i with i = 1,...,5 for inter-layer interactions. The Dirac cone seen in a) for graphene is replaced by parabolic bands in b), while trilayer graphene includes both types of bands.  

In graphene over graphene, the interlayer interaction is due to Van der Waals forces, but the stacking geometry leads to complex electronic phases [11]. For graphene stacked over a metallic substrate as Au, Ag, Cu, etc., strong hybridization is seen. Even impurities can induce bond order as the Kekule patterns [102].

A typical example of a heterostructure hamiltonian is the Slonczewski-Weiss model [103]. Therein, the parameter t 0 is the usual intra-layer graphene’s interaction, as indicated in Fig. 9b). For bilayer graphene, five hoppings γ i with i = 0,...,4 account for the possible overlaps of π-orbitals.

Four on-site energies ϵA1,ϵB1,ϵA2,ϵA2 are used. The parameters are given by: t 0 = 3.16 eV, γ 1 = 0.381 eV, γ 3 = 0.38 eV, γ 4 = 0.14 eV, ϵB1=ϵA2=0.022 eV and ϵA1=ϵB2=0. As Fig. 9c) indicates, the case of trilayer graphene requires two additional parameters (γ 2 and γ 5). The resulting matrix for bilayer graphene is:

H=ϵA1-t0fkγ4fk-γ3f*k-t0f*kϵB1γ1γ4fkγ4f*kγ1ϵA2-t0fk-γ3fkγ4f*k-t0f*kϵB2, (65)

where f(k) is given by Eq. (15). The bands resulting from diagonalization of such matrix are presented in Fig. 9b). Trilayer graphene is treated in a similar way. The energy dispersion seen in Fig. 9c) mixes the monolayer and bilayer case.

5. Superlattices: twisted and displaced graphene bilayers

Now consider the effects of the stacking geometry due to displacements or rotations. Structures known as Moiré pat- terns are produced. At certain angles, periodic superlattice are seen with a unitary cell which usually ranges from 1 to 30 nm [104,105]. Usually, the strain also appears to minimize the elastic free energy. A complete treatment of the problem requires to include not only the electronic energy but the elastic energy as well [13,14,106]. Another interesting way to build superlattices is to consider spatial variations of external electric and magnetic fields [107-110].

Needless to say, the resulting physics is rich. As examples, we have the first experimental observation of the complete Quantum Hall topological phase diagram [111], multiflavored Dirac fermions [31], and quantum phases similar to those seen in high T c superconductors [11]. This last surprising behavior is due to the appearance of flat bands at certain geometrical configurations. Both characteristic behaviors are seen in their most simple form by using a model produced by our group many years ago and that we reproduce in the following subsection [13].

5.1. Quantum Hall effect, flat bands and topology in graphene superlattices

An example of clever use of a Moiré superlattice is the experimental observation of the Hofstadter butterfly [111]. Such fractal spectrum was predicted to occur for electrons in a 2D lattice under a constant magnetic field [112]. From there, the Quantum Hall effect (QHE) was explained in terms of topological phases [113]. Recently, this phase diagram was found using higher-dimensional projections [114], or even more surprisingly, from symbolic sequences, billiards, or Sturmian codings of the magnetic flux, allowing to decode the global fractality of the spectrum and making clear the connection with quasicrystals [115]. Originally, the problem without a lattice was studied by Landau, giving levels with energy E = (n + 1/2)ℏω and n integer. In a lattice, the interatomic distance adds a new scale in the problem that competes with the magnetic length [112]. The spectrum is controlled by the ratio between the magnetic flux and elementary quantum flux. The problem translates into a one dimensional Harper equation, similar to that seen for uniaxial strained graphene [13]. For atomic systems, it was impossible to measure the Hofstadter butterfly as the ratio between fluxes is too low for any real magnetic field. The trick was to build a superlattice using a rotated substrate [111]. This allows increasing the unitary cell size and thus the flux by a dramatic amount while keeping the magnetic field accessible to available laboratory conditions [111].

Let us present a simple model able to capture such behavior. As shown in Fig. 10, suppose graphene over a substrate in which the interaction will induce deformation of graphene. We suppose that the deformation is position-dependent only in one direction, a situation known as uniaxial strain. One can take advantage of the uniaxial property by noting that graphene nanoribbons hamiltonians map into one dimensional effective chains [13] as sketched out in Fig. 10. Consider a nanoribbon in which the atoms’ positions r j are transformed by the deformation into r j +au(r j ), where u(r j ) is any general strain field. Here we demand the field to be uniaxial u(r) = (0,u y (y)). As seen in Fig. 10, the symmetry along the non-strained direction, chosen as the x direction, is not broken. The solution to the Schrödinger equation can be proposed as Ψ(r’) = exp(ik x x)ψ(y). The resulting hamiltonian thus depends upon k x and is labeled as H(k x ).

FIGURE 10 (Color online) A zig-zag nanoribbon modulated by uniaxial strain as described by Eq. (66). The wavy line to the left is the applied strain in the y-direction. The dotted line represents the period in the x direction. A vertical dotted line indicates the limits of the unitary cell in the x direction. Carbon atoms positions in the cells are labeled with A j and B j , corresponding to the bipartite lattice and j is the atom-numbering in the vertical direction. The system maps into a chain joined by effective bonds t j and c(k x )t j . At a special k x value, c(k x ) = 0 and the system decouple into dimers.  

The non-strained case can be thought in terms of cells with four kind of sites (see Fig. 10] and focusing attention to the right of Fig. 10, one notes that all sites within the periodic supercell become inequivalent along the y direction when strain is applied. One labels them as A j and B j , where j is an index for the position in the path, and A or B labels the original bipartite lattice in the absence of strain.

The resulting hamiltonian is that of a one-dimensional modulated chain [13]:

Hkx=j=1N-1tj+1aj+1bj+ckxtjajbj+H.c., (66)

where a j , a j and b j , b j are the annihilation and creation operators in the A and B sublattices, respectively, and N is the number of sites in the A sublattice along the periodic path. For odd j, the effective bonds are defined through:

tj=t0exp-βuyBj-uyAj/a, (67)

and for even j as:

tj=t0exp-βuyAj+1-uyBj/2a, (68)

where the Grunesissen parameter is β ≈ 3. The factor c(k x ) contains the phase in the x−direction:

ckx=2cos3kxa/2. (69)

Figure 11 presents the resulting energy dispersion E(k x ) and the corresponding DOS for the hamiltonian (11) using a periodic substrate using a strain field u(y) = λ(0,cos(2πΛy)). Here Λ controls the spatial wavelength of the deformation and λ its amplitude. Figure 11 panels c) and d) presents the spectrum and DOS for an irrational Λ with respect to the graphene’s unitary cell distance. The spectrum presents a flat band. The model allows to track the origin of the celebrated flat-band to find that they are due to lines of quantum dots which confine states [13]. Also, the Van Hove singularities of unstrained graphene (panels a) and b)) spread along the spectrum producing peaks. For a conmmensurate Λ we can obtain the spectrum and DOS shown in Fig. 11 panels e) and f). The spectrum shows a gap, and the DOS indicates that the system decouples into linear chains, resulting in a Luttinger liquid. Figure 11 indicates the amazing possibilities as we can tune the system to be semimetallic, metallic, or insulator by using strain or a suitable substrate.

FIGURE 11 (Color online) Left columns, energy as a function of k x , and right columns, the corresponding DOS for uniaxial strain in graphene as described by the hamiltonian Eq. (66). Panels a) and b) correspond to unstrained graphene, panels c) and d) contain an incommensurate strain for the unit cell. Panels e) and f) are for conmmensurate strain. On each DOS we include labels to denote the conductivity behavior: semimetallic for graphene, metallic for incommensurate strain and insulator for commensurate strain. The Van Hove singularities and Dirac cones are indicated, as well as the flat band that appears at E = 0 for c) and d). A vertical line going through panels a),c) and d) indicate the special moment k x = 2/√ 3a where the equations are decoupled. Observe how the Van Hove singularity for graphene is a highly-degenerated dimerized system, while the effect of the strain for the incommensurate case is to spread the states forming many peaks and a flat band due to lines of quantum dots. In case f), the commensurate strain produces a gap and Van Hove singularities at band edges. This case results in two effective linear chains forming a Luttinger liquid.  

An interesting situation arises for kx=2/3a, where c(k x ) = 0. In this case, the chain is decoupled into dimers [13]. For unstrained graphene, the dimers produce a massive degeneration leading to the Van Hove singularities observed at E = ±3t 0. They are saddle points of the energy dispersion. For strained graphene, the degeneracy can be completely or partially removed [13], leading to the peaks observed in panels c). When the substrate is another strained and displaced graphene, the resulting hamiltonian is the same as the Harper equation that appears in the Quantum Hall effect problem [112]. Such equation is known by mathematicians as the Almost Mathieu Operator, as is a discretized version of the Mathieu equation which describes the parametric pendulum. The nature of its spectrum for the incommensurate case was one of the famous Simon’s Ten Martini problems, finally solved in 1996 by Jitomiskaya and Last [116].

The previous model contains interesting consequences for electronic localization [13] and non-trivial topological properties. To show the non-trivial topology, consider again a strain field u(y) = λ(0,cos(2πΛy)). For a wavelength such that Λ = 1/(2a), t j takes only two values. By performing a Fourier transform of Eq. (66) using:

aj=1N/2kye-ikyj3/2aky, (70)

and

bj=1N/2kye-ikyj3/2bky, (71)

the hamiltonian becomes:

Hk=hxkσx+hykσy, (72)

where σ x and σ y are the x and y Pauli matrices,

hxk=21-λcos3kx/2+1+λ/2cos3ky/2, (73)

and

hyk=1+λ/2sin3ky/2. (74)

Equation (72) is the celebrated Su-Schrieffer-Heeger model for polyethylene [117]. Therefore, non-trivial topological phases are seen for amplitudes λ > 1/2. For a gapless case (λ < 1/2), the topological modes are flat bands joining Dirac points with opposite topological charges [106,118], as happens in Weyl semimetals [119,120].

Also, this 1D model can be used to study time-dependent strain. The results lead to very complex topological phase diagrams, with new kinds of topological modes which are still in the process of being investigated [121-123].

5.2. Magic angles and superconductivity

As seen in Fig. 12, we can build a Moiré pattern by stacking graphene on top of graphene rotated by an angle θ. Let us derive the corresponding hamiltonian. For simplicity, we only consider a single valley. The hamiltonian is made from two rotated single layer Dirac hamiltonians and an interlayer coupling between them,

H+=HUUHUDHDUHDD. (75)

FIGURE 12 Panel a), twisted graphene over graphene geometry. The stacking points AA, AB, and BA (denoted by r0) are indicated. The Moiré vectors a 1,2 are indicated. Panel b), the Moiré Brillouin zone (usually known as a mini-Brillouin zone). The K point is at the origin. The reciprocal Moiré vectors are b 1,2 and the translation vector b 1,2,3 . The K’ point is given by K’ = K−q1.  

The up (U) and down (D) labels are used to denote each layer, and

HUU=vF-i-KUσθ/2HDD=vF-i-KDσ-θ/2 (76)

are the two rotated versions of the single-layer original Dirac hamiltonian. We have taken into account the shift of the K + Dirac point to produce the high-symmetry points K U and K D on each layer. Also, we defined two rotated versions of the Pauli matrices,

σθ/2=e-iθ4σzσx,σyeiθ4σz. (77)

The Dirac points are now coupled together by an interlayer hopping term such that H DU = H UD with [9],

HUD=T0rσ0+TABrσ++TBArσ- (78)

where T 0(r),T AB (r),T AB (r) are the tunelling matrix elements. To gain a better understanding of how to build these terms, we consider first the case in which electrons tunnel only when the atoms of both layers coincide or when an atom coincides with the center of a hexagon in the other layer. As seen in Fig. 12, there are three cases: both bipartite lattices coincide, as in the stacking points AA or in an equivalent way BB, the A atoms of the upper layer align with the B atoms of the layer (AB stacking) and viceversa (BA stacking). The AA stacking points form a triangular lattice and unfavour tunneling between lattices. The AB and BA points form two sublattices of a big Moiré honeycomb. Therein, the electrons are prone to interlayer tuneling, a fact that will prove to be very relevant for the phenomenology of the system.

Now consider, for example, theP AA stacking point; the tuneling is given by T0r=w0nδr-Rn, where R n is the location of the AA points and w 0 the energetic cost of tunneling there. Using the Moiré reciprocal vectors G m , the delta function is written as a Fourier harmonic superposition δr-Rn=1/AMme-iGmr, where A M is the area of the Moiré unit cell. A similar procedure holds to the other stackings, but due to their displacements from the AA points, the other stacking points acquire a phase,

TABrw1δr-Rn-rA=w1AMme-iGmrAe-iGmrTBArw1δr-Rn-rB=w1AMme-iGmrBe-iGmr, (79)

where w 1 is the energetic cost of tunneling at the AB and BA points.

In the real system, the interlayer hopping t UD (r) is a continuous function of the position, and thus we need to include in the sums over the Moiré reciprocal vectors a modulation factor t~UDGm. A very useful approach to build such hoppings is to identify the coincidence lattice from a higher dimensional cut and projection method [32]. In any case, the interlayer hopping depends on the separation between layers, which is more than twice the intra-layer atomic separations. As the hoppings depend exponentially on the distance [14], the Fourier components decrease rapidly. In the most simple model, only the lowest order harmonics are kept [9]. Also and as the original rotation symmetry is broken, it is customary to move the Dirac points in both layers to the origin by using a unitary transformation with diagonal entries: W = e iKUr ;e iKDr . This leads to the following hamiltonian,

H'UU=vF-iσθ/2H'DD=vF-iσ-θ/2. (80)

Defining the vector q 1 = K U K D = k θ (0,−1) and the vectors related by ϕ = 2π/3 rotations q 2 , q 3, the interaction between layers is now written as,

H'UD=w0U0rσ0+w1UABrσ++w1UBArσ- (81)

with the following definitions,

U0r=e-iq1r+e-iq2r+e-iq3rUABr=e-iq1r+ei2π3e-iq2r+ei4π3e-iq3rUBAr=e-iq1r+e-i2π3e-iq2r+e-i4π3e-iq3r, (82)

where the Moiré lattice vectors a 1 ,2 and reciprocal lattice vectors b 1 ,2 are (see Fig. 12],

a1,2=4π3kθ±32,12,b1,2=q1,2-q1=kθ±32,32. (83)

After all these transformations, we finally obtain the celebrated Bistritzer-Mac Donald (BM) model [9],

H=-ivσθ/2TrTr-ivσ-θ/2 (84)

where,

Tr=w0U0rw1U1rw1U1*-rw0U0rU0r=e-iq1r+e-iq2r+e-iq3rU1r=e-iq1r+eiϕe-iq2r+eiϕe-iq3r (85)

and which has the translational symmetry,

Hr+a1,2=VHrVV=diag1,e-iϕ,1,e-iϕ. (86)

The phase between the top and bottom layers acquired under translations can be understood as an effect of the graphene Brillouin zones twisting, as the momenta in the bottom layer are shifted by q 1 relative to the top layer. This generates the phase shift eiq1a1,2=e-iϕ. The Bloch States inherit this phase and have the following form,

ψkr=ukrei(k-K)rei(k-K')r, (87)

where K’ = K − q1 is the Moiré K’ point wavevector in terms of the Moiré K point wavevector and u k(r) is a periodic function such that u k(r) = u k(r ± a 1,2 ). For zero tunneling, or large angles, the states at K and K’ correspond to the top layer and bottom graphene layer, respectively.

The spectrum of the (BM) model presents flat-bands for several θ; these are known as magical angles. As flat-bands imply localization, Bistritzer-Mac Donald suggested an enhanced electron-electron interaction at these angles and the possibility of having superconductivity [9]. Such hypothesis was experimentally confirmed in 2018, and surprisingly, the quantum phase diagram was found equivalent to those seen in high critical temperature superconductors [11]. There are many open questions in this field, and one of these is the origin of magic angles. To tackle this question, in 2019, Tarnopolsky et al. introduced a simplified version of the BM hamiltonian by making a very useful approximation [12]. They switch off the AA coupling completely by setting w 0 = 0. This results in a chiral version of the original BM hamiltonian,

H=0D*-rDr0 (88)

where the zero-mode operator is,

Dr=-i¯αUrαU-r-i¯, (89)

and

D*-r=-iαU*-rαU*r-i, (90)

with ̿= x + i∂ y , = x i∂ y . The potential is,

Ur=e-iq1r+eiϕe-iq2r+e-iϕe-iq3r. (91)

The corresponding parameters are q1 = k θ (0,−1),

q2=kθ32,12

and

q3=kθ-32,12;

the Moiré modulation vector is k θ = 2k D sin(θ/2) with k D = (4π/3a 0) is the magnitude of the Dirac wave vector and a 0 is the lattice constant of monolayer graphene. The physics of this model depends upon the intensity of the parameter α, defined as α = (w 1 /v 0 k θ ). Here w 1 is the interlayer coupling of stacking AB and BA, and takes the value w 1 = 110 meV, while v 0 is the Fermi velocity, with value v 0 = (19.81 eV/2k D ). At magic angles α =0.586,2.221,3.751,5.276,6.795,8.313,9.829,11.345,... flat-bands appear. These magic α’s have a 3/2 cuantization rule [12] for α > 0.586 .

The chiral hamiltonian captures the “true magic” of the magic angle physics. At point AB the K point wave function has a node [12], and thus flat-bands wave functions are written using Jacobi theta functions, close to the lowest Landau state in the Quantum Hall effect [12]. Interestingly, using this model, one can predict the magic angles for any system with more than two layers [124]. Very recently, this result has been experimentally confirmed as superconductivity was found exactly at the predicted angles for trilayer graphene [125], resulting in the strongest coupled known superconductor with a large Pauli limit violation and reentrant superconductivity phases [10].

Also, recently, Naumis et al. showed that the chiral hamiltonian could be further reduced to a 2 × 2 matrix, bringing to the front the nature of frustration and the topology of the system [53]. To do this, we start with the Schrödinger equation HΦ=EΦ, where Φr=ψ1r,ψ2r,χ1r,χ2rT are the four components of the wave function for twisted graphene bilayer. Here the index 1,2 represent each graphene layer. Then, as we did with graphene, consider the squared hamiltonian H2,

H2=D*-rDr00DrD*-r. (92)

As before, this transformation removes the particle-hole symmetry, which is an anti-unitary anti-commuting symmetry. The resulting 2 × 2 effective hamiltonian H 2 =

D (−r)D(r) is [53],

H2=-2+α2U-r2αArαAr-2+α2Ur2. (93)

The norm of the potential is,

Ur2=3+2cosb1r-ϕ+2cosb2r+ϕ+2cosb3r+2ϕ, (94)

where b1 ,2 = q2 ,3 − q1 are the Brillouin zone Moiré vectors and b3 = q3 − q2. In Fig. 13 we present a contour plot of such effective potential, showing its trigonal symmetry and the structure of minima and maxima that produce a quantum dot effect by confining electrons. Also, notice the similarity of this function with the frustration function F(k) of graphene. The off-diagonal terms are,

Ar=-iμ=13eiqμr2q^μ-kθ, (95)

and

Ar=-iμ=13e-iqμr2q^μ+kθ, (96)

where ∇ = −∇ with ∇ = ( x ,∂ y ) and µ = 1,2,3 as eigenvalues are reals (notice that −A (−r) = A(r)). The unitary vectors q^μ are perpendicular to the set q µ . The term A(r) is a non-Abelian gauge potential as it does not commute with itself at different locations. From here, it is possible to show that at magic angles, the frustration is exactly zero [53]. Moreover, it is possible to obtain in a simple way the spectrum and the topological part of the hamiltonian, showing explicitly the broken time-reversal symmetry as well as a reduction to a friendly rectangular domain by using triangular coordinates [53]. Such renormalization is reminiscent of a phonon problem [53] as the flat-band behaves as a floppy mode band in an equivalent rigidity-phonon problem [126]. These floppy modes occur in flexible systems [127-131]. Finally, we add that some of the electronic properties of Moiré lattices can be obtained through using analogous systems, for example, fluids in a vessel with nested and rotated obstacles [132].

FIGURE 13 Contour plot of the effective potential Eq. (94). Observe the trigonal symmetry. For the first magic angles, the flatband states wave functions tend to confine by tracking the minima of this function, here in blue. The three-light brown maxima around each minimum act as potential barriers which confine the wave function, leading to a quantum dot effect.  

6. Conclusions

A general overview of the electronic properties of 2D materials was provided making contact with the Dirac and Weyl equation approaches. We made special emphasis on the underlying frustration that produces the Dirac cone, a fact usually bypassed in the literature but that our group in Mexico pursued systematically. Moreover, this physical mechanism is in play for twisted graphene bilayers at magic angles, where superconducting phases are observed. We also reviewed a very simple mapping of graphene superlattices into a one-dimensional model that allows to recover the main features of twisted graphene over graphene, like flat-bands, topological states, fractal quantum phase diagrams, etc. Finally, we provided an overview of the physics of twisted bilayer graphene, including a recent derivation that simplifies the problem and brings the physical and topological origin of flat bands to the front [53]. In this review, we emphasized the collaborations made with other groups working in Mexico. Many of these contributions not only helped to understand 2D systems but allowed to develop new techniques to deal with periodic time-driven quantum systems, timedriven topological systems and methods to quantify frustration, build topological phase diagrams and obtain its physical properties. Among these, it is of special importance in all fields of physics the finding of an exact method to build the effective Hamiltonian of quantum time-driven systems, allowing, for example, to solve analytically the quantum harmonic oscillator with time-dependent frequency, i.e., the quantum parametric pendulum, for which solutions were not known [81].

Acknowledgements

I would like to thank my always hard working and bright collaborators: the current students Constanza Cendón, Enrique Aguilar, Saul Herrera, Leonardo Navarro-Labastida, Elías Andrade, Abdiel Champo, and former Ph. D. students Maurice Oliva-Leyva, Pedro Roman-Taboada, Eduardo Barrios Vargas, Hugo Flores, and Francisco López, as well my collegues Salvador Barraza, Ramón Carrillo-Bastos, Alejandro Kunold, Jose Luis Aragón, J.C. Sandoval-Santana and Victor Ibarra for years working together. Finally, I would like to thank Alfredo Raya for the kind invitation to write this review. This work was supported by UNAM-DGAPA project IN102620 and CONACyT project 1564464.

References

1. J. H. C. Scargill. Existence of life in 2 + 1 dimensions. Phys. Rev. Research, 2 (2020) 013217. https://doi.org/10.1103/PhysRevResearch.2.013217 [ Links ]

2. P. W. Anderson. More is different. Science, 177 (1972) 393. 10.1126/science.177.4047.393 [ Links ]

3. K. S. Novoselov, A. K. Geim, S. V. Morozov, D. Jiang, Y. Zhang, S. V. Dubonos, I. V. Grigorieva, and A. A. Firsov. Electric field effect in atomically thin carbon films. Science, 306 (2004) 666. 10.1126/science.1102896 [ Links ]

4. Alexander A. Balandin, Suchismita Ghosh, Wenzhong Bao, Irene Calizo, Desalegne Teweldebrhan, Feng Miao, and Chun Ning Lau. Superior thermal conductivity of single-layer graphene. Nano Lett., 8 (2008) 902. https://doi.org/10.1021/nl0731872 [ Links ]

5. Changgu Lee, Xiaoding Wei, Jeffrey W. Kysar, and James Hone. Measurement of the elastic properties and intrinsic strength of monolayer graphene. Science, 321 (2008) 385. 10.1126/science.1157996 [ Links ]

6. K. S. Novoselov , D. Jiang , F. Schedin, T. J. Booth, V. V. Khotkevich, S. V. Morozov , and A. K. Geim . Two-dimensional gas of massless dirac fermions in graphene. Proc. Natl. Acad. Sci. USA, 102 (2005) 10451. https://doi.org/10.1038/nature04233 [ Links ]

7. Novoselov K. S., Falko V. I., Colombo L., Gellert P. R., Schwab M. G., and Kim K. A roadmap for graphene. Nature, 490 (2012) 192. https://doi.org/10.1038/nature11458 [ Links ]

8. Rafael Roldán, Andrés Castellanos-Gómez, Emmanuele Cappelluti, and Francisco Guinea. Strain engineering in semiconducting two-dimensional crystals. J. of Phys.: Condens. Matter, 27 (2015) 313201. https://doi.org/10.1088/0953-8984/27/31/313201 [ Links ]

9. Rafi Bistritzer and Allan H. MacDonald. Moiré bands in twisted double-layer graphene. Proceedings of the National Academy of Sciences, 108 (2011) 12233. https://doi.org/10.1073/pnas.1108174108 [ Links ]

10. Yuan Cao, Jeong Min Park, Kenji Watanabe, Takashi Taniguchi, and Pablo Jarillo-Herrero. Large pauli limit violation and reentrant superconductivity in magic-angle twisted trilayer graphene, 2021. https://doi.org/10.1038/s41586-021-03685-y [ Links ]

11. Valla Fang Cao, Yuan Fatemi, Shiang, Kenji Watanabe , Takashi Taniguchi , Efthimios Kaxiras, andPablo Jarillo-Herrero . Unconventional superconductivity in magic-angle graphene superlattices. Nature, 556 (2018) 43. https://doi.org/10.1038/nature26160 [ Links ]

12. Grigory Tarnopolsky, Alex Jura Kruchkov, and Ashvin Vishwanath. Origin of magic angles in twisted bilayer graphene. Phys. Rev. Lett., 122 (2019) 106405. https://doi.org/10.1103/PhysRevLett.122.106405 [ Links ]

13. Gerardo G. Naumis and Pedro Roman-Taboada. Mapping of strained graphene into onedimensional Hamiltonians: Quasicrystals and modulated crystals. Phys. Rev. B, 89 (2014) 241404. https://doi.org/10.1103/PhysRevB.89.241404 [ Links ]

14. M. Oliva-Leyva and Gerardo G. Naumis. Understanding electron behavior in strained graphene as a reciprocal space distortion. Phys. Rev. B , 88 (2013) 085430. https://doi.org/10.1103/PhysRevB.88.085430 [ Links ]

15. M. Oliva-Leyva and G. G. Naumis. Anisotropic AC conductivity of strained graphene. J. Phys.: Condens. Matter, 26 (2014) 125302. https://doi.org/10.1088/0953-8984/26/12/125302 [ Links ]

16. M. Oliva-Leyva and Gerardo G. Naumis. Generalizing the fermi velocity of strained graphene from uniform to nonuniform strain. Phys. Lett. A, 379 (2015) 2645. https://doi.org/10.1016/j.physleta.2015.05.039 [ Links ]

17. M. Oliva-Leyva andG. G. Naumis . Tunable dichroism and optical absorption of graphene by strain engineering. 2D Mater., 2 (2015) 025001. https://doi.org/10.1088/2053-1583/2/2/025001 [ Links ]

18. Elias Andrade, Ramon Carrillo-Bastos, and Gerardo G. Naumis. Valley engineering by strain in kekulé-distorted graphene. Phys. Rev. B , 99 (2019) 035411. https://doi.org/10.1103/PhysRevB.99.035411 [ Links ]

19. M.J Katsnelson. Graphene: Carbon in two dimensions. Cambridge University Press, 1st edition, 2012. [ Links ]

20. E. V. Castro et al., Limits on charge carrier mobility in suspended graphene due to exural phonons. Phys. Rev. Lett. , 105 (2010) 266601. https://doi.org/10.1103/PhysRevLett.105.266601 [ Links ]

21. Alonso Contreras-Astorga, V Jakubský , and Alfredo Raya. On the propagation of dirac fermions in graphene with strain-induced inhomogeneous fermi velocity. Journal of Physics: Condensed Matter, 32 (2020) 295301 . https://doi.org/ 10.1088/1361-648X/ab7e5bLinks ]

22. Richard Kerner, Gerardo G. Naumis, and Wilfrido A. Gómez-Arias. Bending and exural phonon scattering: Generalized dirac equation for an electron moving in curved graphene. Physica B: Condensed Matter, 407 (2012) 2002. https://doi.org/10.1016/j.physb.2012.01.129 [ Links ]

23. Anffany Chen, R. Ilan, F. de Juan, D. I. Pikulin, and M. Franz. Quantum holography in a graphene ake with an irregular boundary. Phys. Rev. Lett. , 121 (2018) 036403. https://doi.org/10.1103/PhysRevLett.121.036403 [ Links ]

24. M. I. Katsnelson and K. S. Novoselov . Graphene: New bridge between condensed matter physics and quantum electrodynamics. Solid State Comm., 143 (2007) 3. https://doi.org/10.1016/j.ssc.2007.02.043 [ Links ]

25. Luis E.F. Foa-Torres, Stephan Roche, and Jean-Christophe Charlier. Introduction to Graphene-Based Nanomaterials. Cambridge University Press, 1st edition, 2014. [ Links ]

26. A. H. Castro Neto, F. Guinea, N. M. R. Peres, K. S. Novoselov , andA. K. Geim . The electronic properties of graphene. Rev. Mod. Phys., 81 (2009) 109. https://doi.org/10.1103/RevModPhys.81.109 [ Links ]

27. N.W. Ashcroft and N.D. Mermin. Solid State Physics. Saunders College Publishing, Fort Worth, 1976. [ Links ]

28. Jannik C. Meyer, A. K. Geim , M. I. Katsnelson , K. S. Novoselov , T. J. Booth , and S. Roth. The structure of suspended graphene sheets. Nature, 446 (2007) 60. <DOI: 10.1038/nature05545> [ Links ]

29. Gerardo G Naumis, Salvador Barraza-López, Maurice Oliva-Leyva, and Humberto Terrones. Electronic and optical properties of strained graphene and other strained 2d materials: a review. Reports on Progress in Physics, 80 (2017) 096501. DOI: 10.1088/1361-6633/aa74ef [ Links ]

30. Saúl A. Herrera and Gerardo G. Naumis. Dynamic polarization and plasmons in kekulépatterned graphene: Signatures of broken valley degeneracy. Phys. Rev. B 102 (2020) 205429. DOI: 10.1103/PhysRevB.102.205429 [ Links ]

31. David A. Ruiz-Tijerina, Elias Andrade , Ramon Carrillo-Bastos , Francisco Mireles, and Gerardo G. Naumis. Multi avor dirac fermions in kekulé-distorted graphene bilayers. Phys. Rev. B , 100 (2019) 075431. 10.1103/PhysRevB.100.075431 [ Links ]

32. José L. Aragón, Gerardo G. Naumis, and Alfredo Gómez-Rodríguez. Twisted graphene bilayers and quasicrystals: A cut and projection approach. Crystals, 9 (2019). DOI:10.3390/cryst9100519 [ Links ]

33. S. Das Sarma, Shaffique Adam, E. H. Hwang, and Enrico Rossi. Electronic transport in twodimensional graphene. Rev. Mod. Phys. , 83 (2011) 407. https://doi.org/10.1103/RevModPhys.83.407 [ Links ]

34. V. Meunier, A. G. Souza Filho, E. B. Barros, and M. S. Dresselhaus. Physical properties of low-dimensional sp2-based carbon nanostructures. Rev. Mod. Phys. , 88 (2016) 025005. https://doi.org/10.1103/RevModPhys.88.025005 [ Links ]

35. S. Reich, J. Maultzsch, C. Thomsen, and P. Ordejón. Tight-binding description of graphene. Phys. Rev. B , 66 (2002) 035412. https://doi.org/10.1103/PhysRevB.66.035412 [ Links ]

36. J E Barrios-Vargas andGerardo G Naumis . Doped graphene: the interplay between localization and frustration due to the underlying triangular symmetry. J. Phys.: Condens. Matter , 23 (2011) 375501. https://doi.org/10.1088/0953-8984/23/37/375501 [ Links ]

37. J. E. Barrios-Vargas andG. G. Naumis . Pseudo-gap opening and Dirac point confined states in doped graphene. Solid State Comm ., 162 (2013) 23. <https://doi.org/10.1016/j.ssc.2013.03.006 > [ Links ]

38. J E Barrios-Vargas andGerardo G Naumis . Critical wavefunctions in disordered graphene. J. Phys.: Condens. Matter , 24 (2012) 255305. <DOI: 10.1088/0953-8984/24/25/255305> [ Links ]

39. Andrés R. Botello-Méndez, Juan Carlos Obeso-Jureidini, and Gerardo G. Naumis. Toward an accurate tight-binding model of graphene’s electronic properties under strain. The Journal of Physical Chemistry C, 122 (2018) 15753. https://doi.org/10.1021/acs.jpcc.8b04502 [ Links ]

40. Gerardo G. Naumis, Chumin Wang, and Rafael A. Barrio. Frustration effects on the electronic density of states of a random binary alloy. Phys. Rev. B , 65 (2002) 134203. https://doi.org/10.1103/PhysRevB.65.134203 [ Links ]

41. Y Hatsugai. Topological aspect of graphene physics. J. Phys.: Conf. Ser., 334 (2011) 012004. https://doi.org/10.1088/1742-6596/334/1/012004 [ Links ]

42. G. G. Naumis , M. Terrones, H. Terrones, and L. M. Gaggero-Sager, Design of graphene electronic devices using nanoribbons of different widths. Appl. Phys. Lett., 95 (2009). https://doi.org/10.1063/1.3257731 [ Links ]

43. Gerardo G. Naumis. Internal mobility edge in doped graphene: Frustration in a renormalized lattice. Phys. Rev. B , 76 (2007) 153403. https://doi.org/10.1103/PhysRevB.76.153403 [ Links ]

44. S. Kirkpatrick and Thomas P. Eggarter. Localized states of a binary alloy. Phys. Rev. B , 6 (1972) 3598. https://doi.org/10.1103/PhysRevB.6.3598 [ Links ]

45. D. S. L. Abergel, V. Apalkov, J. Berashevich, K. Ziegler, and Tapash Chakraborty. Properties of graphene: a theoretical perspective. Adv. Phys., 59 (2010) 261. DOI: 10.1080/00018732.2010.487978 [ Links ]

46. A. Bostwick, J. L. McChesney, K. V. Emtsev, T. Seyller, K. Horn, S. D. Kevan, and E. Rotenberg, Quasiparticle Transformation during a Metal- Insulator Transition in Graphene. Phys. Rev. Lett. , 103 (2009) 056404. https://doi.org/10.1103/PhysRevLett.103.056404 [ Links ]

47. J. E. Barrios-Vargas andG. G. Naumis . Electrical conductivity and resonant states of doped graphene considering next-nearest neighbor interaction. Philos. Mag., 91 (2011) 3844. <https://doi.org/10.1080/14786435.2011.594457 > [ Links ]

48. R. Martinazzo, S. Casolo, and G. Franco Tantardini, The effect of atomic-scale defects and dopants on graphene electronic structure. In Dr. Sergey Mikhailov, editor, Physics and Applications of Graphene - Theory, chapter 3. InTech, 2010. [ Links ]

49. Gerardo G. Naumis, Rafael A. Barrio, and Chumin Wang. Effects of frustration and localization of states in the Penrose lattice. Phys. Rev. B , 50 (1994) 9834. https://doi.org/10.1103/PhysRevB.50.9834 [ Links ]

50. J. C. Slonczewski and P. R. Weiss. Band structure of graphite. Phys. Rev., 109 (1958) 272. https://doi.org/10.1103/PhysRev.109.272 [ Links ]

51. Gordon W. Semenoff. Condensed-matter simulation of a three-dimensional anomaly. Phys. Rev. Lett. , 53 (1984) 2449. [ Links ]

52. J. D. Bjorken and S. D. Drell. Relativistic Quantum Mechanics. McGraw-Hill, New York, 1964. [ Links ]

53. Gerardo G. Naumis, Leonardo A. Navarro-Labastida, Enrique Aguilar-Méndez, and Abdiel Espinosa-Champo. Reduction of the twisted bilayer graphene chiral hamiltonian into a 2 × 2 matrix operator and physical origin of at bands at magic angles. Phys. Rev. B , 103 (2021) 245418. https://doi.org/10.1103/PhysRevB.103.245418 [ Links ]

54. V. M. Silkin, J. Zhao, F. Guinea , E. V. Chulkov, P. M. Echenique, and H. Petek. Image potential states in graphene. Phys. Rev. B , 80 (2009) 121408. https://doi.org/10.1103/PhysRevB.80.121408 [ Links ]

55. E. Kogan and V. U. Nazarov. Symmetry classification of energy bands in graphene. Phys. Rev. B , 85 (2012) 115418. <https://doi.org/10.1103/PhysRevB.85.115418 > [ Links ]

56. E. Kogan , V. U. Nazarov , V. M. Silkin , and M. Kaveh. Energy bands in graphene: Comparison between the tight-binding model and ab initio calculations. Phys. Rev. B , 89 (2014) 165430. <https://doi.org/10.1103/PhysRevB.89.165430 > [ Links ]

57. E. Kogan andV. M. Silkin . Electronic structure of graphene: (Nearly) free electron bands versus tight-binding bands. Phys. Status Solidi B, 254 (2017) 1700035. <https://doi.org/10.1002/pssb.201700035 > [ Links ]

58. N. Armbrust, J. Güdde, and U. Höfer. Formation of image-potential states at the graphene/metal interface. New Journal of Physics, 17 (2015) 103043. <https://doi.org/10.1002/pssb.201700035 > [ Links ]

59. D. Niesner and T. Fauster. Image-potential states and work function of graphene. Journal of Physics: Condensed Matter , 26 (2014) 393001. <https://doi.org/10.1088/0953-8984/26/39/393001 > [ Links ]

60. Vyacheslav M. Silkin, Eugene Kogan, and Godfrey Gumbs. Screening in graphene: Response to external static electric field and an image-potential problem. Nanomaterials, 11 (2021). <https://doi.org/10.3390/nano11061561 > [ Links ]

61. M. Hernández, A. Cabo, M. Oliva-Leyva , and Gerardo G. Naumis. How water makes graphene metallic. Physics Letters A, (2019). <https://doi.org/10.1016/j.physleta.2019.125904 > [ Links ]

62. T. E. Weller, M. Ellerby, S. S. Saxena, R. P. Smith, and N. T. Skipper. Superconductivity in the intercalated graphite compounds C6Yb and C6Ca. Nature Physics, 1 (2005) 39. <https://doi.org/10.1038/nphys0010 > [ Links ]

63. B. M. Ludbrook et al., Evidence for superconductivity in Li-decorated monolayer graphene. Proceedings of the National Academy of Sciences , 112 (2015) 11795. <https://doi.org/10.1073/pnas.1510435112 > [ Links ]

64. J. Chapman et al., Superconductivity in Ca-doped graphene laminates. Scientific Reports, 6 (2016) 23254. <https://doi.org/10.1038/srep23254 > [ Links ]

65. T. Scheike, W. Bohlmann, P. Esquinazi, J. Barzola-Quiquia, A. Ballestar, and A. Setzer. Can doping graphite trigger room temperature superconductivity? Evidence for granular high-temperature superconductivity in water-treated graphite powder. Advanced Materials, 24 (2012) 5826. https://doi.org/10.1002/adma.201202219 [ Links ]

66. P. D. Esquinazi, C. E. Precker, M. Stiller, T. R. S. Cordeiro, J. Barzola-Quiquia , A. Setzer , andW. Bohlmann . Evidence for room temperature superconductivity at graphite interfaces. Quantum Studies: Mathematics and Foundations, 5 (2018) 41. <https://doi.org/10.1007/s40509-017-0131-0 > [ Links ]

67. P. M. Echenique and J. B. Pendry. The existence and detection of Rydberg states at surfaces. Journal of Physics C: Solid State Physics, 11 (1978) 2065. <https://doi.org/10.1088/0022-3719/11/10/017 > [ Links ]

68. Kazutoshi Takahashi, Masaki Imamura, Isamu Yamamoto, Junpei Azuma, and Masao Kamada. Image potential states in monolayer, bilayer, and trilayer epitaxial graphene studied with time- and angle-resolved two-photon photoemission spectroscopy. Phys. Rev. B , 89 (2014) 155303. https://doi.org/10.1103/PhysRevB.89.155303 [ Links ]

69. Ferdinand Evers and Alexander D. Mirlin. Anderson transitions. Rev. Mod. Phys. , 80 (2008) 1355. <https://doi.org/10.1103/RevModPhys.80.1355 > [ Links ]

70. Denis Bernard and André LeClair. A classification of 2D random Dirac fermions. J. Phys. A 35 (2002) 2555. <10.1088/0305-4470/35/11/303 > [ Links ]

71. G. Usaj, P. M. Perez-Piskunow, L. E. F. Foa Torres, and C. A. Balseiro. Irradiated graphene as a tunable Floquet topological insulator. Phys. Rev. B , 90 (2014) 115423. <https://doi.org/10.1103/PhysRevB.90.115423 > [ Links ]

72. Horacio Falomir, Marcelo Loewe, Enrique Muñoz, andAlfredo Raya . Optical conductivity and transparency in an effective model for graphene. Phys. Rev. B , 98 (2018) 195430. https://doi.org/10.1103/PhysRevB.98.195430 [ Links ]

73. F. J. López-Rodríguez andG. G. Naumis . Analytic solution for electrons and holes in graphene under electromagnetic waves: Gap appearance and nonlinear effects. Phys. Rev. B , 78 (2008) 201406. <https://doi.org/10.1080/14786431003757794 > [ Links ]

74. F.J. López-Rodríguez and G.G. Naumis. Graphene under perpendicular incidence of electromagnetic waves: Gaps and band structure. Philosophical Magazine, 90 (2010) 2977. <https://doi.org/10.1080/14786431003757794 > [ Links ]

75. D. Valenzuela, S. Hernández-Ortiz, Marcelo Loewe , andAlfredo Raya . Graphene transparency in weak magnetic fields. Journal of Physics A: Mathematical and Theoretical, 48 (2015) 065402. <https://doi.org/10.1088/1751-8113/48/6/065402 > [ Links ]

76. S. Hernández-Ortiz , D. Valenzuela , A. Raya, and S. Sánchez-Madrigal. Light absorption in distorted graphene. International Journal of Modern Physics B, 30 (2016) 1650084. https://doi.org/10.1142/S0217979216500843 [ Links ]

77. M. Oliva-Leyva and Gerardo G. Naumis, Effective dirac hamiltonian for anisotropic honeycomb lattices: Optical properties. Phys. Rev. B , 93 (2016) 035439. <https://doi.org/10.1103/PhysRevB.93.035439 > [ Links ]

78. F. J. Peña andEnrique Muñoz . Magnetostrain-driven quantum engine on a graphene flake. Phys. Rev. E, 91 (2015) 052152. <https://doi.org/10.1103/PhysRevE.91.052152 > [ Links ]

79. E. Muñoz and R. Soto-Garrido, Analytic approach to magneto-strain tuning of electronic transport through a graphene nanobubble: perspectives for a strain sensor. Journal of Physics: Condensed Matter , 29 (2017) 445302. DOI: <10.1088/1361-648X/aa89bc > [ Links ]

80. M. Oliva-Leyva and G. G Naumis. Sound waves induce volkov-like states, band structure and collimation effect in graphene. J. Phys.: Condens. Matter , 28 (2016) 025301. <https://doi.org/10.1088/0953-8984/28/2/025301 > [ Links ]

81. J. Carlos Sandoval-Santana, V. Guadalupe Ibarra-Sierra, J. L. Cardoso, A. Kunold, P. Roman-Taboada, and G. Naumis, Method for finding the exact effective hamiltonian of time-driven quantum systems. Annalen der Physik, 531 (2019) 1900035. <https://doi.org/10.1002/andp.201900035 > [ Links ]

82. Yonatan Betancur-Ocampo, Parisa Majari, Diego Espitia, François Leyvraz, and Thomas Stegmann. Anomalous oquet tunneling in uniaxially strained graphene. Phys. Rev. B , 103 (2021) 155433. <https://doi.org/10.1103/PhysRevB.103.155433 > [ Links ]

83. A. D. Zabolotskiy and Yu. E Lozovik. Strain-induced pseudomagnetic field in the dirac semimetal borophene. Phys. Rev. B , 94 (2016) 165403. <https://doi.org/10.1103/PhysRevB.94.165403 > [ Links ]

84. A. E. Champo andG.G. Naumis , Metal-insulator transition in 8-pmmn borophene under normal incidence of electromagnetic radiation. Phys. Rev. B 99 (2019) 035415. <https://doi.org/10.1103/PhysRevB.99.035415 > [ Links ]

85. A. Kunold , J. C. Sandoval-Santana, V. G. Ibarra-Sierra, andG. G. Naumis . Floquet spectrum and electronic transitions of tilted anisotropic dirac materials under electromagnetic radiation: Monodromy matrix approach. Phys. Rev. B , 102 (2020) 045134. <https://doi.org/10.1103/PhysRevB.102.045134 > [ Links ]

86. J. C. Sandoval-Santana , V. G. Ibarra-Sierra , A. Kunold , andG. G. Naumis . Floquet spectrum for anisotropic and tilted dirac materials under linearly polarized light at all field intensities. Journal of Applied Physics, 127 (2020) 234301. <https://doi.org/10.1063/5.0007576 > [ Links ]

87. T. Horiguchi and C. C. Chen. Lattice green’s function for the diced lattice. Journal of Mathematical Physics, 15 (1974) 659. <https://doi.org/10.1103/PhysRevE.87.032109 > [ Links ]

88. B. Sutherland, Localization of electronic wave functions due to local topology. Phys. Rev. B , 34 (1986) 5208. <https://doi.org/10.1103/PhysRevB.34.5208 > [ Links ]

89. E. V. Gorbar, V. P. Gusynin, and D. O. Oriekhov. Electron states for gapped pseudospin-1 fermions in the field of a charged impurity. Phys. Rev. B , 99 (2019) 155124. <https://doi.org/10.1103/PhysRevB.99.155124 > [ Links ]

90. A. Raoux, M Morigi, J-N Fuchs, F Piéchon, and G Montambaux. From dia-to paramagnetic orbital susceptibility of massless fermions. Physical review letters, 112 (2014) 026402. <https://doi.org/10.1103/PhysRevLett.112.026402 > [ Links ]

91. M. A. Mojarro, V. G. Ibarra-Sierra , J. C. Sandoval-Santana , R. Carrillo-Bastos, andG. G. Naumis . Electron transitions for dirac hamiltonians with at bands under electromagnetic radiation: Application to the α − U3 graphene model. Phys. Rev. B 101 (2020) 165305. <https://doi.org/10.1103/PhysRevB.101.165305 > [ Links ]

92. C. Gutierrez et al., Imaging chiral symmetry breaking from Kekule bond order in graphene. Nat. Phys., 12 (2016) 950. <https://doi.org/10.1038/nphys3776 > [ Links ]

93. O. V. Gamayun, V. P. Ostroukh, N. V. Gnezdilov, I. Adagideli, and C. W. J. Beenakker, Valleymomentum locking in a graphene superlattice with y-shaped kekulé bond texture. New Journal of Physics , 20 (2018) 023016 . <https://doi.org/10.1088/1367-2630/aaa7e5 > [ Links ]

94. S. A. Herrera andG. G. Naumis , Electronic and optical conductivity of kekulépatterned graphene: Intravalley and intervalley transport. Phys. Rev. B , 101 (2020) 205413. <https://doi.org/10.1103/PhysRevB.101.205413 > [ Links ]

95. Qing-Ping Wu, Lu-Lu Chang, Yu-Zeng Li, Zheng-Fang Liu, and Xian-Bo Xiao. Electriccontrolled valley pseudomagnetoresistance in graphene with y-shaped kekul a© lattice distortion. Nanoscale Research Letters, 15 (2020) 46. <https://doi.org/10.1186/s11671-020-3275-5 > [ Links ]

96. S. A. Herrera andG. G. Naumis , Electronic and optical conductivity of kekulépatterned graphene: Intravalley and intervalley transport. Phys. Rev. B , 101 (2020) 205413. <https://doi.org/10.1103/PhysRevB.101.205413 > [ Links ]

97. E. Andrade, G. G. Naumis , andR. Carrillo-Bastos . Electronic spectrum of kekule patterned graphene considering second neighbor-interactions. Journal of Physics Condensed Matter, 2021. https://doi.org/10.1088/1361-648X/abef9a [ Links ]

98. T. Lee et al., Layer-by-Layer Assembly for Graphene-Based Multilayer Nanocomposites: Synthesis and Applications. Chem. Mater. 27 (2015) 3785. <https://doi.org/10.1021/acs.chemmater.5b00491 > [ Links ]

99. C.-J. Kim et al., Chiral atomically thin films. Nat. Nanotech., 11 (2016) 520 . <10.1038/nnano.2016.3 > [ Links ]

100. S. Latil, V. Meunier , and L. Henrard. Massless fermions in multilayer graphitic systems with misoriented layers: Ab initio calculations and experimental fingerprints. Phys. Rev. B , 76 (2007) 201402 . https://doi.org/10.1103/PhysRevB.76.201402 [ Links ]

101. M. Freitag, Graphene: Trilayers unravelled. Nat. Phys., 7 (2011) 596 . DOI: 10.1038/nphys2032 [ Links ]

102. C. Gutiérrez et al., Imaging chiral symmetry breaking from kekulé bond order in graphene. Nature Physics , 12 (2016) 950 . <https://doi.org/10.1038/nphys3776 > [ Links ]

103. E. McCann and Mikito Koshino, The electronic properties of bilayer graphene. Rep. Prog. Phys. 76 (2013) 056503. <https://doi.org/10.1088/0034-4885/76/5/056503 > [ Links ]

104. Ke-Ke Bai et al., Creating One-Dimensional Nanoscale Periodic Ripples in a Continuous Mosaic Graphene Monolayer. Phys. Rev. Lett. , 113 (2014) 086102. <https://doi.org/10.1103/PhysRevLett.113.086102 > [ Links ]

105. A. Artaud et al., Universal classification of twisted, strained and sheared graphene moiré superlattices. Sci. Rep., 6 (2016) 25670 . <https://doi.org/10.1038/srep25670 > [ Links ]

106. Pedro Roman-Taboada and Gerardo G. Naumis. Spectral butter y, mixed Dirac-Schrödinger fermion behavior, and topological states in armchair uniaxial strained graphene. Phys. Rev. B , 90 (2014) 195435. <DOI: 10.1103/PhysRevB.90.195435 > [ Links ]

107. H. García-Cervantes, L. M. Gaggero-Sager , O. Sotolongo-Costa, G. G. Naumis , and I. Rodríguez-Vargas. Angle-dependent bandgap engineering in gated graphene superlattices. AIP Advances, 6 (2016) 035309. <Doi: 10.1063/1.4944495 > [ Links ]

108. D.S. Diaz-Guerrero, I. Rodríguez-Vargas , G. G. Naumis , andL. M. Gaggero-Sager . Self-similar charge transmission in gapped graphene. Fractals, 24 (2016) 1630002. Doi: 10.1142/S0218348X16300026 [ Links ]

109. S. Molina-Valdovinos, J. Martínez-Rivera, N.E. Moreno-Cabrera, and I. Rodríguez-Vargas . Low-dimensional thermoelectricity in graphene: The case of gated graphene superlattices. Physica E: Low-dimensional Systems and Nanostructures, 101 (2018) 188. https://doi.org/10.1016/j.physe.2018.03.005 [ Links ]

110. E. J. Guzmán, O. Navarro, O. Oubram, and I. Rodríguez-Vargas . Transport properties and thermoelectric effects in gated silicene superlattices. Journal of Applied Physics , 124 (2018) 144305. <https://doi.org/10.1063/1.5045479 > [ Links ]

111. C. R. Dean et al., Hofstadter’s butter y and the fractal quantum Hall effect in moiré superlattices. Nature, 497 (2013) 598. https://doi.org/10.1038/nature12186 [ Links ]

112. Douglas R. Hofstadter. Energy levels and wave functions of Bloch electrons in rational and irrational magnetic fields. Phys. Rev. B , 14 (1976) 2239. <https://doi.org/10.1103/PhysRevB.14.2239 > [ Links ]

113. Eduardo Fradkin. Field Theories of Condensed Matter Physics. Cambridge University Press, 2 edition, 2013. [ Links ]

114. Gerardo G. Naumis. Topological map of the Hofstadter butter y: Fine structure of Chern numbers and Van Hove singularities. Phys. Lett. A , 380 (2016 1772). arXiv:1507.08130 [ Links ]

115. Gerardo Naumis. Higher-dimensional quasicrystalline approach to the hofstadter butter y topological-phase band conductances: Symbolic sequences and self-similar rules at all magnetic uxes. Phys. Rev. B 100 (2019) 165101. <https://doi.org/10.1103/PhysRevB.100.165101 > [ Links ]

116. Svetlana Ya. Jitomirskaya and Yoram Last. Dimensional hausdorff properties of singular continuous spectra. Phys. Rev. Lett. , 76 (1996) 1765. <https://doi.org/10.1103/PhysRevLett.76.1765 > [ Links ]

117. W. P. Su, J. R. Schrieffer, and A. J. Heeger. Solitons in polyacetylene. Phys. Rev. Lett. , 42 (1979) 1698. https://doi.org/10.1103/PhysRevLett.42.1698 [ Links ]

118. Sriram Ganeshan, Kai Sun, andS. Das Sarma . Topological zero-energy modes in gapless commensurate aubry-andré-harper models. Phys. Rev. Lett. , 110 (2013) 180403. https://doi.org/10.1103/PhysRevLett.110.180403 [ Links ]

119. T. T. Heikkilä, N. B. Kopnin, and G. E. Volovik. Flat bands in topological media. JETP Letters, 94 (2011) 233. <https://doi.org/10.1134/S0021364011150045 > [ Links ]

120. G. E. Volovik. Flat band in topological matter. J. Supercond. Nov. Magn., 26 (2013) 2887. https://doi.org/10.1007/s10948-013-2221-5 [ Links ]

121. Pedro Roman-Taboada and Gerardo G. Naumis. Topological at bands in time-periodically driven uniaxial strained graphene nanoribbons. Phys. Rev. B , 95 (2017) 115440. <https://doi.org/10.1103/PhysRevB.95.115440 > [ Links ]

122. Pedro Roman-Taboada and Gerardo G. Naumis. Topological edge states on time-periodically strained armchair graphene nanoribbons. Phys. Rev. B , 96 (2017) 155435 . <https://doi.org/10.1103/PhysRevB.96.155435 > [ Links ]

123. Pedro Roman-Taboada andGerardo G Naumis . Topological phase-diagram of time-periodically rippled zigzag graphene nanoribbons. Journal of Physics Communications, 1 (2017) 055023. <https://doi.org/10.1088/2399-6528/aa98fd > [ Links ]

124. Eslam Khalaf, Alex J. Kruchkov, Grigory Tarnopolsky , andAshvin Vishwanath . Magic angle hierarchy in twisted graphene multilayers. Phys. Rev. B , 100 (2019) 085109. <https://doi.org/10.1103/PhysRevB.100.085109 > [ Links ]

125. Jeong Min Park , Yuan Cao , Kenji Watanabe , Takashi Taniguchi , andPablo Jarillo-Herrero . Tunable strongly coupled superconductivity in magic-angle twisted trilayer graphene. Nature, 590 (2021) 249. <https://doi.org/10.1038/s41586-021-03192-0 > [ Links ]

126. Gerardo G. Naumis. Low-frequency vibrational modes anomalies and rigidity: A key to understanding the glass and the electronic properties of exible materials from a Topological Perspective. Front. Mater., 2 (2015) 44. <https://doi.org/10.3389/fmats.2015.00044 > [ Links ]

127. Adrián Huerta and Gerardo G. Naumis. Role of rigidity in the uid-solid transition. Phys. Rev. Lett. , 90 (2003) 145701. <https://doi.org/10.1103/PhysRevLett.90.145701 > [ Links ]

128. Gerardo G. Naumis. Energy landscape and rigidity. Phys. Rev. E, 71 (2005) 026114. <https://doi.org/10.1103/PhysRevE.71.026114 > [ Links ]

129. Adrián Huerta and Gerardo G. Naumis. Evidence of a glass transition induced by rigidity self-organization in a network-forming uid. Phys. Rev. B , 66 (2002) 184204. https://doi.org/10.1103/PhysRevB.66.184204 [ Links ]

130. Hugo M. Flores-Ruiz, Gerardo G. Naumis, and J. C. Phillips. Heating through the glass transition: A rigidity approach to the boson peak. Phys. Rev. B 82 (2010) 214201. <https://doi.org/10.1103/PhysRevB.82.214201 > [ Links ]

131. Cristian F. Moukarzel and Gerardo G. Naumis. Comment on penrose tilings as jammed solids. Phys. Rev. Lett. , 115 (2015) 209801. <https://doi.org/10.1103/PhysRevLett.115.209801 > [ Links ]

132. A. Bazán et al., Quasicrystalline and rational approximant wave patterns in hydrodynamic and quantum nested wells. Phys. Rev. Lett. 97 (2006) 124501. https://doi.org/10.1103/PhysRevLett.97.124501 [ Links ]

Received: July 28, 2021; Accepted: July 28, 2021

Creative Commons License This is an open-access article distributed under the terms of the Creative Commons Attribution License