Bandgap tuning and optical response of twodimensional SiC: A firstprinciples real space study of disordered 2D materials.
Abstract
We present a realspace formulation for calculating the electronic structure and optical conductivity of such random alloys based on the KuboGreenwood formalism interfaced with the augmented space recursion (ASR) [A. Mookerjee, J. Phys. C: Solid State Phys. 6, 1340 (1973)] formulated with the Tightbinding Linear Muffintin Orbitals (TBLMTO) basis with van LeeuwenBaerends corrected exchange (vLB) [Singh et al., Phys. Rev B 93, 085204, (2016)]. This approach has been used to quantitatively analyze the effect of chemical disorder on the configuration averaged electronic properties and optical response of 2D honeycomb siliphene SiC beyond the usual Diraccone approximation. We predicted the quantitative effect of disorder on both the electronicstructure and optical response over a wide energy range, and the results discussed in the light of the available experimental and other theoretical data. Our proposed formalism may open up a facile way for planned band gap engineering in optoelectronic applications.
I Introduction
The discovery of graphene and related two and quasitwo dimensional (2D) materials has accelerated interest towards the fabrication of nanometerscaled twodimensional materials. This is because of the high mobility of their electrons and potential use in both mesoscopic research, materials engineering and nanodevices.geim (); butl (); novo1 (); novo2 () However, the zero band gap of Graphene limits its applications in optoelectronic devices. Thereafter several attempts have been made to increase the band gap, e.g., using bilayer graphene with nonequivalent top and bottom layers,kats (); zhang () or by surface doping and electric gating effect.kaplan () Another effective technique is the hydrogenation of graphene. A maximum band gap of 0.77 eV has been reached with 54% of hydrogen coverage.balog (); elias () Giovannetti et al.giovan () have also reported substrateinduced band gap in graphene on hexagonal boron nitride. However a large band gap satisfying commercial requirements of lightemitting diodes (LEDs) or solar cells is still not available.
In search of new material, recently, focus has been changed from graphene to other graphenelike 2D materials, with the aim of overcoming the shortage of graphene and broadening its range of applications.10 (); 11 () Such 2D layered materials can also be produced using parent 3D bulk materials.12 () Among these silicene, MoS, germanine, etc. have been synthesized successfully and have been found to exhibit new physical properties.liu (); tabert (); qiu () In this work, we focus on 2Dlayered structures of Siliphene, which can be understood in terms of graphene (2DC) or silicene (2DSi) doped with ‘Si’ or ‘C’, respectively, which is another Group IV binary compound displaying interesting properties. Since wurtzite siliconcarbide structure is graphitic in nature and theoretical investigations indicate that a phase transformation from graphinic to graphitelike structure is possible,26 (); 27 (); 28 () and can surely be a good candidate with larger bandgap.
Earlier, firstprinciples calculations reported stable 2DSi buckled honeycomb structures and their charge carriers also behave like massless Dirac fermions and their and * bands cross linearly at the Fermi level. However, unlike 2DSi, the monolayers SiC have stable 2D planar honeycomb structures, which may be attributed to the strong bonding through the perpendicular orbitals,free (); wu (); hsueh () Previously some 2D siliconcarbon binary compounds including hexagonal rings have been proposed, such as hexagonalSiC,lin1 (); lin2 () tetragonalSiC,chao (), gSiC,zhou1 () PtSiC,zhou2 () SiC,ding () which shows interesting properties such as a large direct band gap, improved photoluminescence, and highpower conversion efficiency.free (); wu (); hsueh () Because the honeycomb structure is common to both C and Si, experimentally stable 2DSiC in honeycomb structure has been synthesized successfully. Recent progress in the fabrication of ultra thin layered 2DSiC down to a monolayer (0.51.5nm)hern (); stan (); stan1 (); lin makes it an emerging semiconductor attractive to both fundamental research and practical applications. So, in this work, we emphasize on understanding electronicstructure and optical properties of the semiconducting 2DSiC material both qualitative and quantitative aspect. The choice of the system is based on the existing experiments for direct comparison and validation of our theoretical approach.
We have arranged the paper as follows. In Sec. II, we describe the computational methodology. Results are discussed in Sec. III. When possible, comparisons to previous theoretical and experimental results are made. We summarized our results in Sec. IV.
Ii Computational details
We treat randomdisorder in alloys using augmented space recursion (ASR) technique.asr1 (); asr2 (); asr3 (); asr4 (); asr5 (); asr6 (); asr7 () The ASR allows us to go beyond the local mean field theories like the coherent potential approximation and describe extended disorder like shortranged clustering or ordering as well as longranged disorder like randomly placed StoneWales defect chains forming in a graphene background (left panelFig.1). The ASR can deal both chemical and local structural disorder.SM ()
The ASR method is a realspace recursion technique for treating random disorder of alloys. The recursionviswa (); term () needs a countable basis , which is provided by the tightbinding linearized muffintin orbital (TBLMTO) method.TBLMTO () Here, ‘n’ is the basis and represents atomic positions.heine (); hhk (); hay (); hay2 (); hn (); hn2 () Since the ASR deals with the configuration fluctuations in realspace, e.g. for binary randomness , thus Bloch’s theorem plays no role.
In this work, we use localdensity approximation (LDA) modified with van LeeuwenBaerends (vLB) VLB () correction to exchange part.singh2013 (); singh2016 () The correction give a substantial improvement in the band gaps over the LDA. singh2013 (); singh2016 () We use tightbinding linear muffintin (TBLMTO)vLB approach,TBLMTO (); VLB (); singh2013 (); singh2016 () and combine it with augmented space recursion method asr1 () to carry out configuration averaging in semiconducting random disordered systems. We use LDA correlation energy parameterized by van Barth and Hedin.vBH ()
All calculations were done selfconsistently and nonrelativistically for given experimental geometry till the “averaged relative error” between the converged final and the previous iteration charge density and energy is reached, i.e. 10 and 10 eV/atom respectively. We use the tetrahedron method for kspace integration, and Anderson mixing to facilitate convergence. In TBLMTOvLB core states are treated as atomiclike in a frozencore approximation. Energetically higherlying valence states are addressed in the selfconsistent calculations of the effective crystal potential, which is constructed by overlapping WignerSeitz spheres for each atom in the unit cell.
To assess the optical conductivity of disorder system, we use KuboGreenwood approach (KGA) kubo (); greenwood () with ASR technique.asr1 (); moshi (); singh2016 () The KGAASR technique is used for systematic study of band gap tuning and configuration averaged optical responses in random alloys, which can be a unique formalism for calculating optical response of semiconducting alloys, both bulk and finite size.
The detailed discussion of mathematical approach has been provided in the Appendix.
Iii Results and Discussion.
We performed firstprinciples calculations on SiC ( 0 0.6) and compared our predictions with some existing experimental and theoretical results.azadeh (); lin We chose homogeneously disordered binary 2Dsiliphene (SiC), specifically because of several advantages over graphene.nako (); vogt () One of them is a better tunability of the band gap, very important in designing optoelectronic nanodevices. There has been disagreement between different theoretical and experimental works on band gap engineering in honeycomb 2DSiliphene. A quick review indicates that the main cause for these disagreements, at least among the theoretical studies, has been the dependence of the results on the underlying models and approximations, e.g. the starting structure, the type of exchangecorrelation, the electronic structure method used and how disorder is taken into account, include effects of configuration fluctuations of the immediate neighborhood. These are the focus points of our investigation, as described earlier. In all calculations, we consider the homogenous disorder in SiC.ding (); shi2015 () In their recent work, Ding et al.ding () and Shi et al.shi2015 (), has shown the existence of homogenous disorder mostly for Siconcentration 50 at.% in 2D SiC.
Experimentally, graphene shows zero band gap at a lattice constant of = 2.46 , while siliphene has a 1.90 meV gap with = 3.86 . NJSBKMG2005 () The energy difference of 10.66 eV ( EE) for graphene, and 5.66 eV (EE) for silicene indicate predominantly hybridization due to presence of graphene valence cloud. Thus graphene is more stable as a flat surface, while the valence cloud in siliphene has the possibility of hybridization leading to local corrugations and buckling. In our study we have not included buckling.
The experiments of Lin et al.lin led to lattice constants ranging from 2.4 to 2.8 of 5050 2DSiliphene using high resolution transmission electron microscope. GIXRD analysis confirmed = 2.66 . The two earlier theoretical approaches were by Bekaroglu et al.bekaroglu () and the Azadeh et al.azadeh () Both used space DFT approaches.VASP (); SIESTA (); WIEN2k () This restricted their study to ordered stoichiometric compositions in strict accordance with Bloch’s theorem, ignoring the effects of disorder. For 2DSiC (x=0.5), Bekaroglu et al. bekaroglu () used generalized gradient approximation (GGA) within projector augmented wave approach PAW () and proposed theoretical lattice constant 3.094 nearly equal to that for Wurtzite SiC ( = 3.08 ). On the other hand, Azadeh et al.azadeh () introduce an interesting incrystalsite doping or combination stoichiometric method to generate some ‘semirandom’ compositions and studied them using the WIEN2K code.WIEN2k ()
In Table 1, we compare the results calculated from vLBLMTO with other theoretical methods and measurements. The energy minimized lattice constant of 2DSiC using vLBLMTO is in good agreement with the measured a=2.60.2 . However, we found good agreement between our calculations and the experiments but theoretical results differ almost by 15%. The main source of disagreement of other theoretical approaches with vLBLMTO arise from the use of exchangecorrelation potential. We use vLB+LDA exchangecorrelation,singh2016 () which betters both r0 zero and asymptotic limits. Also, vLB cancels selfinteraction better than other density based functionals.singh2016 ()
In Fig. 2(a)&(c), we plot Nthorder muffintin orbital (NMTO)PRB2K () method calculated (a) bandstructure, and (c) density of states. The maximum DOS contribution comes from the orbital with indirect (direct) band gap of 3.59 (3.90) eV along KM (KK) symmetry direction. However, the vLBcorrected LDAsingh2016 () provides an accurate band structure and density of states, shown in Fig.2 (b)&(d), calculated at equilibrium lattice constant of 2.66 . The vLB+LDA provides indirect band gap of 2.96 eV along K symmetry direction of the Brilluoin zone. Siliphene (SiC, x=0.5), theoretical equivalent of graphene with silicon doping, is also supposed to be gap less. However, the honeycomb lattice of randomly but homogeneously distributed silicon atoms in symmetric semimetallic graphene opens up a gap, and shows semiconducting behavior. This way, the silicon doping allows us to tune the band gap in the SiC.
We show, in Fig. 3(a)&(b), the band gap and TDOS (total density of states) variation with changing lattice constant for given %Si ( see Fig. 4 in appendix, composition dependent lattice parameters). The band gap clearly increases, as shown in Table 2 and Fig. 3(a), with increasing %Si and reaches maximum at 50%, and further increasing %Si reduces the band gap. The large atomic size of Si enhances the interaction strength, and also hybridization between the p orbital for both Si and C. While the band remains unchanged contrary to band.
Bekaroglu et al.bekaroglu () (x=0.5)  Azadeh et al.(x=0.5)  

