Topological phase-diagram of time-periodically rippled zigzag graphene nanoribbons

Graphene, a truly two-dimensional (1D) material, has proven to have very interesting and fascinating properties [1, 2]. Among them, one can mention its extraordinary mechanical features, which can be used to tailor the electronic properties, leading to many novel effects in the static case [3–31]. As a matter of fact, within the tight binding approach and in the absence of interactions between electrons, the effects of a deformation field applied to graphene can be described via a pseudo magnetic field [22, 32–37]. On the other hand, graphene possesses interesting topological properties for both the time-independent [38–52] and the time-dependent cases [53–61]. For instance, in the static case, it has been proven that Dirac cones have a non-vanishing Berry phase, which means that they are robust against perturbations and disorder [62]. In addition, since Dirac cones always come in pairs, each cone has an opposite Berry phase as is companion. Hence, as a consequence of the bulk-edge correspondence, an edge state (such edge state is a flat band for the case of pristine zigzag graphene nanoribbons (ZGNs)) emerges joining two inequivalent Dirac cones (this is, two Dirac cones with opposite Berry phase).

On the other hand, by applying a time-dependent deformation field to graphene, new and novel phenomena appear when compared to the static case [61]. For instance, when a time-dependent in-plane AC electric field is applied to graphene, it is possible to undergo a topological phase transition from a topological semi-metallic phase to a trivial insulator one [63]. Similarly, gaps on the energy spectrum of graphene can be opened by irradiating graphene with a laser by changing its intensity [64, 65]. This gapped phase is also able to host robust topological chiral Floquet edge states, which are highly tunable [56]. These features are similar to the ones observed in topological insulators, which also exhibit robust edge states. However, there is another kind of topological phases akin to gapless systems [66, 67]. Take the kicked Harper model [68] and the kicked SSH model [69], for instance. In the kicked Harper model via periodic driving, one can create many touching band points (i.e. points at where the band edges touch each other following a linear dispersion) that can give rise to highly localized edge states. This occurs because touching band points always come in pairs and each of them have opposite chirality as its companion [68]. These edge states can be flat bands or dispersive edge states. Interestingly enough, one can have the same effect on graphene nanoribbons by applying a time-dependent strain field [61]. The aim of this paper is to show some of these topological properties of gapless systems by studying a periodically driven uniaxial rippled ZGN. To do that, we use a tight-binding Hamiltonian to describe the electronic properties of the periodically driven rippled ZGN within the Floquet formalism. The quasienergy spectrum is then obtained by using an effective Hamiltonian approach.

It is important to remark that the considered deformation field is a corrugation of the graphene ribbon. Here we will restrict ourselves to the case of uniaxial ripples, i.e., the height of carbon atoms with respect to a plane is affected only along one direction (in what follows, we will consider a deformation field applied along the armchair direction). Therefore, it is necessary to take into account the relative change of the orientation between π orbitals [9]. Within such approximation, as will be seen later on, the time-dependent deformation field allows us to create touching band points (touching band points are points at where a band inversion is observed) at zero or $\pm \pi$ quasienergies. Around such points the quasienergy spectrum exhibits a linear dispersion as in the case of Dirac cones. The touching band points originated from the time-dependent deformation field can be of two different kinds: types I and II, each of them giving rise to topologically protected edge states. For the former type, we have found topologically protected flat bands at zero and $\pm \pi$ quasienergy. Such flat bands join two inequivalent touching band points with opposite Berry phase. For the latter, dispersive edge states were found and it was found that they are, at least, topologically weak by obtaining the winding of the effective Hamiltonian.

To finish, it is worthwhile to say that the experimental realization of the deformation pattern here considered can be difficult since it requires very specific hopping parameters values and very fast time scales. In fact, a similar experiment was proposed by us in a previous work [61]. However, this experiment was tailored for in-plane strain [61], and since graphene is almost incompressible, the compressive strain will induce ripples on the nanoribbon. As a result, it is clear that ripple effects are important to be studied. Also, it is possible to have a 1D periodic ripple on graphene. This is done by using thermal enginerring and by growing graphene upon a substrate to induce an anisotropic strain pattern [70]. The time-dependent deformation field can be obtained by applying a time-periodic pressure variation to the whole system [61, 71] (the graphene nanoribbon and the substrate). To observe the results presented below, the pressure needs to be in the frequency range of femto seconds, which can be very challenging in a real experiment. As an alternative, we propose the use of artificial or optical lattices, where the hopping parameters of the graphene nanoribbon lattice can be tailored at will [72–77].

The paper is organized as follows, in section 2 we introduce the model. This is, we briefly discuss how to describe the electronic properties of a rippled ZGN. Then, the time dependence is introduced to the model and the time-evolution operator of the system is defined. In section 3, we analytically obtain the quasienergy spectrum of the system via an effective Hamiltonian approach. Also, the location of both types of touching band points is found and the topological phase diagram of the system is built. The edge states of the system and their topological properties are analyzed in section 4. Some conclusion are given in section 5. Finally, in the appendices some calculations regarding the main text are presented.

We start by considering a ZGN as the one portrayed in figure 1(a), then we apply an out-of-plane uniaxial deformation field (a ripple field) along the y-direction given by,

$$\begin{eqnarray}&&{z}_{j}=z({y}_{j})=\lambda \cos (2\pi \sigma {y}_{j}+\phi ),\end{eqnarray} \tag{ 1 }$$

here, y_j are the positions of the carbon atoms along the y-direction, λ is the amplitude, σ controls the spatial wavelength, and ϕ is a phase. Since such a deformation field modifies the height of the carbon atoms, their positions are also modified and can be written as ${{\bf{r}}}^{{\prime} }=({\bf{r}},z({y}_{j}))$, where ${\bf{r}}$ are the carbon atom positions in unrippled graphene. Within the low energy limit, the electronic properties of a ZGN under a deformation field along the armchair direction, as the one given by equation (1), are well described by the following 1D tight-binding effective Hamiltonian [9],

$$\begin{eqnarray}&&H({k}_{x})=\displaystyle \sum _{j=1}^{N-1}[{\gamma }_{2j}\,{a}_{2j+1}^{\dagger }{b}_{2j}+c({k}_{x}){\gamma }_{2j-1}{a}_{2j-1}^{\dagger }{b}_{2j}]+{\rm{h}}.{\rm{c}}.,\end{eqnarray} \tag{ 2 }$$

where $c({k}_{x})=2\cos (\sqrt{3}{k}_{x}/2)$, the operator a_j (b_j) annihilates an electron at the jth site in the sub lattice A (B), and N is the number of atoms per unit cell (see figure 1, at where the unit cell is indicated by solid red lines). ${\gamma }_{j}$ are the hopping parameters given by [9],

$$\begin{eqnarray}&&{\gamma }_{j}={\gamma }_{0}[1+\alpha (1-{\hat{{\bf{N}}}}_{j}\cdot {\hat{{\bf{N}}}}_{j+1})]\exp (-\beta \delta {l}_{j+1,j}),\end{eqnarray} \tag{ 3 }$$

where ${\gamma }_{0}=2.7\,\mathrm{eV}$ is the hopping parameter for pristine graphene, ${\hat{{\bf{N}}}}_{j}$ is the unit vector normal to the pristine graphene sheet at site j, which has the following form,

