Enhancement of the WS\(_2\) A\(_{1\text{g}}\) Raman Mode in MoS\(_2\)/WS\(_2\) Heterostructures


Abstract

When combined into van der Waals heterostructures, transition metal dichalcogenide monolayers enable the exploration of novel physics beyond their unique individual properties. However, for interesting phenomena such as interlayer charge transfer and interlayer excitons to occur, precise control of the interface and ensuring high-quality interlayer contact is crucial. Here, we investigate bilayer heterostructures fabricated by combining chemical-vapor-deposition-grown MoS\(_2\) and exfoliated WS\(_2\) monolayers, allowing us to form several heterostructures with various twist angles within one preparation step. In case of sufficiently good interfacial contact – evaluated by photoluminescence quenching - we observe a twist-angle-dependent enhancement of the WS\(_2\) A\(_{1g}\) Raman mode. In contrast, other WS\(_2\) and MoS\(_2\) Raman modes (in particular, the MoS\(_2\) A\(_{1g}\) mode) do not show a clear enhancement under the same experimental conditions. We present a systematic study of this mode-selective effect using nonresonant Raman measurements that are complemented with ab-initio calculations of Raman spectra. We find that the selective enhancement of the WS\(_2\) A\(_{1g}\) mode exhibits a strong dependence on interlayer distance. We show that this selectivity is related to the A\(_{1g}\) eigenvectors in the heterolayer: the eigenvectors are predominantly localized on one of the two layers; yet, the intensity of the MoS\(_2\) mode is attenuated because the WS\(_2\) layer is vibrating (albeit with much lower amplitude) out of phase, while the WS\(_2\) mode is amplified because the atoms on the MoS\(_2\) layer are vibrating in phase. To separate this eigenmode effect from resonant Raman enhancement, our study is extended with near-resonant Raman measurements.

1 Introduction↩︎

Since the discovery of graphene [1], two-dimensional (2D) materials have garnered a lot of scientific interest. Besides graphene, transition metal dichalcogenide (TMDC) monolayers, such as MoS\(_2\) and WS\(_2\) are among the most-investigated 2D materials due to their direct bandgap  [2][5] and tightly bound excitons [6], [7]. In addition to the unique properties of individual 2D materials, new emergent properties are present when they are combined into van der Waals (vdW) heterostructures [8]. Material combination and stacking sequence alter the heterostructure’s electronic properties and are therefore parameters for adapting the heterostructure to a desired application. The relative crystallographic orientation of adjacent layers offers an additional new degree of freedom, which is not accessible in epitaxially grown heterostructures. For example, slightly twisted (so-called magic angle) homobilayers of graphene show exciting transport characteristics including superconductivity [9]. In TMDC heterobilayers, inherent phenomena such as (ultrafast) charge transfer [10] and long-lived interlayer excitons [11][15] have been observed, and using twist angle control, their emission quantum yield [16] and energy [17] can be tuned. MoS\(_2\)/WS\(_2\) heterostructures exhibit a type II band alignment [18]. As for other material combinations, charge transfer in MoS\(_2\)/WS\(_2\) heterostructures requires good interlayer coupling between the constituent layers, which can be achieved by annealing of the sample [19]. Thereby, induced photoluminescence (PL) quenching of the monolayer emission in the heterostructure serves as an indicator for efficient charge transfer [10]. Furthermore, Raman spectroscopy represents a non-invasive tool that has been widely used to characterize the quality of interfacial contact. In the low frequency regime, interlayer shear modes have been observed for MoS\(_2\)/WS\(_2\) heterostructures with high-quality interfaces [20][23]. However, due to the enlarged lattice mismatch for larger twist angles, lack of well-defined interlayer atomic registry and hence the missing of a restoring force, the interlayer shear mode is only expected for stacking angles of 0\(^\circ°\) or 60\(^\circ°\) [24], [25]. A more comprehensive picture, that is also applicable to twisted heterostructures, is obtained by evaluating the behavior of the in- (E\(^1_{2g}\)) and out-of-plane (A\(_{1g}\)) high-frequency Raman modes. Indicating good interlayer contact, the stiffening of the A\(_{1g}\) Raman modes [20], [26], the frequency difference A\(_{1g}\)-E\(^1_{2g}\) [22], [27], [28], the A\(_{1g}\) linewidth [28] as well as the intensity ratio A\(_{1g}\)/E\(^1_{2g}\) [22] have been studied extensively for MoS\(_2\)/WS\(_2\) heterostructures. Besides, an enhancement of the WS\(_2\) A\(_{1g}\) mode in the heterostructure compared to the isolated monolayer has been observed and attributed to a strong interlayer coupling [22], [28]. In general, systematic studies of indicators for high-quality interfacial contact, as well as twist-angle-dependent optical signatures, are hampered by the fact that challenges of controlling interlayer twist angle and coupling during sample fabrication typically lead to large variance from sample to sample.

Here, we investigate MoS\(_2\)/WS\(_2\) heterostructures with various interlayer twist angles by means of optical spectroscopy. Using a hybrid fabrication approach combining exfoliated and high-quality chemical-vapor-deposition-grown (CVD) TMDC monolayers [29], [30], we are able to produce a large set of different twist angles in a single preparation process with well-defined parameters. We find a pronounced, selective enhancement of the WS\(_2\) A\(_{1g}\) mode in heterostructures compared to isolated monolayers. This enhancement strongly depends on the interlayer twist. Remarkably, it is observed for both Stokes and anti-Stokes processes, indicating its nonresonant origin. The latter is further examined by contrasting the nonresonant measurements with wavelength- and temperature-dependent Raman experiments, that reveal a typical resonance behavior enhancing all WS\(_2\) modes. The nonresonant, mode-specific Raman enhancement is confirmed using ab-initio calculations of the nonresonant Raman intensities based on density functional theory (DFT).

2 Results and discussion↩︎

The impact of the interfacial quality in MoS\(_2\)/WS\(_2\) heterostructures on PL and Raman features is illustrated in Figure 1. Before the sample was annealed, almost equal PL intensities of the WS\(_2\) A exciton emission in the heterostructure and isolated WS\(_2\) monolayer regions indicate only weak interaction between the heterostructure’s individual layers (Fig. 1a). The small spectral shift observed for the A exciton emission can be attributed to the different dielectric environments provided by the SiO\(_2\) and MoS\(_2\), respectively. In contrast, a pronounced quenching of the WS\(_2\) emission is observed in the heterostructure region after annealing. Interestingly, Raman spectra measured at the same sample position also show significant differences (Fig. 1b). Whereas initially the WS\(_2\) out-of-plane A\(_{1g}\) Raman mode intensity does not differ between the WS\(_2\) monolayer and the heterostructure region, a strong enhancement of this mode occurs after annealing. This suggests that the enhancement of the WS\(_2\) A\(_{1g}\) mode serves as an additional indicator for high-quality MoS\(_2\)/WS\(_2\) heterostructure interfaces. In the following, we further characterize the effect and reveal its microscopic origin.

