Simulating Spin Waves in Entropy Stabilized Oxides

The entropy stabilized oxide Mg$_{0.2}$Co$_{0.2}$Ni$_{0.2}$Cu$_{0.2}$Zn$_{0.2}$O exhibits antiferromagnetic order and magnetic excitations, as revealed by recent neutron scattering experiments. This observation raises the question of the nature of spin wave excitations in such disordered systems. Here, we investigate theoretically the magnetic ground state and the spin-wave excitations using linear spin-wave theory in combination with the supercell approximation to take into account the extreme disorder in this magnetic system. We find that the experimentally observed antiferromagnetic structure can be stabilized by a rhombohedral distortion together with large second nearest neighbor interactions. Our calculations show that the spin-wave spectrum consists of a well-defined low-energy coherent spectrum in the background of an incoherent continuum that extends to higher energies.


I. INTRODUCTION
Entropy stabilization means that a single phase is stabilized by the mixing entropy gained from randomly distributing a number of elements, typically five or more, over a single crystal lattice. [1] In 2015 a new family was added to the class of high entropy materials with the discovery of the high entropy oxides (HEOs). By mixing equimolar amounts of MgO, CoO, NiO, CuO and ZnO Rost et al. [2] found that Mg 0.2 Co 0.2 Ni 0.2 Cu 0.2 Zn 0.2 O (hereafter MgO-HEO) could be stabilized in a simple rock-salt structure in which the oxygen atoms occupy one of the face-centered cubic (FCC) sublattices and the cations are randomly distributed over the other FCC sublattice. The entropy stabilized nature of the phase is apparent when considering that CuO and ZnO do not form in the rock-salt structure that forces Cu and Zn in octahedral coordination. Rost et al. [2] demonstrated that the sample could be switched between a multiphase structure and the aforementioned single phase rock-salt structure by repeated heating and cooling , thereby proving the reversibility of the entropy-driven transition. Furthermore, a combination of X-ray diffraction, X-ray absorption finestructure and scanning transmission electron microscopy with energy dispersive X-ray spectroscopy, showed that the MgO-HEO sample was chemically and structurally homogeneous [2]. The HEOs are not only interesting from a basic scientific point of view, but have also exhibit promising functional properties, such as high Li-ion storage capacity and cycling stability [3], colossal dielectric constants [4], superionic conductivity [5] and a high ratio of elastic modulus to thermal conductivity [6]. In addition to rock-salt HEOs, HEOs with perovskite [7][8][9], fluorite [10,11] and spinel [12] structures have been discovered, together with high entropy carbides [13,14], borides [15] and chalcogenides [16].
The magnetic properties of MgO-HEO had remained tions are gapped, with a spin-gap on the order of ∼ 7 meV [22]. The gap is gradually buried in quasielastic fluctuations upon heating towards T N . Interestingly, the spin-wave excitations persist above T N , up to room temperature [18] . Unlike in binary oxides such as NiO and CoO, no strong lambda-anomaly has been observed in the specific heat [18,19]. Recent muon spin relaxation measurements [22] indicate that MgO-HEO undergoes a broad continuous AFM phase transition starting at 140 K and being fully ordered at 100 K, explaining why the specific heat anomaly is very weak and broad. The discovery of magnetic long-range order in MgO-HEO leads to a number of questions. How can the experimentally observed AFM long range order be stabilized in the presence of the randomly distributed moments of Co, Ni and Cu and the 40% of spin vacancies created by the non magnetic Mg and Zn ions? What is the nature of spin-wave excitations in such a extremely disordered magnetic environment? How do the magnetic properties of MgO-HEO compare with those of Ni 1−x Mg x O [23], Ni 1−x Zn x O [24] and Co 1−x Mg x O [25], which also display antiferromagnetic order in the presence of similarly large spin-vacancy concentrations?
In this paper we theoretically investigate the magnetic ground state and spin-wave excitations in the high entropy oxide MgO-HEO. We model the Co 2+ , Ni 2+ and Cu 2+ cations with spins of size 3/2, 1 and 1/2 respectively and the Mg 2+ and Zn 2+ cations as being spinvacancies. We randomly distribute the cations on a rhombohedrally distorted face-centered cubic lattice and analyze the magnetic properties via linear spin-wave theory in combination with the supercell approximation. We find that large next-nearest neighbor exchange couplings and large rhombohedral distortions are required in order to stabilize the experimentally observed antiferromagnetic ground state with propagation vector q = (1/2, 1/2, 1/2). The spin-wave excitation spectrum of our model for MgO-HEO consists of a coherent component at low energies and an incoherent component at high energies. We also model Co 0.33 Ni 0.33 Cu 0.33 O containing moment-size disorder only and Ni 0.6 Mg 0.4 O containing spin-vacancy disorder only. We find qualitative differences in the spin-wave spectra of these systems compared with that of MgO-HEO. Finally we also investigate the influence of disorder on the spin-wave gap and find that it is proportional to the average moment size per lattice site.