$$\begin{eqnarray}&&{\hat{{\bf{N}}}}_{j}=\displaystyle \frac{{\hat{e}}_{z}-{\rm{\nabla }}{z}_{j}}{\sqrt{1+{({\rm{\nabla }}{z}_{j})}^{2}}},\end{eqnarray} \tag{ 4 }$$

with ${\rm{\nabla }}=({\partial }_{x},{\partial }_{y})$ being the two-dimensional (2D) gradient operator. ${\hat{e}}_{z}$ is a unit vector that is perpendicular to the unrippled graphene sheet, $\alpha \approx 0.4$ is a constant that takes into account the change of the relative orientation between π-orbitals originated from the deformation field, and $\beta \approx 3.37$ is the decay rate (Grüneisen parameter). Finally, the quantity $\delta {l}_{j+1,j}$ is given by,

$$\begin{eqnarray}&&\delta {l}_{j+1,j}=-1+\sqrt{1+{[z({y}_{j+1})-z({y}_{j})]}^{2}}.\end{eqnarray} \tag{ 5 }$$

It is important to say that all distances, here on, will be measured in units of the interatomic distance between carbon atoms (a_c) in pristine graphene. In a similar way, we will set ${\gamma }_{0}$ as the unit of energy. Having said that, it is noteworthy that the energy spectrum of the Hamiltonian equation (2) have been discussed in a previous work for the small amplitude limit and for different ripple's wavelength, see [9]. Also, it is important to say that the deformation field here considered induces a pseudo magnetic field, since such deformation field modifies the relative orientation between π orbitals. In fact, if we assume that $\hat{{\bf{N}}}$ is a smooth function of the position, the magnetic flux through a ripple of lateral dimension l and height z is given by [35],

$$\begin{eqnarray}&&{\rm{\Phi }}\approx \displaystyle \frac{10{\mathring{\rm A} }^{-1}{a}_{{\rm{c}}}^{2}{z}^{2}}{{l}^{3}}.\end{eqnarray} \tag{ 6 }$$

If we introduce all the numerical values, we obtain ${\rm{\Phi }}\approx {10}^{-3}{{\rm{\Phi }}}_{0}$, where ${{\rm{\Phi }}}_{0}=2\pi {\hslash }/c$ and c is the speed of light.

Zoom In Zoom Out Reset image size

Once that the Hamiltonian that describes an uniaxial rippled ZGN has been presented, we proceed to introduce the time-dependence to our model. We will consider a pulse time-driving layout,

