# Strain-induced topological phase transition in phosphorene and in phosphorene nanoribbons

###### Abstract

Using the tight-binding (TB) approximation with inclusion of the spin-orbit interaction, we predict a topological phase transition in the electronic band structure of phosphorene in the presence of axial strains. We derive a low-energy TB Hamiltonian that includes the spin-orbit interaction for bulk phosphorene. Applying a compressive biaxial in-plane strain and perpendicular tensile strain in ranges where the structure is still stable, leads to a topological phase transition. We aslo, examine the influence of strain on zigzag phosphorene nanoribbons (zPNRs) and the formation of the corresponding protected edge states when the system is in the topological phase. For zPNRs up to a width of nm the energy gap is at least three orders of magnitude larger than the thermal energy at room temperature.

## I Introduction

Topological insulators (TIs) with time-reversal symmetry (TRS), have been of increasing interest in condensed matter physics and material science during the last decade. The emergence of robust edge states in two-dimensional (2D) TIs that are protected by TRS, make them promising candidates for potential applications in spintronics and quantum computing Kane and Mele (2005); Hasan and Kane (2010); Moore (2010); Qi and Zhang (2011); Fu et al. (2007); Bernevig et al. (2006). TIs can exist intrinsically or be driven by external factors such as electrical field or by functionalization Ren et al. (2016). Strain engineering is a well known strategy for switching from normal insulator (NI) phase to a TI phase Ren et al. (2016); Ma et al. (2014). Among the wide list of systems that possesses such property, 2D materials with fascinating electronic, mechanical and thermal properties have been in the focus of attention Kane and Mele (2005); Ezawa (2012).

In the past few years, phosphorene, a monolayer of black phosphorus, has emerged as an encouraging 2D semiconducting material for widespread applications. Phosphorene-based field effect transistors (FETs), show a higher ON/OFF ratio in comparison with graphene Koenig et al. (2014); Li et al. (2014a) and has a higher carrier mobility with respect to 2D transition metal dichalcogenides (TMDs) which have recently attracted a lot of attention for FET applications Koenig et al. (2014); Li et al. (2014a); Xia et al. (2014). There exist several works pertinent to the observation of different phases in bulk and multilayer black phosphorous by tuning the lowest energy bands Fei et al. (2015); Kim et al. (2015); Xiang et al. (2015); Liu et al. (2015); Zhang et al. (2015). Using density functional theory (DFT) it was shown that few-layers of phosphorene experiences a NI to TI and then a TI to topological metal (TM) phase transition by applying a perpendicular electric field Liu et al. (2015). In a different DFT study Zhang et al. (2015) such phase transitions for various stacked bilayer phosphorene under in-plane strain has been explored. Owing to the puckered structure of phosphorene, it has a high degree of flexibility. Therefore, it can sustain strain very well specially in the zigzag direction up to about 30% Wei and Peng (2014); Peng et al. (2014). This makes phosphorene promising for possible applications using strain engineering.

In our work, we investigate the effect of strain on the electronic band structure of phosphorene the TB approach. The band gaps of this model Rudenko and Katsnelson (2014) are close to the most reliable DFT and experimental results Liang et al. (2014); Tran et al. (2014) that predict band gaps of eV for phosphorene. In this paper, we propose a model Hamiltonian for the SOC for monolayer phosphorene that can be generalized to few-layers phosphorene. We show that, a model which includes the next-nearest(n-n) neighbors in the upper or lower chains, is sufficient for capturing the main physics. Then, strain engineering of this system is investigated through modifying the hopping parameters of the system. We demonstrate that, by applying particular types of strain, the system can make a phase transition to a TI. Finally, we show numerically that though the topological bulk band gaps induced by SOC is about meV, but the highly anisotropic nature of this material causes the corresponding bulk gaps in large widths zPNRs be at least three order of magnitude larger than room temperature thermal energy ( meV) and makes phosphorene nanoribbons excellent candidates for future applications.

This paper is organized as follows: the effective low-energy TB model Hamiltonian including the SOC terms is obtained in Sec. II. The effect of axial strains on the band structure produced by this model is calculated and our results are compared with DFT results in Sec. III. Demonstration of a topological phase transition in the electronic properties of phosphorene when particular types of strain are applied and the characteristics of corresponding edge states in zPNRs is presented in Sec. IV. The paper is summarized in Sec. V.

## Ii Tight-Binding Model Hamiltonian Including Spin-Orbit Interaction

### ii.1 Structure

The puckered atomic structure of phosphorene and its geometrical parameters are shown in Fig. 1 where the and axes are the armchair and zigzag directions, respectively and the axis is in the normal direction to the plane of phosphorene. With this definition of coordinates, one can indicate the various atom connections which correspond to various hopping parameters that are included in the TB model. The structure parameters have been taken from Qiao et al. (2014) which is very close to experimentally measured parameters Takao and Morita (1981) for its bulk structure. The components of the geometrical parameters as shown in Figs. 1(b) and (c), for bond lengths Å and Å are and , and are simply defined by parameters of and . The two in-plane lattice constants are Å, Å and the thickness of a single layer due to the puckered nature is Å.

### ii.2 Tight-binding model

The phosphorene TB Hamiltonian that has been proposed earlier Rudenko and Katsnelson (2014), without the spin degree of freedom, is given by

(1) |

where the summation is up to fifths neighbors, and are hopping integrals that show the energy transfer between the th and th sites. The hopping terms are shown in Fig. 1(a). and represent the creation and annihilation operators of electrons in sites and , respectively. The numerical values of these hopping parameters are: eV, eV, eV, eV, and eV Rudenko and Katsnelson (2014). Including the spin degree of freedom and SOC the Hamiltonian is modified into

(2) |

where in , the first term is called the usual effective SOC and the second one is the intrinsic Rashba SOC which will be introduced in next subsection. Due to the puckered structure of phosphorene, the Rashba term is rather large as compared to the first term and should be included in our calculations.

### ii.3 Spin-orbit coupling in Phosphorene

The primary goal of this subsection is to introduce a spin-orbit model Hamiltonian for phosphorene which can capture the most important spin-related phenomenon. There exist several studies which showed the anisotropic behaviour in the electronic and optical properties of phosphorene Fei and Yang (2014); Li et al. (2014b); Tran et al. (2014); Sisakht et al. (2015) which are due to the anisotropic nature of the band dispersion of phosphorene. This property is reflected in the effective mass of electrons and holes of phosphorene. As a matter of fact, the corresponding band dispersion of the zigzag direction in real space, is relatively flat near the Fermi energy while it has an approximately linear dispersion in the armchair direction Fei and Yang (2014); Sisakht et al. (2015). One can define two types of n-n neighbors in the phosphorene structure. As shown in Fig. 1(a), each P atom has two intra-chain and four inter-chain n-n neighbors, respectively. The effective mass of electrons in the direction of intra-chain, are at least an order of magnitude larger than the inter-chain direction Fei and Yang (2014). Therefore, electrons usually select the inter-chain path for circular motion, allowing us to ignore the intra-chain neighbors and only consider the four n-n inter-chain P atoms in the SOC model.

In general, the SOC term for a 2D system is given by

(3) |

where , and are Plank’s constant, mass of free electron, and the velocity of light, respectively. is the effective electrostatic force, is the effective momentum and denotes the Pauli matrices. As in the cases of graphene and silicene Liu et al. (2011), the nearest-neighbor SOC is zero in phosphorene, but the SOC terms of the n-n neighbors are nonzero.

As shown in Fig. 1(b), in a honeycomb-like ring of phosphorene, we can define and as vectors that connect the nearest P atoms to each other and the connecting vector of n-n neighbors. Using these vectors, the electrostatic force and momentum can be written as and , with being a prefactor. Rewriting the SOC in terms of the above definitions we obtain

(4) |

Based on experimental and DFT data, and are approximately equal Qiao et al. (2014); Takao and Morita (1981); Castellanos-Gomez et al. (2014); Wei and Peng (2014), therefore and become perpendicular to each other. This leads to

(5) |

where the term will be adjusted to obtain the correct value of SOC as obtained by DFT. Notice that, the above approximations reduce the two parameters of the usual SOC and intrinsic Rashba SOC into a single parameter. Using , where () are the in-plane (out of plane) Pauli matrixes (matrix), we rewrite Eq. (5) as

(6) |

where and is a dimensionless unit vector. The spin-orbit terms in second quantization are given by

(7) |

where and are effective intrinsic SOC and intrinsic Rashba constants, and the summation runs over the inter-chain n-n neighbors. As mentioned before, these two parameters are related to one parameter , which can be estimated by adjusting the TB band structure of phosphorene to the one obtained from DFT. It was shown that in the absence of SOC the energy gap of few-layers phosphorene closes under an external electric field or strain Liu et al. (2015); Zhang et al. (2015). However, including the SOC an energy gap of meV Liu et al. (2015) remains in few-layers phosphorene. This results in the value of meV/ in our TB model.

## Iii Phosphorene Under Strain: Electronic Band Structure

The role of uniaxial and biaxial strain in manipulating the electronic structure of few-layers phosphorene has been investigated via DFT Rodin et al. (2014); Peng et al. (2014); Wang et al. (2015); Zhang et al. (2015); Huang and Xing (2014) and TB approaches Jiang and Park (2015); Mohammadi and Nia (2016); Duan et al. (2016). Applying tensile or compressive strain in different directions results in different modifications of the electronic bands. One can observe a direct to indirect gap transition, or a prior direct band gap closing, depending on the type of applied strain Peng et al. (2014); Wang et al. (2015); Zhang et al. (2015). In this work we consider biaxial compressive strain in the plane of few-layers phosphorene Wang et al. (2015); Zhang et al. (2015), and tensile strain in the normal direction Huang and Xing (2014). This modifies the low energy bands so that the valance and conduction bands approach each other. By further increasing strain, the lower band, coming from orbitals, shifts upward resulting in a semi-metal phase Wang et al. (2015) given that at the band crossing point a mini gap opens due to the SOC. Investigating the local density of states of orbitals Zhang et al. (2015) shows that our one orbital -like TB model is still valid in the low energy limit before the semi-metal phase appears.

In the following, we will first study the bulk band of phosphorene in the presence of axial strains using our TB approach and demonstrate that a band inversion occurs in the energy spectrum of phosphorene in the range where the structure is still stable under strain. It has been shown that the bond lengths and bond angles of phosphorene both change under axial strains Wang et al. (2015); Sa et al. (2014). Therefore, the hopping parameters will change. According to the Harisson rule Harrison (1999); Tang et al. (2009), the hopping parameters for orbitals are related to the bond length as and the angular dependence can be described by the hopping integrals along the and bonds. However, our calculations showed that, though the changes in angles are almost noticeable Wang et al. (2015); Sa et al. (2014), the modification of the hopping parameters due to them is much smaller than the effect of changes of bond lengths. Hence, we consider only changes of the bond lengths in the hopping modulation.

When an axial strain is applied to phosphorene, the rectangle shape of the unit cell with lattice constants of and remains unchanged. Therefore the initial geometrical parameter is deformed as where is the strain in the -direction and is a deformed geometrical parameter. In the linear deformation regime, expanding the norm of to first order of gives

(8) |

where are coefficients related to the structure of phosphorene which are simply calculated via the special geometrical parameters given in previous section. Using the Harrison relation, we obtain the strain effect on the hopping parameters as

(9) |

where is the modified hopping parameter of deformed phosphorene with new lattice constants and .

Let us now study the energy spectrum of strained phosphorene with the modified hopping parameters as given by Eq. (9). The unit cell of monolayer phosphorene is a rectangle containing four atoms as shown in Fig. 1(a). Fourier transform of the strained Hamiltonian of Eq. (2) gives the general Hamiltonian in momentum space as

(10) |

where we have used the basis with being

(11) |

where

(12) |

are matrices

(13) |

whose elements are given by

with , and .

The energy spectrum of pristine phosphorene in the absence of strain has been obtained by numerical diagonalization of the TB Hamiltonian Eq. (10) in different symmetry directions as shown in Fig. 2(a). As we can see in Fig. 2(b), the degeneracies of bands have been removed (black lines) slightly due to the SOC in comparison with the case of zero SOC coupling (red lines) except for the time reversal invariant momentas (TRIMs) which are at least doubly degenerate according to the Kramers theorem.

As seen in Fig. 2 the gap of phosphorene is located at the point which is also a TRIM. At this point, the spin up and spin down valence and conduction bands are degenerate and the change in the gap due to the SOC is very small as compared to the bulk gap. Since axial strain doesn’t break TRS, the bands at this point remain degenerate. Therefore, when the bulk gap is modified by an external factor such as strain, we can safely use the spinless Hamiltonian demonstrating the general trend in changes of the gap. All P atoms in a unit cell have the same on-site energy, so we can project the position of upper and lower chains of phosphorene on a horizontal plane to reduce the spinless Hamiltonian into a two-band TB model Sisakht et al. (2015); Ezawa (2014). Therefore the new -space Hamiltonian of the strained phosphorene in the absence of spin is given by

