Thirteen-band Tight-binding Model for the MoS₂ Monolayer

Abstract

A tight-binding model is fitted to density-functional calculations of the electronic structure of the MoS2 monolayer. The model involves 13 atomic orbitals per unit cell: the 4d orbitals of the molybdenum atom plus the 3s and 3p orbitals of each sulfur atom. The hopping and overlap couplings of each atom with its first nearest neighbors in each crystalline sublattice are considered. Different values are allowed for the intraplane and interplane S-S hopping integrals. A closed-form expression is given for the effective-mass tensor at stationary points. The isotropy of the valence and conduction bands near the edges of the fundamental gap is proven. The role played by the orbital overlapping as well as the crystal-field splitting of the molybdenum 4d level is discussed.

Keywords:
transition-metal dichalcogenide; electronic structure; tight-binding model; density functional theory

1. Introduction

New two-dimensional (2D) materials have attracted a lot of attention since 2004¹. This is due to amazing advances in the investigation of graphene. In this trend, monolayers of transition-metal dichalcogenides (TMDs) such as MoS₂, WS₂, MoSe₂ and WTe₂, have been investigated as well. In the laboratory, the new materials can be obtained by mechanical exfoliation² or molecular vapor deposition³, and can be characterized by different techniques, such as optical and electron microscopy², and Raman and photoluminescence spectroscopy³. The efforts have revealed new physical phenomena, such as Landau quantization in an insulating WTe₂ monolayer⁴, and has led to technological applications, including photodetectors⁵, flexible devices, light-emitting devices, solar cells⁶ and transistors⁷. Biomedical applications in cancer therapy, drug release and tissue-engineering scaffold design have also been considered⁸.

The experimental work supports as well as benefits from the theoretical advances in the description of physical and chemical properties of these materials. Most of those advances rely on ab-initio calculations of the electronic structure, involving the density-functional theory (DFT). Costa-Amaral et al.⁹ reported an investigation of the adsorption properties of $3d$, $4d$, and $5d$ metal adatoms on the MoS₂ monolayer, by applying spin-polarized DFT within the D3 van der Waals correction. Deng and Li¹⁰ predicted a new semiconductor and ferromagnetic 2D crystalline structure of MoS₂. Gabourie et al.¹¹ used molecular dynamics to study the thermal conductivity of monolayer and bilayer structures of MoS₂, comparing suspended, supported, and encased structures. Abdi and Astinchap¹² used the Holstein model to investigate thermoelectric properties of the MoS₂ monolayer.

Further understanding and applications of the electronic structure of TMDs has been gained through semi-empirical tight-binding models (TBMs)^13-15, because of their conceptual, geometrical, and computational simplicity. In short, they approximate the Bloch functions by linear combinations of atomic orbitals (LCAOs), apply variational procedures, explore symmetry, solve eigenvalue problems, and adjust a set of energy and overlap parameters¹⁶.

In 2013, Cappelluti et al.¹³ presented a TBM for the MoS₂ monolayer as well as bulk MoS₂. The $4d$ orbitals of Mo and the $3p$ orbitals of S are considered. This yields 11 energy bands: seven valence bands and four conduction bands. The tight-binding calculations are fitted to DFT energy values, using $11$ ($13$) parameters for the monolayer (bulk) crystal. All bulk parameters are transferred to the monolayer, except for two interlayer hopping integrals between S atoms. In this way, the direct-gap/indirect gap transition is successfully described. Three years later, Silva-Guillén et al.¹⁴ improved their fitting procedures and calculated the electronic structure and the optical conductivity of MoS₂ and other TMDs.

In 2015, Ridolfi et al.¹⁵ used a TBM with $12$ parameters to fit a DFT calculation. Considering the same set of orbitals as Cappelluti et al.¹³, the five energy bands mostly associated to Mo orbitals are fitted very well. This gives satisfactory values for the energy gap, the valence- and conduction-band effective masses, and relevant local extrema of selected bands. However, the six tight-binding valence bands dominated by the $3p$ orbitals of S and their DFT counterparts differ appreciably. The former (latter) set of bands roughly ranges from $-98.5$ eV to $-24.1$ eV ($-6.6$ eV to $-1.1$ eV). As a result, the Mo-S (S-S) hopping integrals may have been overestimated (poorly estimated).

While the above described tight-binding approaches agree with DFT calculations at certain points of the electronic structure, much work remains to be done. Because of simplifications, the fitted parameters can fail to describe the interatomic couplings appropriately. To improve the approach, one can enrich the set of atomic orbitals as well as describe the atomic interactions more realistically. Additionally, the TBM effective-mass values are routinely estimated through quadratic fittings of the energy bands, even though analytical calculations can be performed. To advance in these tasks, the present work investigates a TBM that fits $13$ energy bands of the MoS₂ monolayer obtained by DFT. The main novelties compared with previous works^13,15 rely on: (i) including the $3s$ orbital of S, (ii) considering finite overlaps between the atomic orbitals of nearest neighbors in each crystalline sublattice, (iii) allowing for differences between intraplane and interplane S-S hopping integrals, (iv) adjusting the tight-binding parameters to find good overall agreement between the methods, (v) giving a closed-form expression for the effective-mass tensor, and (vi) proving that the valence and conduction bands are isotropic near the gap edges.

The manuscript is organized as follows. Section 2 describes the crystalline structure of the MoS₂ monolayer, Section 3 describes the DFT procedures, Section 4 presents the tight-binding model, and Section 5 gives the numerical results and discussions. Section 6 summarizes our main findings and perspectives.

2. The MoS₂ monolayer

As depicted in Figure 1, the MoS₂ monolayer has three atoms per unit cell. The structure consists of one Mo atomic plane at $z=0$, sandwiched between two S atomic planes at $z=\pm h$, with¹⁵ $h\approx 1.56$Å. It can be decomposed into three crystalline sublattices. The atomic positions on each plane have the form $R+{\rho }_{\zeta }$. The lattice vectors are given by $R=n_1{a}_1+n_2{a}_2$, where $n_1$ and $n_2$ are integer numbers, $a_1=a\left(+\mathrm{1,}\sqrt{3}\mathrm{,0}\right)/2$ and $a_2=a\left(-\mathrm{1,}\sqrt{3}\mathrm{,0}\right)/2$, with¹⁵ $a\approx 3.16$Å, whereas the relative atomic positions are ${\rho }_0=0$, and ${\rho }_{\pm h}=\left(a_1+a_2\right)/3+\left(\mathrm{0,0,}\pm h\right)$. The Wigner-Seitz cell is the hexagon on the $xy$-plane obtained by requiring that $\left(x,y\mathrm{,0}\right)\cdot R\le R^2/2$ applies for each lattice vector $R$.