II. METHODS
In order to qualitatively investigate the ground state and spin-wave excitations in MgO-HEO we use a simplified model. X-ray absorption spectroscopy [18] and extended x-ray absorption fine structure [26] are indicative of the cations in MgO-HEO being divalent. In addition, electron paramagnetic resonance [19] and Density Functional Theory (DFT) [27] have shown that Co 2+ is in the high spin state. Based on these experimental and theoretical results, we treat the Zn 2+ and Mg 2+ cations as spin-vacancies and the Co 2+ , Ni 2+ and Cu 2+ cations as spins with moments sizes 3/2, 1 and 1/2 respectively, ignoring orbital contributions here. To describe the interactions between the magnetic moments, we take into account the nearest and next-nearest neighbor exchange couplings. In the rock-salt structure, the nearest neighbor cations are connected via a cation-oxygencation group that has a bond angle of 90 degrees. For the second neighboring cations, the mediating oxygen ions are located in between the two cations, such that the cation-oxygen-cation bond angle is 180 degrees. Therefore, according to the Goodenough-Kanamori rules for superexchange [28,29], the nearest neighboring spin interactions are ferromagnetic and the next-nearest neighboring interactions are stronger and antiferromagnetic. This is also what has been concluded from DFT calculations for CoO [30], NiO [31] and MgO-HEO [32] and from some inelastic neutron scattering studies for NiO [33] and CoO [34]. Finally, we incorporate a rhombohedral distortion into our model. On an undistorted FCC lattice the antiferromagnetic order of the second kind will be a mixture of four degenerate ordering wave vectors: q = (1/2, 1/2, 1/2), q = (−1/2, 1/2, 1/2), q = (1/2, −1/2, 1/2) and q = (1/2, 1/2, −1/2). A rhombohedral distortion with the trigonal axis along the [111] direction will naturally stabilize the q = (1/2, 1/2, 1/2) configuration in favor of the other three ordering wave vectors [35]. Putting it all together, we consider to following model: where r labels the sites on the disordered cation FCC lattice, r, r 1p the first nearest neighbors within the ferromagnetic (111) planes, r, r 1a the first nearest neighbors between the ferromagnetic (111) planes and r, r 2 the second nearest neighbors. Overall our simplified model depends on four parameters, J 1 = J 1 (1 − ∆), J 2 and K, where ∆ quantifies the strength of the orthorhombic distortion and K the magnetic anisotropy. This model captures the disorder due to the random distributions of the Co 2+ , Ni 2+ and Cu 2+ moment sizes and the Mg 2+ and Zn 2+ spin-vacancies. For simplicity, disorder in the magnetic exchange couplings [32] has been ignored.
To study the spin wave excitation spectra of our model for MgO-HEO, we use a combination of linear spin-wave theory [36] and the supercell approximation [37]. The linear-spin wave formalism consists of two steps. First the classical ground state configuration is determined. To this end we use the conjugate gradient method as implemented in the YaSpinWave computer program [38]. In the second step, the classical ground state is used to perform a Holstein-Primakoff transformation so that our interacting quantum-spin model (1) can be approximated as a quadratic boson Hamiltonian that can be solved with exact diagonalization. To treat the strong magnetic disorder in our model for MgO-HEO, we use the supercell approximation in which the disordered system is approximated by a large number of supercells, within which the cations are disordered but beyond which they are artificially periodic. For other studies using linear spin-wave theory in combination with the supercell approximation, we refer to Ref. [39,40] and for alternative methods to simulate magnetic excitations in disordered systems we refer to Ref. [41,42].

