Identifying light impurities in transition metal dichalcogenides: the local vibrational modes of S and O in ReSe2 and MoSe2

The transition metal dichalcogenides provide a rich field for the study of two-dimensional materials, with metals, semiconductors, superconductors and charge density wave materials being known. Members of this family are typically hexagonal, but those based on rhenium (ReSe2 and ReS2) and their ternary alloys are attracting attention due to their triclinic structure and their resulting, strong in-plane anisotropy. Here, Raman spectra of dilute ReSe2 - xSx alloys containing low levels of sulfur (x ≤ 0.25) were obtained in order to investigate the distribution of substitutional sulfur atoms over the non-equivalent chalcogen sites of the ReSe2 unit cell. Four different Raman bands arising from the local vibrational modes of sulfur atoms were observed, corresponding to these four sites. One local vibrational mode has a substantially in-plane displacement of the sulfur atom, two are partially out-of-plane and one is completely out-of-plane. The interpretation of the experimental data is based on calculations of the lattice dynamics and non-resonant Raman tensors of a model alloy via density functional theory. For comparison, polarization-dependent Raman spectra of pure ReS2 are also presented; a dramatic increase in the Raman cross-section is found for the out-of-plane modes when the excitation polarization is normal to the layers and the light propagates in the layer plane. A similar increase in cross-section is found experimentally for the local vibrational modes of sulfur in dilute ReSe2 -x Sx alloys and is predicted for dilute sulfur-containing alloys based on MoSe2. The analogous local vibrational modes of substitutional oxygen impurities in ReSe2 were also investigated computationally. The Raman spectra of ReSe2-xSx alloys allows investigation of the distribution of substitutional sulfur atoms within the lattice of ReSe2. A team led by Lewis S. Hart at the University of Bath performed a combination of Raman spectroscopy at visible excitation wavelengths and density functional perturbation theory on ternary ReSe2-xSx alloys obtained by chemical vapor transport. The combination of a traditional normal incidence configuration and an “edge-on” geometry allowed identification of in-plane and out-of-plane Raman features, optically active in the respective experiments owing to their specific selection rules. The frequencies of the observed local vibrational modes arising from sulfur and oxygen atoms were observed to change as a function of the stoichiometric composition of the alloy, thus providing a means of monitoring such impurities. A similar effect was predicted for MoSe2-xSx alloys.


