Free Essay

Multi-valley monolayer materials

By

Fok Hong Ting

2010280804

Under the supervision of

Dr. Wang Yao

Content

Acknowledgement 4

Introduction 5 Motivation 5 Project outline 5

Chapter I Tight binding model 6

Basics 6

The secular equation 8

Conclusion 9

Chapter II Graphene 11

Formulation 11 π energy band of graphene 15 σ energy bands of graphene 17 Conclusion 18

Chapter III Silicene 19

Tight Binding Hamiltonian of silicene without SOC 20 Constructing orthogonal basis 21 1st order SOC inclusion 24

Conclusion 26

Chapter IV Edge state of Group IV elements 27

Formulation 27

Graphene edge 29

Silicene edge 32

Conclusion 42

References 43

Appendix 44

List of Figures 44

Acknowledgement

I hereby would like to express my appreciation and respect to my supervisor Dr. Wang Yao. Although I am not a talented student, Dr. Yao provided me timely support and insight in the field of physics. It is my fortune to take part in this final year project under his guidance.

Moreover, I would like to thank Dr. GuiBin Liu and Mr. We Yue for their support and comments.

Introduction

Motivation

One of the most intriguing phenomena in physics is the edge effect in 2-D systems. With the emergence of 2-D monolayer materials, the study of edge states in such material is of fundamental interest as well as practical interest. A well known example of such material is graphene, the discovery of which has lead to a Nobel Prize 2 years ago. It remains a mystery whether its counterpart, silicene, exhibit a similar property. This gives us great motivation to study the materials.

Project outline

In this project, we try to investigate the band properties of graphene and silicene with tight binding method. The principal of tight binding method is discussed in Chapter I.

In Chapter II, we will formulate a tight binding Hamiltonian for graphene and calculate its energy band structure.

In Chapter III, we will include SOC effect in the calculation of silicene energy bands and look for an analytical solution.

Finally, in Chapter IV, we focus on investigating the edge states of both graphene and silicene.

Chapter I Tight Binding Method

In this section we explain the tight binding method. In our model, we assume the material is a perfect lattice. It has no defects and can be produced by choosing a suitable unit cell.

Basics

Because of the translational symmetry of the unit cells in the direction of the lattice vectors, ai (i = 1,2,3,) , any wave function of the lattice should satisfy Bloch's theorem:

TaiΨ=eik∙aΨ ,(i=1,2,3) where Tai is a translational operation along the lattice vector ai, and k is the wave vector.

There are many possible functional forms of Ψ which satisfy this equation and the most commonly used form for Ψ is a linear combination of plane waves.

Plane waves are used because of the following reasons:

1. The integration of the plane wave wavefunction is easy and can be done analytically,

2. The numerical accuracy only depends on the number of the plane waves used.

However, the plane wave method also has limitations:

1. The scale of the computation is large,

2. It is difficult to relate the plane wave wavefunction to the atomic orbitals in the solid.

Another functional form which satisfies the above equation is based on the jth atomic orbital in a unit cell. This is a tight binding, Bloch function:

Фjk,r=1NRNeik∙Rρj(r-R), (j=1,2,…n)

Where R is the position of the atom and ρj denotes the atomic wavefunction in state j. The number of atomic wavefunctions in the unit cell is denoted by n, and we have n Bloch functions in the solid for a given k. The ρj in the N unit cells are weighted by the phase factor eik∙R and are then summed over the lattice vectors of the whole crystal.

The advantages of using orbitals in Bloch functions are as follows: 1. The number of basis functions can be small compared with the number of plane waves. 2. The formula for many physical properties can be easily derived.

The eigenfunctions in the solid Ψj are expressed by a linear combination of Bloch functions Фj' as follows:

Ψjk,r=j'=1nCjj'(k)Фj' (k, r) where Cjj' are coefficients to be determined.

Then, the j-th eigenvalue as a function of k is given by:

Ejk=ΨjHΨjΨjΨj=j,j'=1nCij*CijФjHФj'j,j'=1nCij*CijФjФj'=j,j'=1nCij*CijHjj'(k)j,j'=1nCij*CijSjj'(k)