$$\begin{eqnarray}H({k}_{x},t)=\left\{\begin{array}{lll}{H}_{0}({k}_{x}) & {\rm{if}} & t\lt {\rm{mod}}(t,T)\lt {t}_{1}\\ {H}_{1}({k}_{x}) & {\rm{if}} & {t}_{1}\lt {\rm{mod}}(t,T)\lt T,\end{array}\right.\end{eqnarray} \tag{ 7 }$$

where T is the driving period and t₁ is a number such that $0\lt {t}_{1}\lt T$. The previous Hamiltonian describes a driving layout in which for times within the interval $({t}_{1},T)$, the deformation field is turned on, whereas it is turned off for times within the interval $(0,{t}_{1})$. For the sake of simplicity, in what follows we will consider the case of short pulses, in other words, we will consider the limit ${t}_{1}\to T$, which resembles the delta driving case. Thus, in the delta driving layout, we turn on the deformation field given by equation (1) at times t = mT , while for $t\ne {mT}$ the deformation field is turned off, here m is an integer number. A graphic representation of this driving layout is shown in figure 1. Within this limit (${t}_{1}\to T$), the time-dependent Hamiltonian (7) takes the following form,

$$\begin{eqnarray}&&H({k}_{x},t)={H}_{0}({k}_{x})+\displaystyle \sum _{m}[{H}_{1}({k}_{x})-{H}_{0}({k}_{x})]\delta (t/T-m),\end{eqnarray} \tag{ 8 }$$

with the Hamiltonians ${H}_{0}({k}_{x})$ and ${H}_{1}({k}_{x})$ given by,

$$\begin{eqnarray}&&{H}_{0}({k}_{x})=\displaystyle \sum _{j=1}^{N-1}{\gamma }_{0}[{a}_{2j+1}^{\dagger }{b}_{2j}+c({k}_{x}){a}_{2j-1}^{\dagger }{b}_{2j}]+{\rm{h}}.{\rm{c}}.,\end{eqnarray} \tag{ 9 }$$

and

$$\begin{eqnarray}&&{H}_{1}({k}_{x})=\displaystyle \sum _{j=1}^{N-1}[{\gamma }_{2j}\,{a}_{2j+1}^{\dagger }{b}_{2j}+c({k}_{x}){\gamma }_{2j-1}{a}_{2j-1}^{\dagger }{b}_{2j}]+{\rm{h}}.{\rm{c}}.\end{eqnarray} \tag{ 10 }$$

Before entering into the details of our model, let us briefly discuss the effect of considering a sinusoidal time perturbation instead of a Dirac delta protocol. The Dirac delta driving is useful because calculations are greatly simplified and because analytical results can be obtained. One can consider a more realistic time perturbation but the system must be treated numerically. Consider for example a cosine-like driving, then the quasienergies of the system are given by the eigenvalues of the so-called Floquet Hamiltonian [78], which is a block diagonal matrix (for our case, each block is N × N matrix with N being the number of atoms per unit cell). By truncating such Hamiltonian (this is, by considering only the first three blocks of such Hamiltonian), one can obtain numerically the quasienergies. By using this kind of driving as we have proven in a previous work [61] for a model quite similar to the one studied here, that the secular gaps are reduced in size when compared with the delta driving. Additionally, the flat bands become dispersive edge states [61]. Summarizing, the emergence of highly localized edge states is not modified if a more realistic driving layout is considered.

To study the time evolution of our system, we define the unitary one-period time evolution operator, $U({k}_{x},T)$, in the usual form,

$$\begin{eqnarray}&&U({k}_{x},T)| {\psi }_{{k}_{x}}(t)\rangle =| {\psi }_{{k}_{x}}(t+T)\rangle ,\end{eqnarray} \tag{ 11 }$$

where $| {\psi }_{{k}_{x}}(t)\rangle$ is the system wave function for a given k_x. The main advantage of using a delta kicking is that the time evolution operator is easy to find. For this case, we have,

$$\begin{eqnarray}U(\tau ) & = & { \mathcal T }{\rm{\exp }}\left[-{\rm{i}}{\displaystyle \int }_{0}^{T}H({k}_{x},t){\rm{d}}{t}/{\hslash }\right]\\ & = & {\rm{\exp }}[-{\rm{i}}\tau ({H}_{1}({k}_{x})-{H}_{0}({k}_{x}))]{\rm{\exp }}[-{\rm{i}}\tau {H}_{0}({k}_{x})],\end{eqnarray} \tag{ 12 }$$

here ${ \mathcal T }$ denotes the time ordering operator and $\tau =T/{\hslash }$. In general Hamiltonians H₁ and H₀ do not commute, therefore, it is a common practice to study the eigenvalue spectrum of the matrix representation of equation (12) via an effective Hamiltonian defined as

$$\begin{eqnarray}&&U({k}_{x},\tau )=\exp [-{\rm{i}}\tau {H}_{\mathrm{eff}}({k}_{x})].\end{eqnarray} \tag{ 13 }$$

Then, the eigenvalues of the time-evolution operator, which we denote by $\tau \omega$, are the eigenvalues of the effective Hamiltonian, $\tau {H}_{\mathrm{eff}}({k}_{x})$. Since $\tau \omega$ are just defined up to integer multiples of $2\pi$, they are called the quasienergies of the system.

Once that the time-dependence have been introduced to our model, we have four free parameters, three owing to the deformation field (λ, σ, and ϕ) and one due to the driving layout (τ). One can study the quasienergy spectrum for a wide range of parameters, however just a few set of parameters allows us to do analytical calculations. Among them, one can mention the case $\sigma =1/3$ and $\phi =0$ for which the system becomes periodic along both the x-direction and the y-direction. This is due to the fact that the hopping parameters, for this particular case, just take two different values, namely,

$$\begin{eqnarray}&&{\gamma }_{j}=\left(1+\alpha -\displaystyle \frac{\alpha }{\sqrt{\tfrac{{\pi }^{2}{\lambda }^{2}}{3}+1}}\right)\exp [\beta (1-\sqrt{{\xi }_{j}{\lambda }^{2}+1})],\end{eqnarray} \tag{ 14 }$$

where ${\xi }_{j}=1/4$ for odd j and ${\xi }_{j}=3/2$ otherwise.

It is noteworthy that for $\sigma =1/3$, our system is quite similar to the system studied in [61], therein a periodically driven uniaxial strained ZGN is studied. The main result of such paper is the emergence of topologically protected flat bands at both zero and $\pm \pi$ quasienergies. The emergence of these flat bands can be understood in terms of a kind of Weyl points that appear each time that the bands are inverted [79]. Therefore, we expect our model to have topological flat bands and Weyl points. This conjecture is confirmed in the next section where the touching band points of the quasienergy spectrum are found.

Our system can be studied numerically for any combination of driving parameters. From an analytical point of view, only few cases are simple enough to carry on calculations. In fact, for inconnmensurate σ, the problem is very complex since quasiperiodicity arises and requires the use of rational approximants and renormalization approaches [80–83]. Here we have chosen to present simple analytical cases and compare it with the numerical results. In particular, we will study the quasienergy touching band points for $\sigma =1/3$, $\phi =0$ and fixed values of λ and τ. For this case, the system becomes periodic along both the x- and y-directions if cyclic boundary conditions are used in the y axis. Nanoribbons are thus studied by changing the boundary conditions. This allows to define the Fourier transformed version of Hamiltonians equations (9) and (10),

$$\begin{eqnarray}{H}_{0}({\bf{k}}) & = & {h}_{0}({\bf{k}}){\hat{{\bf{h}}}}_{{\bf{0}}}({\bf{k}})\cdot {\boldsymbol{\sigma }}\\ {H}_{1}({\bf{k}}) & = & {h}_{1}({\bf{k}}){\hat{{\bf{h}}}}_{{\bf{1}}}({\bf{k}})\cdot {\boldsymbol{\sigma }}\end{eqnarray} \tag{ 15 }$$

by using a vector in reciprocal space ${\bf{k}}=({k}_{x},{k}_{y})$. ${\sigma }_{i}$ ($i=x,y,z$) are the 2 × 2 Pauli matrices and $\hat{{{\bf{h}}}_{{\bf{0}}}}({\bf{k}})={{\bf{h}}}_{0}({\bf{k}})/| {h}_{0}({\bf{k}})|$, $\hat{{{\bf{h}}}_{{\bf{1}}}}({\bf{k}})={{\bf{h}}}_{1}({\bf{k}})/| {h}_{1}({\bf{k}})|$. Here, ${h}_{0}({\bf{k}})$ and ${h}_{1}({\bf{k}})$ denote the norm of ${{\bf{h}}}_{0}({\bf{k}})$ and ${{\bf{h}}}_{1}({\bf{k}})$ respectively. ${{\bf{h}}}_{0}({\bf{k}})$ and ${{\bf{h}}}_{1}({\bf{k}})$ have components which are defined in appendix A. The ${\bf{k}}$-dependent time evolution operator, equation (12), now takes the following form,

$$\begin{eqnarray}&&U({\bf{k}},\tau )=\displaystyle \sum _{{k}_{y}}{ \mathcal U }({\bf{k}},\tau )\otimes | {k}_{y}\rangle \langle {k}_{y}| ,\end{eqnarray} \tag{ 16 }$$

where,

$$\begin{eqnarray}&&{ \mathcal U }({\bf{k}},\tau )={\rm{\exp }}[-{\rm{i}}\tau \delta H({\bf{k}})]{\rm{\exp }}[-{\rm{i}}\tau {H}_{0}({\bf{k}})]\end{eqnarray} \tag{ 17 }$$

and $\delta H({\bf{k}})={H}_{1}({\bf{k}})-{H}_{0}({\bf{k}})$. To obtain the quasienergy spectrum we use an effective Hamiltonian approach. Let us define the effective Hamiltonian as,

$$\begin{eqnarray}&&{ \mathcal U }({\bf{k}},\tau )=\exp [-{\rm{i}}\tau {H}_{\mathrm{eff}}({\bf{k}})].\end{eqnarray} \tag{ 18 }$$

Since the Hamiltonians ${H}_{0}({\bf{k}})$ and ${H}_{1}({\bf{k}})$ are 2 × 2 matrices, it is possible to analytically obtain ${H}_{\mathrm{eff}}({\bf{k}})$ using the addition rule of SU(2) (see appendix A for details). After some calculations and using equations (15) and (17), one gets,

$$\begin{eqnarray}&&{H}_{\mathrm{eff}}({\bf{k}})=\omega ({\bf{k}}){\hat{{\bf{h}}}}_{\mathrm{eff}}({\bf{k}})\cdot {\boldsymbol{\sigma }},\end{eqnarray} \tag{ 19 }$$

and as before, ${\boldsymbol{\sigma }}$ is the Pauli vector. The quasienergies, $\tau \omega ({\bf{k}})$, are given by the following expression,

$$\begin{eqnarray}&&\cos [\tau \omega ({\bf{k}})]=\cos [\tau \,\delta h({\bf{k}})]\cos [\tau {h}_{0}({\bf{k}})]\\ &&\quad -{\hat{{\bf{h}}}}_{{\bf{0}}}({\bf{k}})\cdot \hat{{\boldsymbol{\delta }}{\bf{h}}}({\bf{k}})\sin [\tau \,\delta h({\bf{k}})]\sin [\tau {h}_{0}({\bf{k}})],\end{eqnarray} \tag{ 20 }$$

where ${\boldsymbol{\delta }}{\bf{h}}({\bf{k}})={{\bf{h}}}_{1}({\bf{k}})-{{\bf{h}}}_{0}({\bf{k}})$, and ${\hat{{\bf{h}}}}_{\mathrm{eff}}({\bf{k}})$ is given by,

$$\begin{eqnarray}{\hat{{\bf{h}}}}_{\mathrm{eff}}({\bf{k}}) & = & \displaystyle \frac{-1}{\sin [\tau \omega ({\bf{k}})]}[\hat{{\boldsymbol{\delta }}{\bf{h}}({\bf{k}})}\sin [\tau \,\delta h({\bf{k}})]\cos [\tau {h}_{0}({\bf{k}})]]\\ & & +\displaystyle \frac{-1}{\sin [\tau \omega ({\bf{k}})]}[{\hat{{\bf{h}}}}_{0}({\bf{k}})\sin [\tau {h}_{0}({\bf{k}})]\cos [\tau \delta h({\bf{k}})]]\\ & & +\displaystyle \frac{-1}{\sin [\tau \omega ({\bf{k}})]}[\hat{{\boldsymbol{\delta }}{\bf{h}}({\bf{k}})}\times {\hat{{\bf{h}}}}_{0}({\bf{k}})\sin [\tau \delta h({\bf{k}})]\sin (\tau {h}_{0}[{\bf{k}})]].\end{eqnarray} \tag{ 21 }$$

Since we are looking for touching band points, it is useful to plot the quasienergy spectrum for some characteristic values of λ and τ. In figure 2 we plot the quasienergy band structure for $\sigma =1/3$, $\phi =0$, $\lambda =0.5$, and $\tau =\pi$. Note that apart the Dirac cones (indicated by yellow dots in the figure), there are other touching band points at zero and $\pm \pi$ quasienergies.

Zoom In Zoom Out Reset image size

From figure 2, we can see that touching band points always emerge at zero or $\pm \pi$ quasienergy, then it follows that they can be obtained by imposing $\tau \omega ({{\bf{k}}}^{* })=n\pi$, where n is an integer number and ${{\bf{k}}}^{* }=({k}_{x}^{* },{k}_{y}^{* })$ are the special points where this happens. By substituting ${\bf{k}}={{\bf{k}}}^{* }$ in equation (20), the touching band points are given by the solutions of the following equation,

$$\begin{eqnarray}\pm 1 & = & \cos [\tau \,\delta h({{\bf{k}}}^{* })]\cos [\tau {h}_{0}({{\bf{k}}}^{* })]\\ & & -{\hat{{\bf{h}}}}_{{\bf{0}}}({{\bf{k}}}^{* })\cdot \hat{{\boldsymbol{\delta }}{\bf{h}}}({{\bf{k}}}^{* })\sin [\tau \,\delta h({{\bf{k}}}^{* })]\sin [\tau {h}_{0}({{\bf{k}}}^{* })].\end{eqnarray} \tag{ 22 }$$

A careful analysis of equation (22) shows two possible solutions depending on the value of the dot product ${\hat{{\bf{h}}}}_{{\bf{0}}}({{\bf{k}}}^{* })\cdot \hat{{\boldsymbol{\delta }}{\bf{h}}}({{\bf{k}}}^{* })$. In other words, there are two kinds of touching band points that we have labeled by types I and II. For the type I, it is required that ${\hat{{\bf{h}}}}_{{\bf{0}}}({{\bf{k}}}^{* })\cdot \hat{{\boldsymbol{\delta }}{\bf{h}}}({{\bf{k}}}^{* })=\pm 1$, which is equivalent to ask the commutator $[{H}_{1}({{\bf{k}}}^{* }),{H}_{0}({{\bf{k}}}^{* })]$ to vanish. For type II, it is necessary to impose two simultaneous restrictions, the first one is ${\hat{{\bf{h}}}}_{{\bf{0}}}({{\bf{k}}}^{* })\cdot \hat{{\boldsymbol{\delta }}{\bf{h}}}({{\bf{k}}}^{* })\ne \pm 1$, whereas the second one is given by $\cos [\tau \delta h({{\bf{k}}}^{* })]\cos [\tau {h}_{0}({{\bf{k}}}^{* })]=\pm 1$, this means that type II touching band points never occur at ${k}_{y}^{* }=0,\pm 2\pi /3$. It what follows, we will study the necessary conditions for having these kinds of touching band points. After that, the topological phase diagram of the system is obtained.

3.1. Type I

Although this kind of touching band points have been studied in a previous work for a very particular case of hopping parameters [61], here we obtain the touching band points for the general case of an effective linear chain with two different hopping parameters, say ${\gamma }_{1}$ and ${\gamma }_{2}$. We start our analysis by noticing from equation (A.9), that ${\hat{{\bf{h}}}}_{0}({{\bf{k}}}^{* })\cdot \hat{{\boldsymbol{\delta }}{\bf{h}}}({{\bf{k}}}^{* })=\pm 1$ is fulfilled for ${k}_{y}^{* }=0,\pm 2\pi /3$, needless to say that such values of k_y give the edges of the quasienergy band structure along the y-direction, we stress out the fact that at the edges of the quasienergy band structure, Hamiltonians ${H}_{0}({\bf{k}})$ and ${H}_{1}({\bf{k}})$ commute. By substituting ${k}_{y}^{* }$ into equation (20), one gets,

$$\begin{eqnarray}&&\tau {\omega }_{\pm }({k}_{x})=\tau {\gamma }_{2}\pm 2\tau {\gamma }_{1}\cos (\sqrt{3}{k}_{x}/2),\end{eqnarray} \tag{ 23 }$$

where the plus sign (+) stems for ${k}_{y}^{* }=0$, while the minus sign (−) stems for ${k}_{y}^{* }=\pm 2\pi /3$. Now, in order to have touching band points, two band edges must touch each other. This occurs whenever $\tau \omega ({k}_{x}^{* })=\pm n\pi$ (n being an integer number). By using equation (23), we find that ${k}_{x}^{* }$ has two possible solutions given by,

$$\begin{eqnarray}{k}_{x}^{* (+)} & = & \pm \displaystyle \frac{2}{\sqrt{3}}\arccos \left[\displaystyle \frac{n\pi -\tau {\gamma }_{2}}{2\tau {\gamma }_{1}}\right],\\ {k}_{x}^{* (-)} & = & \pm \displaystyle \frac{2}{\sqrt{3}}\arccos \left[\displaystyle \frac{-n\pi +\tau {\gamma }_{2}}{2\tau {\gamma }_{1}}\right].\end{eqnarray} \tag{ 24 }$$

As before, ${k}_{x}^{* (+)}$ stems for ${k}_{y}^{* }=0$, while ${k}_{x}^{* (-)}$ stems for ${k}_{y}^{* }=\pm 2\pi /3$. From the structure of equation (24), it is easy to see that touching band points always come in pairs, as in the case of Weyl and Dirac points. We have to mention that for n = 0 and for odd n there are two pairs of touching band points, however this is not the case for even n (n different from zero) for which just one pair of touching band points emerge. This can be understood by looking at equation (24). It is readily seen that for even n both ${k}_{x}^{+}$ and ${k}_{x}^{-}$ are the same. On the other hand, the case n = 0 (i.e. the time-independent touching band points) worths special attention, since in this case the touching band points correspond to Dirac cones shifted from their original position due to the deformation field [84]. As is well known, the Dirac cones give rise to flat bands in the time-independent case when the nanoribbon is considered to be finite, this is still true even in the presence of a time-dependent deformation field [61]. As will be seen later on, touching band points for $n\ne 0$ also give rise to topologically protected flat bands.

It is useful to obtain the conditions to have touching band points, since this sheds light about the topological phase diagram of the system. Such information can be readily obtained by observing that in order to have real solutions for equation (24), the following condition must be satisfied,

$$\begin{eqnarray}&&| n\pi -\tau {\gamma }_{2}| \leqslant 2\tau {\gamma }_{1}.\end{eqnarray} \tag{ 25 }$$

In other words, there is a critical treshold for τ, say ${\tau }_{c}$ for having touching band points. Such value depends upon the ripple's amplitude via ${\gamma }_{1}$ and ${\gamma }_{2}$ (see equation (14)). The explicit form of ${\tau }_{c}$ can be obtained from the extremal limits of equation (25), one can prove that is given by,

$$\begin{eqnarray}&&{\tau }_{c}=\displaystyle \frac{\pi }{2{\gamma }_{1}+{\gamma }_{2}}.\end{eqnarray} \tag{ 26 }$$

It is important to say that each time that τ reaches an integer multiple of ${\tau }_{c}$, new touching band points will emerge, in other words, there will be new pairs of touching band points for $\tau =n{\tau }_{c}$, where n is an integer number. Also observe that bands will touch each other at $\pm \pi$ quasienergy if n is odd, whereas they will touch each other at zero quasienergy for even or vanishing n. From equation (25), we can construct the phase diagram of type I touching band points, however, this phase diagram will be incomplete since it will not contain the information of the type II touching band points. Therefore, we leave the construction of the phase diagram to be done after analyzing type II touching band points.

To finish, let us confirm our results. In figure 2 we used $\lambda =0.5$ and $\tau =\pi$, this is, we have $2{\tau }_{c}^{+}\lt \tau \lt 3{\tau }_{c}^{+}$. Therefore, there must be six pairs of touching band points, three pairs at zero quasienergy (two for n = 0 and one for n = 2) and two pairs at $\pm \pi$ quasienergy (n = 1). This is in completely agreement with figure 2.

3.2. Type II

Let us start by determining the location of this kind of touching band points. To do that, we set $\tau \delta h={n}_{1}\pi$ and $\tau {h}_{0}={n}_{2}\pi$ in equation (20), where n₁ and n₂ are integer numbers. After some algebraic manipulations, one obtains,

$$\begin{eqnarray}&&\cos (\sqrt{3}{k}_{x}^{* }/2)\,=\,\pm \,\sqrt{\displaystyle \frac{\tfrac{{n}_{1}^{2}{\pi }^{2}}{{\tau }^{2}}-{({\gamma }_{2}-1)}^{2}+({\gamma }_{1}-1)({\gamma }_{2}-2)\left(1-\tfrac{{n}_{2}^{2}{\pi }^{2}}{{\tau }^{2}}\right)}{4({\gamma }_{1}-1)({\gamma }_{1}-{\gamma }_{2})}},\\ &&\quad \cos (3{k}_{y}^{* }/2)=\displaystyle \frac{{n}_{2}^{2}{\pi }^{2}-4{\cos }^{2}(\sqrt{3}{k}_{x}^{* }/2)-1}{4\cos (\sqrt{3}{k}_{x}^{* }/2)}.\end{eqnarray} \tag{ 27 }$$

Once again, we can obtain the conditions for having these kind of touching band points by noticing that to ensure having real solutions in equation (27), the following conditions need to be held altogether,

$$\begin{eqnarray}&&0\leqslant \displaystyle \frac{\tfrac{{n}_{1}^{2}{\pi }^{2}}{{\tau }^{2}}-{({\gamma }_{2}-1)}^{2}+({\gamma }_{1}-1)({\gamma }_{2}-2)\left(1-\tfrac{{n}_{2}^{2}{\pi }^{2}}{{\tau }^{2}}\right)}{4({\gamma }_{1}-1)({\gamma }_{1}-{\gamma }_{2})}\leqslant 1\\ &&\,\left|\displaystyle \frac{{n}_{2}^{2}{\pi }^{2}-4{\cos }^{2}(\sqrt{3}{k}_{x}^{* }/2)-1}{4\cos (\sqrt{3}{k}_{x}^{* }/2)}\right|\leqslant 1.\end{eqnarray} \tag{ 28 }$$

It is worthwhile to mention that the band edges will touch each other at $\pm \pi$ quasienergy if n₁ is even and n₂ is odd or vice versa, whereas they will touch each other at zero quasienergy for either n₁ and n₂ even or odd.

The conditions given by equation (28) add new phases to the phase diagram of the system. Such diagram will be built in the next section.

3.3. Topological phase diagram

In figure 3, the phase diagram for types I and II touching band points is presented, such diagram was built from the expressions for the critical values of τ obtained from equation (26) and (28). Therein, in figure 3, type I touching band points are labeled by n and each single value of n gives rise to two pairs of this kind of touching band points. For instance, the region label by n = 0, 1 has four pairs of touching band points, two pairs corresponding to n = 0 at zero quasienergy (Dirac cones, as was discussed above) and the others two pairs at $\pm \pi$ quasienergy corresponding to n = 1. Note also that each value of n corresponds to a well defined region in the phase diagram. When it concerns to type II touching band points things become more complicated since each pair of integers (n₁, n₂) results in very intricate regions on the phase diagram, as is clearly seen in figure 3 in the regions labeled by II. Additionally, for having type II touching band points high values of the ripple's amplitude are required, which makes them difficult to be observed experimentally since non-linear effects may appear before reaching this regimen. Finally, note that the fact that both kinds of touching band points always come in pairs suggests that they can give rise to topologically protected edge modes if the system is considered to be finite. In fact, this is the case as is proven below.

Zoom In Zoom Out Reset image size

In this section we discuss the emergence and the topological properties of edge states in a finite ZGN. In the previous section we found touching band points at which the edges of the quasienergy spectrum cross each other, which is a signature for edge states. In order to confirm if edge states emerge, we calculate the quasienergy spectrum for a finite system. To do that, a numerical diagonalization of the matrix representation of the time evolution operator equation (12), as a function of k_x, is done for fixed σ, ϕ, λ and τ. We also study the localization properties of the wave functions of such states. Using the logarithm of the inverse participation ratio (IPR), which is defined as,

$$\begin{eqnarray}&&{\rm{IPR}}(E)=\displaystyle \frac{{\displaystyle \sum }_{j=1}^{N}| \psi (j){| }^{4}}{\mathrm{ln}N},\end{eqnarray} \tag{ 29 }$$

where $\psi (j)$ is the wave function at site j for a given energy (or quasienergy) E. The IPR is a measure of the wave function localization [4]. The closer the IPR to zero the more localized the wave function is. Whereas for the IPR tending to −1, we have completely delocalized wave functions. Having said that, we proceed with the study of the edge states.

4.1. Type I

Let us consider first the case of type I touching band points. We start by obtaining the quasienergy band structure as a function of k_x via the numerical diagonalization of the matrix representation of equation (12). In figure 4 we show the resulting quasienergy band structure for $\sigma =1/3$, $\phi =0$, $\lambda =0.5$, $\tau =\pi$, N = 164 atoms and obtained by using fixed boundary conditions. We used the same condition as in the analytically obtained plot in figure 2(b). Note the excellent agreement between the numerical and the analytical results.

Zoom In Zoom Out Reset image size

In figure 4(a) we also show the winding number of the effective Hamiltonian, which is basically the winding number of the unit vector defined in equation (21) for ${k}_{x}=0.9\pi /\sqrt{3}$, for a phase with flat bands joining two inequivalent Dirac cones. As can be seen, the winding number is one, as expected from the topological properties of a finite ZGN.

The main difference between figures 4 and 2 (apart from the fact that figure 2 is a three-dimensional plot and figure 4 is the projected band structure as a function of k_x) is that, for a finite nanoribbon, highly localized edge modes are clearly seen in figure 4. In addition, we can see more touching band points in figure 2 than in figure 4 since the former is a three-dimensional plot in perspective (we have plotted the front view of the band structure), whereas the latter is a projection of the full band structure. For example, instead of seeing four Dirac cones in figure 4, as happens in figure 2, we just see two Dirac cones because the projection superposes each pair, as happens with other touching band points. The colors used in figure 4 represent the logarithm of the IPR (as defined in equation (29)), blue colors correspond to totally delocalized states and red color represents highly localized wave functions. Also observe how flat bands join two inequivalent touching band points, which suggests that inequivalent touching band points at the same quasienergy have opposite Berry phase. In fact, this is the case for n = 0, which corresponds to Dirac cones, labeled by gray dots in figure 4. This also happens for $n\ne 0$. Before studying the Berry phase of the touching band points and for the sake of clarity, in figure 5 we present the analytical and the numerical band structure of our system for $\sigma =1/3$, $\phi =0$, $\tau =6$, and $\lambda =0.6$. These parameters were chosen in such a way that only type I touching band points appear. In panel figure 5(a) we can observe many touching band points at zero and $\pm \pi$ quasienergies. Each pair produces flat bands as seen in panel (b) of the same figure. It is important to note that the flat bands become more extended as the driving period is increased.

Zoom In Zoom Out Reset image size

To confirm the previous conjecture about the topological nature of the touching band points, we explicitly evaluate the Berry phase for type I touching band points. We start by noticing that near the touching band points the quasienergy spectrum is well described by the one-period time evolution operator, equation (17), expanded up to second order in powers of τ. By using the Baker–Campbell–Hausdorff formula in equation (17), one gets,

$$\begin{eqnarray}&&{ \mathcal U }({\bf{k}},\tau )\approx \exp \{-{\rm{i}}\tau {H}_{1}({\bf{k}})+{\tau }^{2}[{H}_{1}({\bf{k}}),{H}_{0}({\bf{k}})]/2\}.\end{eqnarray} \tag{ 30 }$$

Since we are just interested in what happens in the neighborhood of touching band points, we expand equation (30) around $({k}_{x}^{* },{k}_{y}^{* })$. It is straightforward to show that equation (30) can be written as

$$\begin{eqnarray}&&{ \mathcal U }({q}_{x},{q}_{y},\tau )\approx \exp [-{{\rm{i}}{h}}_{T}\,{\hat{{\bf{h}}}}_{T}\cdot {\boldsymbol{\sigma }}],\end{eqnarray} \tag{ 31 }$$

where ${q}_{x}={k}_{x}-{k}_{x}^{* }$, ${q}_{y}={k}_{y}-{k}_{y}^{* }$, ${h}_{T}=| {{\bf{h}}}_{T}|$, and the vector ${{\bf{h}}}_{T}$ is given by,

$$\begin{eqnarray}&&{{\bf{h}}}_{T}=A(\lambda ,\tau ){q}_{x}{\hat{e}}_{x}+B(\lambda ,\tau ){q}_{y}{\hat{e}}_{y}+C(\lambda ,\tau ){q}_{y}{\hat{e}}_{z},\end{eqnarray} \tag{ 32 }$$

with

$$\begin{eqnarray}A({\gamma }_{1},{\gamma }_{2},\tau ) & = & n\pi +\tau {\gamma }_{1}\sqrt{3+3{\left(\displaystyle \frac{n\pi -\tau {\gamma }_{2}}{2\tau {\gamma }_{1}}\right)}^{2}},\\ B({\gamma }_{1},{\gamma }_{2},\tau ) & = & 3\tau {\gamma }_{2}/2,\\ C({\gamma }_{1},{\gamma }_{2},\tau ) & = & \displaystyle \frac{3\tau }{4{\gamma }_{1}}({\gamma }_{1}-{\gamma }_{2})(n\pi -\tau {\gamma }_{2}).\end{eqnarray} \tag{ 33 }$$

Finally, the Berry phase can be readily obtained from the effective Hamiltonian ${h}_{T}\,{\hat{{\bf{h}}}}_{T}\cdot {\boldsymbol{\sigma }}$. As we prove in appendix B, the Berry phase, ${\gamma }_{C}$, is non-vanishing for touching band points at ${k}_{x}^{+}$, in fact, its value is ${\gamma }_{C}=\pi$. For touching band points at ${k}_{x}^{-}$ the Berry phase takes the opposite value as for ${k}_{x}^{+}$, this is, we have ${\gamma }_{C}=-\pi$. Therefore flat bands joining two touching band points with opposite Berry phase will emerge. Needless to say that these touching band points are topologically protected, so flat bands are topologically non-trivial.

4.2. Type II

Now we analyze the edge states originated from type II touching band points. First, we obtain the quasienergy band structure from the numerical diagonalization of equation (12) for a set of parameters within one of the regions II of the diagram phase figure 3. In figure 6, we show such band structure for $\sigma =1/3$, $\phi =0$, $\lambda =1$, $\tau =5.46$, and N = 164, and obtained using fixed boundary conditions. Observe that in figure 4(b) besides the type I touching band points there is one pair of type II touching band points. As in the case of type I touching band points, edge states emerge from type II touching band points, these edge states seem to be also flat bands. However, as edge states approach k_x = 0, they are no longer flat bands but they become dispersive delocalized states, see the inset in figure 6(b), where a zoom around $\pm \pi$ quasienergy is shown. To get further insight about the edge states that emerge from type II touching band points we plotted, in figure 7, the analytical and numerical quasienergy band structure for $\sigma =1/3$, $\phi =0$, $\lambda =0.9$, and $\tau =7.5$. Observe that the agreement between the numerical (panel b)) and analytic (panel a)) results is quite good. As before, the edge states that appear in panel (b) are dispersive and join two inequivalent touching band points. In addition, the edge states in figure 7 are less localized that the ones in figure 5.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The fact that these edge states start and end at type II touching band points suggest that they have non-trivial topological properties. To study the topological properties of this kind of edge states we cannot proceed as we did with type I touching band points since type II touching band points do not correspond to points at where Hamiltonians (15) commute. Therefore, we analyze the topological properties of a 1D slice of the system, in other words, we study our system for a fixed k_x. Once that we have fixed k_x, the topological properties can be obtained from the winding of the unit vector ${\hat{{\bf{h}}}}_{\mathrm{eff}}$ that appears in the effective Hamiltonian equation (19), since a non-vanishing winding number is a signature of non-trivial topological properties. If ${\hat{{\bf{h}}}}_{\mathrm{eff}}$, for fixed k_x, has a non-vanishing winding number around the origin, then the 1D slice has non-trivial topological properties and the whole 2D system is topologically weak [85–87]. In figure 6(a) we show the winding of the unit vector ${\hat{{\bf{h}}}}_{\mathrm{eff}}$ as a function of k_y obtained from the analytical expression equation (21) for ${k}_{x}=0.9/\sqrt{3}$, $\sigma =1/3$, $\phi =0$, $\tau =5.46$, and $\lambda =1$. As clearly seen in the figure, the winding number is 6, which means that our 1D slice has non-trivial topological properties and that the whole 2D system is topologically weak.

