Eleni Chatzikyriakou*a,
Padeleimon Karafilogloub and
Joseph Kioseogloua
aDepartment of Physics, Aristotle University of Thessaloniki, 54124 Thessaloniki, Greece. E-mail: elchatz@auth.gr; Tel: +30 2310 998109
bLaboratory of Applied Quantum Chemistry, Department of Chemistry, Aristotle University of Thessaloniki, POB 135, 54124 Thessaloniki, Greece
First published on 3rd October 2018
The current carried by a material subject to an electric field is microscopically inhomogeneous and can be modelled using scattering theory, in which electrons undergo collisions with the microscopic objects they encounter. We herein present a methodology for parameter-free calculations of the current density from first-principles using density functional theory, Wannier functions and scattering matrices. The methodology is used on free-standing AB-stacked bilayer penta-silicene. This new Si allotrope has been proposed to have a higher stability than any of its hexagonal bilayer counterparts. Furthermore, its semiconducting properties make it ideal for use in electronic components. We unveil the role of the pz orbitals in the transport through a three-dimensional quantum wire and present current density streamlines that reveal the locations of the highest charge flow. The present methodology can be expanded to accommodate many electron degrees of freedom, the application of electromagnetic fields and many other physical phenomena involved in device operation.
Silicon has historically been the most widely used semiconducting material in the electronics industry both due to its abundance in nature and due to the low-defect interface that it forms with the insulating SiO2 used in transistors and other electronic components. Hexagonal silicene11,12 was one of the first two-dimensional sensations with a plenitude of novel properties, however, when compared to its Dirac cone counterpart, graphene, silicene has a much smaller elastic constant13 and this poses difficulties in its manipulation in the lab.
Other forms of silicon derivatives have been the focus of studies using systematic materials search, revealing novel allotropes with direct band-gaps14 and even high mobility values.15 Penta-silicene is a new form of monolayer silicon that has been observed recently.16,17 Pentagons are rare atomic configurations that are less frequently studied than hexagonal structures mainly due to their difficulties imposed in their fabrication. Using first principles calculations, bilayer penta-silicene, whose layers are stacked in 90° angle between them has been found to be more stable than the most stable form of bilayer hexagonal silicene.18 It is made up of pentagonal rings of Si atoms, similar to penta-graphene19 and other two-dimensional penta-structures.20 Twisting angles between the layers of few-layer materials has been proved to be an efficient method for tuning their properties.21,22 Contrary to the bilayer form with no twisting angle, the AB stacking configuration of penta-silicene has semiconducting properties. This makes it ideal for many different configurations of electronic devices, stemming from the recent surge in heterostructure fabrication.23
In the structures that we examine herein, parts of the system lose their symmetry in favour of stability. Situations such as this, call for the more intuitive chemical description of the system in terms of its valence bonding configuration and this can be very well described with the use of maximally localized Wannier functions (WFs). The latter, acting as alternative representations of Bloch orbitals, have recently gained popularity for examining a plethora of mesoscopic phenomena due to advancements in localization techniques24,25 but also when used as tools for high-throughput screening of topological materials.26 Even more interesting for this study, is the possibility to derive model Tight Binding (TB) Hamiltonians in the WF basis, with specific orbitals involved (i.e. excluding higher conductions bands), arbitrary number of neighbor terms and cut-off values to adjacent cells, and therefore increased accuracy and reduced computational cost, suited to multi-scale modelling of 2D semiconductor transistors.27,28
Therefore, compared to previous work using Kohn–Sham states and mesoscopic transport using Green's functions,29 the use of Wannier functions provides more computational advantages. Finally, the method of scattering matrices is both fast and stable30 and can be expanded to include, among others, spin and other degrees of freedom,31 electromagnetic fields and time-dependent effects.32
Fig. 2 Bandstructure diagram of AB stacked bilayer penta-silicene and orbital projected DoS of atoms 2, 4 and 6. |
Each layer of the structure has tetragonal symmetry and each sublattice is twisted 90° to the other one as shown in Fig. 1a. Aierken et al.18 showed that the structure was stable when the electrons on each of the surfaces of its free-standing form lose their symmetry (Fig. 1b). It is natural that when the material is fabricated on various substrates, depending on the method used and substrate material, bonding and hybridization will be adjusted. Cerdá et al. showed that the pentagons that are formed on Ag(110) host both sp2 and sp3 bonded Si atoms, in which the latter also bond with Ag substrate atoms16 For hexagonal silicene nanoribbons, there has been increased interest in their edge states as doping and hydrogenation has shown to change their electronic and magnetic properties.34,35
Using orbital projected calculations with DFT, the Lowdin charges36 of the atoms in the structure were derived (Table 1). All atoms, except for two at each surface (Si4, 6, 10 and 12) were found to have total charge approximately 4, of which roughly 30% resides in s states and 70% in p states. Then, on each of the surfaces, the Si6 and Si12 atoms lose part their charge to the Si4 and Si10 atoms equivalently, whose hybrid orbitals acquire 35% s and 65% p character.
Atom no. | Total charge (e) | s (e) | p (e) |
---|---|---|---|
Si1–3, 5, 7–9, 11 | 3.95 | 1.16 | 2.78 |
Si4, 10 | 4.12 | 1.45 | 2.67 |
Si6, 12 | 3.74 | 1.29 | 2.44 |
Orbital projected Density of States (DoS) for three representative atoms is plotted next to the band structure in Fig. 2. Most of the atoms have contributions equivalent to that of Si2, except for the two surface atoms that were described previously, for which their pz orbitals show increased contribution around the Fermi level. The bottom of the conduction band is formed by the pz orbitals of Si6 and Si12, while the top of the valence band is formed by the pz orbitals of the Si4 and Si10.
Wannier90 (ref. 37) was used for the Wannierization. In order to derive the effective Hamiltonian from the atomic orbitals, 48 projections were used with atom-centered orbitals, namely, 8 sp3, 4 sp2 and 4 pz, as they were found in the orbital-projected DFT calculations. In Fig. 3, the plots of 5 WFs are given. The results correctly reflect the sp3 nature of the equivalent atoms and the remaining pz orbitals of the distorted symmetry atoms.
Fig. 3 WFs of two atoms. (a)–(d) Show the four orientations of the orbitals in the sp3 hybridized atom (e) shows the pz orbital of an sp2 hybridized atom. |
Transport calculations on a free-standing quantum wire (Fig. 4) with a scattering region of dimensions 9 × 2 × 1 unit cells were performed using Kwant.30 The length expands in the kx direction. The size of the leads is 3 × 2 × 1 with 1D translational symmetry in the directions away from the scattering region. Both the scattering region and the leads are of the same material type, as this allows us to concentrate solely on its properties. A generalization of this can include leads of different material type, or a scattering region with an interface between two different materials. The solution of the Schrodinger equation in the system now corresponds to the plane wave nature of the electrons and not the Bloch waves used in DFT. The Hamiltonian of the system is given by,
(1) |
Fig. 4 Kwant graph showing the sites of a 2 × 9 W × L quantum wire with magenta, green, blue and yellow showing the pz orbitals of Si4, Si6, Si10 and Si12 respectively. |
The on-site energies of the orbitals are given in Fig. 5a. Shown in red are the on-site energies of the pz orbitals of the equivalent atoms, whose energies are higher than the rest of the atomic sp orbitals. Fig. 5b shows the average value of the hopping integrals for each orbital in the unit cell. Both the on-site energies and the hopping integrals show a slight increase in the atoms of the upper layer of the system. The asymmetry of the final TB Hamiltonian is a problem that arises after the disentanglement procedure and which has only recently been addressed.40
Fig. 5 (a) On-site energies of the orbitals of each atom and (b) average value of the hopping integrals showing also a slight increase towards the top orbitals of the system. |
The TB Hamiltonian matrix resulting from the Wannierization is cut-off at three nearest neighbour hoppings. In this context, the nearest neighbours are the hoppings that lie within the home unit cell, second nearest neighbours are the hoppings to adjacent cells in all directions etc. The maximum order of the nearest neighbour interactions is dependent on the k-points used in the Monkhorst–pack mesh during the DFT calculations, as they define the number of periodic images that will appear when the system is translated in real space during Wannierization. As the number of nearest neighbours increases, a large number of neighbour interactions can significantly increase computation time both for creating the system and solving it. However, in most cases, three nearest neighbour interactions provide sufficient accuracy (see ESI†).
The conductance of the wire is calculated from the scattering matrix of the system.30 The dispersion relations within the leads that have 1D translational symmetry are shown together with the conductance plots in Fig. 6. The energies of the modes take on constant values in each direction due to the lack of a 3D translational symmetry. The magnitude of the conductance reflects the probability of transmission of the modes at each electrochemical potential difference between the leads. A value of two shows the there are two modes propagating in the wire at this energy etc.
Local charge and current densities are accessible through solving for individual or pairs of sites in the system, giving access to computations of many properties that exist at surfaces and interfaces. Fig. 7a shows the averaged local density of states from each orbital type which is taken from the expectation value of the local density operator in Kwant using,
(2) |
Fig. 7 (a) Local density of states for each orbital type in the quantum wire (b) orbital-projected DoS derived from DFT calculations using the PBE functional. |
At Ef −0.19 eV, an increased electron density is observed at Si7 sp3 and Si12 pz orbitals. To explain this, the expectation value of the current density at this energy was also extracted for the propagating wavefunctions originating from lead 1 (left lead in Fig. 4). For a hopping from site l to site k, this is calculated from,
(3) |
In the absence of any potential or magnetic field defined explicitly over the 3D region, the distribution of the current density is dictated by the confinement effects in the x, y and z directions,7 that is by the solution of the Schrodinger equation with the translation operator in the leads defining the energies of the propagating modes. Finally, the on-site energies as well as the hopping integral values define the electrostatic landscape for the electrons to find their route through the wire (Fig. 5).‡
Taking into account collective effects from many orbitals, the effects of impurities and screening or the formation of dipoles at interfaces41 can furthermore be examined, opening the road to versatile computational microscopic and topographic studies, while phonon effects can be added as site self-energies using standard procedures.8 Progress in obtaining and engineering model TB Hamiltonians from first principles is currently flourishing42 and is expected to lead to significant advances in predictions of device characteristics.
This methodology presents many advantages for the examination of electron device operation as it allows many microscopic details of the quantum transport to be revealed using realistic values and the addition of disorder, spin, phonon and finite-temperature effects. Bias calculations are also possible by summing the propagating modes around the Fermi level in the scattering region.43
Footnotes |
† Electronic supplementary information (ESI) available. See DOI: 10.1039/c8ra05652h |
‡ While writing this manuscript the authors were made aware of another work that adds to the Wannierization procedure by creating symmetrized tight binding models.40 This is expected to significantly increase the accuracy of the ideal ground state model Hamiltonian. |
§ http://www.quantum-espresso.org/wp-content/uploads/upf_files/Si.blyp-hgh.UPFhttp://www.quantum-espresso.org/wp-content/uploads/upf_files/Si.blyp-hgh.UPF |
¶ http://www.quantum-espresso.org/wp-content/uploads/upf_files/Si.pbe-hgh.UPFhttp://www.quantum-espresso.org/wp-content/uploads/upf_files/Si.pbe-hgh.UPF |
This journal is © The Royal Society of Chemistry 2018 |