Amorphous As2S3 Doped with Transition Metals: An Ab Initio Study of Electronic Structure and Magnetic Properties

Crystalline transition-metal chalcogenides are the focus of solid state research. At the same time, very little is known about amorphous chalcogenides doped with transition metals. To close this gap, we have studied, using first principle simulations, the effect of doping the typical chalcogenide glass As2S3 with transition metals (Mo, W and V). While the undoped glass is a semiconductor with a density functional theory gap of about 1 eV, doping results in the formation of a finite density of states (semiconductor-to-metal transformation) at the Fermi level accompanied by an appearance of magnetic properties, the magnetic character depending on the nature of the dopant. Whilst the magnetic response is mainly associated with d-orbitals of the transition metal dopants, partial densities of spin-up and spin-down states associated with arsenic and sulphur also become slightly asymmetric. Our results demonstrate that chalcogenide glasses doped with transition metals may become a technologically important material.


Introduction
The discovery of semiconducting properties of chalcogenide glasses by Goryunova and Kolomiets has laid the foundation of a new class of semiconductors, viz. amorphous semiconductors [1]. Due to the presence of lone-pair electrons, i.e., paired p-electrons that do not participate in the formation of covalent bonds but occupy the top of the valence band [2], these materials exhibit a number of photo-induced phenomena, such as reversible photostructural changes, photo-induced anisotropy, etc. (for reviews, see [3,4]). In these processes, photo-induced bond modification is caused by bond switching that involves excitation of lone-pair electrons and subsequent formation of a new bonding configuration [5].
On the other hand, transition-metal chalcogenides contain Ch-TM (Ch = chalcogen, TM = transition metal) bonds that are formed by sharing chalcogen lone-pair electrons with empty d-orbitals of a transition metal, which can lead to the appearance of hybridized states [6]. This process results in (i) a decrease in the ability of chalcogenide glasses to undergo photostructural changes and should also result in (ii) an appearance of magnetism in doped glasses. It would also be extremely interesting if one could change the local bonding configuration of TM-doped chalcogenide glasses using electronic excitation in a way that would affect its magnetic properties. This would allow one to design glasses with optically switchable magnetic response. As the first step to addressing this issue, we studied, using first principles simulations, the electronic structure and magnetic properties