PAWLDA  PAWGGA  USLDA  USGGA  LDAGW  VASP/SIESTA  
(Å)  E (eV)  (Å)  E (eV)  (Å)  E (eV)  (Å)  E (eV)  (Å)  E (eV)  (Å)  E (eV) 
3.07  2.51  3.09  2.53  3.05  2.53  3.08  2.54  –  3.90  2.41  2.06 
Lin et al.(x=0.5) Experiment  Our Work (x=0.5)  Our Work (x=0.5)  
TEM  GIXRD  LDALMTO  vLBLMTO  vLBASR  NMTO  
(Å)  E (eV)  (Å)  E (eV)  (Å)  E (eV)  (Å)  E (eV)  (Å)  E (eV)  (Å)  E (eV) 
2.600.2  2.62.9  2.66  2.95  2.66  2.96  2.66  3.49  2.75  3.60  2.66  3.90 
Band Gap (eV)  

%x  Azadeh et al.  This work 
(ordered)  (disordered)  
0.10  –  2.094 
0.17  0.13  2.383 
0.20  –  2.662 
0.25  0.838  3.008 
0.33  1.237  4.387 
0.40  –  5.628 
0.50  2.061  5.835 
0.55  –  5.719 
Repulsive effect between the valence band and the conduction band increases with increasing Si content and reduced CSi bond length. The small bond lengths enhances inplane repulsion, which helps opening up the band gap. The antibonding states moves towards a higher energy level due to increased kinetic energy and repulsion effect of electrons with increasing %Si. Also, delocalization effect due to orbital overlapping between C and Si plays important role in band gap opening with increasing disorder strength. The reduced bondlength enhances delocalization, which broadens both the valence and conduction bands and decreases the band gaps. Assuming silicene as the end point of SiC) with x=1.0, one would expect the downturn of the mentioned trend beyond 50% Si. For graphene the effective mass at Dirac point is at zero, while in other graphenelike materials with finite band gap has larger effective mass and reduced mobility. This establishes that increasing strength of chemical disorder increases the band gap up to certain %Si, in agreement with the previous reports.azadeh ()
Manybody interactions are more significant in 2D materials than in their bulk counterparts.kin1 (); kin2 (); kin3 () This reflects the intrinsic enhancement of the importance of Coulomb interactions in 2D materials and their reduced screening. We note that the theoretical exciton binding energy in bulk 2HSiC is only 0.1 eV, whereas it is 1.17 eV for 2D honeycomb SiC.hsueh () This shows that the reduced dimensionality of a SiC sheet confines the quasiparticles. This significantly enhances the overlap between the electron and the hole wave functions and hence the electronhole (eh) interaction. The presence of the vacuum region reduces the screening and hence provides an extra contribution to the large excitonic effect in the SiC sheet in addition to this quantum confinement,. The role of eh interactions is relevant for the optical response. The repulsive ee interactions shifts the transition peak upwards in energy, and the attractive eh interactions shifts it downwards.AKK2008 (); KOER2010 (); KTB2011 () The optical dielectric function of the 2DSiC sheet is highly anisotropic. Therefore, manybody interaction effects in the band associated with the p orbitals, which extend into the vacuum region from the sheet, would be less screened. Consequently, the band would have a larger quasiparticle corrections. On the other hand, because the inplane bonds are mainly confined to the 2DSiC sheet, screening effects on the bands are more significant, and hence quasiparticle corrections are smaller.
The real and imaginary part of optical conductivity at different Si content are shown in Fig. 3 (c) and (d). We found an optical conductivity peak at 3.96 eV in the lowenergy region, however, at higher energies another peak is found at 6.63 eV. Clearly, increasing %Si increases the band gap (till 50%Si), and optical conductivity peaks exhibit noticeable broadening. The small shift in 55%Si peak (goes below 50%Si) in Fig. 3(b), which is in agreement with band gap change above 50%Si. The optical transition peak implicitly depends on the band gap, and on the factor (E  E + ). The transition peak is shadowed with increasing disorder strength (x) in 2DSiC. The reported transition peak by Hsueh et al.hsueh () from GWBetheSalpeter equation for single 2DSiC sheet have similar magnitude, see Table. 3.
In the lowenergy range (25 eV), the interband optical transitions involve mainly the bands, however, at higher energies, the optical absorption peaks between 511 eV are associated with the interband transitions involving bands. The first prominent peak, located at 3.96 eV arises due to the excitation between the and the states and the pronounced optical peak peak at 6.63 eV is mainly due to the excitation between the and the states. The optical conductivity peak at 3.96 eV is in good agreement with experimental peak at 3.33 eV as reported by Lin et al., lin see Table. 3. The small deviation from their data arises due to multilayers structures. They did their experiment in 2DSiC nanosheets of thickness 10 nm, and the interactions between substrates and samples was taken into consideration in the highresolution transmission electron microscopic measurements.
2D Systems  Theory  Experiment 
(eV)  (eV)  
Monolayer  3.96  – 
SiC  6.63  
[ vLBASR ]  
Single SiC sheet  3.25  – 
5.83  
[ GW+BSE hsueh ()]  
Single SiC sheet  4.42  – 
6.20  
[ GW+RPA hsueh ()]  
Ultra thin SiC nanosheets    3.33 
(10 nm)  lin 
The results shown in Table. 1 followed by subsequent discussions indicates that the calculated electronic structure depends upon the model chosen, the approximations introduced and the calculational methodology followed. The predictions of our proposed model are in good agreement with the experiments. The band gap and optical conductivity calculated at equilibrium lattice parameters are closest to the experiments. The effect of going beyond singlesite is visible in calculated physical properties from ASR used with vLB.
Iv Conclusion.
In summary, we perform DFT calculations with a LDAvLB exchangecorrelation function to calculate the optical conductivity of 2D material, real and imaginary parts, of an isolated singleatomthick layer consisting of a groupIV honeycomb crystal. The LDA and PBE functionals failed to capture important feature in the visible light region, which is important for nanoplasmonic applications. The approach, we introduce reduces the degree of underestimation of the transition energy. Here, we focused our studies on disordered 2DSiC, which shows good agreement with experimental observations. We also found that silicon can be used to tune the bandgap in SiC, which in turn enhances the optical conductivity. The prediction of large optical response in our calculations shows that 2DSiC can be a potential candidate for the solar cell material. The proposed formalism opens up a facile way to band gap engineer material for optoelectronic application.
V Acknowledgment
We are grateful to Prof. Andersen for giving us permission to take apart his LMTOpackage and incorporate our new exchange functional and then coupling it with the Cambridge Recursion Package of Haydock and Nex. BS and AM thank R. Haydock and C.M.M. Nex for permission to use and modify the Cambridge Recursion Package and Yoshiro Nohara, Max Plank Institute, Stuttgart, Germany, for his kind permission to use TBNMTO code developed by our group also. PS and DDJ support is from the U.S. Department of Energy (DOE), Office of Science, Basic Energy Sciences, Materials Science and Engineering Division; our research was performed at the Ames Laboratory, which is operated for the U.S. DOE by Iowa State University under Contract No. DEAC0207CH11358.
*
APPENDIX
In the Appendix we give a concise description of those points which we have introduced in this work. This would give the reader a clear picture of our contribution.
.1 The van LeewenBaerends Exchange.
To address the band gap problem, we have proposed the use of spinpolarized van LeeuwenBaerends (vLB) VLB () corrected exchange potential with a useful addition to satisfy the ionization potential theorem and to make the ionization energy and largest HOMOLUMO difference agree in firstprinciple calculations. The asymptotic behavior of exchange was matched at the atomic sphere boundary (or local interstitial). The vLB exchange has earlier given excellent results in other applications.singh2013 (); singh2016 () The vLB exchangecorrelation can be written as :
(1) 
where and ] are the standard LDA exchange and correlation potentials and is the correction to LDA exchange. The suffix represents the spin degree of freedom. Here, is
(2) 
where parameter was proposed originally by van Leeuwen and Baerends.VLB () Here, signifies the change in mean electronic distance provided density is slowly varying in given region and with strong dependence on gradient of local radius of the atomic sphere .
The effective KongSham potential becomes
where V is the external potential, V is the electronic Hartree potential, V is LDA exchange, V is LDA correlation, and V is the added spinpolarized vLBexchange. The iterative KohnSham scheme now has the effective potential constructed using new electronic density term
(3) 
We solve the KohnSham equation using the tightbinding linear muffintin orbital method with atomic sphere approximation (TBLMTOASA) TBLMTO () to obtain the exchangecorrelation potential in the atomicsphere. The solution is obtained iteratively to selfconsistency.
.2 Lattice constant and Bandgap variation with composition in SiC:
We treat disorder using augmented space formalism, which includes the effect of nearneighbor environment (beyond singlesite unlike CPA). asr1 (); asr2 (); asr3 (); asr4 (); asr5 (); asr6 (); asr7 () The equilibrium lattice constant is calculated using vLB corrected exchange correlation potential within ASR formalism at different Si%. The Vegard’s Law usually works better for elements with similar sizes, while ‘C’ (0.77) and ‘Si’ (1.15) are of different atomicsizes, both locally (in the intraatom sphybridization) and globally (in interatomic bonding) in the solid. This atomicsize difference leads to ‘bowing effect’ in Vegard’s Law for 2DSiC, see Fig. 4.
.3 Real Space and Recursion.
Keeping in mind our aim to develop a formalism without application of Bloch Theorem in any form, we have seamlessly combined the modified LMTO with the Cambridge Recursion Library and the Augmented Space Recursion packages. We propose this new TBLMTOvLBASR to be our technique to tackle disordered semiconductors.
Since we shall be working in real space in Recursion, we express the Kubo response of the valence electrons to optical disturbance in a form so that the necessary correlation functions can also be expressed through a generalized recursion technique.viswa ()
(4) 
where
where is the ground state. In lattices with cubic symmetry, . The FluctuationDissipation Theorem relates the imaginary part of the Laplace transform of the generalized susceptibility to the currentcurrent correlation function :
The correlation function can be obtained most directly using the generalization of the recursion method. This has been described in detail by Viswanath and Müller viswa (); term (). Analogous to the standard Recursion for single particle propagators, we generate a new basis from a three term recurrence :


