First-principles study of the inversion thermodynamics and electronic structure of Fe$M_2X_4$ (thio)spinels ($M=$ Cr, Mn, Co, Ni; $X=$ O, S)

Fe$M_2X_4$ spinels, where $M$ is a transition metal and $X$ is oxygen or sulfur, are candidate materials for spin filters, one of the key devices in spintronics. We present here a computational study of the inversion thermodynamics and the electronic structure of these (thio)spinels for $M=$ Cr, Mn, Co, Ni, using calculations based on the density functional theory with on-site Hubbard corrections (DFT+$U$). The analysis of the configurational free energies shows that different behaviour is expected for the equilibrium cation distributions in these structures: FeCr$_2X_4$ and FeMn$_2$S$_4$ are fully normal, FeNi$_2X_4$ and FeCo$_2$S$_4$ are intermediate, and FeCo$_2$O$_4$ and FeMn$_2$O$_4$ are fully inverted. We have analyzed the role played by the size of the ions and by the crystal field stabilization effects in determining the equilibrium inversion degree. We also discuss how the electronic and magnetic structure of these spinels is modified by the degree of inversion, assuming that this could be varied from the equilibrium value. We have obtained electronic densities of states for the completely normal and completely inverse cation distribution of each compound. FeCr$_2X_4$, FeMn$_2X_4$, FeCo$_2$O$_4$ and FeNi$_2$O$_4$ are half-metals in the ferrimagnetic state when Fe is in tetrahedral positions. When $M$ is filling the tetrahedral positions, the Cr-containing compounds and FeMn$_2$O$_4$ are half-metallic systems, while the Co and Ni spinels are insulators. The Co and Ni sulfide counterparts are metallic for any inversion degree together with the inverse FeMn$_2$S$_4$. Our calculations suggest that the spin filtering properties of the Fe$M_2X_4$ (thio)spinels could be modified via the control of the cation distribution through variations in the synthesis conditions.

antiparallel alignment of the magnetic moments of the ions in the tetrahedral and octahedral sublattices (which is known as collinear Néel configuration), 16 while the hopping of the extra electron in the minority channel of the spins explains the half-metallic properties. 9 Greigite (Fe 3 S 4 ) has been shown to have a similar electronic structure to its oxide counterpart Fe 3 O 4 . 10 Both compounds are sometimes found associated with other transition metals of similar ionic radii 17 and valence as Fe, such as Mn, Co, Ni and Cr, 18 forming spinel compounds of formula FeM 2 X 4 . 19 In these systems, M represents the transition metal and X represents the oxygen or sulfur atom, where the sulfide spinels are usually called thiospinels. [20][21][22] The substituted spinels could retain the type of magnetic behavior of their parent compounds (magnetite or greigite), which is driven by a negative superexchange interaction that is stronger between ions occupying different sublattices than between ions within the same sublattice. 16 The crystal structure of a (thio)spinel is face-centered-cubic and the space group is denote tetrahedral and octahedral sites respectively, while x is the degree of inversion. In normal spinels (x = 0), Fe ions occupy exclusively the A sublattice and M is confined to the B sublattice.
In inverse spinels (x = 1), the A sublattice holds half of the M cations and the B sublattice is equally populated with Fe and M ions. When 0 < x < 1, Fe and M have an intermediate degree of distribution within the A and B sublattices. For all inversion degrees, the 3 Fd m symmetry of the spinel is usually retained, as long as all cations are distributed completely randomly within each sublattice, which makes all the sites within each sublattice effectively equivalent. In such cases the cation distribution is fully characterized by the inversion degree x. The degree of inversion in spinels has been found to be affected by different factors, including the ionic radii of the distributed species, the electronic configuration, the electrostatic energy of the lattice, the short-range Born repulsion energy, crystal field effects, and polarization effects. 18 Structural aspects of FeM 2 X 4 (thio)spinels have been reported extensively in the literature, sometimes also addressing their influence on the magnetic and electronic properties. For example, these studies include: (1) the mixing, non-stoichiometry 24 and magnetic properties as a function of the cation site distribution of the Fe 3 O 4 -FeCr 2 O 4 system; 25 (2) the magnetic ordering of FeCr 2 O 4 26-29 and FeMn 2 O 4 30 upon crystal symmetry lowering; (3) the relevance of the electronic structure to the magnetic properties of FeCr 2 X 4 ; 31 (4) the transport properties based on the half-metallic electronic structure of FeCr 2 S 4 ; 32 (5) 43 as well as (9) the thermodynamic stability 44 and cation distribution of FeNi 2 S 4 . 45,46 Nevertheless, in compounds such as these (thio)spinels, where the other transition metal's atomic number differs only by 1 from Fe, the X-ray diffraction intensities Owing to the experimental limitations for the determination of the cation arrangement in FeM 2 X 4 (thio)spinels, in the present work we have used DFT+U calculations to investigate systematically how modifying the spinel composition affects the equilibrium inversion degree and how this determines the magnetic and electronic properties at a given composition. We study the influence of the nature of the M and X ions (M = Cr, Mn, Co, Ni and X = O, S) on these properties, a type of investigation that has been undertaken previously for other groups of oxide spinels 14 and Heusler alloys 47,48 with potential application in spintronic devices. We will discuss from a thermodynamic point of view the equilibrium cation distribution of these (thio)spinels and the role of the ions' sizes and crystal field stabilization effects. We will also analyze the dependence of the electronic and magnetic structure on the degree of inversion for the normal and completely inverse systems.