Figure 1: (a) Room temperature PL spectra of a MoS_2/WS_2 heterostructure. After annealing, PL quenching of the WS_2 A exciton emission occurs in the heterostructure region, indicating good contact between the constituent layers. Corresponding PL maps are shown in Supplementary Information S1. (b) Initially, Raman spectra of the same spots do not show intensity differences for the WS_2 A_{1g} mode in the monolayer (red line) and heterostructure region (blue area). However, the WS_2 A_{1g} mode is enhanced after annealing. Raman spectra are normalized to the intensity of the Si phonon on the bare Si/SiO_2 substrate. Both Raman and PL spectra were acquired with an excitation wavelength of 532 nm.

Figure 2a shows MoS\(_2\)/WS\(_2\) heterostructures that were fabricated by placing a large exfoliated WS\(_2\) monolayer on top of CVD-grown MoS\(_2\) monolayers (see methods). This yields the advantage of creating several heterostructures with various twist angles within one preparation step. Here, determination of the interlayer twist angle is enabled by optical microscopy. We analyzed the angle between a CVD triangle edge and the long straight WS\(_2\) monolayer edge, both assumed to follow high-symmetry crystal directions. Given the fact that slight variations of preparation parameters can alter a heterostructure’s optoelectronic properties, this method also facilitates a direct comparison between the resulting heterostructures since all structures were produced under the same lab conditions and faced the same annealing conditions. Improvement of the interfacial quality caused by annealing was verified by observation of PL quenching (Supplementary Information S2).

Figure 2: (a) Optical microscope image of a sample consisting of CVD-grown MoS_2 monolayers which were covered with a large exfoliated WS_2 flake. This yields different heterostructures with various twist angles in the same sample. Some CVD-grown monolayers were marked for illustration purposes. (b) Raman scan (532 nm excitation) of the same sample showing the intensity of the WS_2 A_{1g} Raman mode. (c) Averaged Raman spectra of the marked heterostructure with 0^\circ° twist angle. Compared to other Raman modes, the WS_2 A_{1g} is clearly enhanced. Averaged spectra for all heterostructure regions are shown in Supplementary Information S3. (d) The enhancement appears on the Stokes and Anti-Stokes side of the Raman spectrum. The spectrum was measured in the center of the 0^\circ° twist angle heterostructure. (e,f) Maximum and average enhancement factors of the WS_2 A_{1g} mode obtained for heterostructures with different twist angles. For comparison, 2LA(M)/E^1_{2g} average enhancement factors are included, underlining the mode selectivity of the enhancement process. Error bars represent statistical uncertainties determined from all selected individual WS_2 spectra.

The sample was mapped using nonresonant Raman spectroscopy at 532 nm excitation. The corresponding Raman intensity map of the WS\(_2\) A\(_{1g}\) mode reveals a higher intensity in almost all heterostructure regions compared to the WS\(_2\) monolayer, which proves the robustness of the effect (Fig. 2b). Interestingly, a significantly stronger increase is found for stacking angles close to 0\(^\circ°\). Raman spectra belonging to the 0\(^\circ°\) heterostructure marked in (a) were extracted from the scan data, and then summed and divided by the number of selected spectra. The resulting averaged spectrum is depicted in Figure 2c and compared to WS\(_2\) and MoS\(_2\) monolayer spectra that were obtained in the same way. Its characteristic features are the out-of-plane A\(_{1g}\) and in-plane E\(^1_{2g}\) optical modes of both WS\(_2\) and MoS\(_2\) as well as the longitudinal acoustic 2LA(M) mode of WS\(_2\) [31]. For simplicity, we consider the WS\(_2\) 2LA(M) (352 cm\(^{-1}\)) and E\(^1_{2g}\) (356 cm\(^{-1}\)) mode as one combined peak.

In contrast to the WS\(_2\) A\(_{1g}\) mode, the MoS\(_2\) mode intensities as well as the combined WS\(_2\) 2LA(M)/E\(^1_{2g}\) peak do not differ from their monolayers’ counterparts. In the following, we thus primarily focus on the WS\(_2\) A\(_{1g}\) mode behavior. Figure 2d shows that the effect is also apparent on the anti-Stokes side of the Raman spectrum. Together with the absence of enhancement of the WS\(_2\) combined 2LA(M)/E\(^1_{2g}\) peak, this implies that the detected enhancement does not result from a resonant Raman process, as a laser excitation energy matching a WS\(_2\) monolayer optical transition (incoming resonance) would affect all WS\(_2\) modes marked in the spectrum [32].

Figure 3: Calculated Enhancement factors of WS_2 and MoS_2 (blue and red) A_{1g} and E_g (solid and hollow points) Raman modes in WS_2/MoS_2 heterostructures and H_h^h homobilayers as a function of the interlayer distance. The range of equilibrium interlayer distances for each heterostructure is marked by a gray vertical bar. The insets depict the side views of atomic structures of the heterobilayers, where larger (smaller) blue (red) circles represent metal (sulfur) atoms in WS_2 (MoS_2) layer, and the gray dashed lines connect atoms that are vertically aligned in the given stackings.

For further twist angle-dependent analysis, Raman spectra of the different heterostructures and the surrounding WS\(_2\) monolayer were selected from the Raman scan shown in Figure 2b. The WS\(_2\) A\(_{1g}\) mode intensity of every spectrum was then obtained by numerical integration over the energetic region of interest. The maximum enhancement factor (Fig. 2e) is defined as \(I_{HS}/I_{\overline{ML}}\), where \(I_{HS}\) is the highest WS\(_2\) A\(_{1g}\) mode intensity value determined in the respective heterostructure and \(I_{\overline{ML}}\) the average monolayer WS\(_2\) A\(_{1g}\) mode intensity. In Figure 2f, the average enhancement factor \(I_{\overline{HS}}/I_{\overline{ML}}\) for various heterostructures is depicted in dependence of the twist angle. We observe a reduced enhancement for twist angles deviating from 0\(^\circ°\), reaching no observable enhancement close to 30\(^\circ°\), followed by a slight increase up to 60\(^\circ°\). This trend resembles previously reported results [28] for twisted WS\(_2\)/MoS\(_2\) heterostructures where both layers were obtained from CVD growth. Given that the interlayer distance in twisted WS\(_2\)/MoS\(_2\) depends on the twist angle with a smallest layer separation at 0\(^\circ°\) and 60\(^\circ°\) [28], [33], [34], these results suggest that the enhancement is sensitive to the interlayer distance. Moreover, we assume that the different stacking orders of the involved atoms, namely R-type (0\(^\circ°\)) and H-type (60\(^\circ°\)) stacking [33], [35], explain the different enhancement factors at comparable interlayer distance. Average enhancement factors for the WS\(_2\) 2LA(M)/ E\(^1_{2g}\) peaks occurring in the same heterostructures were calculated accordingly. To this end, we integrated over the combined WS\(_2\) 2LA(M)/ E\(^1_{2g}\) peak and the shoulder at 325 cm\(^{-1}\), since the phonon modes overlap in this region. Contrary to the WS\(_2\) A\(_{1g}\) mode, the WS\(_2\) 2LA(M)/E\(^1_{2g}\) peaks do not show any pronounced intensity increase which underlines the mode selectivity of the enhancement process.

