The optical response of artificially twisted MoS\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$_2$$\end{document}2 bilayers

Two-dimensional layered materials offer the possibility to create artificial vertically stacked structures possessing an additional degree of freedom—the interlayer twist. We present a comprehensive optical study of artificially stacked bilayers (BLs) MoS\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$_2$$\end{document}2 encapsulated in hexagonal BN with interlayer twist angle ranging from 0\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{\circ }$$\end{document}∘ to 60\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{\circ }$$\end{document}∘ using Raman scattering and photoluminescence spectroscopies. It is found that the strength of the interlayer coupling in the studied BLs can be estimated using the energy dependence of indirect emission versus the A\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$_\text {1g}$$\end{document}1g–E\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$_\text {2g}^1$$\end{document}2g1 energy separation. Due to the hybridization of electronic states in the valence band, the emission line related to the interlayer exciton is apparent in both the natural (2H) and artificial (62\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^\circ $$\end{document}∘) MoS\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$_2$$\end{document}2 BLs, while it is absent in the structures with other twist angles. The interlayer coupling energy is estimated to be of about 50 meV. The effect of temperature on energies and intensities of the direct and indirect emission lines in MoS\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$_2$$\end{document}2 BLs is also quantified.

R X±,b mn = a 1 m + a 2 n + t X ± ηe z = R mn + t X ± ηe z , respectively. Here we introduced the in-plane primitive lattice vectors of length a 0 the pair of integer numbers (m, n), the short notation for the (m, n)-th lattice vector R mn = a 1 m + a 2 n and unit vectors of Cartesian coordinate system e x , e y , e z . Vectors dene in-plane positions of the metal and chalcogen atoms within a unit cell of S-TMD monolayer, respectively.
Vectors ±ηe z with η > 0 dene the out-of-plane positions of chalcogen atoms in the unit cell. We also introduce the primitive vectors of reciprocal lattice They satisfy the orthogonality property a j b k = 2πδ jk , where δ jk is the Kronecker delta. The top lattice of the bilayer can be obtained from the bottom one as a result of the shift along z direction on some distance l with subsequent rotation around Oz axis on 180 • degree. Then, the positions of the metal and chalcogen atoms of the top lattice become R M,t mn = R mn + le z + t X , (6) R X±,t mn = R mn + le z + t M ± ηe z .
further we will use the notation ±K both for vectors and the positions of the edges of the Brillouin zone (K ± ), for clarity. In the vicinity of ±K points the Bloch states are predominantly made of the d-orbitals of metal atoms. The corresponding valence and conduction bands states of the bottom layer can be presented as Here N is the normalization factor and Y lm (r − R) is the value of the lm-th atomic orbital placed at the point R and calculated at the point r. The operator C 3 , which generates R 2π/3 rotation of the vectors in space, transforms the corresponding Bloch functions as The phases, which appear under the transformation of the basis states of the bottom layer correspond to the notation in Ref. 1.
In addition, the crystal has inversion symmetry I : r ↔ −r + 2R I , with the center of inversion in the point R I = le z /2. This symmetry together with the time-reversal symmetry induces the restriction on the band structure of the crystal. In accordance to the Kramers theorem, all the bands of the crystal become doubly degenerated by spin.
Therefore, it is convenient to dene the second pair of basis states, associated with the top layer, using the abovementioned symmetry operations Ψ t ±K,n (r) = K 0 IΨ b ±K,n (r). Here K 0 and I are complex conjugation and inversion symmetry operators, respectively. Using the tight-binding representation of the basis states of the bottom layer we The states satisfy the following transformation rules under rotation C 3 Ψ t ±K,v (r) = Ψ t ±K,v (r), C 3 Ψ t ±K,c (r) = e ±2πi/3 Ψ t ±K,c (r), with phases which are opposite to the phases of the basis states of the bottom layer . It leads to the fact that bilayer crystal can absorb the light with both circular polarizations in K point as well as in −K one. This feature is a consequence of inversion symmetry of the crystal. Finally, the mirror symmetry operator P acts on the basis states as P Ψ α ±K,n (r) = [Ψ α ±K,n (r)] * = Ψ α ∓K,n (r), where n = c, v and α = b, t.
In further we focus on the states at the K point for brevity. We take into account the spin degree of freedom s =↑, ↓ of electron excitations and introduce the following set of 8 basis states |Ψ t c , s = Ψ t K,c (r)|s , |Ψ t v , s = Ψ t K,v (r)|s .
According to the kp method, developed for S-TMD multilayers 24 the quasiparticles with the momentum k = k x e x + k y e y at the K point are described by the matrix elements Ψ α n , s| H|Ψ α n , s of the one-particle Hamiltonian Here m 0 is electron's mass, c speed of light, Planck's constant, σ = (σ x , σ y , σ z ) are Pauli matrices and p = −i ∇ is the momentum operator. The rst term of the Hamiltonian denes the kinetic energy of an electron which propagates in the crystal eld of the bottom U b (r) and top U t (r) layers of the bilayer. The next two terms describe the spin-orbital interaction in the system, induced by the potentials U b (r) and U t (r), respectively. The last kp term couples valence and conduction bands. This coupling is supposed to be small and we omit its eects for the current study. The detailed analysis of the impact of the kp term can be found in Refs. 3,4.
We consider rst the matrix elements of the states of bottom layer. We present the Hamiltonian as where we introduced the Hamiltonian of the bottom monolayer and the term which aects the motion of quasiparticles of the bottom layer by the crystal eld of the top layer The Hamiltonian H b 0 (r) has the diagonal matrix elements,which are nothing but the position of the conduction and valence bands in monolayer Here E v and E c are positions of the valence and conduction bands without spin splitting, ∆ v and ∆ c are their spin splittings, and σ s = +1(−1) for s =↑ (↓) states respectively. Note that in K point ∆ v is always positive, while ∆ c can be negative (bright type of S-TMD) and positive (darkish type of S-TMD).
The dominant contribution from of H t int (r) is the diagonal matrix elements in spin space We suppose that the corrections to the splitting are small |∆ c | |δ∆ c |, |∆ v | |δ∆ v |, and the type of S-TMD bilayer remains the same as the type of its constituents. H t int (r) term has also non-zero matrix element between valence and conduction bands of the same layer and opposite spins (see Refs. 5,6 for details). However, these matrix elements give the negligibly small contribution to the energies of the bands, proportional to ∼ 1/(E c − E v ). Therefore, we omit them from the study.
The matrix elements between the states of the top layer can be calculated in the same way. We present the total Hamiltonian as where the rst term aects the motion of quasiparticles of the top layer by the crystal eld of the bottom one. With the help of Kramers theorem and the latter result one can immediately get the answer for the matrix elements of the above-mentioned Note that the sign before spin-splitting terms for the states of the top layer is opposite to the sign of the same terms of the bottom layer. This is the manifestation of the double degeneracy by spin of all the bands of the bilayer.
Finally, we calculate the matrix elements of the Hamiltonian between the states of the dierent layers. Namely, we evaluate the following interlayer matrix elements Ψ t n , s|H(r)|Ψ b n , s , where n, n = c, v. The other matrix elements can be obtained by complex conjugation of the considered ones. We present the Hamiltonian in the following way and suppose the orthogonality of the states from the opposite layers Ψ t n , s|Ψ b n , s = 0. Then Ψ t n , s|H(r)|Ψ b n , s = − 1 2m 0 Ψ t n , s| p 2 |Ψ b n , s .
The p 2 operator is a spin singlet, hence the matrix elements are diagonal in spin space Ψ t n , s| p 2 |Ψ b n , s = δ ss Ψ t n | p 2 |Ψ b n , where |Ψ α n = Ψ α K,n (r). Using transformation properties of the basis states under C 3 rotation we get that only one matrix element is non-zero. The P symmetry of the crystal dictates that the parameter t ⊥ is a real number Im Here we introduced he absolute values of the binding energies of corresponding excitons E A , E B , E A , E B , and we the short notations E 2H The sketch of the bands position in MoS 2 bilayer with 2H-stacking is presented in Fig. 1.
For the particular case of MoS 2 the splitting in conduction band is supposed to be much smaller than the splitting in valence band |∆ 2H c | |∆ 2H v |, and the binding energies of A and B excitons are considered to be equal E A = E B . In this approximation we have the following result ∆ 2H Note that the intralayer and interlayer optical transitions in the same K (or −K) point are characterized by opposite circular polarizations and g-factors of corresponding excitons (see Refs. 3,4 for details).