Here Hjj'(k) is the transfer integral matrices and Sjj'(k) is the overlap integral matrices, respectively.

The secular equation

For a given k value, Cij is optimized so that Ej is minimized.

Therefore,

∂Eik∂Cij*=j'=1NHjj'kCij*Cij'j,j'=1NSjj'kCij*Cij'-j,j'=1NHjj'kCij*Cij'j,j'=1NSjj'kCij*Cij'2j'=1NSjj'kCij'=0

Defining a column vector: Ci=Ci1⋮CiN

Then, HCi=EikSCi

H-EikSCi=0

Note that if the matrix H-EikS-1exists, Ci = 0, meaning that there is no wavefunction. Hence, the eigenfunction exists only when H-EikS-1does not exist,

i.e. detH-EikS=0

This is the secular equation. After obtaining the secular equation, the only thing we need to solve are the matrix elements Hαβk and Sαβk.

Note that :

Hαβk=ФαHФβ=1NRaRbeik∙Rb –RaΨαr-RaHΨβr-Rb=1NRaRbeik∙Rb –RaΨαrHΨβr-Rb-RaBy translational symmetry=1NRaRjeik∙Rj ΨαrHΨβr-RjBy letting Rj=Rb-Ra=Rjeik∙Rj ΨαrHΨβr-Rj

Similarly, Sαβk=ФαФβ=Rjeik∙Rj ΨαrΨβr-Rj

Although we now have a simpler expression of Hαβk and Sαβk, they still require the aggregation of integrals with all Rj, which is troublesome. To make it more useful, we consider only the effect of the nearest neighbor. This decision is made because: 1. The number of integrals considered can largely decrease, making calculation feasible. 2. It is believed that the nearest neighbor term contributes a lot more than the higher order terms.

Lastly, for the integrals ΨαrHΨβr-Rj and ΨαrΨβr-Rj, it will be unrealistic for us to calculate them because their exact form may be more complex than the tight binding model itself. In order to overcome this problem, they will be replaced by parameters from estimation or experiments in practice.

Conclusion

In the tight binding method, the one-electron energy eigenvalues E;(l) are obtained by solving the secular equation

H-EikSCi=0

The eigenvalues Eik are a periodic function in the reciprocal lattice, which can be described within the first Brillouin zone. With the equation, we can obtain the dispersion relation easily through the following procedures: 1. Specify the unit cell and the unit vectors, the Brillouin zone and the reciprocal lattice vectors,

2. Calculate the transfer matrix (H) and the overlap matrix (S)

3.Plug in different k values and solve the secular equation to obtain the dispersion relation in the first brillouin zone

To simplify and facilitate the calculation, only the contribution from the nearest neighbor is considered. Finally, the usefulness of tight binding model will be demonstrated in the next few chapters.

Chapter II Graphene

Fomulation

Graphite is a three-dimensional (3D) layered hexagonal lattice of carbon atoms. A single layer of graphite which forms a two-dimensional (2-D) material is called graphene. | a.The unit cell of graphene. | b. Brillouin zone of graphene. |

Here is the primary unit cell (dotted rhombus) and the Brillouin zone (shaded hexagon) of graphene respectively, where a1 and a2, are the unit vectors in real space, and b1 and b2 are the reciprocal lattice vectors.

In the figure, the real space unit vectors a1and a2 of the hexagonal lattice are expressed as a1=32a,a2, a2=32a,-a2

where the lattice constant = a1=a2=2.46Å.

On the other hand, for the reciprocal space unit vectors b1 and b2 : b1=23π,22π, b2=23aπ,-2aπ

Correspondingly, the lattice constant in reciprocal space is 4π3a.

Trivially the direction of the unit vectors b1 and b2 of the reciprocal hexagonal lattice are rotated by 90o from the unit vectors a1 and a2 of the hexagonal lattice in real space. By selecting the first Brillouin zone as shown, the highest symmetry is obtained for the Brillouin zone of graphene. Here three high symmetry points, Γ, K and M, are defined as the center, the corner, and the center of the edge, respectively.