A. Calculation details
We have carried out spin-polarized quantum mechanical calculations using density functional theory (DFT) as implemented in the Vienna Ab-initio Simulation Package (VASP). [49][50][51][52] The Perdew-Burke-Ernzerhof functional revised for solids (PBEsol) 53 was the version of the generalized gradient approximation (GGA) used as exchange-correlation functional for all geometry optimizations and for the calculation of all density of states (DOS), because PBEsol provides a better description of the structure of solids than its parent functional. 54 The semiempirical method of Grimme (D2) was also included in our calculations for modelling the long-range van der Waals interactions. 55 Even when these interactions are not expected to affect significantly the bulk properties of the hard solids investigated here, we have included the D2 correction at this stage because in future work we expect to study the surfaces of these solids and their interactions with adsorbates, where dispersion effects may play a significant role. [56][57][58][59][60][61] The projector augmented wave (PAW) pseudopotential method 62 [68][69][70] In order to improve the description of the highly correlated 3d electrons in the spinels under study, we have included the Dudarev et al. 71 approach for the d orbital correction within the DFT + U method. 72 We report in Table I [73][74][75][76][77][78][79] which provides band gaps of better quality than semi-local functionals. 80 The HSE06 is made up from the Perdew-Burke-Ernzerhof functional (PBE) 81,82 exchange and correlation components, mixed with 25% of short-range Hartree-Fock (HF) exchange. 73 The Coulomb potential exchange is replaced by a screened potential (with screening parameter ω = 0.207 Å −1 ) in order to define the separation between the short-and long-range components of the HF exchange. 79 While the amount of short-range HF is a constant determined by perturbation theory, making HSE06 an adiabatic connection functional in this part of the potential, 83 its screening parameter is a reasonable system-averaged value across a 7 wide variety of systems, giving better agreement with experiments in the case of semiconductors than for metals or insulators. 80 For the fitting, we carried out single-point calculations with both PBEsol + U and HSE06, using unrelaxed structures with normal cation distributions (a 0 and u 0 were taken from experiment for these calculations, values are listed in Table II). In a first step, we determine U eff for Fe, by considering the Fe 3 O 4 and Fe 3 S 4 electronic structures (which have been studied before). 10,44,84 We then keep these values for Fe and perform a set of DFT + U calculations where the effective Hubbard parameter (U eff ) of the M ion was changed in steps of 0.5 eV from 0 to 6.0 eV. In all the HSE06 calculations, we used the same settings as for the PBEsol simulations.  87 The use of the primitive cell ensures that there is a single cation configuration for each of these three degrees of inversion, which simplifies the simulations, allowing us to scan a wide range of FeM 2 X 4 compositions. This approximation follows previous work, where the use of the primitive cell model has been found to adequately describe experimental properties of half-and fully-inverted spinels. [87][88][89][90][91] However, we cannot rule out completely that the use of larger supercells could actually lead to cation configurations with lower energies for the same inversion degree, as found for example in a recent study of CoFe 2 O 4 . 92 Following the collinear Néel model, the initial magnetic moments of the atoms within each sublattice were set parallel among themselves and antiparallel to those of the other sublattice. For each inversion degree, we ran several calculations, specifying different initial magnetic moments, corresponding to different combinations of low-and high-spin states for the transition metal ions in each sublattice, in order to identify the ground state. The magnetic moments were allowed to relax during each of the calculations. It should be noted that the magnetic structure with antiparallel sublattices is not strictly valid in the case of FeCr 2 O 4 , which is known to have a spiral magnetic structure. 33 However, for the sake of comparison with the other spinel systems, we have not considered its different magnetic structure in this study.
The experimental lattice (a 0 ) and anion (u 0 ) parameters (defining the anion position in the crystal) of the spinels are shown on Table II. These were used as the starting structures for our simulations, where a 0 and the internal coordinates were allowed to relax fully for each inversion degree. We kept the cell shape perfectly rhombohedral in such a way that the conventional cell was always cubic. As FeMn 2 S 4 and FeCo 2 S 4 spinels have not been characterized so far, we postulated an initial hypothetical structure for both. For the Mn-and Co-based thiospinels, we kept the same initial anion parameter as in their oxide counterparts and scaled up their initial lattice parameter according to the equation: which gives the estimates shown in Table II.
TABLE II. Summary of the initial unit cell lattice (a 0 ) and anion (u 0 ) parameters of FeM 2 X 4 spinels used in this work. The relaxed a and u are also reported for x = 0, 0.5 and 1. Note that the origin is the center of symmetry. 19 Structure Experimental