II.
BILAYER WITH 0 • -ANGLE ALIGNMENT In order to compare the results of the measurements presented in the main text we describe the optical properties of the bilayer S-TMD with zero-angle alignment (in further bilayer) in the similar way as it was done for 2H-stacked bilayer. Again, we consider the bilayer as a pile of two monolayers (top and bottom), placed in parallel to xy plane.
We assume the positions of metal and chalcogen atoms of the bottom layer are the same as in the previous section The top lattice of the bilayer can be obtained from the bottom one as a result of two consequent shifts: along z direction on distance l (which is not equal to the distance l for the 2H-stacked bilayer) and then along y direction on distance a 0 / √ 3. Then, the position of the metal and chalcogen atoms of the top lattice can be presented as R X±,t mn = R mn + le z + t M ± ηe z .
Note that the half of the chalcogen and half of metal atoms in this bilayer are aligned in z-direction. This type of stacking for hexagonal lattices is called Bernal or AB-stacking. The unit cell of the considering bilayer contains twice more atoms than in monolayer. The positions of metal and chalcogen atoms within the unit cell are dened by vectors {t M , a 1 + le z } and {t X ± ηe z , t M ± ηe z + le z }, respectively.
Note that the considering lattice has neither in-plane mirror symmetry (like AA-stacked case) nor inversion symmetry (like 2H-stacked bilayer). It possesses only C 3 rotation symmetry (with Oz line as a rotational axis) and mirror symmetry P : x ↔ −x (the mirror's plane is yz-plane). Again the crystal has the same hexagonal Brillouin zone as the Brillouin zone of the bottom layer. Hence, we choose the same basis states for valence and conduction bands in the ±K points of the bottom layer as we have in the previous section We dene the Bloch states of the top lattice in the same way They have the corresponding transformation rules Like in previous section we focus on the states at the K point, take into account the spin degrees of freedom and introduce the basis states We use the kp approach and calculate the matrix elements Ψ α n , s|H(r)|Ψ α n , s of the one-particle Hamiltonian All the terms in this Hamiltonian has the same meaning as in previous section. We also omit the kp term, since we are focused on the optical transitions exactly in ±K points for clarity.
Let us consider the states of bottom layer. We present the Hamiltonian in the following form where we introduced the Hamiltonian of the bottom monolayer and the term which aects the motion of quasiparticles of the bottom layer by the crystal eld of the top layer H t int (r) = U t (r) + 4m 2 0 c 2 ∇U t (r), p σ.
to the intralayer A and B excitons in the botom (b) and top (t) layers. |∆v + δ∆ b v | (|∆v + δ∆ t v |) and |∆c + δ∆ b c | (|∆c + δ∆ t c |) denote the splitting in the conduction and valence bands of the bottom(top) layer, respectively.
The Hamiltonian H b 0 (r) has the diagonal matrix elements where all the parameters have the same meaning as in the previous section. The dominating contribution from of H t int (r) has also diagonal matrix elements Ψ b c , s|H t int (r)|Ψ b c , s = δE t c + σ s δ∆ t c /2, This term also has non-zero matrix element between valence and conduction bands of opposite spins. As in previous case they give a negligibly small contribution to the energy of the bands, proportional to ∼ 1/(E c − E v ). Therefore, we omit these matrix element from the current study.
The matrix elements between the states of the top layer can be calculated in the same way. Namely we present the Hamiltonian as with the Hamiltonian of the top monolayer H t 0 (r) = p 2 2m 0 + U t (r) + 4m 2 0 c 2 ∇U t (r), p σ, and the term which aects the motion of quasiparticles of the top layer by the crystal eld of the bottom one H b int (r) = U b (r) + 4m 2 0 c 2 ∇U b (r), p σ.
Again, only diagonal matrix elements of the H t 0 (r) are nonzero Ψ t c , s|H t 0 (r)|Ψ t c , s = E c + σ s ∆ c /2,