We have studied the case of a time-periodically driven rippled ZGN. We obtained the quasienergy spectrum of the time-evolution operator. As a result, two types of touching band points were found for a special value of the corrugation wavelength ($\sigma =1/3$). Each type produces different edge states. For type I edge states, we found that the edge states are flat bands joining two inequivalent touching band points with opposite Berry phase, this was confirmed by the analytical evaluation of the Berry phase. On the other hand, type II edge states were found to have a topological weak nature. This was done by a numerical calculation of the winding number of a 1D slice of the system, in other words, by looking at the topological properties of our system for a fixed k_x. Using this previous information, the phase diagram of the system was built. To finish, we stress out that the experimental realization of our model can be very challenging, however, there are some proposed experiments for similar situations [61, 88, 89]. Experimentally is possible to create a 1D uniaxial ripple of graphene by growing it over a substrate [70]. Then the driving can be achieved by time-periodically applying pressure to the whole system (i.e. to the graphene ribbon and substrate). Time scales of femto seconds are needed to observe the phenomena discussed above, a fact that requires the use of, for example, femto lasers of Ti-Sapphire to induce deformations. As an alternative, optical lattices can be used since the hopping parameters can be tailored at will [88, 89].

Finally, it is important to remark that for observing the edge states studied here, the time driving layout does not need to be a delta driving. Even a cosine-like time perturbation can be used. However, for the case of a cosine-like time-perturbation, the effect could be hard to be observed since the secular gaps are usually smaller [61].

