Electrostatic quantum dots in silicene

We study electrostatic quantum dot confinement for charge carriers in silicene. The confinement is formed by vertical electric field surrounding the quantum dot area. The resulting energy gap in the outside of the quantum dot traps the carriers within, and the difference of electrostatic potentials on the buckled silicene sublattices produces nonzero carrier masses outside the quantum dot. We study the electrostatic confinement defined inside a silicene flake with both the atomistic tight-binding approach as well as with the continuum approximation for a circularly symmetric electrostatic potential. We find localization of the states within the quantum dot and their decoupling from the edge that makes the spectrum of the localized states independent of the crystal termination. For an armchair edge of the flake removal of the intervalley scattering by the electrostatic confinement is found.


B. Szafran, D. Żebrowski & Alina Mreńca-Kolasińska
We study electrostatic quantum dot confinement for charge carriers in silicene. The confinement is formed by vertical electric field surrounding the quantum dot area. The resulting energy gap in the outside of the quantum dot traps the carriers within, and the difference of electrostatic potentials on the buckled silicene sublattices produces nonzero carrier masses outside the quantum dot. We study the electrostatic confinement defined inside a silicene flake with both the atomistic tight-binding approach as well as with the continuum approximation for a circularly symmetric electrostatic potential. We find localization of the states within the quantum dot and their decoupling from the edge that makes the spectrum of the localized states independent of the crystal termination. For an armchair edge of the flake removal of the intervalley scattering by the electrostatic confinement is found.
Silicene 1,2 is a material similar in crystal and electron structure to graphene 3 but with enhanced spin-orbit coupling [4][5][6] that makes this two-dimensional medium attractive for studies of anomalous-5 , spin-6 and valley-quantum Hall effects 7 , giant magnetoresistance 8,9 and construction of spin-active devices 10,11 . The crystal structure of a free-standing silicene is buckled 12 with a relative shift of the triangular A and B sublattices in the vertical direction. The shift allows one to induce and control the energy gap near the charge neutrality point 13,14 . The silicene was first successfully formed on metallic substrates [15][16][17][18][19][20] . For the studies of electron properties of systems based on silicene non-metallic substrates 21 are needed. Theoretical studies have been performed for the silicene on insulating AlN 22 , and semiconducting transition metal dichalcogenides (TMDCs) 1,[23][24][25] . An operating room-temperature field effect transistor was recently realized 26 with silicene layer on Al 2 O 3 insulator. Al 2 O 3 only weakly perturbs the band structure of free-standing silicene near the Dirac points 27 . The silicene islands have also been grown on graphite 28 by van der Waals heteroepitaxy.
In this paper we study formation of an electrostatic quantum dot within the silicene. The electrostatic quantum dots 29 allow for precise studies of the carrier-carrier and spin-orbit interaction. In graphene the electrostatic confinement is excluded since the carriers behave like massless Dirac fermions that evade electrostatic confinement due to the lack of the energy gap in the dispersion relation and chiral Klein tunneling that prevents backscattering 30 . A local electrostatic potential in graphene can only support quasibound states 31,32 of a finite lifetime and cannot permanently trap the charge carriers. Carrier confinement and storage can be realized by finite flakes of graphene [33][34][35][36][37] . However, the electron structure of states confined within the flakes depends strongly on the edge 33,34 that is hard to control at the formation stage and cannot be changed once the structure is grown. The electrostatic confinement 29 is free from these limitations. Finite flakes of silicene as quantum dots were also discussed [38][39][40] . For the graphene, the energy gap 41 due to the lateral confinement or mass modulation by eg. a substrate allows for formation of quantum dots by external potentials [42][43][44][45] . Confinement by inhomogeneity of the magnetic field has also been proposed for graphene 46,47 which removes the edge effects.
The electrostatic quantum dots studied below are formed by an inhomogeneous vertical electric field. We consider a system in which the confinement of the carriers is induced within a region surrounded by strong vertical electric fields [see Fig. 1]. The inhomogeneity of the electric field is translated into position-dependence of the energy gap. Localized states are formed within a region of a small energy gap surrounded by medium of a larger gap. A similar confinement mechanism has previously been demonstrated for bilayer graphene 48 , which also reacts to the vertical electric field by opening the energy gap. The vertical electric field produces potential variation at the A and the B sublattices of the buckled silicene [ Fig. 1(b)]. In this way the system mimics the idea for potential confinement of neutrinos introduced by Berry and Mondragon 49 . A potential of a different sign for the components of spinor wave function was applied 49 that produces a so-called infinite-mass boundary in the limit case of a large potential. The infinite-mass boundary condition is applied for phenomenological modeling of graphene flakes with the Dirac equation 33,34,45,49,50 . The proposed device is a physical realization of this type of the boundary condition. Note that for monolayer TMDCs 51 , materials with hexagonal crystal lattice, the valley degree of freedom and strong spin-orbit coupling, formation of electrostatic quantum dots 52 is straightforward due to the wide energy gap of the system. However, these systems are far from the Dirac physics for massless or light carriers. On the other hand the buckled germanene grown by heteroepitaxy on graphite with the linear Dirac band structure 53 should allow for formation of electrostatic confinement in a manner described below.