We note that the enhancement is not limited to MoS\(_2\)/WS\(_2\) heterostructures consisting of both exfoliated and CVD-grown monolayers, but also occurs in purely exfoliated samples (Supplementary Information S4). Interestingly, we also observe an enhancement of the WS\(_2\) A\(_{1g}\) mode in MoSe\(_2\)/WS\(_2\) heterostructures (Supplementary Information S5).

We complement the nonresonant Raman measurements in MoS\(_2\)/WS\(_2\) heterostructures with DFT calculations. As MoS\(_2\) and WS\(_2\) have nearly identical lattice constants, their heterostructures exhibit significant atomic reconstructions. Instead of moiré patterns, large domains with high-symmetry stackings emerge, as has been observed in MoS\(_2\)/WS\(_2\) and MoSe\(_2\)/WSe\(_2\) heterostructures [24], [35], [36]. In the following, the stacking configurations are labeled according to Yu et al[37], [38]. Samples with twist angles near 0\(^\circ°\) (60\(^\circ°\)) consist mostly of \(R^X_h\) and \(R^M_h\) (\(H^h_h\) and \(H^X_h\)) domains, with characteristic vertical alignment of atoms as depicted by insets in Figure 3. We neglect the \(R^h_h\), \(H^M_h\) and intermediate stacking registries, as they cover a minor part of the sample and their contribution to the Raman signal is negligible. The stackings are modeled within primitive cells and start from finding their equilibrium interlayer distances \(d_{eq}\), which is defined as the vertical distance between Mo and W atomic planes. Then we calculate the Raman intensities \(I=I_{xx}+I_{yy}\) of the MoS\(_2\) and WS\(_2\) A\(_{1g}\) and E\(_g\) modes for a series of interlayer distances \(d\). Other components of the Raman tensor are omitted, since they are not probed in the experimental back-scattering measurements. Figure 3 presents the enhancement factors (EFs) of these modes in function of \(d\), defined as: \[EF=\frac{I(d)}{I(d_{1L})}. \label{EF}\tag{1}\] Here, \(I(d_{1L})\) is the Raman intensity of the respective mode in the monolayer. Distinctly, in all the considered stackings, the enhancement factor of the WS\(_2\) A\(_{1g}\) mode increases when the layers get closer to each other, and reaches \(\approx 2\) at \(d_{eq}\). When \(d\) is further decreased by \(\approx 10\%\) with respect to \(d_{eq}\), the enhancement factor of the WS\(_2\) A\(_{1g}\) mode rapidly grows, exceeding \(\approx 3.5\) in \(R^M_h\). It qualitatively agrees with experimental findings, and supports the hypothesis that annealing leads to a reduction of interlayer distance in MoS\(_2\)/WS\(_2\) samples, which is revealed by the enhancement of the WS\(_2\) A\(_{1g}\) Raman mode intensity. The enhancement factors of other modes exhibit various dependence on \(d\), depending on the stacking. However, their values at \(d_{eq}\) and lower do not exceed 2, which is in line with experimental observations. For comparison, we calculated and plotted the enhancement factors of the discussed modes for MoS\(_2\) and WS\(_2\) homobilayers in \(H^h_h\) stacking, which corresponds to their bulk crystals. Interestingly, intensities of all the modes are enhanced with similar trends, but to a lower extent than in MoS\(_2\)/WS\(_2\) heterostructures.

In order to understand the selective enhancement of the aforementioned Raman modes, let us look at the quantities that contribute to the off-resonant Raman intensity calculated within Placzek’s approximation [39]. For the back-scattering measurement geometry, the Raman intensity (for light scattering between Cartesian directions \(i\) and \(j\)) of the phonon mode \(\nu\) is \[I_{ij}^\nu \sim \left| \sum\limits_{k,n} \alpha'_{ijkn} u_{kn}^\nu \right|^2 \label{raman95formula}\tag{2}\] with \(\alpha'_{ijkn}= \frac{\partial \chi_{ij}}{\partial \tau_{kn}}\) the atom-resolved Raman polarizability defined as derivative of \(\chi_{ij}\), the macroscopic dielectric susceptibility, with respect to \(\tau_{kn}\), the displacement of the \(n\)th atom in the Cartesian direction \(k\). The eigen-displacement of phonon mode \(\nu\) of atom \(n\) in direction \(k\) is denoted by \(u_{kn}\). The prefactor cancels out in the calculation of EFs, and thus is omitted in Equation 2 for brevity.

Let us discuss the \(R^X_h\) heterostructure as a representative case, and WS\(_2\) and MoS\(_2\) \(H^h_h\) homobilayers as reference. The atom-resolved Raman polarizabilities and eigen-displacements of the WS\(_2\) and MoS\(_2\) A\(_{1g}\) and E\(_g\) phonon modes are presented in Figure S6 as a function of interlayer distance. We focus on the non-zero components that contribute to \(I_{xx}\): From symmetry it follows that \(I_{yy}=I_{xx}\) and \(I_{xy}=0\) for the A\(_{1g}\) out-of-plane modes where \(u_{xn}=u_{yn}=0\). For the E\(_{g}\) in-plane modes, \(u_{zn}=0\). As discernible in Figure S6, the eigen-displacements and Raman polarizabilites are more sensitive to the interlayer distance for A\(_{1g}\) modes than for E\(_g\) modes. Furthermore, for the A\(_{1g}\) modes, the interlayer dependence of the \(u_z\) vibration amplitude is more pronounced for the sulfur atoms of the WS\(_2\) layer than for those of the MoS\(_2\) layer.