Note that in graphene, three orbitals 2s, 2px and 2py are hybridized as σ bonds in a sp2 configuration, while the remaining 2pz orbital, which is perpendicular to the horizontal graphene plane, overlaps side by side with that of other atoms and makes π covalent bonds.

The matrix element for the Bloch orbitals between the A and B atoms can be obtained by taking the components of 2px, and 2py in the directions parallel or perpendicular to the bond. Here, we will show how to rotate the 2px atomic orbital and how to obtain the direct overlapping and sideway overlapping components.

For instance, the wavefunction of 2px is decomposed into its σ (direct overlapping) and π (sideway overlapping) components along the bonding direction:

Where 12 comes from cos60o and 32 comes fromsin60o.

By rotating the 2px, and 2py, orbitals in the directions parallel and perpendicular to the desired bonds, the matrix elements appear in only 8 patterns as shown, where shaded and not-shaded regions denote positive and negative amplitudes of the wavefunctions:

The four cases from (1) to (4) correspond to non-vanishing matrix elements and the remaining four cases from (5) to (8) correspond to matrix elements which vanish due to symmetry. The corresponding parameters for both the Hamiltonian and the overlap matrix elements (Hspσ, Hppσ, etc) are shown in the appendix.

Here we demonstrate the calculation of matrix elements 2sAH2pxB and 2pxAH2pyB obtained by the methods mentioned before. Note that in our steps, all vanished integrals are omitted.

2sAH2pxB=-Hspσeikxa3+12Hspσe-ikxa23eikya2+12Hspσe-ikxa23e-ikya2

=Hspσ(-eikxa3+e-ikxa23coskya2)

2pxAH2pyB=34Hppσe-ikxa23eikya2+34Hppπe-ikxa23eikya2-34Hppσe-ikxa23eikya2-34Hppπe-ikxa23eikya2=34Hppσ+Hppπe-ikxa23eikya2-34Hppσ+Hppπe-ikxa23eikya2=3i2Hppσ+Hppπe-ikxa23sinkya2

The resulting matrix element is purely imaginary. Nevertheless, the calculated results for the energy eigenvalues give real numbers.

Recall from our previous analysis, when an orbital is completely vertical, its integral with other orbitals vanishes by symmetry. Since graphene is planar, the overlapping between 2pz and other orbtials (2s, 2px, 2py) will be zero. Therefore we have:

Hgraphene=Hsp2hybridized00H2pz

Where Hsp2hybridized is a 6×6 matrix formed by coupling of 2s, 2px, 2py orbitals, and H2pz is a 2×2 matrix formed by coupling of 2pz orbitals.

π energy band of graphene

Now we consider only 2pz energy bands for graphene, as it is found to have eigenenergy closest to the Fermi energy.

When we consider only nearest-neighbor interactions, then there is only an integration over a single atom in the diagonal elements:

HAA=HBB=ε2p

For off-diagonal matrix element HAB, we consider the three nearest-neighbor B atoms relative to an A atom, which are denoted by the vectors R1, R2 and R3. Then ,

HAB=Hppπeik∙R1+Hppπeik∙R2+Hppπeik∙R3=Hppπ(eikxa3+2e-ikxa23coskya2)=tf(k)

where Hppπ is given by Eq. (2.18)* and f ( k ) is a function of the sum of the phase factors of eik∙Rj (j=1,2,3).

Since f (k) is a complex function, and the Hamiltonian forms a Hermitian matrix, we write HBA=HAB+, in which + denotes the conjugate transpose.

For the overlap integral matrix S , the diagonal matrix elements are given by:

SAA= SBB = 1

And the off diagonal elements are:

SAB = Sppπ f(k) = SBA*

Thus,

H=ε2pHppπf(k)Hppπfk*ε2p, S=1Sppπf(k)Sppπfk*1

Solving the secular equation det(H - ES) = 0 , the eigenvalues Ek are obtained as a function wk, kxand ky:

Egraphenek=ε2p±Hppπwk1±Sppπwk

where the energy with + signs gives the bonding π energy band, and the energy with - signs gives the anti-bonding π* band, and the function wk is given by: wk=fk2=1+4cos3kxa2coskya2+4cos2(kya2) Now we set ε2p=0 ,using Hppπ = -3.033eV and Sppπ = 0.129 , Then we get the following energy band:

The energy dispersion relation for graphene in the First Brillouin zone.

The upper half of the energy dispersion band is the π* energy anti-bonding band, and the lower half curve is the π energy bonding band. We can see that the π* band and the π band are degenerate at K point.

When the overlap integral becomes zero, the π and π * bands become symmetric around ε2p:

Egraphene(kx,ky)=±Hppπ1+4cos3kxa2coskya2+4cos2kya2

In this case, the energies have the values of ±3Hppπ, ± Hppπ and 0, respectively, at Γ, M and K point in the Brillouin zone.

σ energy bands of graphene For the σ bands of graphene, there are three orbitals of sp2 covalent bonding per carbon atom: 2s, 2px and 2py. We have six Bloch orbitals in one unit cell and thus six σ bands. For the eigenvalues obtained, three of the σ bands are bonding σ bands below the Fermi energy, and the other three σ bands are anti-bonding σ* bands above the Fermi energy.

The calculation of Hamitonian elements follows the rules we mentioned before. | The energy dispersion relation for both σ and π energy bands of graphene. |

From the figure, we can see that the π and the some of the σ bands cross each other. However, because of the different group theoretical symmetries between σ and π bands, there is no band separation at the crossing points. [2]

The relative positions of these crossings are known to be important for:

1. photo-transitions from σ to π bands and from π to σ bands

2. charge transfer from alkali metal ions to graphene sheets in graphite intercalation compounds.[2]

Conclusion This tight binding concept for graphene is useful. Since the properties of np and ns orbitals are similar, the method can be generalized to not only carbon nanotubes and fullerenes, but also to other Group IV materials. In the next section we will apply these techniques on silicene.

Chapter III Silicene

Silicene and other monolayer Group IV materials are regarded as a counterpart of graphene. They have similar lattice structure, with one atom surrounded by 3 other atoms and together form a honeycomb structure, but the latter exhibits a totally flat manner while the former one is lower-buckled. planar structure of graphene | lower buckled structure of silicene |

In graphene, due to symmetry, the SOC effect must be zero.

If we assume a force perpendicular to the bonding direction exists, since graphene is completely planar, one can argue the force also exists pointing to the opposite direction. Therefore, the perpendicular force due to SOC effect must be zero.

For silicene, it turns out that the effect of SOC due to the nearest neighbor also vanishes by the same reason. Therefore, only the SOC effect due to self coupling will be included in the following analysis. For higher order terms inclusion such as the intrinsic Rahsba SOC, it may be referred to reference [1].

Tight Binding Hamiltonian of silicene without SOC

In this section we will construct the Hamiltonian of silicene with tight binding model approach. It is known that the silicene bonds are contributed from 4 occupied orbitals with the highest energy level, i.e. 3s, 3px, 3py and 3pz orbitals. In our model we will use this non-orthonormal basis {pZA,pZB,pyA,pxA,sA,pyB,pxB,sB }, where A, B corresponds to the 2 different atoms in a unit cell respectively.

H0= HπHnHn+Hσ Where Hπ= ε3p00ε3p is contributed by atoms’ 3pz orbitals Hσ= ETT+E is contributed by atoms’ 3px, 3py, and 3s orbitals Hn= 000V3'iV3'0V3'-iV3'0000 is contributed by the coupling of H and H T=-V1'-iV1'V2'-iV1'V1'iV2'V2'iV2'0 T=HpyApyB HpxApyB HsApyB HpyApxB HpxApxB HsApyB HpyAsB HpxAsB HsApyB V1'=34sin2θ (Hppπ-Hppσ) V2'=32sinθ (Hspσ) V3'=32sinθcosθ(Hppπ-Hppσ)