Theory
Model system. We consider silicene embedded in a center of a dielectric layer sandwiched symmetrically between metal gates [ Fig. 1(a)]. The distance between the gates is 2h = 2.8 nm outside a circular region of radius 2R = 40 nm, where the spacing is increased to 2H = 28 nm. The model device is a symmetric version of an early electrostatic GaAs quantum dot device 54 . The electrostatic potential near the charge neutrality point can be estimated by solution of the Laplace equation with the Dirichlet boundary conditions at the gates. The solution on the A and B sublattices is shown in Fig. 1(b) for the gate potential V g = 10 V. A potential difference between the sublattices presented in Fig. 1(b) appears as a result of the buckled crystal structure with the vertical distance d = 0.046 nm; between the sublattices [see Fig. 1(a)]. The difference is large outside the central circular area. Beyond this area the potential is U A = eV g d/2h for the A sublattice and U B = −eV g d/2h for the B sublattice [ Fig. 1(b)]. Near the center of the circular area the potential is U A = eV g d/2H, U B = −U A with the gate potential lever arm increased by the larger spacing between the gates. The bottom of the electrostatic potential in the center of the dot in Fig. 1 , with m ∈ (6, 8) i.e. in the Taylor expansion of the potential the parabolic term corresponding to the harmonic oscillator potential is missing. Therefore, in the calculations below we consider a rectangular potential well model A g g and U B = −U A , where r is the in-plane distance from the center of the system. The model potential is plotted with the dotted lines in Fig. 1(b). The results for the exact electrostatic potential are also discussed below. For the discussion of the confinement potential profile depending on the geometry of the gates see ref. 55 . Note that the gate voltage to confinement potential conversion factor depends not only on the spacing between the gates but on the buckling distance 1 which varies for different substrates. The silicene on graphite is characterized by the lattice parameters including the buckling to the free-standing silicene, with the buckling found close to 0.05 nm 28 . On the other hand, for silicene on MoS 2 the buckling can be as large as 0.2 nm 25 . In order to illustrate the dependence on the buckling, the confinement potential at the A sublattices was plotted in Fig. 1  where σ z is the Pauli spin matrix, α † c k is the electron creation operator at ion k with spin α, the symbols 〈 〉 k l , and 〈〈 〉〉 k l , stand for the pairs of nearest neighbors and next nearest neighbors, respectively. The first term of the Hamiltonian accounts for the nearest neighbor hopping with t = 1.6 eV 4,5 . The second term describes the intrinsic spin-orbit interaction 57 with the sign parameter v kl = +1 (v kl = −1) for the counterclockwise (clockwise) next-nearest neighbor hopping and λ so = 3.9 meV 4,5 . The exponents in the first and second sum introduce the Peierls phase, with the vector potential → A . The term with U introduces the model electrostatic potential given by Eq. (1). The last term is the spin Zeeman interaction for perpendicular magnetic field, where μ B is the Bohr magneton and the electron spin factor is g = 2. The applied Hamiltonian is spin diagonal in the basis of σ z eigenstates. We consider the states confined within the confinement potential that is defined inside a finite silicene flake containing up to about 72.5 thousands ions.
where the valley index is η = 1 for the K valley and η = −1 for the K′ valley, τ x , τ y and τ z are the Pauli matrices in the sublattice space, 3 2 is the Fermi velocity with the nearest neighbor distance a = 0.225 nm.
For the isotropic potential U(r) and the symmetric gauge → = − A B y Bx ( /2, /2, 0) the H η Hamiltonian commutes with the total angular momentum operator of the form is the orbital angular momentum operator, and I is the identity matrix. The components of the common H η and J z eigenstates can be put in a separable form where m is an integer. For the K′ valley states we will denote the quantum number by m′. The asymptotic behavior of the radial functions for a given m and η in the center of the potential is 50 . The radial components fulfill the system of equations In the continuum approach we look for the states localized within the confinement potential of radius R within a finite circular flake of radius R′. We are interested in the influence of the type of the flake on the states localized within the electrostatic potential well. In the continuum approximation at the edge of the flake we apply two types of boundary conditions: the infinite-mass boundary condition  Figure 2 shows the energy spectrum and the localization of energy levels obtained with the atomistic tight binding [ Fig. 2(a)] and with the continuum approach [ Fig. 2(b-d)] as functions of the gate voltage. In this plot the spin-orbit interaction was switched off. A vertical magnetic field of 0.5 T is applied, for which splitting of energy levels with respect to the valley but not with respect to the spin is visible on the scale of Fig. 2. One observes the splitting of the energy levels with respect to the orbital angular momentum in the external magnetic field. In the atomistic tight-binding approach [ Fig. 2(a)] a hexagonal flake of side length 43 nm and an armchair boundary were taken. For the continuum approach a circular flake of radius R′ = 2R = 40 nm (b,c) and R′ = 4R = 80 nm (d) were studied. In Fig. 2(b,d) the infinite mass boundary condition is applied at the end of the flake and in Fig. 2(c) the zigzag boundary condition is assumed. The energy levels get localized inside the quantum dot-see the color of the points that indicate the localization of the electron probability density inside the quantum dot. The zigzag edge applied in Fig. 2(c) supports the edge-localized energy levels which correspond to zero energy in the absence of external fields. The edge energy levels for the zigzag flake in Fig. 2(c) are split by the gate potential 38,40 . The energy of the edge states 38 follow the potential energy at the separate sublattices outside the quantum dot. The edge states are missing for the armchair edge of the hexagonal flake adopted for the tight-binding calculations in Fig. 2(a) and for the infinite-mass boundary condition adopted in the continuum model in Fig. 2(b,d). For the larger circular flake [ Fig. 2(d)] the spacing between the energy levels localized outside of the dot is decreased, but the same spectrum of the localized states is found. We can see in all the panels of Fig. 2 that the energy spectra of localized states obtained by the atomistic and continuum approaches with varied boundary conditions become similar for larger V g . The localized states are found in between the two thick gray lines that show potential energy at the A and B sublattices outside the quantum dot U = ±eV g d/2h. A perfect agreement between the energies of the localized states in the tight-binding and Dirac models is obtained for the energy levels that are the closest to the charge neutrality point (E = 0). For the energy levels that are closer to edge states energy [cf. Fig. 2(a) and (b) for E > 100 near V g = 10 V], the wave functions of the localized states penetrate into the region outside of the quantum dot. The external region is different in all the plots of Fig. 2, hence the resolvable difference of the energy levels. The spectrum of the zigzag flake [ Fig. 2(c)] indicates that the confinement of subsequent states within the quantum dot area appears with the crossing of the confined energy levels with the edge states 38 which shift linearly with the external potential. The edge states and thus the crossings are missing in the results obtained with the armchair edge [ Fig. 2(a)], and the infinite-mass boundary conditions [ Fig. 2(b,d)].