But the decisive point for the selective enhancement of the WS\(_2\) A\(_{1g}\) mode is the relative phase of the sulfur atom vibrations on both layers. At large interlayer distance, the vibrations are exclusively localized on the WS\(_2\) or MoS\(_2\) layer, respectively. However, with decreasing interlayer distance, the interlayer coupling leads to a “minority vibration” of the sulfur atoms on the adjacent layer. At equilibrium, the amplitude ratio is \(1:7\) and rapidly increases upon further squeezing of the layers. As can be seen from the orientation of the triangles in Figure S6, for the WS\(_2\) A\(_{1g}\) mode, the vibrations of sulfur atoms (albeit with very different amplitudes) are in phase, while they are out of phase for the MoS\(_2\) A\(_{1g}\) mode. In phase vibration leads to two contributions with the same sign in the summation of Equation 2 and hence to a visible enhancement of the WS\(_2\) mode while out of phase vibration leads to a partial cancellation in the summation and a reduction of the MoS\(_2\) mode intensity. We verified this effect by finite-displacement calculations of the Raman tensor: if only the sulfur atoms on the WS\(_2\) (MoS\(_2\)) layer are displaced, the intensities of both A\(_{1g}\) modes are slightly enhanced (as compared to the isolated single layer). Including the constructive (destructive) interference from the “minority vibration”, one reproduces the experimentally observed and calculated selective mode enhancement of the WS\(_2\) mode. The strong distance dependence of the “minority” amplitude explains then the observed strong distance dependence of the enhancement factor for the WS\(_2\) A\(_{1g}\) mode.

It should be noted that our theoretical reasoning is at the off-resonance limit and neglects resonant effects, which can be significant in layered structures, as shown in [40][42]. However, resonant Raman calculations from first principles are beyond the scope of this study.

In order to gain further experimental insight into the WS\(_2\) A\(_{1g}\) mode’s behavior and to separate the previously described eigenmode effect from resonant effects, we also performed near-resonant and resonant Raman measurements. Resonant conditions were achieved by (i) sweeping the laser excitation energy and (ii) tuning the WS\(_2\) and MoS\(_2\) bandgaps via temperature variation.

For the former, the laser excitation wavelength was tuned below and above the WS\(_2\) A exciton resonance energy, using a Radiant Dye laser (see methods). Spectra at various excitation energies were taken on the WS\(_2\) monolayer and on different heterostructures of the sample shown in Figure 2a. After PL background subtraction, all spectra were normalized to the Si phonon intensity. Intensities of the Raman modes of interest were obtained by numerical integration. The resulting phonon intensities are depicted in Figures 4a and 4b, exemplary Raman spectra above and below the WS\(_2\) A exciton resonance are shown in Figures 4c and 4d.

Figure 4: (a,b) Integrated Raman intensities for the WS_2 out-of-plane and in-plane Raman modes in a WS_2 monolayer (red) and several heterostructures (blue). Raman spectra were obtained from various spots on different heterostructures (symbols) shown in Figure 2a, triangles indicate data originating from the 0^\circ° heterostructure. Data points marked with square and circle represent measurements from the 9^\circ° and 55^\circ° heterostructure, respectively. A laser of tunable wavelength was used for excitation. The lines are fitted resonance curves for the isolated monolayer and the 55^\circ° heterostructure. (c,d) Exemplary Raman spectra excited below (1.98 eV) and above (2.07 eV) the WS_2 A exciton energy. The energy range used for intensity determination via numerical integration is highlighted in gray.

The data points (symbols) constitute the tails of typical Raman resonance curves [32]. Below the A exciton resonance (i.e., at lower excitation energies), larger intensities are detected for both the out-of-plane \(A_{1g}\) Raman mode and the combined 2LA(M)/ E\(^1_{2g}\) peak at all heterostructure spots compared to the isolated monolayer. In contrast, higher intensities occur for the monolayer’s modes above the WS\(_2\) A exciton energy. With increasing deviation from the WS\(_2\) A exciton energy, the intensities of monolayer and heterostructure converge. At the resonance, no data points could be recorded since the available dyes did not allow stable laser emission in the relevant energy range. Resonance profiles (lines) were subsequently generated by fitting the Raman intensities using third-order perturbation theory [40], [43], [44]:

\[I_R(\omega_{ph},E_l) \propto \left| \frac{M_{op}^2 \cdot M_{ep}}{(E_l - E_a +i\gamma_a)(E_l-\hbar\omega_{ph}-E_a+i\gamma_a)} \right|^2\]

We combined the matrix elements \(M_{op}\) (electric dipole transition) and \(M_{ep}\) (electron-phonon interaction) into one fitting parameter. Additional parameters include the energy \(E_a\) of the intermediate resonant state, which is the WS\(_2\) A exciton, and the decay rates \(\gamma_a\) of the incoming and outgoing resonance. As common, the same value \(\gamma_a\) was assumed for both resonance events [43]. \(E_l\) is the laser excitation energy and \(\hbar\omega_{ph}\) the respective phonon energy. In general, the energy of the A exciton in the heterostructure regions is redshifted relative to that in the monolayer [45]. The exemplary resonance profiles (lines) for one heterostructure spot and the isolated WS\(_2\) monolayer clearly reflect this trend. Compared to the monolayer, the resonant condition in the heterostructure regions already sets in at lower excitation energies. The intensities of both out-of-plane and in-plane modes follow the resonance curve for all heterostructure positions. We attribute the rather broad peak for the A\(_{1g}\) intensity profile in the monolayer region to an overlap of incoming and outgoing resonance, which is more distinctly resolved for the in-plane modes. The data were obtained for heterostructures with 0\(^\circ°\), 9\(^\circ°\) and 55\(^\circ°\) twist angle. Interestingly, in none of these heterostructures does the 2LA(M)/ E\(^1_{2g}\) peak show any distinct enhancement at 532 nm excitation (Supplementary Information S3). Therefore, the observed increase of both out-of-plane and in-plane Raman modes - which is typical resonant behavior as previously observed in WS\(_2\) monolayers [32] - provides evidence that different enhancement mechanisms dominate under near-resonant and nonresonant excitation conditions.

Further temperature-dependent measurements were performed on the sample introduced in Figure S4. Figure 5a shows temperature-dependent Raman spectra at 532 nm laser excitation. All spectra were normalized to the Si phonon intensity of each individual spectrum. As the temperature rises, both the WS\(_2\) A\(_{1g}\) and the combined WS\(_2\) 2LA(M)/ E\(^1_{2g}\) mode intensity increase, with the effect starting to become pronounced at about 160 K. Similar temperature-dependent behavior was reported for WS\(_2\) monolayers at 514.5 nm (2.41 eV) excitation and explained by resonance with the WS\(_2\) B exciton [46].