This work was supported by DGAPA-PAPIIT Project 102717. PR-T acknowledges financial support from Consejo Nacional de Ciencia y Tecnología (CONACYT) (México).

In this appendix we analytically obtain the quasienergy spectrum for $\sigma =1/3$, $\phi =0$. As was mentioned in the main text, for $\sigma =1/3$, the system becomes periodic along both the x and y directions. As a result, we can Fourier transform the Hamiltonians (9) and (10) taking advantage of such periodicity. By using the following Fourier transformations,

$$\begin{eqnarray}{a}_{j} & = & \displaystyle \frac{1}{\sqrt{N/2}}\displaystyle \sum _{{k}_{y}}{{\rm{e}}}^{-{\rm{i}}3{k}_{y}j/2}{a}_{{k}_{y}}\\ {b}_{j} & = & \displaystyle \frac{1}{\sqrt{N/2}}\displaystyle \sum _{{k}_{y}}{{\rm{e}}}^{-{\rm{i}}3{k}_{y}j/2}{b}_{{k}_{y}},\end{eqnarray} \tag{ A.1 }$$

and after some algebraic manipulations, one gets the simplified Fourier transformed version of Hamiltonians equations (9) and (10),

$$\begin{eqnarray}{H}_{0}({\bf{k}}) & = & {h}_{0}({\bf{k}}){\hat{{\bf{h}}}}_{{\bf{0}}}({\bf{k}})\cdot {\boldsymbol{\sigma }}\\ {H}_{1}({\bf{k}}) & = & {h}_{1}({\bf{k}}){\hat{{\bf{h}}}}_{{\bf{1}}}({\bf{k}})\cdot {\boldsymbol{\sigma }},\end{eqnarray} \tag{ A.2 }$$