III. GROUND STATE
Before studying the spin-wave excitations in MgO-HEO, we first focus on the stability of the q = (1/2, 1/2, 1/2) AFM ground state depicted in Fig. 1. In the ordered FCC spin lattice this AFM ground state is stabilized in the parameter regime |J 2 /J 1 | > 1 − ∆ [35].
Here J 1 > 0 is the FM nearest and J 2 < 0 the AFM next nearest neighbor exchange couplings respectively, and ∆ parametrizes the strength of the rhombohedral distortion, see Fig. 1. From this inequality we can see that the q = (1/2, 1/2, 1/2) AFM ground state becomes more stable upon increasing the magnitude of the second nearest exchange coupling and the rhombohedral distortion. Figure 1 shows that intuitively this makes sense. By inducing a rhombohedral distortion, the magnitude of the FM exchanges within the (111) planes (J 1 ) will increase relatively to those between the (111) planes (J 1 = J 1 (1 − ∆)). The AFM J 2 exchange only couples moments between the (111) planes and therefore helps to further stabilize the q = (1/2, 1/2, 1/2) AFM ground state.
work the same. In Fig. 2 we consider a 2 × 2 × 2 orthogonal supercell relative to the conventional FCC unit cell. This supercell contains 6 S = 3/2, 6 S = 1 and 6 S = 1/2 moments corresponding to the randomly positioned Co, Ni and Cu cations respectively. That leaves 14 randomly distributed spin vacancies corresponding to the Mg and Zn cation sites. Figure 2(a) illustrates that for J 2 = J 1 and ∆ = 0, the ground state is in some non-collinear configuration. Upon increasing the second nearest exchange to J 2 = −15J 1 and the rhombohedral distortion to ∆ = 0.5, the q = (1/2, 1/2, 1/2) AFM configuration is stabilized despite the disorder (c.f. Fig. 2b). We note that the position of the moments within the supercells in this work have been randomized under the constraint that in the q = (1/2, 1/2, 1/2) AFM configuration the total moment is zero. To accomplish this we divide the FCC lattice into two sets of parallel (111) planes. The planes for which in the q = (1/2, 1/2, 1/2) AFM state the moments are up (e.g. the purple plane in Fig. 2b) contain an equal amount of Co, Ni and Cu cations as the planes for which in the q = (1/2, 1/2, 1/2) AFM state the moments are down (e.g. the blue plane in Fig. 2b). This implies that the number of Co, Ni and Cu cations per supercell needs to be even. In order to satisfy the equiatomic limit as closely as possible under this constraint, we choose the 2 × 2 × 2 FCC supercells to contain 6 magnetic cations of each type and 14 non magnetic cations. We now quantitatively investigate the stability of the q = (1/2, 1/2, 1/2) AFM ground state, by varying the second neighbor exchange coupling and the rhombohedral distortion parameter over a range of values. For each of the 4 × 4 parameter combinations listed in Table I, we simulate 100 orthogonal disordered Co 6 Ni 6 Cu 6 Mg 6 Zn 8 supercells, again under the constraint that in the q=(1/2,1/2,1/2) AFM configuration the total moment is zero. For each supercell we minimize the classical energy via the conjugate gradient method, starting from five different random spin configurations. The configuration with the lowest energy is assigned to be the ground state. The ground state configuration is typically found in three out of five energy minimizations. In total 4 × 4 × 100 × 5 = 8000 energy minimizations have been performed to produce Table I. To deter-mine whether the ground state energy corresponds to the q = (1/2, 1/2, 1/2) AFM configuration, we considered the average of the absolute difference for each of the spin vectors in these configurations. Appendix A contains more details about how we differentiate spin-configurations.
As expected, the percentage of cases for which the q = (1/2, 1/2, 1/2) AFM ground state is realized increases as a function of the next nearest exchange coupling J 2 and the rhombohedral distortion ∆, see Table I. Yet the strength of the parameters required to stabilize the q = (1/2, 1/2, 1/2) AFM ground state is considerably higher than for the ordered case, due to the disorder. For example, for the ordered case with ∆ = 0.2, J 2 = −0.8J 1 is sufficient. Our model for the disordered MgO-HEO with J 2 = −15J 1 realizes the q = (1/2, 1/2, 1/2) AFM ground state in only 40% of the configurations. We have also repeated the same energy minimization in the presence of magnetic anisotropy K = 0.05J 1 , see Table III. The number of configurations, for which the q = (1/2, 1/2, 1/2) AFM ground state is realized, remains essentially the same as for the isotropic case, differing at most by a few percentage points compared to the isotropic case shown in Table III. Having understood the stability regime of the classical q = (1/2, 1/2, 1/2) AFM ground state for our MgO-HEO model, we move on to study the corresponding spinwave excitation spectrum. To that end we will focus on the parameters J 2 = −15J 1 and ∆ = 0.5 for which the q = (1/2, 1/2, 1/2) AFM ground state was found to be most stable within our set of parameters. Before analyzing the influence of disorder in our MgO-HEO model, it is instructive to consider the spin-wave spectra of the ordered CoO, NiO and CuO FCC spin models shown in Fig.  3. We notice two things. First, the energy range of the spin-wave spectra is proportional to the moment magnitude, S=3/2, 1 and 1/2, of the Co, Ni and Cu cations, respectively. Second, we convoluted the spin-wave spectra with a Gaussian of width 0.5J 1 for the purpose of visualization. Other than that, the spin-wave bands in  In Fig. 4 we present the spin-wave spectra for our MgO-HEO model. These spectra are obtained from averaging over 100 configurations each containing 250 FCC lattice sites on average over which the Co, Ni and Cu moments and the Mg and Zn spin-vacancies are randomly distributed. We again use the constraint described in the previous section that in the q = (1/2, 1/2, 1/2) AFM state the total moment vanishes. To better average the disorder configurations, we not only randomly distribute the cations within the supercells, but also randomly change the shapes of the supercells themselves, as App. B explains. Comparing the spin-wave spectrum of the strongly disordered MgO-HEO model in Fig. 4 to those of our models for the ordered CoO, NiO and CuO systems in Fig. 3 reveals notable differences. First, the spin-wave spectrum of the MgO-HEO model displays large broadenings in energy and momentum space corresponding to the short lifetimes and short mean-free-paths of the spin-waves due to the presence of strong disorder. Another striking difference is that for the MgO-HEO model the spin-wave spectrum can be characterized by two regions. In the low energy regime, roughly within [0, 10]J 1 , the spin-wave spectrum exhibits a relatively sharp dispersion. At higher energies, roughly within [10, 100]J 1 , the spin-wave spectrum broadens to such an extent, that there is barely any dispersion left. This featureless incoherent part of the spin-wave spectrum should experimentally be difficult to distinguish from the background. Therefore, most likely the observed spin-wave dispersions in inelastic neutron scattering correspond to the coherent part of the spin-wave spectrum in our theoretical calculations. To further illustrate the coherent and incoherent portions of the spin-wave spec-trum for the MgO-HEO model, Fig. 9 in App. D presents energy distribution curves for eight fixed momenta near the q = (1/2, 1/2, 1/2) AFM ordering wave vector. Note that for the spin-wave spectra of our disordered models, see Fig. 4-6 and 9, a logarithmic scale is used for the spin-wave intensity. Next we dissect the consequences of the two different types of disorder in our MgO-HEO model, the moment size disorder induced by the Co, Ni and Cu cations and the spin-vacancy disorder induced by the Mg and Zn cations. In Fig. 5 we focus on the moment size disorder by presenting the spin-wave spectrum of the Co 0.33 Ni 0.33 Cu 0.33 O model in which there is only moment size disorder, but no spin-vacancy disorder. We note two differences between Fig. 5 and Fig. 4. First, the low energy part of the spectrum of the Co 0.33 Ni 0.33 Cu 0.33 O model is much sharper than for our MgO-HEO model. This illustrates that the scattering of low energy spinwaves against moment size disorder is strongly suppressed as a function of energy. The situation is analogous to case of low energy acoustic phonons scattering against mass disorder, for which in 3D the mass perturbation is proportional to the quadratic power of the energy; see Ref. [43] and references therein. A second difference is that the magnon propagation velocity, that is, the slope of the low-energy spin-wave band close to the q = (1/2, 1/2, 1/2) AFM ordering wave vector, is significantly higher. Similarly, the maximum energy of the spin-wave intensity is also higher. The reason for the reduction of the magnon propagation velocity and the maximum energy of the spin-wave intensity in the MgO-HEO is the presence of the spin vacancies, which make it energetically less costly for the remaining spins to rotate away from the AFM ground state.
In Fig. 6 the spin-wave excitation spectrum for the   Finally, we investigate the influence of magnetic disor-der on the spin-wave gap. Table II lists the spin-gaps induced by an easy-axes anisotropy K = 0.05J 1 for various model systems. For the CoO, NiO and CuO models without disorder, we see that the spin-wave gap is proportional to the size of the cation moments, as expected. To extract the spin-gap for our disordered models, we fit the spin-wave excitation spectrum at the q = (1/2, 1/2, 1/2) AFM ordering wave-vector with a Gaussian. The spinwave gap that we obtain for Co 0.33 Ni 0.33 Cu 0.33 O is nearly equal to that for our NiO model. This suggest that the the spin-wave gap is proportional to the average moment size: 1 3 (3/2 + 1 + 1/2) = 1. This is further confirmed by the fact that the spin-wave gap of the MgO-HEO and the Ni 0.6 Mg 0.4 O model are close in energy. For these two models we also notice that the spin-gap is roughly proportional to (1 − x) with x the vacancy concentration. For example, the spin-wave gap of NiO multiplied by the moment concentration is 0.6 · 4.17 = 2.51, i.e. approximately equal to the spin-wave gap of Ni 0.6 Mg 0.4 O of 2.57, and of Mg-HEO of 2.43. Overall, we find that the spingap in our model is roughly proportional to the average moment size if we count the spin-vacancies as having zero magnetic moment.