B. Configurational free energy of inversion
The calculation of the inversion degree in spinels containing two different cations is based on the thermodynamic considerations of Navrotsky and Keppla, 96 which have proved to agree well with experiments. 88,[97][98][99][100] This methodology is based on the treatment of the spinels' cation distribution as a chemical equilibrium. We calculated the configurational free energy of inversion per formula unit ΔF config as, where ΔE config is the inversion energy per formula unit, T is the temperature and ΔS config is the ideal configurational entropy also per formula unit, which is calculated as, where R is the ideal gas constant. ΔS config = 0 and 11.59 J⋅mol -1 ⋅K -1 for x = 0 and 1 respectively, while it reaches the maximum value 15.88 J⋅mol -1 ⋅K -1 for the complete random distribution at x = 2/3. The above expression means that we have only considered ideal contributions to the configurational entropy, in line with previous work. 88,97,100 We are also ignoring vibrational contributions to ΔF, as their contributions are typically small compared to configurational energies and entropies. 88,100 Table II shows the optimized a and u for the three inversion degrees considered (x = 0, 0.5 and 1.0). In general, the optimized lattice parameter is within 2% from the experimental value, where this is available. However, the relaxed lattice parameter for FeCo 2 S 4 , in the best case (x = 0), is 3.8% away from the initial estimated value, which may be an artefact due to the assumption of linearity between the lattice constants of Fe 3 O 4 , FeCo 2 O 4 and their sulfide counterparts. After relaxation of the structures, u was still different from the value of ¼ that it has in the perfect spinel. This deviation reflects the displacement, in the (111) direction, of the anion from its ideal position in order to accommodate cations of specific volume. The biggest deviation in u in comparison with the experimental value was for the inverse Mn-and Co-based oxide spinels, Table II. We see that in general, u and a values are sensitive to the cation distribution, although no systematic rule can be derived from their dependence.

B. Equilibrium inversion degrees
The configurational inversion energy per formula unit (ΔE config ) was fitted versus the degree of inversion (x) using a quadratic regression curve, see Fig. 2 (a). More details regarding the empirical relationship 101 and theoretical justification 102 of the above fitting in terms of x, a and u can be found elsewhere. In this fitting, we defined the normal spinel as the standard state for a given condition of temperature and pressure.
Using the quadratic equation for ΔE config , it is possible to interpolate the inversion energy for any value of x between 0 and 1. Based on this protocol, we have also estimated the configurational We found the minimum of ΔE config to correspond to a normal distribution of cations, with the exception of Co-based systems and FeMn 2 O 4 , Fig. 2  FeCr 2 X 4 . From Fig. 2 (b), we deduced that at 1000 K, the Cr-based (thio)spinels are normal under equilibrium conditions, as a result of a highly endothermic process of inversion. This normal cation distribution of FeCr 2 X 4 is supported by powder neutron diffraction intensities 33 at room temperature and for the oxide by Mössbauer measurements 25 and DFT calculations, 24 see Table III. FeMn 2 X 4 . At 1000 K, the scenario for FeMn 2 O 4 is unique in this study, as in addition to the global minimum of ΔF config at x = 1, it has a local one at x = 0.1. The local minimum is within a portion of shallow inversion free energy (0 < x < 0.3), which may lead to a metastable inversion degree anywhere within this range, for this spinel's equilibrium structure, see Fig. 2 Table III. The inversion degree has also been found in neutron diffraction experiments to be at x = 0.91 37,38 which is in reasonably good agreement with the global minimum calculated here. We speculate that the two inversion degrees of FeMn 2 O 4 may be hampered by kinetic control, a situation that is outside the scope of this paper, but which explains the different cation arrangements described in the literature.
According to our calculations, samples of FeMn 2 O 4 synthesized above 1150 K can only have inverse cation distribution, as the metastable inversion degree vanishes. FeMn 2 S 4 , on the other hand, is predicted to be mostly normal (x = 0.03) under equilibrium conditions, Fig. 2 (b).  Table III. Its sulfide counterpart, which has not been studied experimentally, shows an equilibrium inversion degree of x = 0.48 in our calculations.  Table III. However, they disagree with the more recent description of synthetic FeNi 2 S 4 samples as completely inverse spinel based on a neutron powder diffraction measurements at temperatures between 100 and 573 K, 45 thermodynamic-based modelling, 46 EXAFS experiment 22 and Mössbauer data. 21 Based on our calculations and the fact that synthetic FeNi 2 S 4 samples cannot be annealed to temperatures higher than 734 K, 104 as they decompose, we propose a rationalization of the different cation arrangements found in natural and synthetic samples of this mineral. We suggest that synthesis produces a kinetic product (with x ~ 1.00) and that these conditions cannot reproduce the hypogene processes occurring in the ores deep below the Earth's surface that lead to the thermodynamic product found in natural samples (x ~ 0).

C. Size of ions and crystal field stabilization effects
We analyze now the effect of cation size and crystal field stabilization energy on the distribution of cations under equilibrium conditions.
Assuming the hard-sphere model, where the ions are spherical, rigid and in contact, the ratio between the tetrahedral (R A ) and octahedral (R B ) bond distances will depend solely on u. The tetrahedral holes are smaller than the octahedral ones for u < 0.2625. 111 Taking into account that for most systems under study here, u is below that value (with a few exceptions in the relaxed structures, see Table II), we can consider that R A < R B is expected for a stable spinel.
According to the Shannon effective radii, 17   reported). These OSPEs clearly show the preference for normal spinels.
The ambiguities in our results are probably not surprising, because earlier attempts to correlate cation distribution with their size and crystal field effects were also not successful, or at least, unable to provide a complete prediction of the degrees of inversion. 114