Simulation Details
The aim of this work was to study theoretically the amorphous magnetism [20][21][22][23] of a prototypical chalcogenide glass a-As 2 S 3 doped with transition metals (Mo, W and V) by means of computer simulations of electronic and magnetic properties within density functional theory (DFT).
The amorphous phase generation procedure used in this study is a "standard" procedure to obtain in silico a melt-quenched amorphous phase of chalcogenides, which has been used in various publications [24,25]. The idea behind this approach is (i) to randomize the structure at a high temperature, (ii) to equilibrate the structure above the melting point of a material and (iii) to quench the melt in order to obtain the amorphous (glassy) phase.
The amorphous phase was generated with the help of ab initio molecular dynamics (AIMD) simulations via the following procedure. As the initial structure we used, the 2 × 2 × 3 supercell of the orpiment [26] (a = 22.95 Å; b = 19.15 Å; c = 12.77 Å; α = γ = 90 • ; β = 90.35 • ) consisting of 240 atoms (48 formula units As 2 S 3 ) with a density of 3.49 g/cm 3 . The initial structure was randomized using an NVT-ensemble with a Nosé thermostat by heating the structure to 3000 K during 20 ps followed by cooling to 900 K (a temperature just above the melting point) over a period of 15 ps. The amorphous phase was generated by quenching from 900 K to 300 K over a duration of 15 ps. The times for each stage we used are similar to the times used by others in the field [24,25,27]. A similar procedure was repeated for Mo doped a-As 2 S 3 . A total of 5 Mo atoms were included in the 240 atom Mo-doped cell corresponding to the Mo concentration of 2%.
The generalized gradient approximation (GGA) for the exchange-correlation functional by Perdew-Burke-Ernzerhof (PBE) [28] along with the Projected Augmented Wave (PAW) [29,30] method to describe the electron-ion interaction as implemented in the plane wave pseudopotential code VASP [31][32][33] were used. A value of 260 eV was chosen for plane wave kinetic energy cutoff E cutoff . Integration over Brillouin zone (BZ) was accomplished with only Γ-point taking into account the relatively large size of the simulation cell.
To address the issue of the effect of different transition metals on the structure and properties of the glass, we used a simplified approach adopting the following strategy. Since the glass formation process is stochastic and the obtained a-As 2 S 3 :Mo structure is just one of the many possible structures, we replaced Mo species by W and V (assuming that such local structures will statistically be formed), after which the structures were relaxed at 0 K. Full geometry optimization (GOpt) at 0 K of all melt-quenched (mq) amorphous structures, both pure and doped, was accomplished via the Broyden-Fletcher-Goldfarb-Shanno (BFGS) [34][35][36][37] optimizer using spin-polarized DFT within the GGA with the exchange-correlation functional PW91 by Perdew and Wang [38] and the van der Waals (vdW) dispersion correction by Tkatchenko and Scheffler (TS) [39], as implemented in the plane-wave pseudopotential code CASTEP [40,41]. The Vanderbilt ultrasoft pseudopotentials (USPs) [42] were chosen to describe the electron-ion interactions. As: 4s 2 4p 3 , S: 3s 2 3p 4 ; Mo: 4s 2 4p 6 4d 5 5s 1 ; V: 3s 2 3p 6 3d 3 4s 2 ; and W: 5s 2 5p 6 5d 4 6s 2 ; these electrons were assigned to the valence space. Geometry optimization (full relaxation) at 0 K was performed at the Γ-point of BZ with E cutoff = 330 eV. Optimization was carried out until the energy difference per atom, the Hellman-Feynman forces on the atoms and all the stress components did not exceed the values of 5 × 10 −7 eV/atom, 1 × 10 −2 eV/Å and 2 × 10 −2 GPa, respectively. The convergence of the SCF-energy was achieved with a tolerance of 5 × 10 −8 eV/atom. While amorphous structures under consideration are not layered and in this sense are not van der Waals (vdW) solids, the presence of strongly polarizable lone-pair electrons requires the use of vdW corrections in calculations. In the CASTEP code, various types of vdW corrections are implemented for different chemical elements and different exchangecorrelation functionals, but not all of them are compatible with all types of GGA-functionals and implemented for all Mo, W and V. The only combination of the GGA functional and the vdW correction implemented in the CASTEP code, available for all three TMs (Mo, W and V), is the combination PW91+TS (Perdew-Wang + Tkatchenko-Scheffler).
Finally, for pure a-As 2 S 3 and doped a-As 2 S 3 :TM structures, both melt-quenched and melt-quenched with subsequent geometry optimization (mq-GOpt) at 0 K, we calculated different electronic and magnetic properties, including band structure, density of states (DOS), partial and local densities of states (PDOS and LDOS), interatomic distances, atomic, bond and spin Mulliken populations [43,44], total and absolute magnetizations. The Monkhorst-Pack [45] 3 × 3 × 5 k-mesh (23 k-points in the irreducible BZ) was used for calculating DOS and PDOS.
The nature of magnetic ordering of a-As 2 S 3 :TM was determined by comparing the calculated values of total and absolute magnetizations. The total magnetization, or in other words the doubled integrated spin density (ISD), gives the doubled total spin of the system. The magnitude of absolute magnetization, or in other words the doubled integrated modulus spin density (IMSD), is a measure of the local unbalanced spin.
Of course, the rather small 240-atom model with just 5 TM atoms may miss some details associated with certain atomic configurations formed in a real glass. Nevertheless, we believe that our model is able to capture the general features of amorphous magnetism in TM-doped chalcogenides glasses.