Figure 1
(a) 3D pictorial view of the MoS₂ monolayer. The blue (yellow) spheres are for the Mo (S) atoms. (b) Top view of panel (a), where the red arrows are for the vectors $a_1$ and $a_2$. (c) Left view of panel (a). (d) Front view of panel (a), where the red arrow are the vectors ${\rho }_h$ and ${\rho }_{-h}$.

The reciprocal lattice is hexagonal as well, with primitive vectors given by $b_1=2\pi \left(\mathrm{1,1}/\sqrt{3}\mathrm{,0}\right)/a$ and $b_2=2\pi \left(-\mathrm{1,1}/\sqrt{3}\mathrm{,0}\right)/a$, and the first Brillouin zone being the hexagon shown in Figure 2. The high-symmetry points $\Gamma =0$, $\text{K}=\left(2b_1+b_2\right)/3=2\pi \left(\mathrm{1,}\sqrt{3}\mathrm{,0}\right)/3a$ and $\text{M}=\left(b_1+b_2\right)/2=2\pi \left(\mathrm{0,1,0}\right)/\sqrt{3}a$ as well as the intermediate points$\Lambda =uK$, $\text{T}=\left(1-u\right)\text{K}+u\text{M}$ and $\Sigma =\left(1-u\right)\text{K}$, where $0\le u\le 1$, are displayed. Our choice for the $\text{M}$ point differs from that in previous works^13,15. It aims to simplify the symmetry analysis of Bloch functions with respect to the plane $x=0$. The triangular path $\Gamma -\text{K}-\text{M}-\Gamma$ is used to partially display the energy bands.

Figure 2
The first Brillouin zone of the MoS₂ monolayer. The points Γ, Κ, M as well as the lines along generic points Λ, Τ, Σ are symmetry elements of the reciprocal lattice.

Besides the translational invariance, the mirror symmetry across the plane $z=0$ is the most apparent symmetry of the MoS₂ monolayer. Hence, the Bloch states can be classified as either symmetric (even) or antisymmetric (odd) under reflection with respect to the $xy$ -plane. Moreover, the monolayer presents reflection and rotation symmetries with respect to planes and axes that are perpendicular to the $xy$-plane, as given by the plane group p3m1^17,18. Together with the time-reversal symmetry, this guarantees that the energy bands present the same symmetries as the Brillouin zone depicted in Figure 2. Then, our calculations can be restricted to the irreducible region limited by the triangular path in the same figure.

The symmetry group has mirror planes at $x=na/2$ and $x\pm \sqrt{3}y=na$, where $n$ takes integer values [see Figure 1b]. They are all equivalent to the plane $x=0$. The reflection symmetry with respect to the latter plane manifests itself in the Bloch functions with wave vector $k$ along the Σ line displayed in Figure 2. They can be classified as either symmetric or antisymmetric across the plane $x=0$. The group also has vertical glide planes given by $x=\left(2n+1\right)a/4$ and $x\pm \sqrt{3}y=\left(2n+1\right)a/2$, where $n$ takes integer values.

Considering vertical rotation axes, the group has third-order axes passing through each Mo atom at $\left(x,y\mathrm{,0}\right)=R$, through each pair of S atoms at $\left(x,y,\pm h\right)=R+{\rho }_{\pm h}$, and through each point given by $\left(x,y\mathrm{,0}\right)=R-\left(a_1+a_2\right)/3$. Here $R$ is a lattice vector [see Figure 1b]. The rotational symmetry around the $z$-axis manifests itself on the Bloch functions of the $\Gamma$ and $\text{K}$ points. They are eigenvectors of the counterclockwise rotation in $120°$, with the eigenvalue being a cubic root of 1, namely $1$ or $\text{exp}\left(\pm 2\pi i/3\right)$.

3. Density-functional-theory

The DFT calculation was run in the PWscf package developed by Quantum ESPRESSO¹⁹ under Linux Ubuntu in a personal computer. Since this code deals with 3D crystalline structures, we consider an infinite stack of MoS₂ monolayers. The spacing between the Mo planes of consecutive monolayers is $5h=7.8$ Å, i.e., we deal with a bulk crystal with a hexagonal lattice generated by the vectors $a_1$, $a_2$, and $\left(\mathrm{0,0,5}h\right)$. The spacing is rather large to ensure that van der Waals interactions between the monolayers is negligible. Perdew-Burke-Ernzerhof functionals for scalar relativistic pseudopotentials with projector-augmented waves and non-linear core corrections of Mo and S atoms are given in files²⁰Mo.pbe-spn-kjpaw_psl.1.0.0.UPF and S.pbe-n-kjpaw_psl.1.0.0.UPF, respectively. Cutoff parameters are given by ecutwfc = 80 and ecutrho = 400. The self-consistent calculation of the crystalline potential is done for $21$ energy bands over a $21\times 21\times 1$ wave-vector mesh. Then, a non-self-consistent calculation of the same bands is performed along the $\Gamma -\text{K}-\text{M}-\Gamma$ path.

4. Tight-binding-model

The TBM approximates each Bloch function of an electron in the MoS₂ monolayer by an LCAO of the constituent elements Mo and S. Here, the spin-orbit interaction is neglected, and the electron Hamiltonian is independent of spin. Considering the $4d$ orbitals of Mo and the $3s$ and $3p$ orbitals of S, the Bloch function reads

${\psi }_k\left(r\right)={\displaystyle \sum }_{\zeta \in \left\{\mathrm{0,}\pm h\right\}}{\displaystyle \sum }_{\phi \in A_{\left|\zeta \right|}}{c}_{\zeta ,\phi ,k}{\varphi }_{\zeta ,\phi ,k}\left(r\right),$ (1)

where

$${\varphi }_{\zeta ,\phi ,k}\left(r\right)={\displaystyle \sum }_R\phi \left(r-R-{\rho }_{\zeta }\right)e^{ik\cdot R},$$ (2)

and $A_0$ and $A_h$ consists of atomic orbitals of Mo and S, namely, $A_0=\left(d_{z^2},d_{xz},d_{yz},d_{x^2-y^2},d_{xy}\right)$ and $A_h=\left(s,p_z,p_x,p_y\right)$. Using the matrix product, the Bloch function can be rewritten as

${\psi }_k\left(r\right)={\displaystyle \sum }_{\zeta \in \left\{\mathrm{0,}\pm h\right\}}{\varphi }_{\zeta ,k}\left(r\right)c_{\zeta ,k}$ (3)

with

$${\varphi }_{\zeta ,k}\left(r\right)={\displaystyle \sum }_R{\text{e}}^{ik\cdot R}{\phi }_{\zeta }\left(r-R-{\rho }_{\zeta }\right).$$ (4)