D. Atomic spin moments and charges
In this and the next section we analyze the electronic and magnetic properties of the spinel materials for the extreme cases of x = 0 and x = 1.
The total magnetization of saturation (M S ) is defined experimentally as the maximum magnetic moment per formula unit of a compound under an increasing magnetic field. This magnitude can also be calculated according to the Néel model as the sum of the atomic spin densities (m s ) in the tetrahedral and octahedral sublattices per formula unit. 16 Table IV  Ni , we also found it to be low spin     In the inverse (thio)spinels, the calculated spin densities for 2 B2 Fe  were slightly overestimated in the oxides compared with the Néel model, while they were more moderately underestimated in the sulfides considering a high-spin electronic distribution of    3 1 2 22 g g g tte for these ions, see Table   IV. The calculated atomic spin densities of the inverse Cr-and Mn-based (thio)spinels agree better, especially in the cations occupying the octahedral (B1) positions, with the high-spin electronic distribution for these ions, as described for the normal spinels. However, in the Coand Ni-based inverse compounds, we found low-spin densities for these atoms in the B1 positions, where the nearly diamagnetic Co B1 in the inverse FeCo 2 O 4 agrees with experiments. 108,109 Notable exceptions are the 3 A Co  in the thiospinel and 3 A Ni  in both oxide and sulfide compounds, where our calculations shift m s by more than 1 μ B /atom with the expected value (in the best case) from a low-spin electronic distribution for these atoms.
The most stable normal cation distribution of FeCr 2 O 4 gave the closest M S to the experimental one, although still overestimated by 1.35 μ B /f.u., as this measurement was carried out at a temperature in which the spins are not collinear anymore. 33 In the case of its normal sulfide counterpart, the difference in spin magnetization of saturation with the experiments is smaller (0.41 μ B /f.u.). 33 The Ni-based (thio)spinels are found experimentally to be paramagnetic. 21,43,103 In the oxide this has been explained as being due to high-spin Fe 3+ ions exclusively localized on the A sublattice whose spins compensate completely the [Ni 2+ Ni 3+ ] occupying the octahedral positions. 43 In the sulfide this has been rationalized on the basis of an A sublattice filled by Ni 3+ and low-spin Fe 2+ occupying octahedral positions. 21 Here, based on our calculated spin magnetization of saturation and assuming intermediate degrees of cation distribution, we present a fresh explanation for the paramagnetism of FeNi 2 X 4 . Considering that M S changes linearly with x, we may postulate that the oxide and sulfide will be paramagnetic when x = 0.33 and 0.30 respectively. Although this suggestion agrees with the oxide and sulfide equilibrium inversion degree calculated in section 3.1, it shows that paramagnetism in these compounds may be due not to the canonical inverse spinel structure with integer oxidation numbers, but to intermediate inversion degrees.
To the best of our knowledge, there is no experimental determination of the saturation magnetization of either FeCo 2 S 4 or FeMn 2 S 4 . Although we found both compounds to be ferrimagnetic, the measurement of M S for FeCo 2 S 4 may be essential to determine the inversion degree of this spinel, as our calculation of the normal and inverse cation distributions show different magnetizations of saturation. Our results agree with the ferrimagnetic behaviour described for FeCo 2 O 4 115 and FeMn 2 O 4 , 30 below the Curie (Néel) temperature. Table V shows the charges (q) gained or lost by an atom with respect to the neutral atom in the FeM 2 X 4 spinels. We clearly appreciate, that charges are systematically underestimated for all the FeM 2 X 4 (thio)spinels. For the Cr-and Mn-based (thio)spinels, q A is clearly smaller than q B for any inversion degree and also for the inverse FeCo 2 O 4 . However, for the Co-and Ni-based systems, the relative charges of the atoms in the tetrahedral and octahedral positions is different in the oxide and sulfide. In the spinel oxides, together with the normal Co-and Ni-thiospinels, q A is bigger than q B, while in the inverse thiospinels, it is the other way round.  The density of states (DOS) in Fig. 3 show that at x = 0, FeCr 2 X 4 is half-metallic, which we confirmed by the integer value of total spin magnetization (M S = 2.00 μ B /f.u.), see Table IV. An integer value of the total spin magnetization discriminates half-metals and insulators from metals.
The total number of electrons of any stoichiometric system is integer (n) and if it has a band gap at least in one spin channel, there is an integer number of electrons ( n ) there too. This makes the difference ( n n n    ), which is the number of electrons on the band that crosses the Fermi level also integer. Therefore, the magnetization of saturation, i.e. the difference of n and n , is also integer. 3,15,116,117 The DOS shows a sharp peak of the partially-occupied e level of Fe A ions in the majority spin channel (α) crossing the Fermi energy, which is weakly hybridized with the empty Cr e g level in the oxide spinel, while the minority spin channel (β) shows a gap near F E . There is a nearly