Results and Discussion
The mass densities of the simulated equilibrated amorphous structures were verified against the experimental data ( Table 1). The presented data demonstrate that they are all slightly lower (3.095-3.321 g/cm 3 ) than that of the crystalline phase, which agrees with the experimental value reported (3.193 g/cm 3 ) [46]. The data presented in Table 1 demonstrate that the discrepancy between the theoretical and experimental a-As 2 S 3 density values is 3%, which is a very reasonable value. Indeed, its magnitude is comparable to the experimental spread of density data reported by different groups for orpiment (crystalline As 2 S 3 ). Obviously for a glass this spread will be even larger.
When doped with TM impurities, the mass densities of amorphous structures increase from 3.095 g/cm 3 for undoped a-As 2 S 3 to 3.159, 3.202 and 3.321 g/cm 3 for a-As 2 S 3 doped with the V, Mo and W impurities, respectively. In this case, we see that the order of increasing mass density upon doping the a-As 2 S 3 with impurity TM-atoms correlates well with an increase in the atomic number of the impurity TM-atom. We also note that undoped a-As 2 S 3 is a semiconductor with a DFT gap of around 1.0 eV as shown in Figure 1(left panel), which is in reasonable agreement with the experimental value of 2.4 eV considering that DFT usually underestimates the gap value by about 50% due to incomplete exclusion of electron self-interaction when using LDA and GGA approximations [47,48]. Figure 2(left) shows the a-As 2 S 3 :Mo melt-quenched structure. At the right of Figure 2, we show fragments of the structure around the Mo atoms. The following observations can be made. One can see that some of the S atoms are three-fold coordinated, i.e., their lonepair electrons are consumed to form covalent (donor-acceptor) bonds with Mo. Analysis of the obtained structure shows that of the 25 S atoms that form covalent bonds with Mo species 13 S are three-fold coordinated. In other words, doping with transition metals results in a decreased concentration of lone-pair electrons. A detailed analysis of the short-range order around molybdenum atoms is difficult as the structure used in the calculations is rather large. That is why to accomplish this analysis fragments of the structure were taken, in which broken bonds were saturated with hydrogen atoms. The differential electron density (charge density difference, CDD), allowing visualization of both covalent bonds (CBs) and lone-pair (LP) electrons, was calculated. CDD is the difference in electron density in the structure under study and the sum of isolated atoms. Consequently the CDD clouds, i.e., an increase in atomic density between atoms, correspond to covalent bonds. Lone pairs are also associated with an increased electron density. An example of electronic distribution visualization is shown in Figure 3. It can be seen that molybdenum atoms are connected with surrounding atoms with covalent bonds. Additionally, lone-pairs formed by s-electrons of arsenic atoms are also visible. At the same time, it is interesting to note that CDD corresponding to p-lonepairs of sulfur atoms in some cases are not observed in agreement with the formation of three-fold coordinated sulfur atoms mentioned above.  The effect of Mo-doping on the electronic structure can be seen from Figure 1, where densities of states for the undoped and Mo-doped mq-GOpt structures are given. One can see that spin-up and spin-down DOSs in pure a-As 2 S 3 are mirror-symmetric in agreement with the experiment (a-As 2 S 3 is diamagnetic). Here it is important to note that spin-up and spin-down DOSs become different in the doped glass, which is an indication that the material became magnetic. Another apparent difference is the formation of a strong tail at the bottom of the conduction band extending down to the valence band and effectively closing the band gap in the Mo-doped glass.
Molybdenum doping also results in the formation of a finite density of states at the Fermi level, i.e., the material becomes a metal. This is in line with the fact that amorphous MoS 2 is metallic [49] while the corresponding crystal in its stable 2H form is a semiconductor. The metallic conductivity of amorphous MoS 2 is due to the formation of a large concentration of metallic Mo-Mo bonds in the amorphous phase. It is remarkable that despite a rather low concentration of Mo atoms in our model, a Mo-Mo dimer is formed, which, again, correlates with the strong tendency of molybdenum to form Mo-Mo bonds in the amorphous phase.
It is not unnatural to assume that isolated Mo atoms and a Mo-Mo dimer contribute differently to the density of states. To address this issue, we removed one of the atoms in the molybdenum dimer, after which the structure was additionally relaxed. In the structure without a dimer, the tail states disappear, while some states in the gap remain. This result allows us to draw a conclusion that the conduction band tail is associated with Mo dimers.
Since Mo atoms all possess different local structures, it is interesting to see which of them (and how) contribute to the magnetic properties. This was done via DFT-method calculating the contributions from individual molybdenum atoms, as well as vanadium and tungsten atoms, to the total spin of the a-As 2 S 3 :TM structures using the well-known Mulliken approximation for the analysis of atomic, bond and spin populations. The results of the Mulliken analysis of spin populations of individual atoms, along with the atomic and bond populations, are presented in Appendix A. The performed DFT calculations showed that due to the different local structure of impurity TM-atoms, not all of them contribute to the total spin of the system. In particular, as can be seen from Tables A1 and A2, only Mo1 and Mo3 atoms contribute noticeably to the total spin of a-As 2 S 3 :Mo. Similar conclusions can also be drawn for vanadium atoms, where four out of five V atoms make a significant contribution to the total spin for mq a-As 2 S 3 :V (see Table A3) and only two out of five-for mq-GOpt a-As 2 S 3 :V (see Table A4). As regards the antiferromagnetic structure mq a-As 2 S 3 :W with zero total spin, the main noticeable contributions to the total spin come from the projections of spins of opposite signs of atoms W1 and W3 (see Table A5). For equilibrated paramagnetic mq-GOpt a-As 2 S 3 :W structure, all tungsten atoms have zero spins (see Table A6). It should be noted that, in addition to the Mulliken atomic and spin populations, Appendix A also lists the bond populations along with the bond lengths between impurity TM atoms and ligand atoms (see Tables A7-A12). An analysis of these data allows us to conclude that in the glassy (amorphous) structures a-As 2 S 3 :TM there are mainly two types of hybridization of the orbitals of impurity TM atoms, namely, trigonalpyramidal (d 1 sp 3 ) and octahedral (d 2 sp 3 ). The first hybridization type involves s, p x , p y , p z and d z 2 orbitals in the combination (s+p x +p y )+(p z +d z 2 ). The other one contains another d x 2 −y 2 in addition to the same orbitals in the combination (s+p x +p y +d x 2 −y 2 )+(p z +d z 2 ). Of course, in glass hybridization is distorted due to the spread in bond lengths and angles.
While the calculation of the contributions of individual atoms to the total spin of the system using the Mulliken approach is very approximate, the total spin can also be calculated more accurately by integrating the total spin density. Indeed, the CASTEP code enables one to calculate two quantities, which qualitatively determine the type of magnetic ordering, namely the integrated spin density (total spin of the system) and the integrated modulus spin density (a measure of the local unbalanced spin), which allow one to characterize magnetism in materials. Depending on the obtained values, a material can be characterized as being paramagnetic, ferromagnetic, ferrimagnetic or antiferromagnetic. Thus, if both quantities are equal to zero, the material is paramagnetic. Both non-zero quantities and the same magnitude determine ferromagnetic ordering. In the case of non-zero value for the first number and a bigger value for the second, the material is ferrimagnetic. Finally, a zero value for the first number and non-zero for the second characterizes antiferromagnetic ordering of material (the total spin is zero, but the local spin density varies).
In order to verify the appearance of magnetism, we performed spin-up and spin-down PDOS (alpha and beta) calculations. As can be seen from Figure 4, spin-up and spin-down contributions from Mo d-electrons are clearly different. At the same time, PDOS of spin-up and spin-down states associated with As and S atoms also become slightly asymmetric which indicates the presence of some contribution to the magnetic properties of the a-As 2 S 3 :TM from As and S atoms when doping with a TM impurity. It is interesting to note that a similar result, namely, an appearance of magnetism on normally diamagnetic S atoms was also observed for edge states of WS 2 nanosheets [50].
Of significant interest is the observation that the contribution of d-electrons substantially changes upon relaxation of the structure. Since the relaxation is done at 0 K when bond breaking is very unlikely to happen and only the existing bond lengths and bond angles slightly change, the result suggests that magnetism in doped a-As 2 S 3 :TM is rather "fragile". In Figure 5, we compare spin-up and spin-down contributions to PDOS from d-electrons of Mo, W and V atoms before and after relaxation. One can see that-as in the case of Mo doping-PDOS shape changes after relaxation. Of special interest is the fact that in the relaxed W-doped a-As 2 S 3 spin-up and spin-down contributions become identical (paramagnetic state). We believe that the observed fragility of magnetism in doped a-As 2 S 3 is related to the absence of dark ESR in chalcogenide glasses, which is due to the absence of dangling bonds with unpaired spins caused by the fact that a melt-quenched relaxed glassy structure adjusts to have all bonds saturated [51]. It should be noted that magnetism is also not observed in single crystal transition-metal chalcogenides but appears at sample edges and grain boundaries [6]. Similarly, magnetism in the studied glasses may be associated with the presence of soft modes [52], where some of the bonds may be electron deficient. The latter may be easily affected by slight changes in the local environment.
Analysis of the integrated spin density and integrated modulus spin density shows that the nature of magnetic ordering in a-As 2 S 3 :TM is different for different dopants, namely, the melt-quenched a-As 2 S 3 :Mo and melt-quenched a-As 2 S 3 :V are ferrimagnetic, while the melt-quenched a-As 2 S 3 :W structure is antiferromagnetic. After subsequent relaxation at 0 K of the melt-quenched structure geometry, the Mo-and V-doped structures remain ferrimagnetic, while the W-doped structure becomes paramagnetic. A summary of the results is given in Table 2, which shows the integrated spin densities (total spin of the system) and integrated modulus spin densities (a measure of the local unbalanced spin) for the materials studied. From Table 2, one can see that despite the observed differences in PDOS of as-quenched and 0K-optimized structures the type of magnetism is preserved (ferrimagnetic) in Mo-and V-doped a-As 2 S 3 , while in W-doped a-As 2 S 3 structural optimization (relaxation) changes the type of magnetic ordering from antiferromagnetic to paramagnetic.
We attribute this unusual behavior of the W-doped amorphous structure to the fact that when the Mo atom is substituted by a bigger W atom, the lattice becomes locally stressed and when this stress is removed during the subsequent optimization of the geometry at 0 K, the nature of the magnetic ordering changes significantly, as a result of which the optimized W-doped a-As 2 S 3 structure loses magnetic properties. This opens up the interesting possibility of controlling the magnetic properties of the TM-doped amorphous structure by applying external pressure. The study of pressure effect on the magnetic properties of a-As 2 S 3 is currently underway. Table 2. Type of magnetic ordering of a-As 2 S 3 :TM structure.

Conclusions
The electronic configuration of sulphur, s 2 p 1 x p 1 y p 2 z , means that two of the four pelectrons are used to form covalent bonds and two are left as a lone-pair. The accomplished DFT simulations indeed showed that for both structures of a pure glass and a-As 2 S 3 :Mo the majority of As and S atoms satisfy the 8-N rule, i.e., As atoms are three-fold coordinated and S atoms are two-fold coordinated. In the doped glass, some lone-pair electrons are consumed to form Ch-TM bonds.
We demonstrate a strong effect of TM-dopants on the electronic structure of the a-As 2 S 3 as well as an appearance of a magnetic response. Of special interest is the fact that doping with different transition metals results in the formation of materials with different magnetic properties. Our results suggest that chalcogenide glasses doped with transition metals may become a technologically important material. Data Availability Statement: Additional data can obtained from the authors upon reasonable request.

Acknowledgments:
The authors are grateful to Yuta Saito and Paul Fons for their involvement during the initial stage of this project.

Conflicts of Interest:
The authors declare no conflict of interest.

Abbreviations
The following abbreviations are used in this manuscript: