Data-driven design of B20 alloys with targeted magnetic properties guided by machine learning and density functional theory

Chiral magnets in the B20 crystal structure host a peculiar spin texture in the form of a topologically stable skyrmion lattice. However, the helical transition temperature ( T C ) of these compounds is below room temperature, which limits their potential in spintronics applications. Here, a data-driven approach is demonstrated, which integrates density functional theory (DFT) calculations with machine learning (ML) in search of alloying elements that will enhance the T C of known B20 compounds. Initial DFT screening led to the identi ﬁ cation of chromium (Cr) and tin (Sn) as potential substituents for alloy design. Then, trained ML models predict Sn substitution to be more promising than Cr-substitution for tuning the T C of FeGe. The magnetic exchange energy calculated from DFT validates the promise of Sn as an effective alloying element for enhancing the T C in Fe(Ge,Sn) compounds. New B20 chiral magnets are recommended for experimental investigation.


Introduction
Transition metal compounds that crystallize in the cubic B20type structure belong to the space group P2 1 3 (# 198). One of the intriguing characteristics of a P2 1 3 space group is that it lacks the spatial parity operation of inversion (i.e., noncentrosymmetric). Many interesting properties emerge due to the inversion asymmetry, such as magnetoelectricity, second harmonic generation, pyroelectricity, ferroelectricity, and optical activity [1]. The main focus of this article is on chiral magnets in the bulk B20 crystal structure, which have been shown to host a peculiar spin texture in the form of a topologically stable skyrmion lattice. In this materials class, the skyrmion phase appears in a small pocket of the temperaturefield phase diagram [2]. The origin of chiral spin modulation in noncentrosymmetric crystals is attributed to the energetic competition between the stronger magnetic exchange (J) interactions and the relatively weaker Dzyaloshinskii-Moriya (DM) interactions. An additional contributing factor to DM Prasanna V. Balachandran Prasanna V. Balachandran is an assistant professor at the University of Virginia, with a joint appointment in the Department of Materials Science and Engineering and the Department of Mechanical and Aerospace Engineering. He was a postdoctoral researcher in the Theoretical Division at Los Alamos National Laboratory. He received his PhD degree in materials science and engineering from Iowa State University in 2011. He was a recipient of the Research Excellence Award that recognizes "the best of the best" graduating thesis work. His research focuses on computational materials science and information science applied to accelerate the search and discovery of novel functional and structural materials. Balachandran can be reached at pvb5e@virginia.edu.
ª Materials Research Society 2020. This is an Open Access article, distributed under the terms of the Creative Commons Attribution licence (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted re-use, distribution, and reproduction in any medium, provided the original work is properly cited. interaction is induced by the relativistic spin-orbit coupling that scales roughly with Z 2 , where Z is the atomic number of the chemical element [3].
Skyrmions are regarded as promising for spintronics applications for at least three reasons [4]: (i) nanometer size, (ii) topological protection, and (iii) appearance of spin transfer torque, where magnetic structures and textures are manipulated by electric currents. Typical length of a magnetic domain wall that is used for hard disk recording is of the order of 30 nm [5]; however, sub-30 nm size skyrmions are experimentally observed [6,7], thus showing the exciting potential for encoding vast amount of information relative to the magnetic domain wall technologies [4]. The topological protection contributes to superior thermal stability [8,9]. The spin transfer torque appears when the current density exceeds ;10 6 A/m 2 , which is 4-5 orders of magnitude smaller than the current density required to move the ferroelectric domain walls [10]. However, not all chiral magnets that form stable skyrmions satisfy the stringent property requirements for use in spintronics applications.
In Fig. 1, a plot is shown that summarizes the experimentally measured helical period (denoting skyrmion size) and helical transition temperature (T C ) of several bulk B20 alloys [11]. It is interesting to note that none of the known B20 alloys satisfy the combined requirement of smaller helical period (,50 nm) and higher T C (.300 K) for practical spintronics applications. The best-known material in the B20 family is FeGe, which has a helical period and T C of 70 nm and 278 K, respectively. To date, only ;20 compositions have been experimentally explored in the B20 alloy family. This indicates the challenge and the low-throughput nature of the synthesis and characterization efforts in the materials discovery cycle. The low-throughput nature of data generation from experiments is justified by the need for performing expensive and time-consuming synthesis (sometimes requiring nonequilibrium processing routes [12,13]). Moreover, nontrivial characterization studies are also needed to reliably report the crystal structure and magnetic properties of these alloys, adding to the challenge. However, the periodic table and the principles of alloy theory offer more avenues for composition design than what has been explored so far in the literature. We ask the following question-Are there novel B20 alloy compositions in the unexplored search space with T C better than that of known materials? Alloy design using computational methods have not attracted a lot of attention. Most of the computational efforts that involve density functional theory (DFT) calculations have focused on either providing an explanation of the observed properties from the viewpoint of the band theory or extracting properties such as the strength of the DM term as a function of alloy substitutions and epitaxial strain [14,15,16,17,18,19,20,21,22]. In this work, a novel computational approach, built on the foundations of machine learning (ML) and DFT, is developed to accelerate the design of B20-based chiral magnets with improved T C . Although ML methods have been used in the past to predict the ferromagnetic Curie temperature of alloys [23,24], properties of hard permanent magnets [25], two-dimensional materials [26] and magnetic properties of single-molecule magnets [27,28], no a priori rules exist that link alloy compositions to T C for the B20 alloys.
The objectives are accomplished by developing a datadriven approach that synergistically integrates insights from DFT calculations and published experimental literature with regression-based ML methods. The role of DFT is 2-fold. (i) To perform initial screening and identify promising alloying elements for ML. (ii) To validate the predictions from ML and provide a physical basis for meaningful interpretation of the ML outcome. The task of ML methods is to establish a quantitative relationship between the alloy compositions and the experimentally measured T C of known B20 alloys. The trained ML models are then used to rapidly predict promising new alloy compositions. The final outcome is the recommendation of new B20 chiral magnets for experimental validation.

Results and discussion
The overarching data-driven strategy is schematically shown in Fig. 2. The data-driven search was initiated by identifying a total of 36 binary compounds in the AB stoichiometry, where A is a transition metal atom and B 5 Si, Ge, or Sn. All compounds were constrained in the B20 bulk crystal structure (although B20 is not the lowest energy crystal structure for many of the AB compounds explored in this work). Spinpolarized DFT calculations were performed to screen for AB compounds, such that the converged structure is ferromagnetic Figure 1: A summary of experimentally determined helical period in nanometer (x-axis) and helical transition temperature in K (T C , y-axis) of bulk B20 alloys taken from the published literature [11]. None of the known B20 alloys satisfy the combined requirement of smaller helical period (,50 nm) and higher T C (.300 K) for practical spintronics applications, which is denoted as "Target" in the figure. In this study, the focus is on the computational design of new B20 alloys that have T C . 300 K.
(see "Methods" for more details), which is a necessary condition for stabilizing a skyrmion phase in the temperature-field phase diagram. The purpose of these initial DFT calculations is not to accurately calculate the atomic magnetic moment of the transition metals in the B20 structure, because DFT is known to overestimate the magnetic moment for some of the B20 compounds (e.g., MnSi [18]). Instead, this step serves as an accelerated screening criterion to downselect promising alloy substituents. The key parameter of interest is the atomic magnetic moment and only those B20 compounds with total magnetization .0.1 Bohr magneton per cell were downselected as candidates from DFT. Out of the 36 compounds that were explored during this initial screening stage, only seven compounds (MnSi, MnGe, MnSn, FeGe, FeSn, CrGe, and CrSn) satisfied this criterion. Three of the seven screened compounds correspond to MnSi, MnGe, and FeGe, which represent the most widely studied chiral magnets in the B20 materials class that also host a skyrmion phase. This gives confidence in the initial screening step. The DFT optimized lattice constant and atomic magnetic moment data of the A-site that contains the transition metal atoms are given in Table I. One of the important outcomes from the initial DFT screening step is the identification of Cr and Sn as potential substituents for tuning the magnetic properties of B20 alloys. The impact of Cr and Sn alloying elements on the T C of B20 alloys have not been explored in the literature. The next key question is then the following: Can we rapidly assess the potential of Cr and Sn when added as a substitutional impurity to known B20 alloys in impacting the T C of the alloy? We choose ML methods to address this question.
Starting from a training dataset of 18 alloy compositions, we built descriptors based on orbital radius [29], effective principal quantum number [30], ratio of orbital radius to effective principal quantum number, and valence electron number [31], which are known to capture the chemical trends of solid compounds [29,30,31,32,33,34,35,36]. The full list of input descriptors is given in Table II. The alloy compositions used for training and testing the ML models are given in Table III. Each B20 chemical composition is represented in terms of the linear combination of the weighted contribution of these elemental descriptors. For example, the A-site valence electron number for Fe 0.2 Mn 0.8 Ge is described as 0.2 Â VEN Fe 1 0.8 Â VEN Mn . Separate sets of descriptors were constructed to independently represent A-and B-sites (see Table II). Pairwise statistical correlation analysis indicated a strong Figure 2: The overarching strategy for accelerating the search for new B20 alloy compositions with improved T C is shown. The initial screening step involves running spin-polarized DFT calculations on the chosen chemical space of 36 AB compounds in the constrained bulk B20 crystal structure, where A is the transition metal atom and B is Si, Ge, or Sn. Out of 36 AB compounds, only 7 converged in the magnetic structure, which is a necessary condition for stabilizing a skyrmion phase in B20 compounds. From the seven downselected compounds, two promising alloying elements were identified, namely Cr and Sn, which had not been explored in the literature. ML models were used to establish a relationship between descriptors that represent the compositions of experimentally known B20 and their measured T C values. The trained models were, in turn, used to predict the potential of Cr-and Sn-substituted B20 alloys in the affecting the T C property. ML models identified Sn as a candidate to improve the T C of B20 alloys. The data-driven predictions were validated using DFT calculations, where the magnetic exchange energy was calculated using a supercell approach. Promising and previously unexplored B20 alloy compositions were identified and are recommended for experimental validation.  Table II. It is important to note that we did not include either J or DM interactions as input descriptors for building the ML models. This is mainly because calculation of J and DM interactions from DFT requires supercells and non-collinear spins, which are computationally challenging even for a simple compound. Furthermore, our composition space also includes solid solutions of B20 alloys (see Table III), and DFT calculation of J and DM interactions for solid solutions is a nontrivial task. An ensemble of 100 ML models was trained based on the support vector regression (SVR) method for predicting the T C as a function of three input descriptors (see "Methods" for more details). Balachandran and co-workers have demonstrated the potential of ensemble-based SVR methods for building reliable ML models from small data [27,38,39,40,41,42]. The trained models were validated using Mn 0.75 Rh 0.25 Ge B20 alloy composition as a test case [37], which was recently synthesized by Sidorov et al. using a nonequilibrium high-pressure synthesis method. In the family of Mn 1Àx Rh x Ge alloys prepared by Sidorov et al., Mn 0.75 Rh 0.25 Ge had the highest T C (125 K). The trained ML models predicted the T C value for Mn 0.75 Rh 0.25 Ge compound as 105.6 6 40.2 K. As shown in Fig. 3, the measured T C falls within the error bar of the predicted values, which is encouraging. Note that Mn 0.75 Rh 0.25 Ge was not part of the dataset, that was used for training the ML models. This simple, yet powerful, exercise showed that the trained ML models have predictive capabilities and can be used to rapidly predict the T C for previously unexplored B20 alloy compositions.
We specifically focused on the effect of Cr and Sn substitutions on some of the experimentally known B20 alloys, namely MnSi, MnGe and FeGe, to maximize the impact of our predictions and motivate new experiments. On the basis of chemical similarity, we assumed that Cr randomly substitutes the Mn-and Fe-site, whereas Sn randomly substitutes the Si-or Ge-sites. The trained ML models consistently predicted Cr to decrease the T C of FeGe (best in the training set) at a rapid rate, i.e., the ML predicted T C values for Fe 0.9 Cr 0.1 Ge and Fe 0.8 Cr 0.2 Ge with increasing Cr content were 220 6 14.4 and 178.4 6 19.9 K, respectively. On comparing this result with the ML predicted T C for unsubstituted FeGe (261.1 6 13.2 K), it can be inferred that Cr substitution is unfavorable for improving the T C of FeGe compound. However, Cr substitution increased the T C of MnSi, i.e., the ML predicted T C values for Mn 0.9 Cr 0.1 Si and Mn 0.8 Cr 0.2 Si were 32.9 6 9.1 and 44.9 6 16.4 K, respectively. Furthermore, Cr substitution only had a marginal or no effect on the T C of MnGe (the ML predicted T C values for Mn 0.9 Cr 0.1 Ge and Mn 0.8 Cr 0.2 Ge were 141.5 6 36.6 and 141.54.9 6 45.3 K, respectively). Large ML prediction uncertainties indicate lack of knowledge on the (Mn,Cr)Ge material system.
In contrast, Sn-substituted FeGe alloys showed promise. The ML predicted T C values for FeGe 0.9 Sn 0.1 , FeGe 0.8 Sn 0.2 , and FeGe 0.75 Sn 0.25 were 259 6 13.2, 257 6 13.2, and 255 6 13.3 K, respectively. Intriguingly, Sn substitution did not have any discernible impact on the T C of FeGe alloys. Unlike the Cr substitution of FeGe that resulted in a sharp T C decrease, the addition of Sn is predicted to only marginally decrease the T C  of FeGe compound. This outcome warrants an explanation. It is important to recognize that extrapolation is considered a challenging problem for ML methods, especially when these methods are used without an adaptive learning scheme or iterative feedback loop [43]. Balfer and Bajorath [44] showed that when SVR methods are used for building quantitative structure-activity relationships, these methods systematically underestimate the true value of high-potency compounds due to the regularization term in the loss function that balances the trade-off between model complexity and its ability to generalize to an unseen data point. Therefore, our ML prediction for Snsubstituted FeGe needs further attention and cannot be discarded solely on the basis of ML predictions. Similar to Cr substitution, the addition of Sn increased and only had a marginal effect on the T C for MnSi and MnGe, respectively. The ML predicted T C values for MnSi 0.9 Sn 0.1 , MnSi 0.8 Sn 0.2 , MnGe 0.9 Sn 0.1 , and MnGe 0.8 Sn 0.2 were 32.5 6 8.7, 41.6 6 9.9, 141.9 6 28.6, and 140.6 6 28.4, respectively. The next step involved validating the data-driven predictions of Snsubstituted B20 alloys using DFT calculations.
A new set of DFT calculations were then performed to further understand the predictions from the data-driven ML models (see "Methods" for more details). A Heisenberg Hamiltonian was used to describe the magnetic exchange energy difference in the ferromagnetic ground state at T 5 0 K, which can be written as follows [45]: where DE is the total energy difference between ferromagnetic and antiferromagnetic spin states, z is the coordination number for (oppositely) ordered moments, S is the local magnetic moment at the atomic sites, and J is the exchange interaction energy between two spins. Thus, Eq. (1) relates DE and J. When the value of J is known, one can use the well-known mean field approximation formula to estimate the paramagnetic Curie temperature [45,46]: where k B is the Boltzmann's constant. From Eqs. (1) and (2) [44]. On the other hand, the promising prediction of Sn-substituted FeGe from DFT indicates that the FeGe 0.75 Sn 0.25 composition is worthy of experimental investigation. In the literature, Seow and Ziegler [47] suggested a remedial measure to overcome the SVR underprediction problem. This involves synthetically increasing the proportion of the high value points using bootstrapping or oversampling, so that the relative proportion of the high value data points becomes large. Seow and Ziegler demonstrated the efficacy of this approach on an environmental engineering problem. However, this method was not considered for the current study because of the smaller size of the training dataset. The skyrmion dataset is at least thirty times smaller than the environmental engineering problem, which had 638 data points to train the SVR models. Therefore, experimental validation remains one of the promising avenues to test the underprediction problem and provide feedback for ML model improvement.
In addition to the known B20 compounds, Table IV also lists a prediction for one more compound, namely FeSn, in the hypothetical B20 crystal structure. The novelty of this result stems from the fact that this compound has never been synthesized in the B20 crystal structures. Its ML predicted T C of 234.8 6 14.2 is greater than that of MnGe (which has an experimentally determined T C of 170 K), but is lower than that for FeGe 0.75 Sn 0.25 . Another encouraging aspect of bulk FeSn compound is that it has been experimentally synthesized, although in a centrosymmetric B35 crystal structure (space group, P6/mmm). To determine if the B20 structure is energetically close to B35, the total energies for FeSn in both B20 and B35 crystal structures were calculated. The B20 structure was found to be 1211.8 meV/atom higher in energy than the B35 crystal structure. Thus, FeSn in B20 structure is predicted to be highly metastable. Similar to FeSn, FeGe has also been synthesized in the B35 crystal structure [48]. Nevertheless, the B20 FeGe is only 117.6 meV/atom higher in energy than the B35 crystal structure.
As discussed earlier, nonequilibrium synthesis routes [12] were able to surmount the energy barrier and stabilize FeGe in the B20 crystal structure.

Conclusions
A data-driven computational approach is demonstrated to accelerate the design of novel chiral magnets in the B20 crystal structure with targeted T C property. This work highlights how ML methods can be used to efficiently guide the more computationally intensive physics-based DFT calculations toward promising directions in the composition space. This is especially beneficial when the costs of running the high-fidelity calculations are high. New alloy compositions are identified as promising with likely improved T C values, namely (Cr, Mn)Si, Mn(Si, Sn), and Fe(Ge, Sn), relative to the parent compounds, MnSi and FeGe. These results motivate new experiments to validate the computational predictions and provide feedback for model improvement. One of the major outcomes of this work is the prediction of Fe(Ge, Sn) with improved T C values, which can have key implications in the design of skyrmions for room temperature spintronics applications.

Density functional theory
DFT calculations were performed within the generalized gradient approximation, as implemented in Quantum ESPRESSO [49]. The PBEsol exchange-correlation functional [50] was used and the core and valence electrons were treated with ultrasoft pseudopotentials [51]. The Brillouin zone integration used Marzari-Vanderbilt smearing [52], with a smearing width of 0.27 eV and a 12 Â 12 Â 12 Monkhorst-Pack [53] k-point mesh centered at C. We used 60 Ry plane-wave cutoff for wave functions and 600 Ry kinetic energy cutoff for charge density and potential. The scalar relativistic pseudopotentials were taken from the PSLibrary [54]. The atomic positions and the cell volume were allowed to relax until an energy convergence threshold of 10 À8 eV and Hellmann-Feynman forces less than 2 meV/Å were achieved. For the initial accelerated screening step, collinear ferromagnetic spin configurations were imposed on the A-site transition metal atoms, where all spins are constrained to point in the same direction.
The magnetic exchange energy difference (DE) was calculated as the total energy difference between the ferromagnetic and the lowest energy antiferromagnetic spin states. Several antiferromagnetic spin configurations were explored and the one with the lowest energy was taken as the most favorable antiferromagnetic spin configuration for the calculation of DE. A 2 Â 2 Â 1 supercell with 32 atoms was considered for the DE calculation. In the case of calculations that involved Sn substitution, four Si or Ge atoms were substituted by the Sn atoms. The Sn-Sn bond distances were maintained at 4.639 and 4.715 Å in the supercells of MnSi 0.75 Sn 0.25 and FeGe 0.75 Sn 0.25 , respectively.

Machine learning
All ML analyses were formed using the R statistical environment [55]. The SVR was performed using the method implemented in the e1071 package [56]. More details about the SVR method can be found in the literature [57,58]. We used the nonlinear radial basis function (RBF) kernel to establish a quantitative relationship between the input descriptors (see Table II) and experimentally determined T C . Error bars for each prediction were estimated using the bootstrap resampling method [59,60]. The hyperparameters for the SVR were optimized using leave-one-out cross-validation method for each bootstrap sample. The trained SVR models are used to predict the T C for unexplored alloy compositions. From 100 SVR models, we get 100 predicted T C values for each composition. The mean (l) and standard deviations (error bar, r) are estimated from the 100 SVR RBF models.