(15) |

Diagonalizing this Hamiltonian at the point gives the band gap as

(16) | |||||

where denotes the summation over , , components. The first bracket is the unstrained band gap i.e. eV and the second one indicates the structural dependent values of changes in the band gap due to the axial strains. Inserting the numerical values of the structural parameters in Eq. (16) we obtain a compact form for the gap equation

(17) |

where eV, eV and eV. Eq. (17) shows that by applying in-plane compressive biaxial strain and perpendicular tensile strain, the band gap decreases which is consistent with DFT calculations Rodin et al. (2014); Peng et al. (2014); Wang et al. (2015); Zhang et al. (2015); Huang and Xing (2014). It is shown that DFT calculations using the PBE exchange correlation functional anticipate properly the general trends of the band structure when applying axial strains on phosphorene Peng et al. (2014); Wang et al. (2015). A comparison between the band gaps as function of axial strains using available DFT data Peng et al. (2014); Wang et al. (2015); Huang and Xing (2014) and TB model demonstrate that the modification of the hopping parameters in the linear regime are valid for rather large strains and show that the modified TB model predicts correctly the variation of the low energy spectrum. Figure 3 shows the band gap values evaluated at the point in the presence of (a) uniaxial perpendicular tensile strain (b) uniaxial compressive strain in armchair direction, and (c) biaxial compressive in-plane strain, respectively. In both DFT and TB approaches the band gaps exhibit linear dependence with applied strain. The discrepancy between the values of the band gaps originate from the specific calculation method.

As a particular case we consider the modification of energy the spectrum under a perpendicular tensile strain. By increasing the tensile strain, a band inversion occurs at the critical value of . This is a signal of a topological phase transition. Figs. 2(c), (d) show the low energy bands just before and after band closing at 11.5% and 12.5% tensile strain, respectively. As shown in the inset of Fig. 2(d), the SOC opens a small gap of about 5 meV after band closing preventing the formation of a Dirac like-cone.

## Iv Topological Phase Transition of Phosphorene Under Strain

The classification provides a very strong distinction between two different time reversal topological and trivial phases. Pristine phosphorene as a trivial insulator when the intrinsic SOC effect is included preserves the TRS and can exhibit a quantum spin Hall (QSH) phase when its electronic properties is influenced by external factors e.g. electric field or strain. In the following, we first briefly describe our approach for calculating the invariant. This approach, when working in the frame of the TB model Soluyanov (2012) is quite efficient for 2D materials such as phosphorene. Then, we will demonstrate numerically a topological phase transition in strained phosphorene and calculate the phase diagrams accordingly. Finally we will show the existence of protected edge states in zPNRs and discuss their fascinating properties.

### iv.1 Calculation of invariant

Fu and Kane Fu and Kane (2006) showed that an equivalent way to calculate the invariant is as an integral over half the Brillouin zone given by

(18) |

where HBZ denotes half the Brillouin zone. is the Berry gauge potential and the Berry field strength is written as where is the periodic part of the Bloch state with band index and the summation runs over all occupied states. According to Stoke’s theorem, it is obvious that if and have the same gauge which is smooth over HBZ, the result will vanish. Therefore, one needs to fix the gauge with some additional constraints Soluyanov and Vanderbilt (2011). By choosing a gauge, in which the corresponding states fulfills the TRS constraints in addition to the periodicity of the points, that are related by a reciprocal lattice , the gauge fixing procedure is complete and the returned results of or represents the trivial and topological phases, respectively. In the case of phosphorene, where bands cross or degeneracies are present in the energy spectrum, the Berry potential and Berry field strength must be extended to non-Abelian gauge field analogies Fukui et al. (2005) associated with a ground state multiplet in the equation .

Based on the above extension, the discretized Brillouin zone version Fukui and Hatsugai (2007) of Eq. (18) for numerical computing the invariant, is written as

(19) |

where each site in the square lattice of the Brillouin zone of phosphorene is labeled by and l specifies a plaquette with so-called unimodular link variable

(20) |

where denotes a unit vector in - plane. Such a link variable allows us to define the Berry potential and Berry field as

(21) | |||||

(22) |

Berry potential and Berry field strength are both defined within the branch of and .

Figure 4 shows the results of corresponding to the energy bands in Fig. 2. As can be seen, at the critical strain of , which is consistent with the condition of for band inversion, the invariant jumps from to . This, demonstrates a topological phase transition in the electronic properties of phosphorene.

According to Eq. (17), another way to observe a topological phase transition in phosphorene, is by applying in-plane compressive biaxial strain at a fixed value of tensile strain in the direction. Figs. 5 show the numerically computed phase diagrams as a function of and at a fixed value of . As can be seen, there is a linear border between two distinct topological phases that corresponds to the regimes before and after the gap closing condition of , where is a fixed value of strain in the direction of .

It is worth mentioning that, the relatively large bulk band gap of monolayer phosphorene necessitates a rather large value of strain in order to observe band inversion. As mentioned before, according to DFT calculations, this is accompanied by an upward shift of a new VBM. After a critical percentage of strain, a direct band touching occurs, which is characterized by a TI phase. However, further increase of strain leads to a metal phase and because the topological nature does not change, the system may fall into the TM phase. Our model can not predict the VBM upward shift, hence, in spite of demonstrating the change of the topological phase, it can not distinguish between the TI and TM phases.

Note that our approach can be simply extended to the case of few-layers phosphorene in which we expect to observe the topological phase transition at lower strain values, due to the fact that the inter-layers hoppings result in a smaller gap Cai et al. (2014).

### iv.2 Electronic properties of phosphorene nanoribbons under strain

In this subsection, we investigate the evolution of the band structure of phosphorene nanoribbons in the presence of in-plane and perpendicular strain. In the following, we refer to the width of zPNRs as -zPNR with being the number of zigzag chains across the ribbon width. As we showed in the previous section, a topological phase transition occurs in the band spectrum of phosphorene. This should lead to the formation of topologically protected edge states in the band structure of the corresponding nanoribbons. We obtain the eigenvalues and eigenvectors using the following matrix

(23) |

where are the 1D Bloch wave functions. , denote super-cells; , are the basis sites in a super-cell and , denote the spin degree of freedom. is the wave vector, and represents a Bravais lattice vector. are the hopping integrals with usual SOC or intrinsic Rashba coupling that are conveniently defined between the basis site with spin of super cell and the basis site with spin of unit cell .

Note that, Eq. (23) is related to the energy spectrum of nanoribbons that are not edge passivated. The experimental realization of such nanoribbons with pristine edges in low dimensional materials as graphene is well known Zhang et al. (2013) and may be extended to the case of phosphorene nanoribbons. However, the stability of such ribbons is important from the experimental point of view. Formation energy studies Carvalho et al. (2014) showed that pristine phosphorene nanoribbons are stable specially for ribbon widths which we have considered in this paper.

The emergence of quasi-flat bands which are detached completely from the bulk bands due to the special structure of phosphorene are well known Ezawa (2014); Sisakht et al. (2015); Grujić et al. (2016). As shown in Fig. 6(a), there are topologically non-protected edge modes in the 1D bands of a typical zPNR (the results are for ). These quasi-flat bands have been used to propose a field-effect transistor driven by an in-plane electric field Ezawa (2014); Sisakht et al. (2015). However, since pristine bulk phosphorene is a trivial insulator, the existence of topologically non-protected edge modes in the corresponding nanoribbons which can be affected by environmental conditions such as disorder or impurities, may not be a good candidate for practical use. As an example, we consider the zigzag nanoribbon in the presence of perpendicular strain. The behaviour in the presence of other types of strain is similar to this case. As can be seen in Figs. 6(b) and (c), by increasing strain the bulk gap of the nanoribbon gradually decreases and after a critical strain, where a band inversion occurs in the bulk spectrum, the corresponding edge states in the ribbon cross the gap which demonstrates a topological insulator phase. Owing to the dependence of the nanoribbon gap on the ribbon width, the critical strain for driving it to a topological insulator phase depends on the width as well. If we consider ribbons with very large widths, the critical value approaches the critical strain value of bulk 11.8% that we have calculated in previous section.

The anisotropic structure of phosphorene results in a large bulk gap for zigzag nanoribbons with experimentally accessible widths. This makes strained zPNRs ideal systems for observing topological states even at room temperature. As shown in Fig. 6(c) for a zigzag nanoribbon of width nm this gap is about 200 meV which is much larger than room temperature thermal energy. We have calculated numerically these bulk gaps for relatively large ribbons up to a width of nm and found that the mentioned gaps are at least three orders of magnitude larger than the thermal energy at room temperature.It is worth mentioning that, such a typical ribbon width is wide enough to prevent from overlapping of edge states living on opposite sides of the ribbon. The corresponding amplitude probability of the topological edge modes of Fig. 6(c) which have amplitude on opposite edges are shown in Fig. 6(d) for a definite point. The amplitude of the wave functions drop very quickly along the width of the ribbon demonstrating that the nanoribbon width is wide enough to prevent quantum tunneling. Such excellent properties can pave the way for utilizing it in device applications.

## V Conclusions

In summary, we derived a spin-orbit model Hamiltonian based on the structural and electronic properties of phosphorene that captures the main physical properties of spin-orbit related subjects. Then we showed in the frame of this TB model that gap engineering of phosphorene by axial strains can lead to a topological phase transition in the electronic properties of phosphorene. In spite of the relatively small gap induced by SOC in bulk monolayer phosphorene, we predict that due to the special puckered structure of phosphorene, zigzag nanoribbons in the regime of TI have topologically protected edge states with rather large bulk band gaps of about meV for a typical ribbon of width nm. Such gaps are larger that the thermal energy at room temperature and are therefore sufficiently large for practical device engineering at room temperature.

## Acknowledgement

This work was supported by Iran’s ministry of science. M.Z. is a postdoc fellow of the Felamish Research Foundation (FWO-Vl).

## References

- Kane and Mele (2005) C. L. Kane and E. J. Mele, Phys. Rev. Lett. 95, 226801 (2005).
- Hasan and Kane (2010) M. Z. Hasan and C. L. Kane, Rev. Mod. Phys. 82, 3045 (2010).
- Moore (2010) J. E. Moore, Nature (London) 464, 194 (2010).
- Qi and Zhang (2011) X.-L. Qi and S.-C. Zhang, Rev. Mod. Phys. 83, 1057 (2011).
- Fu et al. (2007) L. Fu, C. L. Kane, and E. J. Mele, Rev. Mod. Phys. 98, 106803 (2007).
- Bernevig et al. (2006) B. A. Bernevig, T. L. Hughes, and S.-C. Zhang, Science 314, 1757 (2006).
- Ren et al. (2016) Y. Ren, Z. Qiao, and Q. Niu, Reports on Progress in Physics 79, 066501 (2016).
- Ma et al. (2014) Y. Ma, Y. Dai, W. Wei, B. Huang, and M.-H. Whangbo, Sci. Rep. 4, 7297 (2014).
- Ezawa (2012) M. Ezawa, Phys. Rev. Lett. 109, 055502 (2012).
- Koenig et al. (2014) S. P. Koenig, R. A. Doganov, H. Schmidt, A. C. Neto, and B. Oezyilmaz, Appl. Phys. Lett. 104, 103106 (2014).
- Li et al. (2014a) L. Li, Y. Yu, G. J. Ye, Q. Ge, X. Ou, H. Wu, D. Feng, X. H. Chen, and Y. Zhang, Nat. Nanotechnol. 9, 372 (2014a).
- Xia et al. (2014) F. Xia, H. Wang, and Y. Jia, Nat. Commun. 5, 4458 (2014).
- Fei et al. (2015) R. Fei, V. Tran, and L. Yang, Phys. Rev. B 91, 195319 (2015).
- Kim et al. (2015) J. Kim, S. S. Baik, S. H. Ryu, Y. Sohn, S. Park, B.-G. Park, J. Denlinger, Y. Yi, H. J. Choi, and K. S. Kim, Science 349, 723 (2015).
- Xiang et al. (2015) Z. J. Xiang, G. J. Ye, C. Shang, B. Lei, N. Z. Wang, K. S. Yang, D. Liu, F. B. Meng, X. G. Luo, L. J. Zou, Z. Sun, Y. Zhang, and X. H. Chen, Phys. Rev. Lett. 115, 186403 (2015).
- Liu et al. (2015) Q. Liu, X. Zhang, L. Abdalla, A. Fazzio, and A. Zunger, Nano Lett. 15, 1222 (2015).
- Zhang et al. (2015) T. Zhang, J.-H. Lin, Y.-M. Yu, X.-R. Chen, and W.-M. Liu, Sci. Rep. 5, 13927 (2015).
- Wei and Peng (2014) Q. Wei and X. Peng, Appl. Phys. Lett. 104, 251915 (2014).
- Peng et al. (2014) X. Peng, Q. Wei, and A. Copple, Phys. Rev. B 90, 085402 (2014).
- Rudenko and Katsnelson (2014) A. N. Rudenko and M. I. Katsnelson, Phys. Rev. B 89, 201408 (2014).
- Liang et al. (2014) L. Liang, J. Wang, W. Lin, B. G. Sumpter, V. Meunier, and M. Pan, Nano Lett. 14, 6400 (2014).
- Tran et al. (2014) V. Tran, R. Soklaski, Y. Liang, and L. Yang, Phys. Rev. B 89, 235319 (2014).
- Qiao et al. (2014) J. Qiao, X. Kong, Z.-X. Hu, F. Yang, and W. Ji, Nat. Commun. 5, 4475 (2014).
- Takao and Morita (1981) Y. Takao and A. Morita, Physica B+C. 105, 93 (1981).
- Fei and Yang (2014) R. Fei and L. Yang, Nano Lett. 14, 2884 (2014).
- Li et al. (2014b) Y. Li, S. Yang, and J. Li, J. Phys. Chem. C 118, 23970 (2014b).
- Sisakht et al. (2015) E. T. Sisakht, M. H. Zare, and F. Fazileh, Phys. Rev. B 91, 085409 (2015).
- Liu et al. (2011) C.-C. Liu, H. Jiang, and Y. Yao, Phys. Rev. B 84, 195430 (2011).
- Castellanos-Gomez et al. (2014) A. Castellanos-Gomez, L. Vicarelli, E. Prada, J. O. Island, K. Narasimha-Acharya, S. I. Blanter, D. J. Groenendijk, M. Buscema, G. A. Steele, J. Alvarez, H. W. Zandbergen, J. J. Palacios, and H. S. J. van der Zant, 2D Materials 1, 025001 (2014).
- Rodin et al. (2014) A. S. Rodin, A. Carvalho, and A. H. Castro Neto, Phys. Rev. Lett. 112, 176801 (2014).
- Wang et al. (2015) C. Wang, Q. Xia, Y. Nie, and G. Guo, Appl. Phys. Lett. 117, 124302 (2015).
- Huang and Xing (2014) G. Huang and Z. Xing, arXiv:1409.7284 (2014).
- Jiang and Park (2015) J.-W. Jiang and H. S. Park, Phys. Rev. B 91, 235118 (2015).
- Mohammadi and Nia (2016) Y. Mohammadi and B. A. Nia, Superlatt. and Microstruct. 89, 204 (2016).
- Duan et al. (2016) H. Duan, M. Yang, and R. Wang, Physica E 81, 177 (2016).
- Sa et al. (2014) B. Sa, Y.-L. Li, J. Qi, R. Ahuja, and Z. Sun, J. Phys. Chem. C 118, 26560 (2014).
- Harrison (1999) W. A. Harrison, Elementary electronic structure (World Scientific, Singapore 1999).
- Tang et al. (2009) H. Tang, J.-W. Jiang, B.-S. Wang, and Z.-B. Su, Solid State Commun. 149, 82 (2009).
- Ezawa (2014) M. Ezawa, New. J. Phys. 16, 115004 (2014).
- Soluyanov (2012) A. A. Soluyanov, Topological aspects of band theory, Ph.D. thesis, Rutgers University-Graduate School-New Brunswick (2012).
- Fu and Kane (2006) L. Fu and C. L. Kane, Phys. Rev. B 74, 195312 (2006).
- Soluyanov and Vanderbilt (2011) A. A. Soluyanov and D. Vanderbilt, Phys. Rev. B 83, 035108 (2011).
- Fukui et al. (2005) T. Fukui, Y. Hatsugai, and H. Suzuki, J. Phys. Soc. Jpn. 74, 1674 (2005).
- Fukui and Hatsugai (2007) T. Fukui and Y. Hatsugai, J. Phys. Soc. Jpn. 76, 053702 (2007).
- Cai et al. (2014) Y. Cai, G. Zhang, and Y.-W. Zhang, Sci. Rep. 4, 6677 (2014).
- Zhang et al. (2013) X. Zhang, J. Xin, and F. Ding, Nanoscale 5, 2556 (2013).
- Carvalho et al. (2014) A. Carvalho, A. S. Rodin, and A. H. Castro Neto, Europhys.Lett 108, 47005 (2014).
- Grujić et al. (2016) M. M. Grujić, M. Ezawa, M. Ž. Tadić, and F. M. Peeters, Phys. Rev. B 93, 245413 (2016).