where ${\bf{k}}=({k}_{x},{k}_{y})$, ${\sigma }_{i}$ ($i=x,y,z$) are the 2 × 2 Pauli matrices, $\hat{{{\bf{h}}}_{{\bf{0}}}}({\bf{k}})={{\bf{h}}}_{0}({\bf{k}})/| {h}_{0}({\bf{k}})|$, $\hat{{{\bf{h}}}_{{\bf{1}}}}({\bf{k}})={{\bf{h}}}_{1}/| {h}_{1}({\bf{k}})|$ $[{h}_{0}({\bf{k}})$ (${h}_{1}({\bf{k}})$) being the norm of ${{\bf{h}}}_{0}({\bf{k}})$ (${{\bf{h}}}_{1}({\bf{k}})$)]. ${{\bf{h}}}_{0}({\bf{k}})$ and ${{\bf{h}}}_{1}({\bf{k}})$ have components given by

$$\begin{eqnarray}{h}_{0}^{(x)}({\bf{k}}) & = & 2\cos (\sqrt{3}{k}_{x}/2)+\cos (3{k}_{y}/2),\\ {h}_{0}^{(y)}({\bf{k}}) & = & \sin (3{k}_{y}/2),\\ {h}_{1}^{(x)}({\bf{k}}) & = & 2{\gamma }_{1}\cos (\sqrt{3}{k}_{x}/2)+{\gamma }_{2}\cos (3{k}_{y}/2),\\ {h}_{1}^{(y)}({\bf{k}}) & = & {\gamma }_{2}\sin (3{k}_{y}/2).\end{eqnarray} \tag{ A.3 }$$