INTRODUCTION
Rhenium disulfide (ReS 2 ) and rhenium diselenide (ReSe 2 ) have attracted considerable recent attention because of their unusual structure. They belong to the transition metal dichalcogenide (TMD) family reported by Wilson and Yoffe 1 but, unlike more wellknown TMDs, their structure is highly anisotropic in the layer plane and this is reflected in all their physical properties. [2][3][4] The class of van der Waals layered semiconductors with in-plane anisotropy is a topic of great current interest and, besides ReS 2 and ReSe 2 , it includes black phosphorus, [5][6][7] GeS, 8 transition metal trichalcogenides 9 and Sb 2 Se 3 . 10 Proposed applications of ReX 2 (X = S, Se) include plasmonic materials, 11 polarization-sensitive photodetectors with high sensitivity, [12][13][14][15] inverters, 16 catalytic devices, 17,18 and few-layer field effect 19,20 or heterojunction 21 transistor structures.
The anisotropy of the ReX 2 compounds arises because, instead of adopting the octahedral (2H) or trigonal (1T) polymorphs typical of TMDs, they undergo a distortion to form a trigonal 1T′ polymorph. 22 The rhenium atoms are bonded to each other to form diamond-like chains, with the direction of the chains defined here as the a direction 23 (Supplementary Fig. S1). The inter-layer coupling is unusually weak compared to other TMDs, [24][25][26] though recent angle-resolved photoemission spectroscopy (ARPES) studies show that the band structure of the bulk ReX 2 family still has a three-dimensional character. [27][28][29][30] ReSe 2 is an indirect band gap semiconductor for all numbers of layers [31][32][33][34] with a band gap of 1.36 eV 34 whilst there is disagreement as to the nature of the band gap for bulk ReS 2 , with some groups reporting an indirect band gap of 1.41 eV and others claiming a direct band gap of 1.5 eV. [31][32][33][34][35][36] In this paper, we focus on the alloys of composition ReSe 2 -x S x . These alloys are of interest for two major reasons. Firstly, by changing the composition of sulfur within these alloys it is possible to tune the band gap, 32,37 as found for other TMD alloys 38 and alloys can be produced over the whole composition range, with the anisotropic trigonal structure being preserved throughout. Secondly, the more well-known TMDs exhibit strong spinorbit splittings (SOS) in the valence band, with the magnitude of this SOS being 180 meV at the K and K' points of the Brillouin zone for bulk MoS 2 . 39 ReSe 2 has space group p1; with inversion symmetry even for monolayers and thus its spin-orbit splitting is expected to be zero throughout the Brillouin zone. However, DFT calculations with and without the inclusion of spin-orbit coupling (SOC) show that the SOC of the rhenium d bands still perturbs the band structure of ReSe 2 significantly. 40 By introducing a sulfur impurity, the center of inversion is removed; as a result of this symmetry breaking, the host crystal field acts as an asymmetric potential and will lift the spin degeneracy at least locally. 41 Rhenium has a high atomic number (higher than molybdenum and tungsten) and thus the spin-orbit splitting is expected to develop rapidly with this symmetry breaking. However, a random distribution of impurity atoms may effectively restore inversion symmetry if the crystal field averages to zero. It is therefore essential to understand whether particular chalcogen sites are preferred and whether non-equilibrium growth techniques or post-growth methods can be used to influence the distribution of sulfur over the possible sites.
Here, we address the key question of how to identify sulfur impurities on the different chalcogen sites in ReSe 2 -x S x (Supplementary Fig. S1) and we show that this is possible via Raman spectroscopy. In addition, we consider whether substitutional oxygen impurities on chalcogen sites may be detected by the same method. Because all ternary alloys in the TMD family are of considerable potential interest, we consider also the case of sulfur impurities in a more representative and conventional TMD, MoSe 2 , for which experimental Raman data already exist, identifying similarities in the structure of the Raman tensors for the impurity local vibrational modes (LVMs) in each case. ReSe 2 has 4 formula units per unit cell (12 atoms) and 36 vibrational normal modes are therefore expected of which 33 have non-zero frequencies at the Brillouin zone center and 18 are Raman-active A g species. Once the center of inversion is eliminated in a random alloy, all 33 modes become Ramanactive in principle. There have been very few Raman studies of sulfur-doped ReSe 2 32,42 though Liu et al. reported the Raman spectrum of the ReSeS ternary alloy. However, with 33 modes spread over a range of 100-500 cm −1 , the analysis of this spectrum is not trivial. 25 Wen et al. reported Raman spectroscopy of ReSe 2 - x S x alloy monolayers suggesting, for low sulfur compositions, that the higher frequency modes disappear. 43 It is shown here that these modes can be observed clearly in dilute alloy samples, and appropriate experimental geometries for this are identified.