Results
The effects of the spin-orbit coupling and the results for the exact confinement potential are given in Fig. 3. The plot of Fig. 3(b)-with the external field 0.5 T and without the spin-orbit interaction is the zoom of Fig. 2(b). The Zeeman spin splitting is still not resolved at this energy scale, but the splitting of the energy states with respect to the valley is evident. The K (K′) valley states with the indicated angular momentum quantum number m (m′) for the A sublattice is given in the Figure. In Fig. 3(b) all the energy levels are nearly degenerate with respect to the spin. For comparison the result for 0 T is plotted in Fig. 3(a), where all the energy levels are strictly degenerate. The degeneracy is fourfold: with respect to both the spin and the valley. The results with the intrinsic spin orbit coupling are displayed in Fig. 3(c) for B = 0. The intrinsic spin-orbit interaction introduces an effective valley-dependent magnetic field which forms spin-valley doublets. The energy effects of the splitting is comparable to the external magnetic field of 0.5 T given in Fig. 3(b).
The results of the present manuscript were obtained with the rectangular quantum well potential [Eq. (1)] approximation to the actual electrostatic potential [see Fig. 1(b)]. The results for the rectangular potential [ Fig. 3(b)] can be compared to the ones with the exact potential [ Fig. 3(d)]. The energy levels for the exact potential are shifted up on the energy scale-since the rectangular potential well is a lower bound to the exact potential [cf. Fig. 1(b)]. However, the order of the energy levels and the relative spacings obtained with the exact potential are close to the ones obtained for the quantum well ansatz. Figure 2 demonstrates that the dot-localized states are insensitive to the type of the edge and the size of the flake, which results from their decoupling from the edge. In particular, the valley mixing by the armchair edge is removed. The removal of the valley mixing has distinct consequences for the electron wave functions as described within the atomistic approach. Figure 4 shows the absolute value of the probability amplitude at the A (left column) and B (right column) sublattices for varied values of the gate voltage and the lowest-energy conduction-band state for B = 0.5 T. For V g = 0 the electron density at both the sublattices undergoes rapid oscillations which result . The probability density is then (r) ( r) (r) Due to the large distance between K and K′ in the reciprocal space the exponent term induces rapid variation of the density from one atom to the other even when φ | | A 2 and φ | | ′ A 2 densities are smooth. A smooth |Ψ | A amplitude can only be obtained provided that one of the valley components φ A or φ ′ A is zero. Figure 4 shows that indeed as the gate voltage is increased the rapid oscillations of the density disappear. The valley mixing disappears along with the coupling to the edge.
In Fig. 4 a circular symmetry of the confinement potential is reproduced by the electron density for larger V g . In the lowest-energy conduction band state that we plot in Fig. 4 the density is locally maximal in the center of the quantum dot in the A sublattice. For the B sublattice a zero of the density is found. In the continuum approach the ground state at B > 0 corresponds to K′ valley with the total angular momentum j′ = −1/2 or m′ = 0 in Eq. (4). The A component of the wave function corresponds to an s state and the B component to a p state, which agrees with the results of Fig. 4(g,h). For V g = 10 V the electron density far from the dot disappears. A penetration of the electron density outside the nominal radius of the dot is still present, but short range.
In Fig. 5(a) we plot the absolute value of the wave function for the same state as obtained with the continuum approach with the infinite-mass boundary condition at the end of the flake R′ = 80 nm. The vanishing derivative of the probability amplitude at r = 0 is found for the A sublattice and a linear behavior for the radial function on the B sublattice. Equations (5, 6) translate the potential step into a jump in the derivative of the radial functions. V g shifts most of the probability amplitude of the lowest-energy conduction band states to the A sublattice (see also Fig. 4). However, for large V g the radial functions for both sublattices tend to the same amplitude (see Fig. 5(a)), since in the limit of infinite V g the variation of the electrostatic potential at the outside of the dot induces an infinite-mass boundary at r = R, which implies equal amplitudes of the wave functions therein 33,34,49 . The results obtained with the exact potential are given in Fig. 5(b). The derivatives of the wave function are continuous for the smooth potential variation. The maxima of the wave function amplitude on the B sublattice no longer exactly coincide with R for the exact potential. Also, Fig. 1(b) shows that the energy difference between the sublattices is larger for the exact potential than in the rectangular quantum well approximation.

Summary and Conclusions
We found formation of states localized by external electrostatic potential within a silicene flake. The potential used for this purpose results from the inhomogeneity of the vertical electric field that induces an energy gap outside the quantum dot and the buckling of the silicene surface. The energy spectrum for a finite flake can be separated into quantum-dot localized states and the states delocalized over the rest of the flake. The localized and delocalized states appear in separate parts of the energy spectrum limited by the electrostatic potential on the separate sublattices of silicene. We have demonstrated that the states localized within the quantum dot are separated from the edge and independent of the boundary condition applied therein. A very good agreement between the atomistic tight-binding and continuum model results have been obtained. The electrostatic confinement opens perspectives for studies of localized states in the anomalous, spin and valley quantum Hall effects conditions. Data availability. All data generated or analysed during this study are included in this published article. . The infinite mass boundary conditions are applied at the end of the flake at r = R′ = 80 nm. Panel (a) shows the results for the rectangular potential well [Eq. (1)] and (b) for the smooth electrostatic potential of Fig. 1(b).