Two-dimensional oxygen functionalized honeycomb and zigzag dumbbell silicene with robust Dirac cones

The optimized geometry of PDB-h silicene monolayer is shown in figure 1(a). In this configuration, the Si atoms can be divided into two types, fourfold coordinated Si(α) atoms, similar to Si atoms in a bulk crystal, and threefold coordinated Si(β) atoms, similar to atoms in LB silicene. Every three Si(α) and two Si(β) atoms form one cage-like unit, and the cage-like units are arranged in a hexagonal honeycomb pattern. Our calculated lattice constant is 7.43 Å, which is in good agreement with reference [26]. The thickness of the configuration given by the distance between up-plane and down-plane Si(β) atoms is 2.71 Å.

Zoom In Zoom Out Reset image size

With all three-fold silicon atoms bonded with oxygen atoms, we get the ODB-h structure, as shown in figure 1(c). In ODB-h, not only Si(α) but also Si(β) atoms are fully coordinated. The Si–O bond length is 1.55 Å, roughly the same as that in β-cristobalite [42], and smaller than that of Si–O bond in α-quartz (1.61 Å). The lattice constant is increased to 7.49 Å, with a longer distance of 4.32 Å between two same-plane Si(β). The thickness of the monolayer characterized by the distance between up-plane and down-plane O atoms increases to 5.92 Å, in which the distance between up-plane and down-plane Si(β) atoms increases to 2.82 Å, implying stronger sp³ hybridization in Si(β).

By ab initio structure-prediction calculations using a constrained crystal-structure prediction algorithm, Borlido et al. predicted that a zigzag like DB silicene has the lowest energy among all DB silicene structures [29]. As shown in figure 1(b), PDB-z silicene is also realized with cage-like units, but the units are arranged in a zigzag pattern. According to our PBE functional calculations, the primitive cell is a rectangle, with its lattice parameter being 6.40 Å and 7.42 Å, and the area being 4.75 Ų/atom, almost the same as in PDB-h, i.e. 4.78 Ų/atom.

The oxidized form, ODB-z silicene is shown in figure 1(d). With a similar rectangular unit cell, the lattice constants are enlarged to 6.47 Å and 7.48 Å. The distance between two Si(β) layers is 2.83 Å, and the Si–O bond length is 1.55 Å, both of which are very close to that in ODB-h, implying similar hybridization mechanism in ODB-z and ODB-h silicene.

The Si-Si bond lengths are illustrated by the background colors in figures 1(a)–(d). In all the PDB and ODB structures, the bonds between two Si(α) atoms (labelled as r_α,α ) are shorter than 2.365 Å, showing red colors, and the bonds between Si(α) and Si(β) atoms (labelled as r_α,β ) are longer than r_α,α , showing white or blue colors. Comparing the Si–Si bond lengths in pristine DB silicene and oxidized DB silicene, oxidization significantly changes the Si–Si bonds. In PDB structures, r_α,α and r_α,β are quite close to the Si–Si bond length in bulk silicon (r_bulk). In contrast, in ODB structures, the two lengths are very different. By oxidization, r_α,α is shortened to roughly 2.28 Å, close to the Si–Si bond length in LB silicene (r_LB), while r_α,β are increased to about 2.45 Å.

To compare the bonding nature in PDB-h, PDB-z, ODB-h, and ODB-z, electron localization functions (ELF) are calculated and analyzed. To highlight the bonds between Si(α), Si(β) atoms, and between Si(β) and O atoms, we plot the sliced ELF contour map in the top panel of figure 2. The region with an ELF value of 1 indicates fully localized electrons, the region close to 0 is typical of a low electron density area, and the region of 0.5 is an area with a homogeneous electron gas. Interestingly, the ELF contour plots of ODB-h and ODB-z are quite similar, which implies that the bonding mechanism of ODB-h and ODB-z have no significant difference. Comparing PDB and ODB monolayers, it is found that with O atoms adsorbed, the localization of electrons inside the DB silicene layers is weakened. Compared to PDB, electrons between Si(α) and Si(β) in ODB silicene are less localised. Furthermore, in ODB, the electrons around O atoms behave mainly as a homogeneous electron gas. The binding energy of O atoms in PDB-h and PDB-z silicene is, respectively, −5.19 eV/O atom and −5.17 eV/O atom, implying that O atoms are chemisorbed on pristine DB-h and DB-z.

Zoom In Zoom Out Reset image size

Contrary to carbon atoms, silicon atoms prefer to be fourfold co-ordinated compared to threefold co-ordination. This is seen from the relatively higher stability of pristine DB-h and DB-z compared to LB silicene. With dangling bonds filled by O atoms, the stability of ODB-h and ODB-z is quite promising. We will investigate the stability by their energetic preference, which is measured by their formation energy, and through their dynamical stability, which is confirmed by their phonon spectra.