Despite the lower excitation energy in our experiment (2.33 eV), we still observe the onset of a resonant process due to temperature-dependent modifications in the WS\(_2\) bandstructure in both monolayer and heterostructure. To illustrate this effect, we performed white-light reflectance measurements on the MoS\(_2\)/WS\(_2\) heterostructure and the isolated WS\(_2\) monolayer (Fig. 5b,c). Here, several relevant aspects are apparent: First, the characteristic absorption features of WS\(_2\) A and B exciton are broader and redshifted in the heterostructure compared to the monolayer; and second, in both monolayer and heterostructure a broadening and redshift of the absorption dips occur with increasing temperature. Thus, resonance with the WS\(_2\) B exciton is established in the heterostructure region close to room temperature. Consequently, we observe a pronounced enhancement of both in- and out-of-plane Raman modes compared to the isolated monolayer (Supplementary Information S4). Remarkably, the WS\(_2\) A\(_{1g}\) mode intensity is more affected than the combined 2LA(M)/ E\(^1_{2g}\) mode intensity. Analogous to the the calculations for the heterostructures presented in Figure 2, we extract from the room temperature Raman scan (Fig. S4) an average enhancement factor of 4.55 \(\pm\) 0.46 for the WS\(_2\) A\(_{1g}\) mode and an average enhancement factor of 1.487 \(\pm\) 0.045 for the combined WS\(_2\) 2LA(M)/ E\(^1_{2g}\) peak. This is in contrast to previously reported resonant Raman profiles for WS\(_2\) monolayers, where the enhancement at the WS\(_2\) B exciton energy is slightly smaller for the WS\(_2\) A\(_{1g}\) than for the WS\(_2\) 2LA(M) and E\(^1_{2g}\) modes [32], which indicates that we observe a superposition between resonant enhancement and the mode-selective effect occurring in MoS\(_2\)/WS\(_2\) heterostructures. We note that, with both WS\(_2\) peaks enhanced at room temperature compared to the isolated monolayer, the sample’s spectral behavior differs from that of the zero-degree twist angle heterostructure spectrum shown in Figure 2c. However, the comparability of various samples with regard to their resonance energy is complex due to multiple factors. The exciton binding energy - and hence the A exciton resonance - is sensitive to dielectric screening which is governed by the interlayer distance [47]. The latter in turn is highly affected by contaminants introduced during sample preparation such as hydrocarbons [48] and PDMS residues [49]. Furthermore, band gaps can be modified by strain occurring in the heterostructure region [50][53]. Since the energetic difference between A and B exciton is expected to remain almost constant for distinct samples [54], these factors have a direct impact on the resonance condition discussed above. The WS\(_2\) A\(_{1g}\) mode enhancement in MoS\(_2\)/WS\(_2\) heterostructures is thus driven by a complex interplay of interlayer distance between the constituent layers and resonant conditions that depend on the individual heterostructure’s properties.

Figure 5: (a) Temperature-dependent Raman spectra of a MoS_2/WS_2 heterostructure at 532 nm excitation. All spectra are normalized to the Si phonon. Raman modes of interest are the WS_2 A_{1g} mode at 420 cm^{-1} and the combined 2LA(M)/ E^1_{2g} peak at 355 cm^{-1}. (b,c) Low-temperature white-light reflectance contrast (RC) measurements of the same sample. RC=(R_{Sample}-R_{Ref})/R_{Ref} was used for normalization, with R_{Ref} being the reflectance of the bare Si/SiO_2 substrate. Furthermore, all spectra were normalized to the minimum of the WS{_2} A exciton for better illustration of energetic shifts. The dotted line indicates the WS_2 monolayer’s A and B exciton energies at 4 K.

In summary, we have fabricated vdW heterostructures by combining CVD grown MoS\(_2\) monolayers with large exfoliated WS\(_2\) monolayers. This approach allows us to achieve comparable preparation conditions and facilitates comparison between various samples. Under nonresonant conditions, we observe a twist-angle-dependent mode-selective enhancement of the WS\(_2\) A\(_{1g}\) mode which serves as an easily accessible indicator for high-quality interfacial contact. DFT calculations reveal that it originates from an in-phase oscillation of the MoS\(_2\) sulfur atoms with those in WS\(_2\) for the A\(_{1g}\) mode displacement. Near-resonant Raman measurements, realized by varying the excitation wavelength or the temperature, demonstrate that the mode-selective enhancement is independent of resonant effects. Our results highlight the complex interplay of phonon modes across the van der Waals gap.

3 Methods↩︎

3.1 Sample fabrication↩︎

CVD-grown [30] MoS\(_2\) monolayers were picked up from the growth substrate by capillary-force assisted transfer [55]. For this, we either used deionized water vapor to wet a PDMS stamp or a deionized water droplet that was carefully placed on the growth substrate prior to pick-up. MoS\(_2\) monolayers were then transferred to a Si/SiO\(_2\) substrate by viscoelastic stamping [56]. Several MoS\(_2\) monolayer triangles were covered with a large WS\(_2\) monolayer, that was previously mechanically exfoliated onto PDMS, in a second transfer step. All exfoliated flakes were obtained from bulk crystals by HQ graphene. To increase interlayer contact, samples were annealed in high vacuum at 420 K for 3 hours.

3.2 Optical measurements↩︎

3.2.1 Raman and Photoluminescence Measurements↩︎

For PL and nonresonant Raman measurements, a continuous-wave 532 nm diode-pumped solid state laser was used for excitation. The laser beam was focused on the sample through a 100x (room temperature) or a 80x (low temperature) microscope objective. The scattered/emitted light was collected through the same objective. Mounted on a motorized xy-stage, the sample was moved with \(\mu\)m precision under the fixed laser spot position. A set of three successive Bragg filters suppressed Rayleigh-scattered light in front of a spectrometer equipped with a CCD sensor operated at -70\(^\circ°\)C. For low-temperature Raman measurements, the sample was placed in a flow cryostat cooled with liquid helium.

For (near-)resonant Raman measurements the sample was excited with a wavelength-tunable Radiant Dye Laser, using two dyes (DCM and R6G). With an excitation power of 1.1 mW the laser beam was focused through an 100x microscope objective under which the sample was placed. The scattered light was collected using the same objective and guided into a Horiba T64000 micro-spectrometer with triple-grating configuration. The Raman signal was dispersed by a grating with 900 grooves per mm, and a CCD sensor was used for detection. All (near-)resonant Raman measurements were performed at room temperature.

3.2.2 White Light Reflectance Measurements↩︎