${\gamma }_{1}$ and ${\gamma }_{2}$ have been defined in equation (14). By using equation (A.2), the time evolution operator equation (12) can be written as

$$\begin{eqnarray}&&U({\bf{k}},\tau )=\displaystyle \sum _{{k}_{y}}{ \mathcal U }({\bf{k}},\tau )\otimes | {k}_{y}\rangle \langle {k}_{y}| .\end{eqnarray} \tag{ A.4 }$$

Here $\delta H({\bf{k}})={H}_{1}({\bf{k}})-{H}_{0}({\bf{k}})$, and

$$\begin{eqnarray}&&{ \mathcal U }({\bf{k}},\tau )={\rm{\exp }}[-{\rm{i}}\tau \delta H({\bf{k}})]{\rm{\exp }}[-{\rm{i}}\tau {H}_{0}({\bf{k}})].\end{eqnarray} \tag{ A.5 }$$

Even though, H₁ and H₀ generally do not commute, one can rewrite equation (A.5) as follows,

$$\begin{eqnarray}&&{ \mathcal U }({\bf{k}},\tau )=\exp [-{\rm{i}}\tau {H}_{\mathrm{eff}}({\bf{k}})],\end{eqnarray} \tag{ A.6 }$$

where the effective Hamiltonian is given by

$$\begin{eqnarray}&&{H}_{\mathrm{eff}}({\bf{k}})=\omega ({\bf{k}}){\hat{{\bf{h}}}}_{\mathrm{eff}}({\bf{k}})\cdot {\boldsymbol{\sigma }},\end{eqnarray} \tag{ A.7 }$$

the quasienergies $\tau \omega ({\bf{k}})$ are given by the next relation,

$$\begin{eqnarray}\cos [\tau \omega ({\bf{k}})] & = & \cos [\tau \,\delta h({\bf{k}})]\cos [\tau {h}_{0}({\bf{k}})]\\ & & -{\hat{{\bf{h}}}}_{{\bf{0}}}({\bf{k}})\cdot \hat{{\boldsymbol{\delta }}{\bf{h}}}({\bf{k}})\sin [\tau \,\delta h({\bf{k}})]\sin [\tau {h}_{0}({\bf{k}})],\end{eqnarray} \tag{ A.8 }$$

where ${\boldsymbol{\delta }}{\bf{h}}({\bf{k}})={{\bf{h}}}_{1}({\bf{k}})-{{\bf{h}}}_{0}({\bf{k}})$, and

$$\begin{eqnarray}&&{\hat{{\bf{h}}}}_{{\bf{0}}}({\bf{k}})\cdot \hat{{\boldsymbol{\delta }}{\bf{h}}}({\bf{k}})=\displaystyle \frac{1}{{h}_{0}({\bf{k}})\delta h({\bf{k}})}\left[4({\gamma }_{1}-1){\cos }^{2}\left(\displaystyle \frac{\sqrt{3}}{2}{k}_{x}\right)\right]\\ &&\quad +\,\displaystyle \frac{1}{{h}_{0}({\bf{k}})\delta h({\bf{k}})}\left[2({\gamma }_{1}+{\gamma }_{2}-2)\cos \left(\displaystyle \frac{\sqrt{3}}{2}{k}_{x}\right)\cos \left(\displaystyle \frac{3{k}_{y}}{2}\right)\right]\\ &&\quad +\,\displaystyle \frac{{\gamma }_{2}-1}{{h}_{0}({\bf{k}})\delta h({\bf{k}})}.\end{eqnarray} \tag{ A.9 }$$

Finally, the unit vector ${\hat{{\bf{h}}}}_{\mathrm{eff}}({\bf{k}})$ is given by

$$\begin{eqnarray}{\hat{{\bf{h}}}}_{\mathrm{eff}}({\bf{k}}) & = & \displaystyle \frac{-1}{\sin [\tau \omega ({\bf{k}})]}[\hat{{\boldsymbol{\delta }}{\bf{h}}({\bf{k}})}\sin [\tau \,\delta h({\bf{k}})]\cos [\tau {h}_{0}({\bf{k}})]]\\ & & +\displaystyle \frac{-1}{\sin [\tau \omega ({\bf{k}})]}[{\hat{{\bf{h}}}}_{0}({\bf{k}})\sin [\tau {h}_{0}({\bf{k}})]\cos [\tau \delta h({\bf{k}})]]\\ & & +\displaystyle \frac{-1}{\sin [\tau \omega ({\bf{k}})]}[\hat{{\boldsymbol{\delta }}{\bf{h}}({\bf{k}})}\times {\hat{{\bf{h}}}}_{0}({\bf{k}})\sin [\tau \delta h({\bf{k}})]\sin (\tau {h}_{0}[{\bf{k}})]].\end{eqnarray} \tag{ A.10 }$$

In this appendix, the explicit evaluation of the Berry phase for type I touching band points is done. The Berry phase is defined as

$$\begin{eqnarray}&&{\gamma }_{C}={\oint }_{C}{\bf{A}}\cdot \,{\rm{d}}{\bf{k}},\end{eqnarray} \tag{ B.1 }$$

where ${\bf{A}}=-{\rm{i}}\langle {\psi }_{k}| {{\boldsymbol{\nabla }}}_{k}| {\psi }_{k}\rangle$ is the so-called Berry connection (a gauge invariant quantity), and ${{\rm{\nabla }}}_{k}=({\partial }_{{k}_{x}},{\partial }_{{k}_{y}})$ is the gradient operator in the momentum space. Since we are interested in what happens in the neighborhood of type touching band points, it is enough to calculate the Berry phase of ${\hat{{\bf{h}}}}_{T}\cdot {\boldsymbol{\sigma }}$, which is the effective Hamiltonian in the neighborhood of type I touching band points and that is defined in equation (32).

To obtain the Berry phase, we first need to calculate the eigenvectors of Hamiltonian equation (32), it can be proven that such eigenvectors are given by the following spinors,

$$\begin{eqnarray}| {\psi }_{q^{\prime} }^{\uparrow }\rangle & = & \displaystyle \frac{1}{\sqrt{2}}\left(\begin{array}{c}\sqrt{1+\tfrac{C\,q{{\prime} }_{y}}{B\,{h}_{T}}}\\ {{\rm{e}}}^{{\rm{i}}\xi {\alpha }_{q^{\prime} }}\sqrt{1-\tfrac{C\,q{{\prime} }_{y}}{B\,{h}_{T}}}\end{array}\right)\\ | {\psi }_{q^{\prime} }^{\downarrow }\rangle & = & -\displaystyle \frac{1}{\sqrt{2}}\left(\begin{array}{c}{{\rm{e}}}^{-{\rm{i}}\xi {\alpha }_{q^{\prime} }}\sqrt{1-\tfrac{C\,q{{\prime} }_{y}}{B\,{h}_{T}}}\\ -\sqrt{1+\tfrac{C\,q{{\prime} }_{y}}{B\,{h}_{T}}}\end{array}\right),\end{eqnarray} \tag{ B.2 }$$

where

$$\begin{eqnarray}q{{\prime} }_{x} & = & {q}_{x}/A\\ q{{\prime} }_{y} & = & {q}_{y}/B\end{eqnarray} \tag{ B.3 }$$

and ${\alpha }_{{q}^{{\prime} }}$ is given by,

$$\begin{eqnarray}&&{\alpha }_{q^{\prime} }={\tan }^{-1}\left(\displaystyle \frac{q{{\prime} }_{y}}{q{{\prime} }_{x}}\right).\end{eqnarray} \tag{ B.4 }$$

ξ can take the values $\xi =+1$ which corresponds to $+{k}_{x}^{* (+)}$ and $\xi =-1$ to $-{k}_{x}^{* (+)}$. Now, the Berry connection can be calculated using such spinors, for simplicity we set $\xi =1$, however the result does not depend upon ξ. After some calculations, one obtains that the Berry connection is,

$$\begin{eqnarray}&&{\bf{A}}=\displaystyle \frac{1}{2}\left(1-\displaystyle \frac{C}{B\,{h}_{T}}q{{\prime} }_{y}\right){{\rm{\nabla }}}_{q^{\prime} }{\alpha }_{q^{\prime} },\end{eqnarray} \tag{ B.5 }$$

where

$$\begin{eqnarray}&&{{\rm{\nabla }}}_{q^{\prime} }{\alpha }_{q^{\prime} }=\displaystyle \frac{-q{{\prime} }_{y}\,{\hat{e}}_{x}+q{{\prime} }_{x}\,{\hat{e}}_{y}}{{(q{{\prime} }_{x})}^{2}+{(q{{\prime} }_{y})}^{2}}.\end{eqnarray} \tag{ B.6 }$$

Finally, we calculate the Berry phase along a circumference centered at $q{{\prime} }_{x}=q{{\prime} }_{y}=0$. By using polar coordinates, defined as, $q{{\prime} }_{x}=q^{\prime} \cos \theta$ and $q{{\prime} }_{y}=q^{\prime} \sin \theta$ where ${(q^{\prime} )}^{2}={(q{{\prime} }_{x})}^{2}+{(q{{\prime} }_{y})}^{2}$, the Berry connection is readily obtained,

$$\begin{eqnarray}{\gamma }_{C} & = & {\displaystyle \int }_{0}^{2\pi }\,{\bf{A}}\cdot {\rm{d}}{\bf{q}}^{\prime} \\ & = & \displaystyle \frac{1}{2}{\displaystyle \int }_{0}^{2\pi }\left(1-\displaystyle \frac{\tfrac{C}{B}{\rm{\sin }}\theta }{\sqrt{1+\tfrac{{C}^{2}}{{B}^{2}}{{\rm{\sin }}}^{2}\theta }}\right){\rm{d}}\theta =\pi .\end{eqnarray} \tag{ B.7 }$$

A similar calculation can be done for ${k}_{x}^{* (-)}$, which gives ${\gamma }_{C}=-\pi$.