RESULTS
Raman scattering of high-frequency vibrational modes of ReS 2 Before considering the alloys, we discuss a useful but uncommon experimental geometry in which the laser excitation propagates in the layer plane. Such "edge-on" spectroscopic experiments have identified out-of-plane vibrational modes in the Raman spectra of other 2D materials 44,45 and the technique can give a better understanding of photoluminescence selection rules. 46 This geometry is particularly valuable for identifying the highestfrequency vibrational mode of ReSe 2 -x S x ; to show this, we consider first the case of pure ReS 2 . The Brillouin zone-center phonon modes of ReS 2 include 18 Raman-active modes, 24,47 several of which can be seen in the experimental data of Fig. 1a.
Their intensities I are determined by Raman tensors R according to I / je i :R:e s j 2 ; (1) where e i and e s are the incident and scattered light polarization vectors respectively. As in the case of ReSe 2 , all Raman-active modes in ReS 2 belong to the non-degenerate symmetry species A g with a Raman tensor R which has only non-zero elements and is given in Cartesian axes in the laboratory frame by We used a polarization-insensitive spectrometer with no analyzer before the detector. 48 Then, with propagation of light along the normal to the layer plane (i.e, along z in the coordinate system of equation 2), the Raman scattering of intensityIa given mode is θ is the polarization angle of the excitation with respect to the inplane x direction, here parallel to the rhenium chains and the crystallographic direction a as shown in Fig. 1b. The use of equation 3 to model the angle-dependent Raman scattering of bulk ReSe 2 and ReS 2 was tested previously. 48,49 For ReS 2 , the highest-frequency predicted mode (438 cm −1 ) has a relatively large component t of the Raman tensor R as defined in equation 2 (atomic displacement patterns and numerical values of R are given in Supplementary Fig. S3). Therefore, when the incident light is polarized perpendicular to the layer plane, the 438 cm −1 Raman mode will be significantly more intense than for polarization in the layer plane. To test this, we use the "edge-on" geometry with the propagation of light along the crystallographic a direction, Fig. 1c. Analogous to equation 3, the intensity of each Raman band with unpolarized detection now follows the relationship for excitation polarized parallel and perpendicular to the layer normal (φ = 0°and φ = 90°respectively). Figure 1a shows the resulting Raman spectra of ReS 2 . The intensities of the lowerfrequency modes (from 100 to 400 cm −1 ) change by less than a factor of two between these two polarizations; however, the intensity of the 438 cm −1 mode increases very strikingly at φ = 0°a nd fitting the experimental data gives an enhancement of Fig. 1 Enhancement of the Raman signals of ReS 2 in the "edge-on" geometry. a Raman spectra of ReS 2 for light polarized parallel (φ = 0°, black line) and perpendicular (φ = 90°, red line) to the normal to the layers z with the direction of propagation along the crystallographic a direction (the direction of the rhenium chains). Inset: experimentally-determined (points) and fitted (solid line) dependence of the intensity of the 438 cm −1 Raman mode of ReS 2 on the angle φ in this geometry. b and c Schematic diagrams of the conventional and "edge-on" Raman experiments respectively; the solid black arrows represent the direction of the incident laser beam (back scattered light is detected in both cases) and the dashed arrows define the angles used to specify the laser polarization I para =I perp 6. The true enhancement in this geometry may be even larger, given that the thickness of the flake is comparable to the laser spot size, so that the effective volume of sample is lower in the "edge-on" geometry. The inset to Fig. 1a shows the measured dependence of the intensity of the 438 cm −1 Raman signal on the incident polarization angle φ. A fit using equation 4 has been added, showing the predicted two-fold symmetry. Supplementary Fig. S3 shows the simulated angle dependence based on the Raman tensor obtained from DFT calculations, which give a predicted value of I para =I perp 14, in reasonable agreement with experiment. Finally, the inset of Fig. 1 shows that the polarization direction at maximum signal is inclined at 6 ± 1°to the normal to the layers; this is also reproduced in the simulated angle-dependence ( Supplementary Fig. S3). These results demonstrate, firstly, that the calculated Raman tensors for the ReS 2 system (and the atomic displacements on which they depend) are reliable and, secondly, that the highestfrequency Raman active phonon of ReS 2 involves the out-of-plane motion of the sulfur atoms, with a strong Raman cross-section for excitation polarized normal to the layers. ReSe 2 shows very similar polarization-dependence for its highest frequency Raman-active mode with predicted and measured I para =I perp 12 and 4, respectively ( Supplementary Fig. S4).
Lattice dynamics of ReSe 1.75 S 0.25 : computational results We now turn to alloys dilute enough to display a sulfur LVM in a ReSe 2 lattice; this facilitates experimental recognition of this mode, and provides an input for modified random element isodisplacement (MREI) models, 50,51 as applied recently to a range of TMDs. [52][53][54][55] Figure 2 (inset) shows the eight chalcogen sites that substitutional sulfur may occupy in the ReSe 2 unit cell. The number of arrangements of r sulfur atoms on the n chalcogen sites of the primitive unit cell (considered as non-equivalent) is given by n C r with n = 8. This number should be halved since each configuration has an experimentally-equivalent one related by inversion (Fig. 2, inset). Therefore, for one sulfur atom per unit cell, composition ReSe 1.75 S 0.25 , the sulfur atom can occupy one of only four distinct chalcogen positions; for two sulfur atoms per unit cell, there are 14 non-equivalent arrangements and for three there are 28. It is feasible but computationally expensive to calculate the phonon modes for all these configurations but at low sulfur concentrations this is unnecessary, as follows. Ho et al. concluded from X-ray studies that sulfur in bulk crystals grown by vapor transport is distributed homogenously. If we assume sulfur is randomly distributed with a constant occupation probability of x/2 for any chalcogen site, then by the binomial theorem, the probability of r sulfur atoms being within a single unit cell is P(r) = 8 C r (1 − x/2) 8-r (x/2) r . Thus, for our experimental concentration of one sulfur per eight formula units, x = 0.125, P(0) = 0.60, P(1) = 0.32, P(2) = 0.07, and P(m) < 0.01 for all m > 2. Thus, the most probable unit cell configuration is pure ReSe 2 and only~8% of unit cells contain more than one sulfur atom.
The lattice modes of ReSe 1.75 S 0.25 were calculated via density functional perturbation theory (DFPT) and Fig. 2 summarizes the results: the frequencies up to the 29th mode are similar to those of ReSe 2 49 and are insensitive to the site occupied by sulfur (changes in frequencies vary by less than typical DFPT accuracy,~2-3%). From the 30th mode upwards, the frequencies diverge from those of ReSe 2 . The most dramatic change is seen at mode 36, which occurs at four remarkably different frequencies: 350, 390, 410 and 433 cm −1 , depending on the site occupied by the sulfur (Fig. 2, right). These frequencies cover a large enough range that Raman spectroscopy should be capable of identifying sulfur in all four sites, in contrast to the case of the binary compound ReSeS, where there are many possible unit cell configurations, so that its Raman bands overlap strongly. The Raman tensors associated with these modes are given in Supplementary Fig. S5, and we compare these to experiment later Importantly, our LVM calculations are not significantly affected by the use of a small unit cell. We have calculated the frequencies of the phonon modes for a 2 × 1 × 1 and a 2 × 2 × 1 supercell, each containing one sulfur atom (compositions ReSe 1.825 S 0.125 and ReSe 1.9375 S 0.0625 ) and find the spread of high-frequency modes is similar to those calculated for ReSe 1.75 S 0.25 (frequencies are tabulated in Supplementary Fig. S5, and compared graphically, Fig. S9). We do see small systematic differences in the frequencies of the spatially-extended phonons as the unit cell size changes ( Supplementary Fig. S9) as expected since the net chalcogen mass in the supercell is changing.
An important observation is that the predicted frequency of the sulfur LVM at site D (433 cm −1 ) is extremely close to the highest- . This is not a coincidence; the atomic displacements are related for both modes and predominantly involve the same chalcogen site D (compare Fig. 2 and Fig S1). This phonon branch in ReS 2 is essentially dispersionless; the calculated phonon dispersion for ReS 2 24 shows that the motions of sulfur atoms in neighboring unit cells are not strongly correlated, so that a large supercell is not required. Another way of demonstrating this is that the displacements of neighboring atoms tend to zero rapidly with distance from the sulfur atom, particularly for the D site ( Supplementary Fig. S8) so that the highfrequency sulfur vibrational modes are localized and independent of concentration.
We conclude that the highest-frequency LVM of sulfur in ReSe 2 is related to the highest-frequency phonon mode of ReS 2 and its predicted frequency is the same to within the typical accuracy of DFPT calculations. Thus, the appearance of this mode in the Raman spectra of ReSe 2 -x S x alloys is not evidence of ReS 2 phase segregation, but does indicate the presence of substitutional sulfur. Similar behavior is expected for light substitutional S (and O) impurities in other TMDs, as exemplified next by the case of MoSe 2 .
Lattice dynamics of MoSe 2 -x S x : computational results We now consider the related question of the LVM arising from sulfur as a light impurity in a more conventional, hexagonal TMD, choosing 2H-MoS x Se 2 -x , for which experimental Raman studies have been reported. 38,52 Here, only one chalcogen site exists in the primitive unit cell, and localized modes can be of both A or B (out-of-plane, non-degenerate) and E (in-plane, two-fold degenerate) symmetry types. 56 Our findings are summarized in Fig. 3; as in the case of ReSe 2 , the highest-frequency mode involves an outof-plane motion of the S impurity whilst other atomic displacements are small. Again, the largest deviation in frequency from the bulk modes of MoSe 2 is seen for the highest-frequency modes, though there are now four modes that are strongly distinguished from the host phonons. The highest-frequency sulfur-related mode is non-degenerate (type A), there is a nearly-degenerate pair of in-plane modes, and there is one lower-frequency nondegenerate mode. The highest-frequency sulfur mode is once more above the top of the MoSe 2 phonon dispersion, forming an out-of-plane LVM.
Experimental studies of 2H-MoS x Se 2 -x showed that the observed phonon frequencies can be described well by an MREI model from which the impurity (local) modes of sulfur in MoSe 2 Fig. 3 Local vibrational modes of substitutional S in MoSe 2 .Top: atomic displacements of the four highest-frequency vibrational modes of monolayer MoSe 2 containing one S impurity atom per 12-atom 2 × 2 × 1 supercell (MoSe 1.75 S 0.25 ), calculated via DFPT (vectors are scaled proportionately to the atomic displacement and are drawn for all atoms but some are too small to be seen). Bottom: vibrational frequencies calculated by DFPT for a single MoSe 2 unit cell, a 2 × 2 × 1 MoSe 2 supercell, and a 2 × 2 × 1 supercell with one Se atom replaced by S. The horizontal dashed lines indicate the symmetry species and the frequencies of the bulk MoSe 2 zone-center modes (adjustable parameters in MREI) were inferred to have frequencies of 343 cm-1 (for the in-plane E 2g mode) and 391 cm-1 (A 1g mode). These values are in good agreement with those we obtain from DFPT as shown on Fig. 3, which are 367 and 417 cm −1 for the Etype modes (34,35) and A-type mode (36), respectively. The agreement is particularly good if we take into account the fact that our calculations for MoSe 2 systematically overestimate the phonon frequencies by about 5%, or~16 cm −1 in this frequency range.
These modes were not observed unambiguously at low concentrations in the previous study, although there is some weak structure in the experimental spectra for x = 0.2 at around 350-360 cm −1 which strengthens and evolves towards the MoS 2 E 2g mode with increasing sulfur content. 52 A difficulty in backscattering experiments is that numerous overtone modes fall within the frequency range of the sulfur LVM, for example E 1g + A 1g and E 2g + LA(M). 52 An experiment in the "edge-on" geometry may be able to resolve this. The calculated Raman tensors for modes 33-36 are given in Supplementary Fig. S14 and, for the highest-frequency LVM, t is~17 times larger than the next largest component of the Raman tensor, so that its "edge-on" signal should be enhanced by up to~17 2 = 270. This agrees qualitatively with data on vertically-aligned MoS 2 layers (for which the backscattering geometry is now "edge-on"), where the intensity of the A 1g mode at 408 cm −1 increased relative to that of the E 2g mode (383 cm −1 ) in the "edge-on" geometry. 45 Experimental Raman spectra of ReSe 1.875 S 0.125 We now test the above predictions experimentally for a dilute alloy and compare to the pure binary materials, taking ReSe 1.875 S 0.125 as a composition that is dilute enough for a true LVM to exist but concentrated enough that signals are easily detected. Measured Raman spectra of ReSe 2 , ReSe 1.875 S 0.125 and ReS 2 are shown in Fig. 4. We note first that we observe bands in the ReSe 2 spectrum which have not previously been reported (seen on the left of Fig. 4b). These lie between 320 and 360 cm −1 and do not correspond to first-order Raman-active modes of ReSe 2 . 48 These bands are most probably due to overtones or combinations of phonons, as are known to be observed in some TMDs 57-59 though, owing to the large number of bands (36) in the phonon dispersion of ReSe 2 , it is difficult to identify the exact branches giving rise to these overtones.
By comparing the alloy spectra to that of ReSe 2 , we observe four peaks which can be attributed to the substitutional sulfur LVM, which we label A-D; the frequencies of these are shown in Table 1. Peak A is predicted from our DFPT calculations to be at 350 cm −1 and there is indeed a peak in the Raman spectra of the alloy near this position (345 cm −1 ); however, this overlaps with a broad band in the ReSe 2 spectrum, and so it is not possible to attribute this peak unambiguously to a LVM. The other three Raman bands of the alloy, B-D, are expected from calculations of their Raman tensors ( Supplementary Fig. S5) to be more intense than band A, and are observed in a region of the Raman spectra where there are no modes of the pure material (Fig. 4a, b), so their Fig. 4 Experimental Raman spectra of ReS 2 , ReSe 2 , and ReSe 1.875 S 0.125 . a Raman spectra for excitation normal to the layers and the incident polarization in-plane and perpendicular to the a direction; b Raman spectra as in a but showing the high frequency region on an expanded scale; c Raman spectra recorded using the "edge-on" geometry with the incident polarisation normal to the layers; d Raman spectra in the "edge-on" geometry as in c in the frequency region of the LVMs. Spectra are displaced vertically for clarity. Vertical arrows in a and c indicate the highest-frequency phonon modes of pure ReSe 2 (296 cm −1 ) and ReS 2 (438 cm −1 ) identification is clear. There is also a strong peak at 322 cm −1 , close to the frequencies of the 34th and 35th alloy modes, which is particularly strong in this geometry. This mode almost certainly arises from the incorporation of sulfur since it is not seen in ReSe 2 but, because the 34th and 35th modes are so close-lying for the different sulfur sites A-D, it is difficult to assign this peak to a specific site or mode with certainty ( Fig. 2 shows that these modes are closely grouped between 310 to 350 cm −1 ). However, the agreement of these observations with DFPT results (shown in Table 1) to well within the limits of error (here, a~3% overestimate) gives us confidence in the assignment of at least three of the four highest-frequency Raman bands to sulfur LVMs. Table 1 includes calculated results with three supercells, giving sulfur concentrations above, below and equal to the experimental concentration; the results do not vary significantly, as expected of an LVM. Supplementary Fig. S9 gives further comparison of the results for different supercells.
The highest-frequency vibration was predicted above to have a large out-of-plane displacement and a large element t in its Raman tensor (Eq. 2). Therefore this mode should give a strong signal for incident light polarized perpendicular to the layer. To test this, Raman spectra were obtained in the "edge-on" geometry as for ReS 2 above; results are shown in Fig. 4c, d). Using this geometry, these modes are clearly identifiable, demonstrating that this experimental technique allows for the identification of the highest vibrational mode of lightly sulfur doped ReSe 2 . The Raman tensors, normal incidence θ-dependence and "edge-on" φ-dependence of bands A-D are presented in Supplementary  Figs. S5, S6 and S7 and, as found experimentally, our calculations predict that the near in-plane mode B is stronger than band D for polarization in the layer plane and normal to a (see Fig. 4b, S7) whilst D is stronger than B for polarization normal to the layers (see Fig. 4d, S5).
One comment concerns the frequencies of the spatially extended phonon modes (below 300 cm −1 ) for the alloy, Fig. 4a, c. Fitting of these spectra reveals that all bands shift to higher frequency by about 1% with respect to the corresponding modes of ReSe 2 , confirming the incorporation of sulfur. This shift is known 43 and, over the whole composition range, is approximately linear. It is reproduced in our calculations and is principally accounted for by the mass difference of the chalcogen atoms 47 though, interestingly, the sulfur site plays a role; the calculated shifts of representative lattice modes (Supplementary Fig. S10) are dependent on the location of sulfur atom (Supplementary Figs. S11 and S12). We chose modes #9 and #10 (124.8 and 159.7 cm −1 in ReSe 2 ) for this, and we used their calculated atomic displacement patterns to identify the analogous modes in the alloy (Supplementary Fig. S10). Thus, interpolation between the binary compounds can only be used to infer composition from alloy mode frequencies if some assumption is made about the distribution of sulfur over sites A-D, such as the random distribution normally assumed in MREI. 52,55 Oxygen impurities in ReSe 2 : computational results It is interesting to ask whether oxygen can form a similar substitutional impurity to sulfur, either via contamination during bulk crystal growth in silica or via atmospheric oxidation. Previous computational studies have investigated oxygen incorporated into ReX 2 , either as substitutional atomic O 60 or molecular O 2 , and have studied the resulting electronic band structure. Here, we investigate what vibrational modes would be expected from substitutional oxygen impurities; this is important for our understanding of the dilute alloys with sulfur in order to establish that the LVMs we observe are not due, for example, to native oxide at the surface. Figure 5 shows that the highest vibrational mode again has four frequencies corresponding to oxygen occupying sites A-D and these lie in a different range to those of sulfur. The displacement patterns are, however, similar to those of sulfur and the Raman tensors again have large components for polarization normal to the layers.
We annealed ReSe 2 in air with times and temperatures following a previous work, 61 where X-ray photoelectron spectroscopy (XPS) revealed surface oxidation. However, we did not observe any oxygen-related modes, presumably because Raman is not as surface-sensitive as XPS and the total quantity of oxide is small. There are reports of the Raman spectra of Re oxides, including partially oxidized ReS 2 containing ReO 3 , 62 pure ReO 3 63,64 and Re 2 O 7 , 50,65,66 but neither of these oxides possesses a configuration analogous to the O atom bridging three Re atoms (site D) that provides the highest frequency LVM here, and neither shows strong Raman bands in the region of 450 to 650 cm −1 . It is not our aim here to investigate all configurations of oxygen in ReSe 2 and, certainly, more substantial structural modifications and also the role of chalcogen vacancies 67,68 should be considered.
In summary, the LVMs observed in dilute sulfur alloys cannot be attributed to oxygen, and we do not detect significant oxidation of micron-scale flakes even after annealing. However, our results provide a means of identifying the specific configuration of a substitutional oxygen atom on a ReSe 2 chalcogen site.