White-light reflectance contrast measurements were performed using a quartz tungsten halogen lamp. A collimated beam was focused on the sample by a 80x objective. The reflected light was collected through the same objective and was guided into a spectrometer equipped with a CCD detector. Enhanced spatial resolution was achieved by a spatial-filtering module in the detection path (see [57] for details). The sample was mounted in a flow cryostat cooled with liquid helium.

3.3 Theoretical calculations↩︎

Density functional theory calculations were performed with the Quantum Espresso package [58], [59] within local density approximation (LDA) to exchange-correlation functional and norm-conserving scalar-relativistic pseudopotentials of version 1.2 [60]. The energy cut-off for the wave functions was 90 Ry. The Brillouin zone was sampled with a \(\Gamma\)-centered Monkhorst-Pack grid of 12\(\times\)​12\(\times\)​1 k-points. The in-plane and out-of-plane lattice constants of 3.1270 Å and 40 Å were used for all the structures. The former represents the average of experimental values for WS\(_2\) and MoS\(_2\). The latter yields an effective vacuum distance between repeated hetero-layers of 20 Å, which is sufficient to reproduce the phonons of isolated hetero-layers, as verified by convergence tests. The atomic positions were optimized until all the forces acting on the atoms were less than \(10^{-6}\) Ry/Bohr. When varying the interlayer distance, positions of metal atoms were kept fixed during the optimization. The phonons at \(\Gamma\) were calculated using density functional perturbation theory (DFPT). The nonresonant Raman intensities were evaluated with DFPT, using the implementation of Lazzeri and Mauri[61]. For the proper treatment of 2D boundary conditions, the Coulomb cut-off technique was used in all calculations [62].

4 Acknowledgements↩︎

T.K. acknowledges financial support by the DFG via the following grants: SFB1477 (project No. 441234705) and SPP2244 (project No. 443361515). L.W. acknowledges funding by the FNR Luxembourg through project C22/MS/17415967/ExcPhon. The DFT calculations were carried out using the HPC facilities of the University of Luxembourg. PK and OG acknowledge the DFG for funding (KU4034 2-1). Z.G., A.G. and AT acknowledge financial support of the DFG through SPP2244 "2DMP" (Projects TU149/13-1 and TU149/21-1) and CRC 1375 "NOA" (Project B2).

5 Supplementary experimental data↩︎

Figure 6: (a) Optical Microscope image of the MoS_2/WS_2 sample on which the spectra shown in Figure 1 (main text) were measured. CVD-grown MoS_2 monolayers were covered with a large exfoliated WS_2 monolayer, yielding several heterostructure regions. The dots indicate the measurement positions of the heterostructure (black) and WS_2 monolayer (yellow) spectra presented in the main text. (b) PL maps (acquired with an excitation wavelength of 532 nm) show the spatial distribution of the WS_2 A exciton emission. PL quenching is observed after annealing, indicating an improvement of the heterostructure’s interface.
Figure 7: (a) Optical microscope image and (b) PL maps of the MoS_2/WS_2 sample introduced in Figure 2 in the main text. Good interfacial contact between the heterostructures’ layers was established by annealing and verified by the observation of quenching of the WS_2 A-Exciton emission. A 532 nm laser was used for excitation.
Figure 8: All averaged Raman spectra (532 nm excitation wavelength) for the heterostructures in the sample presented in the main text (Figure 2).
Figure 9: (a) Optical microscope image of a MoS_2/WS_2 heterostructure where both layers were obtained by exfoliation from a bulk crystal. (b) The fluorescence microscope shows quenching in the heterostructure region where an enhancement of the WS_2 A_{1g} mode is observed in the corresponding Raman map (c). (d) Room temperature averaged Raman spectra of MoS_2 monolayer, WS_2 monolayer and heterostructure. Raman spectra were measured using an excitation wavelength of 532 nm.
Figure 10: Characterization of a MoSe_2/WS_2 heterostructure. (a) Optical microscope image of the sample, consisting of an exfoliated MoSe_2 monolayer (bottom) and an exfoliated WS_2 monolayer (top). (b) PL map showing the WS_2 A exciton emission after preparation. While some heterostructure regions ((1),(2)) lack good interfacial contact, others already show PL quenching ((3), (4)). Thus, the sample was not annealed. (c) Raman map displaying the spatial distribution of the WS_2 A_{1g} mode intensity. Enhancement occurs in regions with pronounced PL quenching. (d) Averaged Raman spectra of heterostructure regions. Average enhancement factors of 3.34 \pm 0.53 (3) and 2.26 \pm 0.30 (4) are obtained. Maximum detected enhancement factors are 4.73 (3) and 4.26 (4). PL and Raman spectra were measured using an excitation wavelength of 532 nm.

6 Supplementary theoretical data↩︎

Figure 11: Calculated non-zero contributions to I_{xx} components of Raman intensity tensor for A_{1g} and E_g phonon modes of MoS_2 and WS_2 layers (top row) in R^X_h MoS_2/WS_2 heterostructure and (bottom row) in H^h_h MoS_2 and WS_2 bilayers: u_z(y), the out-of(in)-plane atoms eigendisplacements and \alpha'_{xxz(y)}, the atomic-resolved Raman polarizabilities. Blue (red) symbols and lines correspond to WS_2 (MoS_2) layers: dots (triangles) are used for metal (sulfur) atoms with upper (lower) triangles for upper (lower) sulfur atoms in a layer. Horizontal dashed lines mark some of the separated layers values for reference. Insets present the atomic displacements for the considered phonon modes.

References↩︎