V. DISCUSSION
In our study we have qualitatively modeled the magnetic ground state and spin-wave excitations in MgO-HEO with a rhombohedrally distorted FCC spin model that depends on four global parameters and fully takes into account the moment size and spin-vacancy disorder. One aspect not modeled here and to be considered in future studies is the influence of disorder in the exchange couplings on the spin-wave excitations. Recently, the exchange couplings in MgO-HEO have been computed from first principles [32]. From fitting the total energy of 3 magnetic configurations in 6 ordered cation distributions the ratio of the second to first nearest neighbor exchange couplings was found to vary from J 2 = −11.4J 1 for Ni-Ni pairs to J 2 = −2.6J 1 for Cu-Cu pairs, with an average ratio of J 2 = −6.6J 1 . This is reasonably consistent with our conclusion that J 2 = −5J 1 is needed to stabilize a significant portion of the disordered configurations in the q = (1/2, 1/2, 1/2) AFM ground state. In this same study [32] an average fitting of the exchange coupling was performed via calculations in a large disordered supercell and led to the conclusion that J 2 = −15.1J 1 , very close to the value used in our study to study spin-wave excitations. Another aspect that would be significant to consider is the influence of local distortions on the exchange couplings. Local Jahn-Teller distortions of the oxygen anions surrounding the Cu 2+ cations have been observed in MgO-HEO via extended X-ray absorption fine structure spectroscopy [26], electron paramagnetic resonance and X-ray diffraction studies [44]. From DFT [45] it is concluded that these Jahn-Teller distortions are randomly orientated. Such distortions could induce a large spread in the exchange couplings even for first or second nearest neighbor couplings of a fixed inter-species pair. For example, for the high entropy alloy CrFeCoNi, DFT calculations show that the nearest neighbor Fe-Fe interactions can fluctuate roughly between +10 and −10 meV depending on the local chemical environment [46]. A realistic spin model can be bench-marked against inelastic neutron scattering. Features that can be used to constrain the model are the experimentally observed spinwave gap of 7 meV [22] and the magnon propagation velocity [18]. In addition, our current study has found that the spin-wave excitations in MgO-HEO consist of a coherent part and a featureless incoherent part. Inelastic neutron scattering experiments on MgO-HEO are most likely only able to resolve the coherent portion of the spin-wave spectrum. The maximum of the coherent spectra (in our current simulations roughly equal to 10J 1 ) is an additional quantitative aspect of the spin-wave spectrum that can be benchmarked against inelastic neutron scattering.
Another interesting aspect to consider in future investigations is the magnetic frustration induced by competition between the first and second nearest neighbor exchange couplings. In this study we have worked under the assumption that below the Néel temperature MgO-HEO orders in the q = (1/2, 1/2, 1/2) AFM ground state. This is what has been reported based on neutron diffraction studies on MgO-HEO and it is also the ground state that has been experimentally observed for the end members NiO and CoO [20,21]. However, it is possible that no global rhombohedral distortion takes place. Indeed, as of yet, no global distortion has been observed experimentally, although this could also be due to limited resolution. In the absence of a rhombohedral distortion, the AFM ground state of MgO-HEO could be a mixture of the four ordering wave vectors q = (1/2, 1/2, 1/2), q = (−1/2, 1/2, 1/2), q = (1/2, −1/2, 1/2) and q = (1/2, 1/2, −1/2) as is the case for the ordered FCC model with weak nearest FM and strong next nearest AFM couplings [35]. Note that such a ground state is not a mixture of single q domains; even within the domains the four q vectors are mixed. The AFM ground state has been resolved from powder diffraction for which it is not possible to distinguish between the single or the mixed q state. Interestingly, in the absence of a rhombohedral distortion, the first nearest neighbor interaction is frustrated. As can be seen in Fig. 1, half of the nearest FM exchange couplings, labeled J 1 , are frustrated in the q = (1/2, 1/2, 1/2) configuration. A moderate frustration, expressed as the ratio of the Curie-Weiss temperature over the Neel temperature, has in fact been reported for MgO-HEO [19]. Frustration of the nearest neighbor interaction could possibly also explain some of the other curious magnetic properties of MgO-HEO, such as for example the broad phase transition ranging from 140 to 100 K [22], the persistence of spin-wave excitations up to room temperature [18], and the bifurcation of the field-cooled and zero field cooled magnetic susceptibil-ity [18,19].