The coordinates of both the row matrix ${\phi }_{\zeta }\left(r\right)$ and the column matrix $c_{\zeta ,k}$ are arranged as the orbital $\phi$ runs over the sequence $A_{\left|\zeta \right|}$. One may also write ${\psi }_k\left(r\right)={\varphi }_k\left(r\right)c_k$, with a row matrix ${\varphi }_k\left(r\right)=\left({\varphi }_{\mathrm{0,}k}\left(r\right),{\varphi }_{h,k}\left(r\right),{\varphi }_{-h,k}\left(r\right)\right)$ and a column matrix $c_k=\left(c_{\mathrm{0,}k},c_{h,k},c_{-h,k}\right)$.

The TBM requires finding the best values for the coefficients $c_{\zeta ,\chi }\left(k\right)$ in Equation 1. This is done by minimizing the expected value of the electron energy in the Bloch state, which is given by the Rayleigh quotient

$E_k=\frac{{{\displaystyle \int }}_{-\infty }^{+\infty }{{\displaystyle \iint }}_{\text{WS}}^{}{\psi }_k^{*}\left(r\right)\widehat{H}{\psi }_k\left(r\right)dxdydz}{{{\displaystyle \int }}_{-\infty }^{+\infty }{{\displaystyle \iint }}_{\text{WS}}^{}{\psi }_k^{*}\left(r\right){\psi }_k\left(r\right)dxdydz}=\frac{c_k^{†}{ℍ}_k{c}_k}{c_k^{†}{\mathbb{S}}_k{c}_k}.$ (5)

Here WS denotes the Wigner-Seitz cell of the MoS₂ monolayer,

${ℍ}_k={{\displaystyle \int }}_{-\infty }^{+\infty }{{\displaystyle \iint }}_{\text{WS}}^{}{\varphi }_k^{†}\left(r\right)\widehat{H}{\varphi }_k\left(r\right)dxdydz$ (6)

and

${\mathbb{S}}_k={{\displaystyle \int }}_{-\infty }^{+\infty }{{\displaystyle \iint }}_{\text{WS}}^{}{\varphi }_k^{†}\left(r\right){\varphi }_k\left(r\right)dxdydz.$ (7)

These 13×13 matrices are Hermitian, and can be given in terms of their blocks