The formation energy of oxidized DB silicene can be hinted at by the reaction: Si₁₀ + O₂ → Si₁₀O₄, in which Si₁₀ denotes PDB silicene, O₂ is oxygen gas, and Si₁₀O₄ is ODB silicene. The negative formation energy (−0.675 eV/atom for ODB-h and −0.667 eV/atom for ODB-z) implies that the reaction is exothermic. We further confirmed their dynamic stability by the phonon spectra calculated along the high symmetry directions in the first BZ, as shown in figure 2 of the bottom panel. It is observed that the phonon spectra of ODB-h and ODB-z are quite similar in the low energy part. In both cases, the optical phonon modes have a gap at around 6 THz, and most of their optical phonon modes below 7 THz are concentrated around 3.5 THz and 5 THz. There are several high energy modes (of about 34 THz) in both structures, whose eigenvectors are shown in the supplementary files (https://stacks.iop.org/NJP/23/023007/mmedia) ODBh.mp4 and ODBz.mp4, respectively [43]. The high energy modes are corresponding to the strong chemical bonds between Si and O atoms [44]. Notably, it can be seen that there are no or only tiny imaginary frequencies near the Γ point, and the structures are therefore expected to be dynamically stable.

To examine the electronic structure of oxygen functionalized DB silicene, we first calculate the local DOS (LDOS) for Si(α), Si(β), and O layers. As shown in figure 3, in both ODB-h and ODB-z, p_z orbital of O atoms contributes to the energy states from −3 eV to −1 eV and hybridizes strongly with p_z orbitals of Si(β). Near the Fermi level p_x,y orbitals of O atoms significantly contribute to the electronic structure and hybridize with the p orbitals of Si(α). The strong hybridization in ODB can decrease the binding energy of O on DB silicene, and thus stabilize the structures. Interestingly, in both cases, the DOS shows a gapless electronic structure with zero DOS at the Fermi level.

Zoom In Zoom Out Reset image size

Here, we analyze the electronic band structures of honeycomb and zigzag PDB and ODB. As shown in figure 4(a), PDB-h monolayer shows an indirect PBE bandgap of ∼0.24 eV with valence band maximum (VBM) at K point and conduction band minimum at Γ point, in agreement with the results of previous works [27, 28]. By oxidization, it is found that ODB-h is a gap-less semiconductor where VBM and CBM meeting at the K point. Near Γ point, the local valence band maximum (LVBM) and the local conduction band minimum (LCBM) are at the same place, Γ point, and LVBM degenerates with second local valence band maximum (LVBM-1), which features are also found in PDB-h. The energy difference between LVBM and LCBM at Γ point is about 0.32 eV, 0.27 eV smaller than that in PDB-h. With two bands crossing linearly at the Fermi level, the charge carriers can be characterized by massless Dirac fermions. By a linear fitting of the band structure, following

$$\begin{equation}{\nu }_{F}=\frac{1}{\hslash }{\nabla }_{k}E\left(\mathbf{k}\right),\end{equation} \tag{ 2 }$$

we obtain the Fermi velocity of ODB-h to be 3.1 × 10⁵ m s⁻¹, which is comparable to those of the group IV elemental monolayers (4.7–8.5 × 10⁵ m s⁻¹), according to the definition, the effective mass of charge carriers

$$\begin{equation}{m}^{{\ast}}=\hslash \left[\frac{{\mathrm{d}}^{2}E\left(k\right)}{\mathrm{d}{k}^{2}}\right]\end{equation} \tag{ 3 }$$

equal to zero near the Fermi level. Therefore, ODB-h silicene is expected to have ultra–high carrier mobility.

Zoom In Zoom Out Reset image size

As shown in figure 4(b), PDB-z monolayer shows an indirect bandgap of roughly 0.62 eV with CBM both at Γ point and VBM near S point. At the Γ point, the energy level of LVBM is only 0.04 eV lower than VBM near S, which is in the correct range for applications in photovoltaics [29]. Interestingly, similar to PDB-h silicene, its oxide, ODB-h, is also a gapless semiconductor with two bands crossing linearly at the Fermi level at a point on the path between Y and Γ points. Since the Dirac cone is not situated at high symmetry points, similar to 1T-YN₂ and 1T-LaN₂ [45, 46], H_4,4,4-graphyne [47], PAI-graphene [48], and Tungsten boride [49], the Fermi velocities are anisotropic. By a linear fitting of the two bands crossing the Dirac point, we get two different Fermi velocities for the two bands (with slopes S > 0 and S < 0) ν_F(S > 0) = 3.4 × 10⁵ m s⁻¹, and ν_F(S < 0) = 2.9 × 10⁵ m s⁻¹. The high Fermi velocity and massless carrier character make the ODB-h and ODB-z ideal materials for building high–speed electronic devices.

As illustrated in figure 3, for ODB-h and ODB-z, the states in the vicinity of the Fermi level are mainly from p_x,y orbitals of O atoms, p_x,y orbitals of Si(α) atoms and p_z orbitals of Si(α) atoms. To reveal the relationship between electronic band structures and atomic orbitals, we have shown the projected orbitals to bands in figures 5(a) and (c), where the size of a dot demonstrates the contribution of the respective atomic orbital to the electronic bands. It is found that in ODB-h, the LVBM and LVBM-1 at Γ point are dominated by p_x,y orbitals of Si(α) with a little contribution by p_x,y orbitals of O atoms, while the VBM/CBM at K points is mainly contributed by p_x,y orbitals of O atoms, and p_z orbitals of Si(α) atoms. Moreover, in ODB-z, p_x,y orbitals of Si(α) atoms mainly contribute to the bands near Γ point, while the Dirac cones are made up by p_x,y orbitals of O atoms and p_z orbitals of Si(β) atoms. These features of the bands at K point are also demonstrated in the density profiles of VBM/CBM electronic wave functions shown in figures 5(b) and (d).

Zoom In Zoom Out Reset image size

To reveal the origin of the Dirac cone explicitly, tight-binding (TB) models are helpful. As illustrated in figures 5(c) and (d), both Si(α) and O atoms contribute to the Dirac cone, the electronic wave functions of O atoms are much larger than those of Si(α) atoms. By proposing simple TB models, which only include the p_x,y orbitals of O atoms, we are able to fit the bands near the Fermi level.

To further illustrate the origin of the Dirac points in the two monolayers, we propose two TB models involving the p_x and p_y atomic orbitals. The effective Hamiltonian is taken as in reference [50]

$$\begin{equation}{H}_{\mathrm{T}\mathrm{B}}=\sum\limits _{i,\alpha }{\varepsilon }_{i}^{\alpha }{c}_{i}^{\alpha +}{c}_{i}^{\alpha }+\sum\limits _{\left\langle i,j\right\rangle ,\alpha ,\beta }{t}_{ij}^{\alpha \beta }\left({c}_{i}^{\alpha +}{c}_{j}^{\beta }+\mathrm{h}.\mathrm{c}.\right),\end{equation} \tag{ 4 }$$

where ${\varepsilon }_{i}^{\alpha }$, ${c}_{i}^{\alpha +}$, and ${c}_{i}^{\alpha }$ represent the on-site energy, creation and annihilation operators of an electron at the α-orbital of the ith atom, respectively. The ${t}_{ij}^{\alpha \beta }$ parameter is the nearest-neighbour hopping energy of an electron between an α-orbital of ith atom and β-orbital of jth, and α, β ∈ (p_x , p_y ). Considering the hexagonal lattice first, we neglect the interactions between the O atoms in different z planes, because the nearest-neighbor distance in the z-direction is more than 20 Å, which is much less than that in the xy plane. We only consider the O atoms in the up/down plane, which is a graphene-like lattice. This graphene-like lattice with p_x,y atomic orbitals leads to a Dirac point in the K point, which is similar to in previous investigated systems of BiX/SbX (X = H, F, Cl, Br, I) [51]. According to the TB theory, only one ${t}_{ij}^{\alpha \beta }$ parameter should be considered, and the hopping energies can be evaluated on the basis of two parameters (V_ppσ and V_ppπ ). Since all atoms are oxygen atoms, we can neglect the on-site energy difference of these atoms and set it as 0. Taking V_ppσ = 0.4 eV and V_ppπ = −0.2 eV, the band structure from TB calculation can reproduce the DFT band structures close to the Fermi level and the Dirac point in the K point (see figure S1 in the supplementary information (SI)). For the rectangular lattice, we still chose to neglect the interactions between the O atoms in different planes and set the on-site energy to zero. Different from the hexagonal lattice, there are two nearest-neighbour hopping energies, leading to four parameters (V_ppσ , V_ppπ , V'_ppσ , and V'_ppπ ). Taking V_ppσ = 0.7 eV, V_ppπ = −0.3 eV, V'_ppσ = 0.0 eV and V'_ppπ = −0.2 eV, we obtained the band structure from TB calculations. The location of the Dirac point on the line from Y point to Γ is in good agreement with the result of DFT. From the above TB calculations for the two lattices, we prove that the Dirac points in the two lattices mainly originate from the p_x and p_y atomic orbitals of O atoms.

External compressing/tensile strain is reported as a feasible method to tune the electronic structures of 2D materials and transform them between metal, semiconductor, and semimetal [28, 52, 53]. Here we investigate the mechanical and electronic properties of ODB-h and ODB-z under strain, which is applied by changing the lattice constant from a₀ to (1 + τ)a₀. Under biaxial compressing/tensile strains of −5% < τ < 10%, as shown in figures 6(a) and (b), the energy changes continuously without any abrupt changes, which indicates that no structural transition occurs. Notably, it does not mean that the two structures are dynamically stable under the strains. Due to the freedom in out-of-plane direction, two-dimensional materials are not dynamically stable under compressive strain [54]. Under tensile strain, as shown in figure S7 in the SI, ODB-h is demonstrated to be still stable when τ = 10%. In contrast, ODB-z is not dynamically stable when τ = 10%, but it is stable when τ = 8%. The structural robustness against intense strain is related to its buckling configuration, similar to the low-bucking structure of silicene.

Zoom In Zoom Out Reset image size

Under such strains, the Dirac cones near the Fermi level in ODB-h and ODB-z are always conserved. The band structure evolution is presented in the SI (figures S2 and S3). In ODB-h, under tensile strain, the Dirac cone is exactly at the Fermi level, and the bandgap at the Γ point is maintained until τ > 9%. In contrast, under compressing strain, the Dirac point moves below the Fermi level, and the bandgap at Γ point vanishes for τ < −2%, making the nanosheets behave as n-type doped. The self-doping behavior arises from the buckling structures, which is also found in LB silicene and germanene [55]. This is because the Dirac cone at K point is dominated by p_z orbitals of Si(α) and p_x,y orbitals of O, which are not sensitive to the in-plane strains. In contrast, the LVBM and LVBM-1 at the Γ point made up of p_x,y orbitals of Si(α), which are affected by the change of Si-Si bond length. With similar structural features, ODB-z behaves quite similarly. Under tensile strain, ODB-z preserves the semimetallic property until 9%, while under compressing strain, the self-doping behavior makes it metallic. Overall, for τ < 9%, the semimetallic features of ODB-h and ODB-z are preserved, with the Fermi velocities changing less than 2.5%, and the insensitivity to tensile strain is advantageous for applications in complicated mechanical environments.

We have also performed band structure calculations with SOC taken into account. Considering that the Si and O are quite light, the SOC effect is generally not strong. For both ODB-h and ODB-z, with or without SOC, the band structures do not change significantly, except that SOC can introduce a bandgap of roughly 10 meV. This bandgap is much more significant than that of LB silicene, i.e. 1.55 meV [5]. The difference is due to the Dirac cone formation mechanism. In LB silicene, the Dirac cone at K point is only formed by p_z orbital of Si atoms, while in ODB-h and ODB-z silicene, p_x,y orbitals of O atoms contribute much more than other atomic orbitals.

Since SOC opens a band gap at the Dirac cones in ODB-h and ODB-z, the topological properties are investigated. It is noteworthy that the chemical adsorption of atoms or molecules also leads to realizing topological insulators [56, 57]. As both ODB-h and ODB-z are inversion symmetric, the topological invariant Z₂ can be calculated using the parity criteria method [58], which claims that the topological invariant can be determined by the parities of the time-reversal and parity invariant points in the BZ, as shown in figure 7. For 2N occupied bands, the Z₂ index ν is defined by

$$\begin{equation}{\left(-1\right)}^{\nu }=\prod\limits _{i}{\delta }_{i},\quad \mathrm{w}\mathrm{i}\mathrm{t}\mathrm{h}\quad {\delta }_{i}=\prod\limits _{n=1}^{N}{\xi }_{2n}\left({{\Gamma}}_{i}\right),\end{equation} \tag{ 5 }$$

where ${\xi }_{2n}\left({{\Gamma}}_{i}\right)$ is the parity eigenvalue of the 2nth occupied energy band at the time-reversal and parity invariant points Γ_i . For ODB-h, δ_i values of the four time-reversal invariant momenta are −1 at Γ point, and +1 at three M points. For ODB-z, the δ_i values are −1 at Γ point, and +1 at X, S, and Y points. Thus the Z₂ indexes of ODB-h and ODB-z are nontrivial, e.g. ν = 1, demonstrating that they are both 2D topological insulators.

Zoom In Zoom Out Reset image size

Considering the underestimation of the bandgap in PBE calculations, we perform HSE06 calculations to examine the results. By HSE06 calculation, in both ODB-h and ODB-z, the electronic bands near the Fermi level are changed a bit: in most regions, the distance between the lowest conduction bands and highest valence bands is enlarged, while the Dirac cones are preserved. Considering the SOC effect, bandgaps of roughly 20 meV open up at the Dirac cones. These details are shown in figures S5 and S6 in the SI.