VI. CONCLUSION
In this study we have investigated the spinwave excitations in the high entropy oxide Mg 0.2 Co 0.2 Ni 0.2 Cu 0.2 Zn 0.2 O (MgO-HEO). To model the magnetic interactions in MgO-HEO, we utilize a FCC spin model with FM nearest, and AFM next nearest neighbor interactions in the presence of a rhombohedral distortion and magnetic anisotropy. The Co, Ni and Cu cations are modeled as spins with a moment size of 3/2, 1 and 1/2 respectively, whereas the Mg and Zn cations are taken to be spin-vacancies. We employ linear spin-wave theory in combination with the supercell approximation to treat the disorder in order to theoretically investigate the ground state and the spin-wave excitations. Stabilizing the experimentally reported q = (1/2, 1/2, 1/2) AFM ground state in our disordered model for MgO-HEO requires significantly larger second nearest neighbor interactions and rhombohedral distortions compared to the case without disorder. The calculated spin-wave spectrum of MgO-HEO consists of a coherent dispersive part at low energies, whereas at higher energies the spectrum is incoherent due to the disorder induced broadening. To differentiate the influence of moment size disorder and the spin-vacancy disorder we compute the spin-wave spectrum of Co 0.33 Ni 0.33 Cu 0.33 O and Ni 0.6 Mg 0.4 O and find these to be qualitatively different from each other and from the one of MgO-HEO. We also investigated the spin-wave gap for various cation compositions and find that for the considered model it is proportional to the average moment size per lattice site. States Government retains a non-exclusive, paid-up, irrevocable, world-wide license to publish or reproduce the published form of this manuscript, or allow others to do so, for United States Government purposes. The Department of Energy will provide public access to these results of federally sponsored research in accordance with the DOE Public Access Plan (http://energy.gov/downloads/doe-public-access-plan).

Appendix A: Differentiating spin configurations
To quantify how close the ground state of a disordered supercell is to the q = (1/2, 1/2, 1/2) AFM configuration one could compute the q = (1/2, 1/2, 1/2) AFM energy relative to the ground state energy. However, for comparing among different parameter sets the energy isn't a good measure, because for a given configuration the energy varies as a function of the parameters. Therefore instead we use the spin-configuration of the ground state to determine its closeness to q = (1/2, 1/2, 1/2) AFM configuration. To this end we define the difference between two spin configurations as the absolute difference of the spin-vectors averaged over the spins in the configurations. While doing this, we have to be mindful of the rotational invariance when computing this difference. For example, in the isotropic case the difference between two spin-configurations should not depend on an overall rotation of either of the spin-configurations. To this end we define the difference between spin-configurations A and B as follows. We first perform an overall rotation of configuration B such that its first spin aligns as much as possible with the first spin of configuration A and then compute the average absolute difference of the spin-vectors. To not bias our measure of difference toward the first spin, we also repeat the procedure with aligning each of the other spins as much as possible. Following the above, the difference of spin-configurations B and A is given by: where N is the total number of spins, S A r is the spinvector of spin r in configuration A and is the spin-vector of spin r in configuration B that is rotated to align its spin r as much as possible with the corresponding spin r in configuration A. For isotropic interactions, the rotation matrix is an element of the 3D rotation group: U A r ∈ SO(3). In that case, the spin r of configuration A and B can be perfectly aligned. In the presence of an easy axis anisotropy K > 0, the rotation matrix is an element of the product of the group of 2D rotations around the easy axis and the group of reflections in the plane perpendicular to the easy axis : U A r ∈ SO(2) × Z 2 . For example, suppose we have an easy axis along the z-direction and we have two spins: S A = (cos(φ A ) sin(θ A ), sin(φ A ) sin(θ A ), cos(θ A )) and S B = (cos(φ B ) sin(θ B ), sin(φ B ) sin(θ B ), cos(θ B )).
Having defined a measure to differentiate two spinconfigurations, we have to choose a cut-off below which we consider two spin-configurations to be equal. In this work, we consider a spin-configuration A to have the q = (1/2, 1/2, 1/2) AFM ground state if ∆S < 10 −4 . For example, the configuration in Fig. 2(b) (Fig. 2(a)) has ∆S equal to 0.00002(2.12636) and therefore is considered to be (not to be) in the q = (1/2, 1/2, 1/2) AFM ground state. To better average the disorder configurations we don't only randomly distribute the cations within the supercells, but also randomly change the shapes of the supercells themselves [37]. To construct the non-orthogonal supercells, we first consider the primitive cell of the antiferromagnetic order of the second kind with ordering wave vector q = (1/2, 1/2, 1/2) illustrated in Fig. 7.

Co
Ni Cu Mg Zn more stable in these non-orthogonal supercells compared to the orthogonal ones used to determine the ground state configuration in section III. We select the sizes of these supercells to be within 115 and 135 q = (1/2, 1/2, 1/2) AFM primitive cells, corresponding to 230 and 270 cation sites. We allow the unit cell angles to vary between 80 and 100 degrees. In Fig. 8 we illustrate one of the nonorthogonal supercells used for computing the disorderconfiguration averaged spin-wave spectra. It is spanned by A 1 = −a 1 + 5a 2 − 5a 3 , A 2 = −4a 1 + 4a 2 + 2a 3 and A 3 = 5a 1 + a 2 − 3a 3 . It contains 248 cations and its unit cell angles are α = 83.1354, β = 90 and γ = 85.904.