${\mathbb{H}}_{\zeta ,{\zeta }^{\mathrm{\text{'}}},\mathit{k}}={\mathrm{}\mathrm{e}}^{i\mathrm{}\mathit{k}\bullet ({\mathit{\rho }}_{\zeta }-{\mathit{\rho }}_{\zeta }\mathrm{\text{'}})}\mathrm{}\sum _{\mathit{d}\in D_{\zeta ,{\zeta }^{\mathrm{\text{'}}}}}{\mathrm{}\mathrm{e}}^{i\mathrm{}\mathit{k}\bullet \mathit{d}}\mathrm{}{\mathbb{𝕙}}_{\zeta ,{\zeta }^{\mathrm{\text{'}}},\mathit{d}},\mathrm{}\mathrm{}$ (8)

and

${\mathbb{S}}_{\zeta ,{\zeta }^{\prime },k}={\text{e}}^{ik\cdot \left({\rho }_{\zeta }-{\rho }_{\zeta }\text{'}\right)}{\displaystyle \sum }_{d\in D_{\zeta ,{\zeta }^{\prime }}}{\text{e}}^{ik\cdot d}{\mathbb{S}}_{\zeta ,{\zeta }^{\prime },d},$ (9)

with $D_{\zeta ,{\zeta }^{\prime }}$ being the set of vectors $d=R+{\rho }_{{\zeta }^{\prime }}-{\rho }_{\zeta }$ joining an atom of sublattice $\zeta$ at ${\rho }_{\zeta }$ to each atom of sublattice ${\zeta }^{\prime }$ at $R+{\rho }_{{\zeta }^{\prime }}$,

${\mathbb{𝕙}}_{\zeta ,{\zeta }^{\mathrm{\text{'}}},\mathit{d}}=\underset{-\mathrm{\infty }}{\overset{+\mathrm{\infty }}{\iiint }}{\mathit{\phi }}_{\zeta }^{†}\left(\mathit{r}-{\mathit{\rho }}_{\zeta }\right)\mathrm{}\widehat{H}\mathrm{}{\mathit{\phi }}_{{\zeta }^{\mathrm{\text{'}}}}\left(\mathit{r}-{\mathit{\rho }}_{\zeta }-\mathit{d}\right)\mathrm{}dV$ (10)

and

$${\mathbb{S}}_{\zeta ,{\zeta }^{\prime },d}=\underset{-\infty }{\overset{+\infty }{{\displaystyle {\displaystyle {\iiint }_{}}}}}{\phi }_{\zeta }^{†}\left(r\right){\phi }_{{\zeta }^{\prime }}\left(r-d\right)dV.$$ (11)

The kinetic-energy operator is isotropic, whereas the interaction energy can be approximated by a superposition of isotropic interactions with the atomic cores. In Equation 10, the atomic cores at ${\rho }_{\zeta }$ and ${\rho }_{\zeta }+d$ dominate the contribution of the interaction energy.

For $d=0$, the atomic energy levels are affected by the nearest neighbors. Then, the $4d$ energy level of Mo splits into three sublevels $E_{d0}$, $E_{d1}$ and $E_{d2}$, associated to the orbital subspaces $\left\{d_{z^2}\right\}$, {$d_{xz},d_{yz}\}$ and {$d_{x^2-y^2},d_{xy}$}, respectively. In the same way, the $3p$ energy level of S splits into a couple of sublevels $E_{p0}$, $E_{p1}$ for the orbital subspaces $\left\{p_z\right\}$ and {$p_x,p_y$}, respectively. Additionally, the S atom has the energy level $E_s$ of the $3s$ orbital. This yields six on-site energies for the TBM.

For $d\ne 0$, the integral can be calculated under the assumption that $\widehat{H}$ presents cylindrical symmetry around the line through the atomic positions ${\rho }_{\zeta }$ and ${\rho }_{\zeta }+d$. In this way, one can apply the Slater-Koster procedure to obtain ${\mathbb{𝕙}}_{\zeta ,{\zeta }^{\text{'}},\mathit{d}}$ and ${\mathbb{S}}_{\zeta ,{\zeta }^{\prime },d}$ in terms of $\sigma$, $\pi$ and $\delta$ integrals^14,16. We consider the coupling of each atom to its first nearest neighbors in each crystalline sublattice of the MoS₂ monolayer. Only independent hopping and overlap integrals are listed below.

Each Mo atom couples to the six nearest Mo atoms through three hopping integrals ($V_{dd\sigma }$, $V_{dd\pi }$, $V_{dd\delta }$) and three overlap integrals ($S_{dd\sigma }$, $S_{dd\pi }$, $S_{dd\delta }$). It also couples to the three nearest S atoms at $z=h$ and the three nearest S atoms at $z=-h$. The coupling is given by three hopping integrals ($V_{ds\sigma }$, $V_{dp\sigma }$, $V_{dp\pi }$) and three overlap integrals ($S_{ds\sigma }$, $S_{dp\sigma }$, $S_{dp\pi }$).

The coupling of each S atom to the three nearest Mo atoms is given by three hopping integrals ($V_{sd\sigma }$, $V_{pd\sigma }$, $V_{pd\pi }$) and three overlap integrals ($S_{sd\sigma }$, $S_{pd\sigma }$, $S_{pd\pi }$). Each S atom also couples to its six nearest S atoms on the same plane through four hopping integrals ($V_{ss\sigma },V_{sp\sigma }$, $V_{pp\sigma }$, $V_{pp\pi }$) and four overlap integrals ($S_{ss\sigma },S_{sp\sigma }$, $S_{pp\sigma }$, $S_{pp\pi }$), as well as to the nearest S atom on the other S plane, with four overlap integrals (${\overline{V}}_{ss\sigma },$ ${\overline{V}}_{sp\sigma }$, ${\overline{V}}_{pp\sigma }$, ${\overline{V}}_{pp\pi }$) and five overlap integrals (${\overline{S}}_{ss\sigma },$ ${\overline{S}}_{sp\sigma }$, ${\overline{S}}_{pp\sigma }$, ${\overline{S}}_{pp\pi }$). Since the distance between coplanar S atoms is $a\approx 3.16$ Å, while that between the S planes is $2h\approx 3.12$ Å, we assume that ${\overline{S}}_{sp\sigma }=S_{sp\sigma }$, ${\overline{S}}_{pp\sigma }=S_{pp\sigma }$ and ${\overline{S}}_{pp\pi }=S_{pp\pi }$. This gives $14$ independent hopping integrals: $V_{dd\sigma }$, $V_{dd\pi }$, $V_{dd\delta }$, $V_{ss\sigma },V_{sp\sigma }$, $V_{pp\sigma }$, $V_{pp\pi }$, ${\overline{V}}_{ss\sigma },$ ${\overline{V}}_{sp\sigma }$, ${\overline{V}}_{pp\sigma }$, ${\overline{V}}_{pp\pi }$, $V_{sd\sigma }$, $V_{pd\sigma }$, $V_{pd\pi }$, and $10$ overlap integrals $S_{dd\sigma }$, $S_{dd\pi }$, $S_{dd\delta }$, $S_{ss\sigma },S_{sp\sigma }$, $S_{pp\sigma }$, $S_{pp\pi }$, $S_{sd\sigma }$, $S_{pd\sigma }$, $S_{pd\pi }$. Collecting the on-site energies, the hopping integrals, and the overlap integrals, our TBM deals with $30$ adjustable parameters. The dependence of ${ℍ}_k$ and ${\mathbb{S}}_k$ on those parameters is explained in Appendix 1.

To ensure the existence of a global minimum of $E_k$ as a function of $c_k$ we require ${\mathbb{S}}_k$ to be non-singular²¹. Since neglecting the orbital overlapping would reduce ${\mathbb{S}}_k$ to the identity matrix, the feasible region of the overlap parameters is such that ${\mathbb{S}}_k$ is a positive definite matrix for each $k$. To find the optimal $c_k$ we derive the second line of Equation 5 with respect to $c_k^{†}$. In fact, instead of the real and imaginary parts of $c_k$, the complex variables $c_k$ and $c_k^{†}$ can be treated as independent variables, and the corresponding derivatives of the real-valued function $E_k$ are the complex conjugates of each other. The derivative vanishes whenever¹⁶

${ℍ}_k{c}_k=E_k{\mathbb{S}}_k{c}_k$.(12)

Hence, the energy bands and the coefficients in the LCAO are the solutions of the generalized eigenvalue problem for the pair $\left({ℍ}_k,{\mathbb{S}}_k\right)$. This is solved by using the Eigensystem command in Wolfram Mathematica²², giving 13 energy bands.

The mirror symmetry across the plane $z=0$ is related to the reflection operator $\widehat{\sigma }$ given by $\widehat{\sigma }\psi \left(x,y,z\right)=\psi \left(x,y,-z\right)$. Even and odd Bloch states are eigenvectors of $\widehat{\sigma }$ with eigenvalue $+1$ and $-1$, respectively. The even states will be given by ${\psi }_k^{\mathcal{E}}\left(r\right)={\varphi }_k^{\mathcal{E}}\left(r\right)c_k^{\mathcal{E}}$, with^13,15

$$\begin{array}{l}{\varphi }_k^{\mathcal{E}}=({\varphi }_{\mathrm{0,}{d}_{z^2},k},{\varphi }_{\mathrm{0,}{d}_{x^2-y^2},k},{\varphi }_{\mathrm{0,}{d}_{xy},k},\frac{{\varphi }_{h,s,k}+{\varphi }_{-h,s,k}}{\sqrt{2}},\frac{{\varphi }_{h,p_z,k}-{\varphi }_{-h,p_z,k}}{\sqrt{2}},\\ \frac{{\varphi }_{h,p_x,k}+{\varphi }_{-h,p_x,k}}{\sqrt{2}},\frac{{\varphi }_{h,p_y,k}+{\varphi }_{-h,p_y,k}}{\sqrt{2}})={\varphi }_k\left(r\right){\mathbb{M}}_{\mathcal{E}}\end{array}$$,(13)

with ${\mathbb{M}}_{\mathcal{E}}$ being a $13\times 7$ matrix, whereas the odd states read ${\psi }_k^{\mathcal{O}}\left(r\right)={\varphi }_k^{\mathcal{O}}\left(r\right)c_k^{\mathcal{O}}$, with^13,15

$$\begin{array}{l}{\varphi }_k^{\mathcal{O}}=({\varphi }_{\mathrm{0,}{d}_{xz},k},{\varphi }_{\mathrm{0,}{d}_{yz},k},\frac{{\varphi }_{h,s,k}-{\varphi }_{-h,s,k}}{\sqrt{2}},\frac{{\varphi }_{h,p_z,k}+{\varphi }_{-h,p_z,k}}{\sqrt{2}},\\ \frac{{\varphi }_{h,p_x,k}-{\varphi }_{-h,p_x,k}}{\sqrt{2}},\frac{{\varphi }_{h,p_y,k}-{\varphi }_{-h,p_y,k}}{\sqrt{2}})={\varphi }_k\left(r\right){\mathbb{M}}_{\mathcal{O}},\end{array}$$ (14)

where ${\mathbb{M}}_{\mathcal{O}}$ is a $13\times 6$ matrix. They satisfy the generalized eigenvalue problem

${ℍ}_k^{\mathcal{E},\mathcal{O}}{c}_k^{\mathcal{E},\mathcal{O}}=E_k^{\mathcal{E},\mathcal{O}}{\mathbb{S}}_k^{\mathcal{E},\mathcal{O}}{c}_k^{\mathcal{E},\mathcal{O}},$ (15)

where ${ℍ}_k^{\mathcal{E},\mathcal{O}}={\mathbb{M}}_{\mathcal{E},\mathcal{O}}^{†}{ℍ}_k{\mathbb{M}}_{\mathcal{E},\mathcal{O}}$ and ${\mathbb{S}}_k^{\mathcal{E},\mathcal{O}}={\mathbb{M}}_{\mathcal{E},\mathcal{O}}^{†}{\mathbb{S}}_k{\mathbb{M}}_{\mathcal{E},\mathcal{O}}$. Here ${ℍ}_k^{\mathcal{E}}$ and ${\mathbb{S}}_k^{\mathcal{E}}$ are $7\times 7$ matrices, whereas ${ℍ}_k^O$ and ${\mathbb{S}}_k^O$ are $6\times 6$ matrices. They lead to seven (six) energy bands with even (odd) Bloch functions, provided ${\mathbb{S}}_k^{\mathcal{E}}$ and ${\mathbb{S}}_k^{\mathcal{O}}$ are positive definite matrices²¹ for each $k$.

The tight-binding parameters should fit the DFT results. They depend on the set of DFT $\left(k,E_k\right)$ points involved. We first use the $13$ band energies of Γ, Κ and Μ, i.e., $39$ $\left(k,E_k\right)$ points, and minimize an objective function denoted by $f_{\Gamma KM}$. Then we double the amount of information by including the $13$ band energies at the midpoints of the line segments Λ, Τ and Σ in Figure 2. The new objective function is denoted by $f_{\Gamma \Lambda K{\rm T}M\Sigma }$. Each objective function is the root mean square of the energy differences between the TBM and DFT at the corresponding points. The optimal set of parameters obtained for $f_{\Gamma KM}$ can be used as a starting point for the minimization of $f_{\Gamma \Lambda K{\rm T}M\Sigma }$.

Regarding the tight-binding parameters, we run two configurations: (TBM-20), where the overlap integrals between neighbor atoms are neglected, and (TBM-30) where the full set of parameters is considered. In both cases, additional conditions reduce the number of parameters by two. The energies at the $\text{K}$ point of both the fifth and sixth TBM bands with even Bloch functions are pinned the DFT values. This gives the on-site energies $E_{d0}$ and $E_{d2}$ in terms of the remaining parameters.

The numerical fitting is performed by combining random multi-start and local minimization procedures²³. Only starting points with all the on-site energies between $-20$ eV and $10$ eV, and $E_s$ being the lowest among them, are considered. The local minimization uses the FindMinimum command in Wolfram Mathematica²². In the TBM-30 configuration, each overlap integral $S$ between neighbor atoms should satisfy $-1<S<1$. To simplify this, we map each $S$ into an auxiliary variable $Z$ according to the one-to-one relation $Z=S/\left(1-\left|S\right|\right)$. Nevertheless, the feasibility region for the overlap parameters is given by the restriction that both ${\mathbb{S}}_k^{\mathcal{E}}$ and ${\mathbb{S}}_k^{\mathcal{O}}$ should be positive definite matrices²¹ for each $k$.

The effective-mass values at the selected stationary points are among the main information usually extracted from the band structure. They affect thermal, transport and optical properties of the investigated material. The corresponding tensor is given by

$\frac{{\hslash }^2}{m^{*}}=\nabla \nabla E_k,$ (16)

where the Hessian matrix is the result of applying the dyadic product of $\nabla$ with itself. It can be obtained by implicit derivation of the secular equation

$\det \left({\mathbb{D}}_k\right)=0$,(17)

where ${\mathbb{D}}_k={ℍ}_k-E_k{\mathbb{S}}_k$. Since the partial derivative of $\det \left({\mathbb{D}}_k\right)$ with respect to ${\left({\mathbb{D}}_k\right)}_{j,j^{\prime }}$ equals²⁴ the cofactor ${\left({\mathbb{F}}_k\right)}_{j,j^{\prime }}$ of ${\mathbb{D}}_k$, the gradient of Equation 17 with respect to $k$ reads

${\displaystyle \sum }_{j,{\dot{j}}^{\text{'}}=1}^2{\left({\mathbb{F}}_k\right)}_{j,j^{\prime }}{\left(\nabla {\mathbb{D}}_k\right)}_{j,j\text{'}}=\mathrm{0,}$ (18)

i.e., $\text{Tr}\left({\tilde{F}}_k\nabla {\mathbb{D}}_k\right)=0$, where ${\mathbb{F}}_k$ is the cofactor or adjugate matrix of ${\mathbb{D}}_k$. This leads to

$\nabla E_k=\frac{\text{Tr}\left({\tilde{\mathbb{F}}}_k\left[\nabla {ℍ}_k-E_k\nabla {\mathbb{S}}_k\right]\right)}{\text{Tr}\left({\tilde{\mathbb{F}}}_k{\mathbb{S}}_k\right)}$.(19)

Differentiating again with respect to $k$, and considering a stationary point, Equation 16 becomes

$\frac{{\hslash }^2}{m^{*}}=\frac{\text{Tr}\left(\nabla {\tilde{\mathbb{F}}}_k\left[\nabla {ℍ}_k-E_k\nabla {\mathbb{S}}_k\right]+{\tilde{\mathbb{F}}}_k\left[\nabla \nabla {ℍ}_k-E_k\nabla \nabla {\mathbb{S}}_k\right]\right)}{\text{Tr}\left({\tilde{\mathbb{F}}}_k{\mathbb{S}}_k\right)}.$ (20)

Since ${ℍ}_k$ and ${\mathbb{S}}_k$ are linear combinations of exponentials like $\text{exp}\left(ik\cdot R\right)$ in Equations 8 and 9, their gradients ($\nabla {ℍ}_k$ and $\nabla {\mathbb{S}}_k$) and their Hessians ($\nabla \nabla {ℍ}_k$ and $\nabla \nabla {\mathbb{S}}_k$) can be obtained by substituting each exponential by $iR\text{exp}\left(ik\cdot R\right)$ and $-RR\text{exp}\left(ik\cdot R\right)$, respectively. Moreover, $\nabla {\left({\mathbb{F}}_k\right)}_{j,j^{\prime }}={\left(-1\right)}^{j+j^{\prime }}\text{Tr}\left({\tilde{\mathbb{G}}}_k^{\left(j,j^{\prime }\right)}\nabla {\mathbb{D}}_k^{\left(j,j^{\prime }\right)}\right)$, where ${\mathbb{D}}_k^{\left(j,j^{\prime }\right)}$ is the submatrix of ${\mathbb{D}}_k$ obtained by crossing out row $j$ and column $j\text{'}$, and ${\mathbb{G}}_k^{\left(j,j^{\prime }\right)}$ is the cofactor matrix of ${\mathbb{D}}_k^{\left(j,j^{\prime }\right)}$.

The high-symmetry points $\Gamma$ , $\text{K}$ and $\text{M}$ are stationary points. Because of the three-fold rotational symmetry of the energy bands around $\Gamma$ and $\text{K}$, the effective masses at such points are scalars. In other words, each effective-mass tensor equals a scalar multiple of the identity matrix, i.e., the effective mass is isotropic at these points. Unexpectedly, claims of the effective-mass isotropy at $\text{K}$ being only approximate are available in the literature^15,25. This can be an artifact a sparse $k$ mesh, truncation in data storage, or limited precision of fitting procedures.

At the $\text{M}$ point, the bands present mirror symmetry with respect to the $\Gamma -\text{M}$ line (see Figure 2). There, the principal axes of each energy band are parallel to the Cartesian axes. If the corresponding effective-mass values have the same sign (opposite signs), then $\text{M}$ gives a local extremum (saddle point). The effective mass anisotropy can also be found at other stationary points, such as the $\text{Q}$ point considered by Ridolfi et al.¹⁵.

5. Numerical results and discussions

Figure 3 shows $13$ bands of the electronic structure of the MoS₂ monolayer along the $\Gamma -\text{K}-\text{M}-\Gamma$ path. They have been produced by PWscf¹⁹, as described in Section 3. We have discarded the lowest four bands, which are associated to the $4s$ and $4p$ orbitals of Mo atoms, as well as the highest four bands. The latter ones can be dominated by either the $5s$ and $5p$ orbitals of Mo atoms or the $3d$ orbitals of S atoms. Colors are used to distinguish even from odd Bloch functions with respect to the Mo plane. This classification is done by means of a code developed in Wolfram Mathematica²². It considers the weights of atomic orbitals in the output of the code projwfc.x code provided by Quantum ESPRESSO¹⁹. Even (odd) Bloch functions have no contribution from $d_{xz}$ and $d_{yz}$ ($d_{z^2},d_{x^2-y^2},$ and $d_{xy}$) orbitals of the Mo atoms.

Figure 3
The electronic structure of the MoS₂ as given by the DFT calculations. Blue (red) dots are for Bloch states that are symmetric (antisymmetric) with respect to the Mo plane.

The lowest two bands in Figure 3 are dominated by the $3s$ orbitals of S atoms. The first six bands above them are mainly associated to the $3p$ orbitals of S. The remaining five bands are dominated by the $4d$ orbitals of Mo. A fundamental energy gap of $1.76$eV splits this latter group of bands into one branch of the valence band and four branches of the conduction band. The gap width is near experimental photoluminescence results^15,26, namely close to $1.9$eV. The whole band structure is shifted in energy to put the valence-band top at $E=0$ eV.

The valence and conduction effective-mass values at the $\text{K}$ point, denoted by $m_v$ and $m_c$, respectively, have been extracted from our DFT calculations. This is done by using the same quadratic term but different cubic terms along the $\text{K}-\Gamma$ and $\text{K}-\text{M}$ directions, giving $m_v\approx -0.5986m_0$ and $m_c\approx 0.4691m_0$, where $m_0$ is the bare electron mass.

The DFT results at the $\Gamma$, $\Lambda$ , $\text{M}$, $\text{T}$, $\Sigma$ and $\text{K}$ points associated to the vertical lines in Figure 3 have been fitted by two TBM configurations. Tables 1, 2 and 3 display the parameters that were found for the configurations TBM-20 (neglecting overlaps) and TBM-30 (with finite overlaps). For 20 (30) parameters, the approximate root-mean-square difference with the DFT results at the $\Gamma$, $\Lambda$, $\text{M}$, $\text{T}$, $\Sigma$ and $\text{K}$ points equals $0.34$eV ($0.14$eV).

Table 1
Optimized on-site energies of the tight-binding model. The TBM-20 (TBM-30) columns are for the calculation that neglects (considers) overlap integrals between neighbor atoms.

On-site energy (eV)	TBM-20	TBM-30	On-site energy (eV)	TBM-20	TBM-30
$E_{d0}$	$-0.004$	$-0.392$	$E_s$	$-10.455$	$-11.515$
$E_{d1}$	$0.165$	$-1.740$	$E_{p0}$	$-1.966$	$-2.031$
$E_{d2}$	$0.100$	$-0.536$	$E_{p1}$	$-2.125$	$-2.287$

Table 2
As of Table 1, but for the optimized hopping integrals.

Hopping integral (eV)	TBM-20	TBM-30	Hopping integral (eV)	TBM-20	TBM-30
$V_{dd\sigma }$	$-0.739$	$-0.834$	$V_{sd\sigma }$	$2.405$	$-2.015$
$V_{dd\pi }$	$0.716$	$0.375$	$V_{pd\sigma }$	$2.105$	$2.062$
$V_{dd\delta }$	$-0.065$	$0.011$	$V_{pd\pi }$	$-1.014$	$-1.097$
$V_{ss\sigma }$	$-0.463$	$-0.500$	${\overline{V}}_{ss\sigma }$	$-0.733$	$-0.892$
$V_{pp\sigma }$	$0.768$	$1.206$	${\overline{V}}_{pp\sigma }$	$1.488$	$1.945$
$V_{pp\pi }$	$-0.228$	$-0.217$	${\overline{V}}_{pp\pi }$	$-0.419$	$-0.605$
$V_{sp\sigma }$	$-0.423$	$-0.498$	${\overline{V}}_{sp\sigma }$	$1.451$	$-0.585$

Table 3
As of Table 1, but for the optimized overlap integrals.

Overlap integral	TBM-20	TBM-30	Overlap integral	TBM-20	TBM-30
$S_{dd\sigma }$	0	$0.110$	$S_{sd\sigma }$	0	$0.023$
$S_{dd\pi }$	0	$0.011$	$S_{pd\sigma }$	0	$-0.166$
$S_{dd\delta }$	0	$0.031$	$S_{pd\pi }$	0	$0.104$
$S_{ss\sigma }$	0	$0.013$	$S_{pp\pi }$	0	$-0.019$
$S_{pp\sigma }$	0	$-0.064$	$S_{sp\sigma }$	0	$0.008$

As shown in Table 1, the splitting of the $4d$ Mo energy level given by the TBM-20 configuration fulfills $E_{d0}<E_{d2}<E_{d1}$. This agrees with a simple model of crystal-field theory where the nearest S atoms are substituted by negative point charges^27,28. In the TBM-30 case, the sublevels are ordered as $E_{d1}<E_{d2}<E_{d0}$. Here we note that the orbital overlap between neighbor atoms affects the splitting appreciably. Pasquier and Yazyev²⁹ have reported still another ordering, namely $E_{d2}<E_{d0}<E_{d1}$. This should be investigated in more detail, by improving the optimization procedures, increasing the amount of DFT information to be fitted, and considering more realistic tight-binding models.

Our optimal tight-binding parameters in Tables 1 and 2 do not correlate well with the values found in previous works^13,15. On one hand, our values have the same order of magnitude as those reported by Cappelluti et al.¹³. However, the ordering of on-site energies and the signs of various hopping integrals differ. On the other hand, comparing with the values reported by Ridolfi et al.¹⁵, our S on-site energies are more than $35$ eV higher, and several of our hopping integrals are roughly $10$ times smaller. The cited report prioritizes the five bands dominated by the $4d$ Mo orbitals over the six bands associated to the $3p$ S bands. We also note that the largest overlap integrals in Table 3 are $S_{dd\sigma }$, $S_{pp\sigma }$, $S_{pd\sigma }$ and $S_{pd\pi }$.

Figure 4 displays the TBM and DFT energy bands of the MoS₂ monolayer. In Figure 4a, the 20 optimal parameters roughly fit the DFT energy bands. Neglecting the overlap integrals, the direct gap at the $\text{K}$ point is not well reproduced. Much better agreement is apparent in Figure 4b, where $30$ parameters include finite overlap integrals between nearest neighbors. The two lowest bands, dominated by the $3s$ orbitals of S are very well fitted. The next six bands, associated to the $3p$ orbitals of S, are reasonably well reproduced. The main differences are near the $\text{K}$ point, in the first band with odd Bloch functions. The top five bands, dominated by the $4d$ orbitals of Mo, display good agreement. Noticeable differences are near the $\text{M}$ point and/or in the highest two bands.

Figure 4
Calculated tight-binding (curves) and DFT (dots) energy bands of the MoS₂ monolayer. Blue (red) curves and dots are for even (odd) Bloch states with respect to the Mo plane. Panel (a) [(b)] is for the computational configuration TBM-20 [TBM-30] which neglects [considers] orbital overlap between neighbor atoms.

The curvatures of the valence and conduction bands near their edges seem to be in good agreement. However, the valence and conduction effective masses at the $\text{K}$ point, as given by Equation 20, are $m_v=-0.7557m_0$ and $m_c=0.5637m_0$, respectively. This can be improved, for instance, by fitting the TBM with additional DFT data¹⁵. Moreover, the $5s$ and $5p$ orbitals of Mo atoms, which are not included in our TBM, can interact with the $4d$ orbitals of Mo. This can affect the shape of the conduction bands with a strength that should be assessed in further investigations.

6. Conclusions

A $13$-band tight-binding model has been fitted to a DFT electronic structure of the MoS₂ monolayer. The model was run in two configurations: (i) neglecting the overlap between neighbor atoms and dealing with $20$ parameters, and (ii) considering the overlap, thus increasing the number of parameters to $30$. The first one is more realistic than previous works^13,15 because the $3s$ S orbitals are included, and different hopping integrals are considered for intraplane and interplane couplings between S atoms. It gave a crystal-field splitting of the $4d$ Mo level in agreement with crystal-field theory^27,28 but did not reproduce the DFT bands appropriately. The second configuration fitted the DFT results very well but reversed the sublevel ordering in the crystal-field splitting.

A closed-form expression for the effective-mass tensor at stationary points was derived. Besides giving relevant information accurately, it can be used to adjust the tight-binding parameters whenever effective-mass values appear explicitly in the objective function. In contrast with previous works^15,25, the valence and conduction effective masses at the $\text{K}$ point were shown to be exactly isotropic.

Our tight-binding parameters are not well correlated to the values found in previous reports^13-15. We have used more adjustable parameters, and chose different DFT information, another objective function, and alternate optimization procedures. Our DFT results are in good qualitative agreement with the literature and yields an energy gap close to the experimental values. They have served the purpose of providing enough information to fit the tight-binding energy bands. The present approach is expected to successfully fit more accurate DFT calculations as well. The tight-binding method can be applied to multilayer structures with different stacking patterns³⁰ (1T, 2H, 3R), and to twisted bilayer structures³¹. This requires additional discussions on symmetry as well as the inclusion of new interlayer hopping and overlap matrix elements¹³.

More realistic tight-binding models should be investigated in detail. A small correction can be done by considering different values for the intraplane and interplane overlap integrals between S atoms. Crystal-field theory needs further consideration to better understand the splitting of the Mo atomic level as well as to better describe the effects on the $3s$ and $3p$ levels of S. For instance, the three neighbor Mo atoms of each S atom can induce an on-site $3s-3p_z$ coupling. Considering the lack of rotation symmetry around either the intraplane S-S lines or the Mo-S lines should introduce new hopping integrals between atomic orbitals having different projections of the angular momentum along the line. Including the $5s$ and $5p$ Mo orbitals, and possibly the $3d$ S orbitals in both the DFT and the tight-binding descriptions should give better account of the conduction band. This claim is motivated by a previous study on the electronic structure of the Mo atom that places the $5s$ level below the $4d$ level³². The interatomic couplings can also be extended to the second nearest neighbors. Moreover, there is plenty of room for improvements in the computational methods^13,14, in the selection of what experimental or ab-initio information will be fitted, and in choosing an appropriate weight for each information in the objective function. Efforts along these lines will contribute to a better understanding of the physical properties of the new 2D materials under consideration and should impact the development of theoretical and computational methods for Materials Science.

Supplementary material

The following online material is available for this article:

Appendix 1

7. Acknowledgments

This study was financed in part by the Coordenação de Aperfeiçoamento de Pessoal de Nível Superior - Brasil (CAPES) - Finance Code 001. L.A.M.J. thanks POSMAT for all opportunities. A.B-A. is grateful to the Brazilian Agency CNPq for financial support (Grant # 312770/2017-0). He is also indebted to Prof. Dr. Humberto Rodríguez Gutierrez and late Prof. Dr. Edméa Cássia Baptista for encouraging discussions on transition-metal dichalcogenides and on optimization techniques, respectively.

8. References

• ¹
  Cahangirov S, Sahin H, Le Lay G, Rubio A. Introduction to the physics of silicene and other 2D materials. Berlin: Springer; 2017.
• ²
  Wypych F. Dissulfeto de molibdênio, um material multifuncional e surpreendente. Quim Nova. 2002;25(1):83-8.
• ³
  Sahoo PK, Memaran S, Nugera FA, Xin Y, Díaz Márquez T, Lu Z, et al. Bilayer lateral heterostructures of transition metal dichalcogenides and their optoelectronic response. ACS Nano. 2019;13(11):12372-84.
• ⁴
  Wang P, Yu G, Jia Y, Onyszczak M, Cevallos FA, Lei S, et al. Landau quantization and highly mobile fermions in an insulator. Nature. 2021;589(7841):225-9.
• ⁵
  Thakar K, Lodha S. Optoelectronic and photonic devices based on transition metal dichalcogenides. Mater Res Express. 2020;7(1):014002.
• ⁶
  Vikraman D, Patil SA, Hussain S, Mengal N, Kim HS, Jeong SH, et al. Facile and cost-effective methodology to fabricate MoS₂ counter electrode for efficient dye-sensitized solar cells. Dyes Pigments. 2018;151:7-14.
• ⁷
  Radisavljevic B, Kis A. Mobility engineering and a metal-insulator transition in monolayer MoS₂ Nat Mater. 2013;12(9):815-20.
• ⁸
  Anju S, Mohanan PV. Biomedical applications of transition metal dichalcogenides (TMDCs). Synth Met. 2021;271:116610.
• ⁹
  Costa-Amaral R, Forhat A, Caturello NAMS, Silva JLFD. Unveiling the adsorption properties of 3d, 4d, and 5d metal adatoms on the MoS₂ monolayer: A DFT-D3 investigation. Surf Sci. 2020;701:121700.
• ¹⁰
  Deng X, Li Z. Intrinsic ferromagnetic semiconductivity realized in a new MoS₂ monolayer. Phys Chem Chem Phys. 2020;22(24):13363-7.
• ¹¹
  Gabourie AJ, Suryavanshi SV, Farimani AB, Pop E. Reduced thermal conductivity of supported and encased monolayer and bilayer MoS2. 2D Materials. 2021;8:011001.
• ¹²
  Abdi M, Astinchap B. Thermoelectric properties of monolayer MoS₂ in the presence of magnetic field and electron/hole doping by using the Holstein model. ECS J Solid State Sci Technol. 2020;9(7):073004.
• ¹³
  Cappelluti E, Roldán R, Silva-Guillén JA, Ordejón P, Guinea F. Tight-binding model and direct-gap/indirect-gap transition in single-layer and multilayer MoS₂ Phys Rev B Condens Matter Mater Phys. 2013;88(7):075409.
• ¹⁴
  Silva-Guillén JA, San-Jose P, Roldán R. Electronic band structure of transition metal dichalcogenides from ab initio and Slater–Koster tight-binding model. Appl Sci. 2016;6(10):284.
• ¹⁵
  Ridolfi E, Le D, Rahman TS, Mucciolo ER, Lewenkopf CH. A tight-binding model for MoS₂ monolayers. J Phys Condens Matter. 2015;27(36):365501.
• ¹⁶
  Martin RM. Electronic structure: basic theory and practical methods. Cambridge: Cambridge University Press; 2012.
• ¹⁷
  Landwehr K. Visual discrimination of the 17 plane symmetry groups. Symmetry. 2011;3(2):207-19.
• ¹⁸
  Schattschneider D. The plane symmetry groups: their recognition and notation. Am Math Mon. 1978;85(6):439-50.
• ¹⁹
  Giannozzi P, Baseggio O, Bonfà P, Brunato D, Car R, Carnimeo I, et al. Quantum ESPRESSO toward the exascale. J Chem Phys. 2020;152(15):154105.
• ²⁰
  Quantum Espresso [homepage on the Internet]. 2020 [cited 2020 Oct 5]. Available from: https://www.quantum-espresso.org/pseudopotentials
  » https://www.quantum-espresso.org/pseudopotentials
• ²¹
  Parlett BN. The symmetric eigenvalue problem. Philadelphia: Society for Industrial and Applied Mathematics; 1998.
• ²²
  Wolfram Research, Inc. Mathematica, Version 11.3. Champaign, IL; 2018.
• ²³
  Zhigljavsky A, Zilinskas A. Bayesian and high-dimensional global optimization. Cham: Springer International Publishing; 2021.
• ²⁴
  Horn RA, Johnson CR. Matrix analysis (2nd ed.). Cambridge: Cambridge University Press; 2012.
• ²⁵
  Peelaers H, Van de Walle CG. Effects of strain on band structure and effective masses in MoS₂ Phys Rev B Condens Matter Mater Phys. 2012;86(24):241401.
• ²⁶
  Mak KF, Lee C, Hone J, Shan J, Heinz TF. Atomically thin MoS₂: a new direct-gap semiconductor. Phys Rev Lett. 2010;105(13):136805.
• ²⁷
  Mulak J, Gajek Z. The effective crystal field potential. Amsterdam: Elsevier; 2000.
• ²⁸
  Yan S, Qiao W, He X, Guo X, Xi L, Zhong W, et al. Enhancement of magnetism by structural phase transition in MoS₂ Appl Phys Lett. 2015;106(1):012408.
• ²⁹
  Pasquier D, Yazyev OV. Crystal field, ligand field, and interorbital effects in two-dimensional transition metal dichalcogenides across the periodic table. 2D Materials. 2019;6:025015.
• ³⁰
  Wypych F. Dissulfeto de molibdênio, um material multifuncional e surpreendente: doze anos depois. Quim Nova. 2014;37(7):1220-6.
• ³¹
  Venkateswarlu S, Honecker A, Laissardière GT. Electronic localization in twisted bilayer MoS₂ with small rotation angle. Phys Rev B. 2020;102(8):081103.
• ³²
  Orofino H, Machado SP, Faria RB. The use of Rich and Suter diagrams to explain the electron configurations of transition elements. Quim Nova. 2013;36(6):894-6.

Publication Dates

• Publication in this collection
  01 Nov 2021
• Date of issue
  2021

History

• Received
  20 Jan 2021
• Reviewed
  10 Aug 2021
• Accepted
  23 Sept 2021