[1]
K. S. Novoselov, A. K. Geim, S. V. Morozov, D. Jiang, Y. Zhang, S. V. Dubonos, I. V. Grigorieva, and A. A. Firsov, https://doi.org/10.1126/science.1102896.
[2]
K. F. Mak, C. Lee, J. Hone, J. Shan, and T. F. Heinz, https://doi.org/10.1103/physrevlett.105.136805.
[3]
A. Splendiani, L. Sun, Y. Zhang, T. Li, J. Kim, C.-Y. Chim, G. Galli, and F. Wang, https://doi.org/10.1021/nl903868w.
[4]
W. Zhao, Z. Ghorannevis, L. Chu, M. Toh, C. Kloc, P.-H. Tan, and G. Eda, https://doi.org/10.1021/nn305275h.
[5]
P. Tonndorf, R. Schmidt, P. Böttger, X. Zhang, J. Börner, A. Liebig, M. Albrecht, C. Kloc, O. Gordan, D. R. T. Zahn, S. Michaelis de Vasconcellos, and R. Bratschitsch, https://doi.org/10.1364/oe.21.004908.
[6]
A. Molina-Sánchez, K. Hummer, and L. Wirtz, https://doi.org/10.1016/j.surfrep.2015.10.001.
[7]
A. Chernikov, T. C. Berkelbach, H. M. Hill, A. Rigosi, Y. Li, B. Aslan, D. R. Reichman, M. S. Hybertsen, and T. F. Heinz, https://doi.org/10.1103/physrevlett.113.076802.
[8]
A. K. Geim and I. V. Grigorieva, https://doi.org/10.1038/nature12385.
[9]
Y. Cao, V. Fatemi, S. Fang, K. Watanabe, T. Taniguchi, E. Kaxiras, and P. Jarillo-Herrero, https://doi.org/10.1038/nature26160.
[10]
X. Hong, J. Kim, S.-F. Shi, Y. Zhang, C. Jin, Y. Sun, S. Tongay, J. Wu, Y. Zhang, and F. Wang, https://doi.org/10.1038/nnano.2014.167.
[11]
H. Fang, C. Battaglia, C. Carraro, S. Nemsak, B. Ozdol, J. S. Kang, H. A. Bechtel, S. B. Desai, F. Kronast, A. A. Unal, G. Conti, C. Conlon, G. K. Palsson, M. C. Martin, A. M. Minor, C. S. Fadley, E. Yablonovitch, R. Maboudian, and A. Javey, https://doi.org/10.1073/pnas.1405435111.
[12]
M. Palummo, M. Bernardi, and J. C. Grossman, https://doi.org/10.1021/nl503799t.
[13]
S. Gao, L. Yang, and C. D. Spataru, https://doi.org/10.1021/acs.nanolett.7b04021.
[14]
R. Gillen and J. Maultzsch, https://doi.org/10.1103/PhysRevB.97.165306.
[15]
E. Torun, H. P. C. Miranda, A. Molina-Sánchez, and L. Wirtz, https://doi.org/10.1103/PhysRevB.97.245427.
[16]
P. K. Nayak, Y. Horbatenko, S. Ahn, G. Kim, J.-U. Lee, K. Y. Ma, A.-R. Jang, H. Lim, D. Kim, S. Ryu, H. Cheong, N. Park, and H. S. Shin, https://doi.org/10.1021/acsnano.7b00640.
[17]
J. Kunstmann, F. Mooshammer, P. Nagler, A. Chaves, F. Stein, N. Paradiso, G. Plechinger, C. Strunk, C. Schüller, G. Seifert, D. R. Reichman, and T. Korn, https://doi.org/10.1038/s41567-018-0123-y.
[18]
K. Kośmider and J. Fernández-Rossier, https://doi.org/10.1103/physrevb.87.075451.
[19]
S. Tongay, W. Fan, J. Kang, J. Park, U. Koldemir, J. Suh, D. S. Narang, K. Liu, J. Ji, J. Li, R. Sinclair, and J. Wu, https://doi.org/10.1021/nl500515q.
[20]
J. Zhang, J. Wang, P. Chen, Y. Sun, S. Wu, Z. Jia, X. Lu, H. Yu, W. Chen, J. Zhu, G. Xie, R. Yang, D. Shi, X. Xu, J. Xiang, K. Liu, and G. Zhang, https://doi.org/10.1002/adma.201504631.
[21]
M. Okada, A. Kutana, Y. Kureishi, Y. Kobayashi, Y. Saito, T. Saito, K. Watanabe, T. Taniguchi, S. Gupta, Y. Miyata, B. I. Yakobson, H. Shinohara, and R. Kitaura, https://doi.org/10.1021/acsnano.7b08253.
[22]
Y. Saito, T. Kondo, H. Ito, M. Okada, T. Shimizu, T. Kubo, and R. Kitaura, https://doi.org/10.35848/1347-4065/ab9400.
[23]
K. H. Shin, M.-K. Seo, S. Pak, A.-R. Jang, and J. I. Sohn, https://doi.org/10.3390/nano12091393.
[24]
J. Holler, S. Meier, M. Kempf, P. Nagler, K. Watanabe, T. Taniguchi, T. Korn, and C. Schüller, https://doi.org/10.1063/5.0012249.
[25]
C. H. Lui, Z. Ye, C. Ji, K.-C. Chiu, C.-T. Chou, T. I. Andersen, C. Means-Shively, H. Anderson, J.-M. Wu, T. Kidd, Y.-H. Lee, and R. He, https://doi.org/10.1103/physrevb.91.165403.
[26]
K.-G. Zhou, F. Withers, Y. Cao, S. Hu, G. Yu, and C. Casiraghi, https://doi.org/10.1021/nn5042703.
[27]
L. Liang and V. Meunier, https://doi.org/10.1039/c3nr06906k.
[28]
L. Wu, C. Cong, J. Shang, W. Yang, Y. Chen, J. Zhou, W. Ai, Y. Wang, S. Feng, H. Zhang, Z. Liu, and T. Yu, https://doi.org/10.1007/s12274-020-3193-y.
[29]
S. Shree, A. George, T. Lehnert, C. Neumann, M. Benelajla, C. Robert, X. Marie, K. Watanabe, T. Taniguchi, U. Kaiser, B. Urbaszek, and A. Turchanin, https://doi.org/10.1088/2053-1583/ab4f1f.
[30]
A. George, C. Neumann, D. Kaiser, R. Mupparapu, T. Lehnert, U. Hübner, Z. Tang, A. Winter, U. Kaiser, I. Staude, and A. Turchanin, https://doi.org/10.1088/2515-7639/aaf982.
[31]
A. Berkdemir, H. R. Gutiérrez, A. R. Botello-Méndez, N. Perea-López, A. L. Elías, C.-I. Chia, B. Wang, V. H. Crespi, F. López-Urías, J.-C. Charlier, H. Terrones, and M. Terrones, https://doi.org/10.1038/srep01755.
[32]
E. del Corro, A. Botello-Méndez, Y. Gillet, A. L. Elias, H. Terrones, S. Feng, C. Fantini, D. Rhodes, N. Pradhan, L. Balicas, X. Gonze, J.-C. Charlier, M. Terrones, and M. A. Pimenta, https://doi.org/10.1021/acs.nanolett.5b05096.
[33]
W. Yang, H. Kawai, M. Bosman, B. Tang, J. Chai, W. L. Tay, J. Yang, H. L. Seng, H. Zhu, H. Gong, H. Liu, K. E. J. Goh, S. Wang, and D. Chi, https://doi.org/10.1039/c8nr07498d.
[34]
M. Tebyetekerwa, J. Zhang, S. E. Saji, A. A. Wibowo, S. Rahman, T. N. Truong, Y. Lu, Z. Yin, D. Macdonald, and H. T. Nguyen, https://doi.org/10.1016/j.xcrp.2021.100509.
[35]
M. R. Rosenberger, H.-J. Chuang, M. Phillips, V. P. Oleshko, K. M. McCreary, S. V. Sivaram, C. S. Hellberg, and B. T. Jonker, https://doi.org/10.1021/acsnano.0c00088.
[36]
X. Wang, K. Yasuda, Y. Zhang, S. Liu, K. Watanabe, T. Taniguchi, J. Hone, L. Fu, and P. Jarillo-Herrero, https://doi.org/10.1038/s41565-021-01059-z.
[37]
H. Yu, G.-B. Liu, J. Tang, X. Xu, and W. Yao, https://doi.org/10.1126/sciadv.1701696.
[38]
H. Yu, G.-B. Liu, and W. Yao, https://doi.org/10.1088/2053-1583/aac065.
[39]
G. Placzek, in Handbuch der Radiologie. Band 6, Teil 2, edited by E. Marx(Akademische Verlagsgesellschaft, Leipzig, 1934) pp. 209–374.
[40]
O. Garrity, T. Brumme, A. Bergmann, T. Korn, P. Kusch, and S. Reich, https://doi.org/10.1021/acs.nanolett.4c02757.
[41]
S. Reichardt and L. Wirtz, https://doi.org/10.1126/sciadv.abb5915.
[42]
M. Nalabothula, L. Wirtz, and S. Reichardt, https://doi.org/10.1021/acs.nanolett.5c00355.
[43]
A. Jorio, R. Saito, G. Dresselhaus, and M. S. Dresselhaus, in Raman Spectroscopy in Graphene Related Systems(Wiley, 2011) pp. 103–119.
[44]
C. Hyeonsik and J.-U. Lee, in Raman Spectroscopy of Two-Dimensional Materials, edited by P.-H. Tan(Springer, 2019) pp. 185–201.
[45]
A. F. Rigosi, H. M. Hill, Y. Li, A. Chernikov, and T. F. Heinz, https://doi.org/10.1021/acs.nanolett.5b01055.
[46]
X. Huang, Y. Gao, T. Yang, W. Ren, H.-M. Cheng, and T. Lai, https://doi.org/10.1038/srep32236.
[47]
M. Florian, M. Hartmann, A. Steinhoff, J. Klein, A. W. Holleitner, J. J. Finley, T. O. Wehling, M. Kaniber, and C. Gies, https://doi.org/10.1021/acs.nanolett.8b00840.
[48]
S. J. Haigh, A. Gholinia, R. Jalil, S. Romani, L. Britnell, D. C. Elias, K. S. Novoselov, L. A. Ponomarenko, A. K. Geim, and R. Gorbachev, https://doi.org/10.1038/nmat3386.
[49]
A. Jain, P. Bharadwaj, S. Heeg, M. Parzefall, T. Taniguchi, K. Watanabe, and L. Novotny, https://doi.org/10.1088/1361-6528/aabd90.
[50]
K. He, C. Poole, K. F. Mak, and J. Shan, https://doi.org/10.1021/nl4013166.
[51]
H. J. Conley, B. Wang, J. I. Ziegler, R. F. Haglund, S. T. Pantelides, and K. I. Bolotin, https://doi.org/10.1021/nl4014748.
[52]
G. Plechinger, A. Castellanos-Gomez, M. Buscema, H. S. J. van der Zant, G. A. Steele, A. Kuc, T. Heine, C. Schüller, and T. Korn, https://doi.org/10.1088/2053-1583/2/1/015006.
[53]
X. He, H. Li, Z. Zhu, Z. Dai, Y. Yang, P. Yang, Q. Zhang, P. Li, U. Schwingenschlogl, and X. Zhang, https://doi.org/10.1063/1.4966218.
[54]
B. Zhu, X. Chen, and X. Cui, https://doi.org/10.1038/srep09218.
[55]
X. Ma, Q. Liu, D. Xu, Y. Zhu, S. Kim, Y. Cui, L. Zhong, and M. Liu, https://doi.org/10.1021/acs.nanolett.7b03449.
[56]
A. Castellanos-Gomez, M. Buscema, R. Molenaar, V. Singh, L. Janssen, H. S. J. van der Zant, and G. A. Steele, https://doi.org/10.1088/2053-1583/1/1/011002.
[57]
S. Deb, J. Krause, P. E. Faria Junior, M. A. Kempf, R. Schwartz, K. Watanabe, T. Taniguchi, J. Fabian, and T. Korn, https://doi.org/10.1038/s41467-024-52011-3.
[58]
P. Giannozzi, S. Baroni, N. Bonini, M. Calandra, R. Car, C. Cavazzoni, D. Ceresoli, G. L. Chiarotti, M. Cococcioni, I. Dabo, A. Dal Corso, S. de Gironcoli, S. Fabris, G. Fratesi, R. Gebauer, U. Gerstmann, C. Gougoussis, A. Kokalj, M. Lazzeri, L. Martin-Samos, N. Marzari, F. Mauri, R. Mazzarello, S. Paolini, A. Pasquarello, L. Paulatto, C. Sbraccia, S. Scandolo, G. Sclauzero, A. P. Seitsonen, A. Smogunov, P. Umari, and R. M. Wentzcovitch, http://www.quantum-espresso.org.
[59]
P. Giannozzi, O. Andreussi, T. Brumme, O. Bunau, M. B. Nardelli, M. Calandra, R. Car, C. Cavazzoni, D. Ceresoli, M. Cococcioni, N. Colonna, I. Carnimeo, A. D. Corso, S. de Gironcoli, P. Delugas, R. A. D. Jr, A. Ferretti, A. Floris, G. Fratesi, G. Fugallo, R. Gebauer, U. Gerstmann, F. Giustino, T. Gorni, J. Jia, M. Kawamura, H.-Y. Ko, A. Kokalj, E. Küçükbenli, M. Lazzeri, M. Marsili, N. Marzari, F. Mauri, N. L. Nguyen, H.-V. Nguyen, A. O. de-la Roza, L. Paulatto, S. Poncé, D. Rocca, R. Sabatini, B. Santra, M. Schlipf, A. P. Seitsonen, A. Smogunov, I. Timrov, T. Thonhauser, P. Umari, N. Vast, X. Wu, and S. Baroni, http://stacks.iop.org/0953-8984/29/i=46/a=465901.
[60]
D. R. Hamann, https://doi.org/10.1103/PhysRevB.88.085117.
[61]
M. Lazzeri and F. Mauri, https://doi.org/10.1103/PhysRevLett.90.036401.
[62]
T. Sohier, M. Calandra, and F. Mauri, https://doi.org/10.1103/PhysRevB.96.075448.