We generate the rest of the basis recursively :

Mutual orthonormality leads to :
We can now expand any wave ket in terms of the orthogonal basis just generated :
(5) 
with =0 and . Taking Laplace transform we get,
(6) 
Recursive steps now give us :
So that :
(7) 
.4 Disorder and the Augmented Space approach.
We model local disorder with associating with the tightbinding Hamiltonian elements a set of random parameters . For binary alloying, for example, takes the values 0 and 1 with probabilities and .
(8)  
with the projection and transfer operators and . Following the ideas of measurement theory, the prescription is to replace the random parameters by the operator asr2 () whose spectrum are the possible results of measurement and whose spectral density is the probability of obtaining these results. For example, if takes the values 0 and 1 with probabilities and then is of rank 2 with a representation :
The basis used here is and
The full configuration space is and the augmented Hamiltonian in is :
(9)  
The Augmented Space Theorem asr2 () tells us that the configuration average of any function of is :
(10) 
With and
Now disorder can be combined with Recursion using the ideas of Augmented Space with all operators and states in It follows that :
We calculate this using Generalized Recursion in the full augmented space . The dielectric function and optical conductivity are related :
.5 The Optical Current.
The calculations begin with a thorough electronic structure calculation using the tightbinding linear muffintin orbitals method. The technique, like many others, is based on the density functional theory, where the energetics depends entirely upon the charge and spin densities. However, the current operator which is characteristic of an excited electron : excited optically, electronically or magnetically, is described by transition probabilities which depend upon the wavefunction. The wavefunctions are expressed as linear combinations of the linearized basis functions : the muffintin orbitals.
TBLMTOvLBASR calculations were done selfconsistently and nonrelativistically for a given geometry till the “averaged relative error” between the converged final and the previous iteration charge density and energy is reached. Nowhere do we use lattice translation symmetry so that Bloch’s theorem plays no role.
The calculations begin with a thorough electronic structure calculation using the tightbinding linear muffintin orbitals method. The technique, like many others, is based on the density functional theory, where the energetics depends entirely upon the charge and spin densities. However, the current operator which is characteristic of an excited electron : excited optically, electronically or magnetically, is described by transition probabilities which depend upon the wavefunction. The wavefunctions are expressed as linear combinations of the linearized basis functions of the LMTO. The wavefunction representation in the LMTO basis is :
where
The coefficients are available from our TBLMTO secular equation. To obtain the currents we follow the procedure of Hobbs et al. hobs (). The optical current operator is given by :
(15) 
where,
We have followed the prescription of Hobbs et al.hobs () to evaluate the integrals above. For details we again refer to reader to the above reference. We should note that we have so far not introduced the idea of reciprocal space. We expect our methodology to be applicable to situations w The ASR method is a realspace technique for treating random disorder of alloys. Recursion viswa (); term () needs a countable basis which is provided by the TBLMTOASAvLB basis with labelling the atomic positions .
References
 (1) A.K. Geim and K.S. Novoselov, Nat. Mater 6, 183 (2007).
 (2) S.Z. Butler et al., ACS Nano 7, 2898 (2013).
 (3) K.S. Novoselov, A.K. Geim, S.V. Morozov, Z. Jiang, Y. Zhang, S.V. Dubonos, I.V. Grigorieva, and A.A. Firsov, Science 306, 666 (2004).
 (4) K.S. Novoselov, A.K. Geim, S.V. Morozov, Z. Jiang, M.I. Katsnelson, I.V. Grigorieva, S.V. Dubonos, and A.A. Firsov, Nature 438, 197 (2005).
 (5) M.I. Katsnelson, K.S. Novoselov and A.K. Geim , Nat. Phys. 2, 620 (2006).
 (6) Y. Zhang, T.T. Tang, C. Girit, Z. Hao, M.C. Martin, A. Zettl, M.F. Crommie, Y.R. Shen, and F. Wang, Nature 459, 820 (2009)
 (7) D. Kaplan, V. Swaminathan, G. Recine, R. Balu, and S. Karna, J. Appl. Phys. 113, 183701 (2013).
 (8) D.C. Elias, R.R. Nair, T.M.G. Mohiuddin, S.V. Morozov, P. Blake, M.P. Halsall, A.C. Ferrai, D.W. Boukhvalov, M.I. Katsneslon, A.K. Geim, and K.S. Novoselov, Science 323, 610 (2009)
 (9) R. Balog, B. Jorgensen, L. Nilsson, M. Andersen, E. Rienks, M. Bianchi, M. Fanetti, E. Egsgaard, A. Baraldi, S. Lizzit,Z. Sljivancanin, F. Besenbacher, B. Hammer, T.G. Pedersen, P. Hofmann, and L. Hornekaer, Nat. Mater. 9, 315 (2010).
 (10) G. Giovannetti, P.A. Khomyakov, G. Brocks, P.J. Kelly, and J. van den Brink, Phys. Rev. B 76, 073103 ( 2007).
 (11) R. Ganatra, and Q. Zhang, ACS Nano 8, 4074 (2014).
 (12) M. Chhowalla, H.S. Shin, G. Eda, L.J. Li, K.P. Loh, H. Zhang, Nat. Chem. 5, 263 (2014).
 (13) A.K. Geim, and I.V. Grigorieva, Nature, 499, 419 (2013).
 (14) H. Liu, J. Gao and J. Zhao, ACS Nano 7, 10353 (2013).
 (15) C.J. Tabert, and E.J. Nicol, Phys. Rev. Lett. 110, 197402 (2013).
 (16) D.Y. Qiu, F.H. da Jornada, and S.G. Louie, Phys. Rev. Lett. 111, 216805 (2013).
 (17) E. Hernndez, C. Goze, P. Bernier, and A. Rubio, Phys. Rev. Lett. 80, 4502 (1998).
 (18) S.M. Lee, Y.H. Lee, Y.G. Hwang, J. Elsner, D. Porezag, and T. Frauenheim, Phys. Rev. B 60, 7788 (1999).
 (19) M. Zhao, Y. Xia, D. Zhang, and L. Mei, Phys. Rev. B 68, 235415 (2003).
 (20) C.L. Freeman, F. Claeyssens, N.L. Allan, and J.H. Harding, Phys. Rev. Lett. 96, 066102 (2006).
 (21) I. J. Wu and G. Y. Guo, Phys. Rev. B 76, 035343 (2007).
 (22) H.C. Hsueh, G.Y. Guo and Steven G. Louie, Phys. Rev. B 84, 085404 (2011).
 (23) S. Lin, J. Phys. Chem. C 116, 3951 (2012).
 (24) X. Lin, S. Lin, Y. Xu, A.A. Hakro, T. Hasan, B. Zhang, B. Yu, J. Luo, E. Li, and H. Chen, J. Mater. Chem. C 1, 2131 (2013).
 (25) C. Yang, Y. Xie, L.M. Liu, and Y. Chen, Phys. Chem. Chem. Phys. 17, 11211 (2015).
 (26) L. Zhou, Y. Zhang and L. Wu, Nano Lett. 13, 5431 (2013).
 (27) Y. Li, F. Li, Z. Zhou and Z. Chen, J. Am. Chem. Soc. 133, 900 (2010).
 (28) Y. Ding and Y. Wang, J. Phys. Chem. C 118, 4509 (2014).
 (29) Y. Hernandez et al., Nat. Nanotechnol. 3, 563 (2008).
 (30) S. Stankovich, D.A. Dikin, G.H.B. Dommett, K.M. Kohlhaas, E.J. Zimney, E.A. Stach, R.D. Piner, S.T. Nguyen, and R.S. Ruoff, Nature 442, 282 (2006).
 (31) X. Li, X.R. Wang, L. Zhang, S.W. Lee, and H.J. Dai, Science 319, 1229 (2008).
 (32) A. Mookerjee, in: A. Mookerjee, D.D. Sarma (Eds.), Electronic Structure of Clusters, Surfaces and Disordered Solids, Taylor Francis, (2003).
 (33) A Mookerjee J. Phys. C: Solid State Phys 6 L205 (1973).
 (34) A. Mookerjee J. Phys. C 6, 1340 (1973).
 (35) A. Mookerjee J. Phys. C 8, 1524 (1974).
 (36) A. Mookerjee J. Phys. C 8, 2688 (1976).
 (37) A. Mookerjee J. Phys. C 9, 1225 (1976).
 (38) S.Chowdhury, D.Jana, B.Sadhukhan, D.Nafday, S.Baidya, T.SahaDasgupta, A.Mookerjee, Indian J. Phys 90, 649 (2016).
 (39) T. Saha Dasgupta and A. Mookerjee, J. Phys Condens Matter 8, 2915 (1995)
 (40) V. Viswanath and G. Müller, sl The Recursion Method, Applications to ManyBody Dynamics Germany : SpringerVerlag (1994)
 (41) V. S. Viswanath, and G Müller, J. Appl. Phys. 67, 5486 (1990).
 (42) O. Jepsen and O.K. Andersen, The Stuttgart TBLMTOASA program, version 4.7, MaxPlanckInstitut für Festkörperforschung, Stuttgart, Germany (2000).
 (43) V. Heine, in Solid State Physics (Academic Press) Vol 35 1 (1980)
 (44) R. Haydock, V. Heine and M.J. Kelly, J. Phys. C:Solid State 5, 2845 (1972).
 (45) R. Haydock, in Solid State Physics (Academic Press) Vol 35 216 (1980).
 (46) R. Haydock, Philos. Mag. B43, 203 (1981).
 (47) R. Haydock and C.M.M. Nex, J. Phys. C:Solid State 17, 4783 (1984).
 (48) R. Haydock and C.M.M. Nex, J. Phys. C:Solid State 18, 2285 (1985).
 (49) R. van Leeuwen, and E.J. Baerends, Phys. Rev. A 49, 2421 (1994).
 (50) P. Singh, M. K. Harbola, B. Sanyal, and A. Mookerjee, Phys. Rev. B 87, 235110 (2013).
 (51) P. Singh, M. K. Harbola, M. Hemanadhan, A. Mookerjee, and D. D. Johnson, Phys. Rev. B 93, 085204 (2016).
 (52) U. von Barth and L. Hedin, J. Phys. C 5, 1629 (1972).
 (53) R. Kubo, J. Phys. Soc. Jpn. 12, 570 (1957).
 (54) D.A. Greenwood, Proc. Phys Soc 71, 585 (1958).
 (55) M. Rahaman, K. Tarafder, B. Sanyal, and A. Mookerjee, Physica B: Condensed Matter 406, 11 (2011).
 (56) M.S. Sharif Azadeh, A. Kokabi, M. Hosseini, and M. Fardmanesh, Micro & Nano Letters 6, 582 (2011).
 (57) H. Nakano, T. Mitsuoka, M. Harada, K. Horibuchi, H. Nozaki, N. Takahashi, T. Nonaka, Y. Seno, and H. Nakamura, Angew. Chem 118 6451 (2006).
 (58) P. Vogt, P. DePadova, C. Quaresima, J. Avila, E. Frantzeskakis, M. C. Asensio, A. Resta, B. Ealet, and G. LeLay, Phys. Rev. Lett. 108, 155501 (2012).
 (59) Z. Shi, Z. Zhang, A. Kutana, and B. I. Yakobson, ACS Nano 9, 9802 (2015).
 (60) K.S. Novoselov, D. Jiang, F. Schedin, T.J. Booth, V.V. Khotkevich, S.V. Morozov, A.K. Geim, Proc. Natl. Acad. Sci. 102, 10451 (2005).
 (61) E. Bekaroglu, M. Topsakal, S. Cahangirov and S. Ciraci, Phys. Rev. B 81, 075433 (2010).
 (62) J. Kunes, R. Arita, P. Wissgott, A. Toschi, H. Ikeda, and K. Held, Comp. Phys. Commun. 181, 1888 (2010).
 (63) G. Kresse, and J. Hafner, Phys. Rev. B 47, RC558 (1993).
 (64) J.M. Soler, E. Artacho, J.D. Gale, A. Garca, J. Junquera, P. Ordejn, and D. SnchezPortal, Journal of Physics: Condensed Matter 14, 2745 (2002).
 (65) P. E. Blöchl, Phys. Rev. B 50, 17953 (1994).
 (66) O.K. Andersen and T. SahaDasgupta, Phys. Rev. B 62, R16219 (2000).
 (67) K.F. Mak, J. Shan, and T.F. Heinz, Phys. Rev. Lett. 106, 046401 (2011).
 (68) B.E. Feldman, J. Martin, and A. Yacoby, Nat. Phys. 5, 889 (2009).
 (69) Y. Zhao, P. CaddenZimansky, Z. Jiang, and P. Kim, Phys. Rev. Lett. 104, 066801 (2010).
 (70) R. Armiento, S. Kümmel, and T. Körzdörfer, Phys. Rev. B 77, 165106 (2008).
 (71) M. Kuisma, J. Ojanen, J. Enkovaara, and T.T. Rantala, Phys. Rev. B 82, 115106 (2010).
 (72) D. Koller, F. Tran, and P. Blaha, Phys. Rev. B 83, 195134 (2011).
 (73) D. Hobbs, E. Piparo, R. Girlanda, and M. Monaca, J. Phys. Condens. Matter 7 2541 (1995).