FeMn 2 X 4
When the FeMn 2 X 4 (thio)spinel is normal (x = 0), the Fe A e and t 2 levels, which are very close or hybridized, appear in the maximum of the valence and minimum of the conduction bands in the minority and majority channel respectively, see Fig. 4 left panels. The half-metallic character of the normal FeMn 2 X 4 (thio)spinels is also confirmed by the spin density analysis, showing a spin magnetization per formula unit of M S = 4.00 μ B , see Table IV. At the Fermi energy, the spin up partially-occupied e g band of Mn B appears highly hybridized with the X p orbitals and with the Fe A e and t 2 levels in the oxide and sulfide respectively, in agreement with the bigger atomic volume, enhancing orbital overlapping. The rest of the density of states is essentially the same as in the Cr-material, while in the Mn-based spinels the valence band is slightly shifted towards the Fermi energy and the Fe A t 2 and e bands in the β channel of the spins appear more prominently.
With the normal cation distribution, the bands of FeMn 2 S 4 ( Fig. 4 bottom-left panel), as in the case of FeCr 2 S 4 , are shifted towards the Fermi energy, except for the Mn B t 2g and e g bands in the α channel of the spins.
When FeMn 2 O 4 is a completely inverse spinel (x = 1), the system becomes half-semiconductor with a negligible band-gap, which has also been found experimentally, 35 see Fig. 4 top-right panel and also note in Table IV the

FeCo 2 X 4
When FeCo 2 O 4 has a normal distribution, all the bands are pushed slightly towards the Fermi energy and especially those due to Co B , see Fig. 5 top-left panel. As a result, the partiallyoccupied Co B t 2g level that crosses the Fermi energy has a minimal band gap in the minority channel of the spins, making the normal FeCo 2 O 4 spinel almost a half-metal, see also the integer value of the spin magnetization of saturation in Table IV. On the other hand, the sulfide counterpart has all the bands closer to the Fermi energy with symmetrical Co B bands in the minority and majority spin channels (Fig. 5 bottom-left panel), due to the fully-occupied t 2g level which indicates non-magnetic behaviour for this atom (see Table IV). There is a peak at the Fermi energy, in the majority spin channel, with contributions from the partially-occupied Fe A e and Co B e g levels. In the minority spin channel, the normal FeCo 2 S 4 spinel is weakly conducting, as there is a small Co B fully-occupied t 2g band strongly hybridized with S p orbitals that ends shortly after the Fermi energy in the conduction band side. Overall, the sulfide counterpart is metallic which is confirmed by a decimal spin magnetization of saturation,  Table IV, typical of metals). These properties are due to the hybridized partially-occupied Co A t 2 and fully-occupied Co B1 t 2g levels and to the merged Co A t 2 , Co B1 t 2g and Fe B2 t 2g levels in the majority and minority channel of spins, respectively.

FeNi 2 X 4
For both Ni-based (thio)spinels, when x = 0, the bands' pattern is similar and follows the same distribution as described in previous cases, see Fig. 6 left panels. The oxide is half-metal due to a strong hybridization of the partially-occupied Ni B e g band with the O p orbitals that cross the Fermi energy in the majority channel of the spins, see also the integer value of M S in Table IV. The main difference between oxide and sulfide lies in the fact that bands in the latter are closer to the Fermi energy in both spin channels, thus becoming a metallic system. The thiospinel's metallic character is given by the S p orbitals with a small hybridization (only in the majority channel of the spins) with the partially-occupied Fe A e and t 2g levels, see the decimal M S in Table   IV.
When Ni is filling the A sublattice, the Ni-based spinel becomes half-semiconductor with a band gap of 0.20 and 2.05 eV in the majority and minority channels of the spins, respectively, see Fig.   6

IV. SUMMARY AND CONCLUSIONS
We have performed systematic electronic structure calculations for a series of (thio)spinels, which elucidate the cation distribution as well as the magnetic and electronic properties of these materials.
We have determined the thermodynamic inversion degree for the FeM 2 X 4 (thio)spinels at temperatures used typically in their synthesis, which agrees reasonably well with available No single factor among those analyzed, i.e. neither crystal field stabilization effects nor the size of the cations, can account by themselves for the equilibrium inversion degree.
For the two extreme scenarios, namely the completely normal and inverse spinels, we have calculated the electronic and magnetic properties of the metal atoms as well as the electronic properties of the bulk phase. We found that the majority of the spinels for any extreme inversion degree are half-metals in the ferrimagnetic state. Notable exceptions are the inverse Co and Ni oxide spinels, which are insulators, and their sulfide counterparts that are metallic for any inversion degree, together with the inverse FeMn 2 S 4 . Notably, we found that hard anions stretch the band structure, giving the biggest band gaps and therefore the best half-metallic properties.
Finally, we have proposed a theoretical structure for FeMn 2 S 4 and FeCo 2 S 4 and have predicted their electronic and magnetic properties and equilibrium inversion degree.