By setting ε3p = 0, we have: H0=00000V3'-iV3'000V3'iV3'00000V3'000-V1'-iV1'V2'0-iV3'000-iV1'V1'-iV2'0000Δ-V2'iV2'0V3'0-V1'iV1'-V2'000iV3'0iV1'V1'-iV2'00000V2'iV2'000Δ Constructing orthogonal basis In order to get a diagonal matrix from Ho, we define a unitary matrix U1: U1 = 100000000001000000000-i√2 i2 -i2 00000-1√2 -12 12 0100000000-i√2 000-i2 -i2 001√2 000-12 -12 00001000 Then we can have H1=U1+H0U1 : H1= 00-i2V3'000000Δi2V2'00000i2V3'-i2V2'00000000000-i2V3'000000Δ-i2V2'00000i2V3'i2V2'0000000002V1'00000000-2V1' Under the transformation, we have new orthonormal basis: {pZA,sA,℘2B,pZB,sB,℘1A,℘3,℘4 } Where: ℘1A=-12pxA+ipyA=p+A ℘2B=12pxB-ipyB=p-B ℘3=12-12pxA-ipyA-12pxB+ipyB ℘4=1212pxA-ipyA-12pxB+ipyB Under the new basis, we can see H1 can be represented by 3 block matrix: H1= HA000HB000HC Where

HA= 00-i2V3'0Δi2V2'i2V3'-i2V2'0 HB= 00-i2V3'0Δ-i2V2'i2V3'i2V2'0 HC= 2V1'00-2V1' For HA, the characteristic equation is: λ3- ∆∙λ2-2λV2'2+V3'2+2∆∙V3'2=0 Then the 3 roots are: ε1=∆3-13312[-2∆3-18∆V2'2+V3'2+54∆∙V3'2+-2∆3-18∆V2'2+V3'2+54∆∙V3'22-4∆2+6λV2'2+V3'23-13312[-2∆3-18∆V2'2+V3'2+54∆∙V3'2--2∆3-18∆V2'2+V3'2+54∆∙V3'22-4∆2+6λV2'2+V3'23 ε2=∆3+1+i36312[-2∆3-18∆V2'2+V3'2+54∆∙V3'2+-2∆3-18∆V2'2+V3'2+54∆∙V3'22-4∆2+6λV2'2+V3'23+1-i36312[-2∆3-18∆V2'2+V3'2+54∆∙V3'2--2∆3-18∆V2'2+V3'2+54∆∙V3'22-4∆2+6λV2'2+V3'23 ε3=∆3+1-i36312[-2∆3-18∆V2'2+V3'2+54∆∙V3'2+-2∆3-18∆V2'2+V3'2+54∆∙V3'22-4∆2+6λV2'2+V3'23+1+i36312[-2∆3-18∆V2'2+V3'2+54∆∙V3'2--2∆3-18∆V2'2+V3'2+54∆∙V3'22-4∆2+6λV2'2+V3'23 And the corresponding eigenvector is 1αiV2'εiαiΔ-εiV3'iεi2αiV3' , where the normalization constant is αi=1+V2'εiΔ-εiV3'2+εi2V3'2 Define the unitary matrix UA=1α1V2'εiα1Δ-εiV3'iεi2α1V3'1α2V2'εiα2Δ-εiV3'iεi2α2V3'1α3V2'εiα3Δ-εiV3'iεi2α3V3' For UB, since we observe that HB can be obtained by replacing V2’ by -V2’ in HA, we define UB=1α1-V2'ε1α1Δ-ε1V3'iε12α1V3'1α2-V2'ε2α2Δ-ε2V3'iε22α2V3'1α3-V2'ε3α3Δ-ε3V3'iε32α3V3' Therefore, we construct a unitary transformation matrix U2 by aggregating UA, UB, and a 2x2 identity matrix as diagonal blocks and rearranging the columns and rows. The overall transformation matrix is obtained as =U1U2 . H'0= ε100000000ε100000000ε200000000ε200000000ε300000000ε3000000002V1'00000000-2V1' And the final basis are: {ϕ1,ϕ4,ϕ2,ϕ5,ϕ3,ϕ6,ϕ7,ϕ8 }={pZA,pZB,pyA,pxA,sA,pyB,pxB,sB }∙U For small θ, ε1≈∆V32V22, ε2≈∆+∆2+4V222,ε3≈∆-∆2+4V222. And it is found that the basis {ϕ1,ϕ4} gives the eigenvalue nearest to the Fermi energy.[2] By performing 1sr order k expansion around K points under the basis {ϕ1,ϕ4}, we have: HK=ε1I+0vF(kx+iky)vF(kx-iky)0 Where vF=3a2[u112Hppπsin2θ+Hppσcos2θ-u212Hssσ+2u11u21cosθHspσ-12u312sin2θHppσ-Hppπ][2] 1st order SOC inclusion Unlike graphene, SOC is in effect in silicone. In this section we will show how the effect of SOC can be amended into the model. Hso is calculated with basis:

{pZA,pZB,pyA,pxA,sA,pyB,pxB,sB }⊗{↑,↓} After calculating all elements, we obtain the following SOC value table for different coupling: | pzi | pzi | pzi | pzi | pzi | 0 | iσx | -iσy | 0 | pzi | -iσx | 0 | iσz | 0 | pzi | iσy | -iσz | 0 | 0 | pzi | 0 | 0 | 0 | 0 | Where i=A,B corresponds to a particular atom in a unit cell, σx, σy, σz are Pauli matrices and 0 is zero matrix. All cross atom SOC terms are zero. Under the basis ϕ1,ϕ4⊗↑,↓=ϕ1↑,ϕ1↓,ϕ4↑,ϕ4↓ The diagonalized Hamiltonian Hso'=U+HsoU is of the form: H'so=-λso1st0000λso1st0000λso1st0000-λso1st And it is found that λso1st=ξ02u312 [1] Therefore the whole Hamiltonian becomes: Htotal=-λso1st+ε10vF(kx+iky)00λso1st+ε10vF(kx+iky)vF(kx+iky)0λso1st+ε100vF(kx+iky)0-λso1st+ε1

Conclusion

From the analysis we can see that silicene is quite different from grpahene. Due to the low buckled structure of silicene, the basis corresponding to the eigenenergy near the Fermi energy is no longer only the pz orbital, but a hybrid of {pZ,py,px,s}. Nevertheless, with a small k expansion, we can still observe that the energy band has a linear dependence on k near K points. Also, after including SOC effect, energy shows spin splitting, which can be seen in the diagonal of the Hamiltonian under basisϕ1↑,ϕ1↓,ϕ4↑,ϕ4↓. For higher order treatment to the SOC effect, one can look into reference [1].

Chapter IV Edge state of Group IV elements Formulation The most common and simple edges are zigzag edge and armchair edge. Zigzag edge | Armchair edge | | |

Note that for the ribbons, although the unit cell (and hence the Hamiltonian) is dependent on the ribbon width, we can still break it down in to several cases.

For zigzag edge, consider the following unit cell inside the blue rectangle, with each atoms numbered from the top to the bottom:

The bonding between atom 1 and atom 2, is the same as the bonding between atom 3 and 4. Other bondings which are similar are painted in purple. Recall previous results, we know that off diagonal tight binding Hamiltonian elements depends on the bonding with the nearest neighbours. With identical bonding arrangement, the matrix elements would be the same.

From the graph, obviously there are 2 types of off diagonal elements. Denote it by H12 (corresponding to the purple bond), and H23 (corresponding to the green bond), we have:

H0=H11H1200⋯00H12*H11H2300000H23*H11H12⋱0⋮00H12*⋱⋱00⋮0⋱⋱H11H2300000H23*H11H1200⋯00H12*H11

Where H11 is the self-coupling matrix, and H12* and H23* can be obtained by the conjugate transpose of the upper triangular Hamiltonian.

In the case of armchair edge,

The unit cell contains more atoms, yielding more variations in terms of the matrix elements. After categorizing, we have: H12 (corresponding to the pink bond), H13 (corresponding to the purple bond) and H24 (corresponding to the green bond). Therefore we get:

H0=H11H12H130⋯00H12*H110H24000H13*0H11H12⋱0⋮0H24*H12*⋱⋱H130⋮0⋱⋱H110H24000H13*0H11H1200⋯0H24*H12*H11

Nevertheless, the above equations consider only orbitals as basis without taking SOC effect into account. For materials that exhibit significant SOC effect, modification will be required.

Finally, for simplicity, in the coming up calculations we use the approximation S=I, where S is the overlapping matrix and I is the identity matrix.

Graphene edge

As we concern on the energy near Fermi surface, only 2p orbital is included in calculation.

For armchair ribbon, we have:

Harmchair=ε2pⅇ-ⅈak3Hppπⅇⅈak23Hppπ0⋯00ⅇⅈak3Hppπε2p0ⅇ-ⅈak23Hppπ000ⅇ-ⅈak23Hppπ0ε2pⅇⅈak3Hppπ⋱0⋮0ⅇⅈak23Hppπⅇ-ⅈak3Hppπ⋱⋱ⅇⅈak23Hppπ0⋮0⋱⋱ε2p0ⅇ-ⅈak23Hppπ000ⅇ-ⅈak23Hppπ0ε2pⅇⅈak3Hppπ00⋯0ⅇⅈak23Hppπⅇ-ⅈak3Hppπε2p

Hence we have the following dispersion relation:

We can see that the energy gap is nearly zero at k=0. Furthermore, there is no band gap for every width that equals 3M-1 (with M an integer). It is consistent with the analytical result obtained in reference [3].

For zigzag ribbon, Hzigzag=ε2p-2Hppπcoska00⋯00-2Hppπcoskaε2pHppπ00000Hppπε2p-2Hppπcoska⋱0⋮00-2Hppπcoska⋱⋱00⋮0⋱⋱ε2pHppπ00000Hppπε2p-2Hppπcoska00⋯00-2Hppπcoskaε2p

And the following energy dispersion is obtained:

Here we observe a horizontal line at the 2 sides of the Brillounin zone. To further confirm electron are localized at the edge under these k values, we have also plotted the extension length:

For the extension length, we can see that it is almost zero at the boundaries of Brillounin zone, indicating that at that k values the edge contributes the most to the edge state band. In other words, the electron is localized at the edge. However, as k is closer to 0 the electron is more likely to appear within the ribbon, so proportion of inner atoms contributes more to the edge energy band.

Silicene edge Recall from previous result of silicene, we can decompose the total Hamiltonian H into tight binding Hamiltonian and SOC Hamiltonian. In this section we will use similar approach to calculate the edge state Hamiltonian for silicene. Now we let: HsoAA=0iσx-iσy0-iσx0iσz0iσy-iσz000000 Then the total SOC Hamiltonian Hso is simply: Hso=HsoAA000⋱000HsoAA By combining H=H0⊗1001+Hso for either edge, we have a 8n×8n matrix for a silicene ribbon with width n atoms. After calculations, we obtain the following energy dispersion relation for silicene edge state, with –πa≤k≤πa (i.e. in the Brillouin zone) and width 20: To find the energy band near the Fermi energy, we consider the following reasons: The energy dispersion relation for silicene. | The energy dispersion relation for both σ and π energy bands of graphene. |

Due to silicene’s half filling, there should be 4 eigenvalues below the Fermi energy, and 4 above it. It is very likely that the middle 2 bands are those in our interest. Comparing the silicene bulk band with graphene’s, we also see these bands are located in a similar manner like the π and π* band of grpahene.

To consolidate the argument, we also compare their edge state. Here we alter the original low-buckled model to a planar one and compare it with graphene’s edge state. original silicene model (θ=101.7o) | | silicene model with θ=90o | | graphing model | |

We can observe that the zigzag edge state for θ=90o is similar with that of graphene, which is consistent with our analysis before.

From the comparison, an interesting observation is that as θ increases, the bands cross at a narrower k separation. Also, the edge state is no more localized at the peripherals of the Brillouin zone. Rather, they emerge from the bulk band and form a distinctive energy state when k is near zero.

Now we check whether Kramers degeneracy[4] exists in silicene. As SOC is included in our analysis and the basis is doubled, the number of edge state should also be doubled compared to the pure tight binding model. Here we show plots of squared wavefunction vs x of 2 lower edge states (band A and B) for a ribbon with width 50. Band | k value | spin up contribution | spin down contribution | A | -4π3a | | | B | 4π3a | | | A | 0 | | | B | 0 | | |

From the graphs, considering the state with major contribution, we can see the squared wave function of a certain spin in one band at k is very similar to that of another spin in the other band at –k. Although the previous observation is not true for the squared wavefunction with minor contribution, we believe it is due to rounding errors as they are of the order ~10-4. Therefore, it is verified that Kramers degeneracy exists in silicene.

Now we add a staggered potential V to the previous zigzag edge with width 20. V (eV) | Dispersion relation | -1.2 | | -0.9 | | -0.6 | | -0.3 | | 0 | | 0.3 | | 0.6 | | 0.9 | | 1.2 | |

As the potential magnitude increases, the edge states which are originally very close to each other splits in the middle of the Brillouin zone. Furthermore, it is interesting to see that the seemingly overlapping bands at V=0eV, are actually 2 bands at V=12eV.

Now we add a B-field to the edge without applying staggered potential. Note that we found it is pretty meaningless to plot the dispersion relation and squared wavefunction as there is no observable difference between the graphs. Hence we choose to plot the energy difference of the lower edge state vs k in the Brillouin zone.

B-field (meV) | energy difference vs k | -2 | | -1 | | 0 | | 1 | | 2 | |

Here we observe that no splitting occurs at k=0. The splitting is also close to zero at k≈±9π60a. Finally, there is nearly no splitting for B=0.

For armchair edge, we obtain the dispersion relation below:

With similar argument, we believe the 2 bands inside the blue rectangle are near Fermi energy. original silicene model (θ=101.7o) | | silicene model with θ=90o | |

Here shows the armchair edge energy bands for original and planar model with width 32. The latter model obeys the band property of graphene armchair edge, which exhibits metallic property for width=3M-2[5]. However, the low-buckled geometry opens a gap of 0.1897eV at k=0.

Conclusion

Although silicene and graphene are often viewed as counterparts, their difference in structure has led to distinctive properties in silicene and graphene. We believe as a semi-conductor, silicene can play crucial role in making nano-devices.

References

1. L. Brey and H. A. Fertig, Physical Review B 73, 235411 (2006)

2. R. Saito, G. Dresselhaus, and M. S. Dresselhaus, Physical Properties of Carbon Nanotubes , Imperial College, London, 1998

3. Kyoko Nakada and Mitsutaka Fujita, Physical Review B 54.17954 (1996)

4. Y. Yafet, g Factors and Spin-Lattice Relaxation of Conduction Electrons, 1963

5. Wang Yao, Sheungyuan A. Yang, and Qian Niu, Physical Review Letters 102, 096801 (2009)

6. Cheng-Cheng Liu, Hua Jiang and Yugui Yao, 1108.2933v1 (2011)

7. M. Fujita, M. Yoshida, and K. Nakada, Fullerene Sci. Technol. 4, 565 , 1996

8. Iijima, Sumio, Nature (354), 56 (1991)

9. M. Fujita, K. Wakabayashi, K. Nakada and K. Kusakabe, JPSJ (65), 7 (1996)

Appendix

List of figures

Values for coupling parameters in Hamiltonian matrix Material | Hppσ | Hppπ | Hssσ | Hspπ | ∆ | ξ0 | Graphene | 5.0372 | -3.0332 | -6.7692 | -5.5802 | -8.8682 | N/A | Silicene | 4.4471 | -1.121 | -1.931 | 4.471 | -7.031 | 34×10-31 |

Values for coupling parameters in overlap integral matrix Material | Sppσ | Sppπ | Sssσ | Sspπ | Graphene | 0.1462 | 0.1292 | 0.2122 | 0.1022 |