DISCUSSION
The frequencies of the LVMs of substitutional oxygen and sulfur atoms on the four non-equivalent chalcogen sites A-D of ReSe 2 have been identified computationally and those of sulfur have been verified experimentally. The four LVMs of sulfur have remarkably different frequencies, providing a means of monitoring the occupation of these sites in ReSe 2 -x S x alloys. Analogous behavior has been predicted for the cases of oxygen impurities in ReSe 2 and sulfur impurities in MoSe 2 . The relationship between the ReSe 2 -x S x LVMs and the highest-frequency phonon mode of ReS 2 has also been demonstrated. The agreement between experimental and computational results extends also to the enhanced Raman scattering cross-sections when the experimental geometry is chosen so that the light propagates in the layer plane, In a dilute alloy sample, all four LVMs were observed, though the identification of the band arising from site A was tentative since it overlaps an overtone band of ReSe 2 . The observation of all four predicted bands suggests a random distribution of substitutional sulfur impurities in bulk alloys grown by vapor transport.

METHODS
The ReS 2 crystals were supplied by 2D Semiconductors, CA, and have a purity of 99.9995% characterized by secondary ion mass spectroscopy. ReSe 2 samples were supplied by HQgraphene, Netherlands. A sulfurcontaining ReSe 2 sample with a nominal composition of ReSe 1.9 S 0.1 was grown via chemical vapor transport following the technique described by Huang et al. 69 and energy-dispersive X-ray analysis (EDX) indicated a sulfur concentration of~1 sulfur atom per 8 formula units (ReSe 1.875 S 0.125 ), consistent with the nominal composition. EDX results also suggested a small chalcogen deficiency in this sample, as found by Ho et al. 42 From EDX in imaging mode, a uniform distribution of sulfur was found on a macroscopic scale ( Supplementary Fig. S13). Synchrotron-based ARPES studies of the binary compounds and associated XPS measurements 28,29 revealed no significant impurities except for native surface oxide for samples cleaved in air. Samples cleaved in UHV did not show the presence of oxygen. Thick (micron-scale) flakes of the rhenium dichalcogenides were exfoliated using blue Nitto tape and placed onto Gel-Film® WF 6.0 mil in order to perform out-of-plane Raman measurements. To explore the effects of oxygen incorporation, samples were annealed in air following the procedure reported for ReS 2 , 61 but no changes were observed in the Raman spectra. Small changes in the lower-frequency ReS 2 phonons were reported after exposure to O 2 plasma, 62 but this was coupled with the appearance of ReO 3 in XPS spectra, rather than substitutional oxygen. Raman spectra were recorded using a Renishaw inVia Raman microscope with a 50× lens with a spatial resolution of~1 µm. A 532 nm laser was used for all the results presented here and the intensity was kept below 300 µW. No polarization analyzer was used; the detector was previously demonstrated to be polarization insensitive for Raman shifts below 500 cm −1 . 48 For "edge-on" measurements, thick flakes were exfoliated onto Gel-Pak film using micromechanical exfoliation. A single crystal flake partially suspended from the edge of the film was then identified (Supplementary Fig. S2) and the sample was placed upon a support to orient it easily. It would be challenging to observe the signals we discuss below in few-layer structures but suitable structures could be prepared by exfoliation onto graphite, supporting the resulting heterostructure in the same way as above.
Data presented here were obtained with a photon energy of 2.33 eV, above the excitonic band gaps of ReSe 2 and ReS 2 (1.3-1.6 eV, discussed above), so we do not expect strong resonant or optical interference effects in the Raman scattering. Spectra obtained using 1.58 eV excitation (above the band gap of ReSe 2 but close to that of ReS 2 ) did not show significant changes at room temperature though we observe clear resonant effects in the Raman scattering of bulk ReS 2 at low temperatures in the energy range of its excitonic transitions (results will be presented elsewhere).
Phonon modes of all materials were calculated using DFPT using a plane wave basis as implemented in the Quantum Espresso code. 70,71 Calculation of the non-resonant Raman tensors necessitated the use of normconserving pseudopotentials; the valence of Re was taken as 7 with atomic configuration 5d 5 6 s 2 and a local density approximation exchangecorrelation functional was used of the PZ form. 72 Convergence of the total energy was checked with respect to the kinetic energy cutoff (60 Rydberg) and the Monkhorst-Pack k-point grid 73 for a 12-atom unit cell of (8 × 8 × 8). Bulk, three-dimensional material was simulated with three different cell sizes; 12-atoms, 24-atoms and 48-atoms, corresponding to the primitive unit cell and cells doubled along one or both in-plane dimensions. In each case, one S atom was introduced, with separate calculations being made for each of the four non-equivalent chalcogen positions. This gave effective concentrations of twice, equal to, and half of the experimental concentration of S. Initial unit cell parameters and atomic coordinates were taken from Lamfers et al. 23 or Ho et al. 32 with similar results; these were relaxed to obtain atomic forces less than 10 −3 eV/Å. Similar criteria gave good agreement for the phonons of other TMDs. 74,75 The atomic configuration in the vicinity of the S impurity undergoes significant relaxation in order to shorten the Re-S bond distance (e.g, to 2.33 Å compared to a distance of 2.45 Å for a Se atom on the same site). The resulting structure is stable, with no imaginary phonon frequencies. Calculated frequencies were overestimated by up to 2-3% compared to experiment, typical of the above methodology, and confirmed by our experimental results for ReS 2 and ReSe 2 . The crystal structures were drawn using the Jmol or XCrysden packages. 76 It is beyond our present scope to consider the dependence of the phonons on layer thickness, but studies of ReS 2 77,78 and ReSe 2 48,79 as well as Re(S,Se) 2 43 show that the phonons shift only weakly in frequency with thickness down to the monolayer level. For the LVMs that are our focus, we find no significant dependence on unit cell size.

Data availability
Data that supports the findings of this study is available from the University of Bath archive (DOI: 10.15125/BATH-00373)