Synthetic accessibility and stability rules of NASICONs
Bin Ouyang, Jingyang Wang, Tanjin He, Christopher J. Bartel, Haoyan Huo, Yan Wang, Valentina Lacivita, Haegyeom Kim, Gerbrand Ceder
1 Synthetic accessibility and stability rules of NASICONs
Bin Ouyang † , Jingyang Wang † , Tanjin He , Christopher J. Bartel , Haoyan Huo , Yan Wang , Valentina Lacivita , Haegyeom Kim and Gerbrand Ceder * Materials Sciences Division, Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA Department of Materials Science and Engineering, University of California, Berkeley, CA 94720, USA Advanced Material Lab, Samsung Research America, 10 Wilson Rd., Cambridge, MA 02138, USA † These authors contributed equally to this work * corresponding author: Prof. Gerbrand Ceder (Email: [email protected]) Abstract
In this paper we develop the stability rules for NASICON-structured materials, as an example of compounds with complex bond topology and composition. By applying machine learning to the ab-initio computed phase stability of 3881 potential NASICONs we can extract a simple two-dimensional descriptor that is extremely good at separating stable from unstable NASICONS. This machine-learned “tolerance factor” contains information on the Na content, the radii and electronegativities of the elements, and the Madelung energy. We test the predictive capability of this approach by selecting six predicted NASICON compositions. Five out of the six resulted in a phase pure NASICON while the sixth composition led to a NASICON that coexisted with other phases, validating the efficacy of this approach. This work not only provide tools to understand synthetic accessibility of NASICON-type materials, but also demonstrate an efficient paradigm for discovering new materials with complicate composition and atomic structure. Introduction
The advancement of high-throughput computation and data mining techniques (
1, 2 ) now enables us to rapidly search for new materials with interesting properties (
1, 3-6 ). This success has now become limited by our lack of understanding which predicted new materials are experimental accessible. A trial-and-error approach to synthesis is tedious, if not intractable, due to the exponential scaling of compositional possibilities in a multi-component spaces ( ). Therefore, the development of a rapid approach to evaluate the stability of compounds with complicated composition and structure would be a crucial step to accelerate materials discovery. Stability is a non-local property in the composition space of a system as it is a competition between the energy of a compound, and a large number of possible combinations of competing phases, usually at different compositions (
8, 9 ). For this reason, models that solely focus on reproducing the formation energy of a compound have limited ability to predict its stability ( ), and phenomenological models that focus on physical factors seem to be better at capturing phase stability. One of the best known such models, the Goldschmidt tolerance factor, proposed by Victor Goldschmidt in 1926 ( ), uses the ratio of ionic radii to estimate the stability of a perovskite. Similar models have been formulated for other materials with relatively simple bonding topology and composition ( ). For materials with more complicated bonding topology and composition, however, it is more challenging to predict the stability with one simple physical factor such as ionic radius, electronegativity, or by empirical rules such as the Pauling’s rules ( ). Sodium (Na) Super Ionic Conductor (NASICON) compounds are a good example of materials with complex composition and bonding topology. These materials have the general chemical formula Na x M ( A O ) and consist of a 3D corner-sharing framework made of A O tetrahedra and M O octahedra ( Fig. 1A ). Within the polyanion framework, there are four different cation sites, as illustrated in
Fig. 1B . All the oxygens sit in the 36f Wyckoff positions. In such materials, the Na content ( x ) can vary from 0 to 4, and a variety of cations can coexist in both the M site and A site ( Fig. 1C ) ( ). Given such compositional and bonding complexity, multiple physical factors can affect the phase stability, including an element’s stability in the cation site environment, and the bond compatibility among different coordination environments. As a result, any stability rules for NASICON materials are likely to be nontrivial. We focus on NASICONs as an example because their high ionic conductivity makes them of interest for many applications including selective ion membranes, gas sensing devices, and electrodes or solid-state electrolytes for Na-ion batteries ( ). Since the first report of the prototype NASICON Na Zr (SiO ) (PO ) in the 1960s ( ), the search for new NASICONs has been hindered by the large search space, and most of the successfully synthesized new materials are compositionally similar to the original Na Zr (SiO ) (PO ), as illustrated by the data in Fig. 1C . In the current work, we developed and applied a suite of materials discovery, physical interpretation, and machine-learning tools to uncover the stability rules of NASICON-type materials across a very broad compositional space. We first performed high-throughput ground state diagram calculations to sample the energy of NASICON compounds in the chemical space Na- M - M - A - B -O. The high-throughput screening covered 21 metal elements and evaluated 3881 compositions, only 32 of which have been previously explored. Guided by the high-throughput phase stability predictions, we successfully synthesized 5 out of 6 predicted new NASICONs in the Na-rich silicon phosphate sub-space. We also rationalize the stability trend of NASICON in terms of bond compatibility and site miscibility. Finally, to enable the prediction of stability from basic physical properties, we combined Sure Independent Screening (SIS) ( ) and machine-learned ranking (MLR) (
38, 39 ) to identify the best 2D descriptors that can separate NASICONs that are likely to be synthesized from unstable NASICONs. With the constructed features, 𝑡1 = %𝑁 !" ! + 𝑋 and 𝑡2 = 𝐸 ’(")&* ∙ 𝑋 +!" ∙ 𝑅 +$%& , where 𝑁 !" being Na content, 𝑋 being standard deviation of A site electronegativity, 𝐸 ’(")&* being Ewald summation of the electrostatic energy, 𝑋 +!" being electronegativity difference between M site and Na, 𝑅 +$%& being the standard deviation of the radius in M site, the synthetic accessibility with solid-state synthesis at 1000K can be estimated using the simple relationship as . Fig. 1. Structure and chemical space of NASICON. (A)
A Rhombohedral unit cell of NASICON-structured Na M ( A O ) with blue units represent M O octahedra, grey A O tetrahedra, and yellow/orange spheres Na ions. (B) Local structure showing distinct cation sites. (C)
The number of reports of NASICON containing a specific element. Compounds are either grouped by the central cation A in the polyanion group (left), or by the redox-inactive metal M (right). The color bar shows the maximum observed site occupancy of each element in the literature. The tolerance factors we developed for NASICONs is one of the first machine-learned phenomenological model to estimate the stability of complex oxides within a large combinatorial space. Our work will thus not only facilitate the search for NASICON-type materials, but from a broader perspective, it will also provide insight to help apply the state-of-the-art machine-learning techniques and data-mining tools to accelerate the discovery of new inorganic materials with complicated compositions and bond topologies.
Results Stability map of NASICON based on high-throughput DFT calculations
To computationally explore NASICONs with the chemical formula Na x M y M’ (AO ) z (BO ) , we sampled two compositional variables: the Na content x from 0.0 to 4.0 (in steps of 0.5) and the average M/M’ charge state from +2 to +5 (in steps of 0.5). For each sampled Na content and M/M’ charge state, the possible combinations of (SiO ) , (PO ) and (SO ) polyanions resulting in a charge-balanced composition were enumerated, leading to 3881 possible NASICON compounds. The phase stability of all compounds was evaluated by computing the phase diagram in their respective compositional spaces, taking the ground states in the Materials Project ( ) as possible competing phases. Relative stability is quantified by the energy above the convex hull ( 𝐸 ,-)) ) ( ). Fig. 2. The thermodynamic stability of 3881 calculated NASICONs. (A)
Stability map. The color of each block in the lower triangle represents the median E hull (bottom-left) of NASICONS that contain the metals on the x and y axis. The upper triangle shows the fraction of compounds ( f accessible ) with E hull – S ideal × ≤ (B) The distribution of E hull (top) and f accessible (bottom) as a function of Na content and polyanion composition. Polyanions are represented by their central cation, e.g., Si P indicates (SiO ) (PO ). Green (red) indicates good (poor) NASICON stability in (A) and (B). In Fig. 2 , the effect of chemical composition on the phase stability of the NASICON compounds is described with two factors: the median of the E hull values in the relevant compositional space, and the frequency ( 𝑓 accessible ) by which E hull of the NASICON compound is less than S Ideal × S ideal is the ideal entropy of mixing as explained in Supplementary Note 2 . The second factor gauges which NASICONs are close enough to the 0 K hull so that they may be stabilized by configurational entropy at the synthesis temperature. These two factors are presented in
Fig. 2 . In
Fig. 2A , the lower triangle corresponds to the median E hull values for all the NASICONs that contain both the metal on the x -axis and the metal on the y -axis. The upper triangle shows the corresponding 𝑓 "../0012)/ values. Fig. 2B present the median E hull values and 𝑓 "../0012)/ for all NASICONs with a specific Na content and polyanion composition. The results in Fig. 2A indicate that NASICONs containing metals such as Hf , Zr , Ta and Sc have higher stability and compatibility with other metals, whereas Ca , Zn , Ge and La are metal species that are difficult to stabilize in a NASICON structure. Moreover, as indicated by the green (dark) region in the upper (bottom) panel of Fig. 2B , NASICONs are relatively more stable with low Na content or with a single type of polyanion. In particular, mixing different polyanions appears to have a negative effect on stability.
Table 1. Statistics of calculated NASICONs.
All calculated NASICONs are separated into three categories: thermodynamic ground states (GS) ( 𝐸 )*++ ≤ 0 ), likely synthesizable (LS) ( )*++ ≤ 𝑆 ,-./+ × 1000𝐾 ) and unlikely synthesizable (US) ( 𝐸 )*++ > 𝑆 ,-./+ × 1000𝐾 ). Materials category Number/All Percentage
Thermodynamic ground state (GS) ( E hull ≤ E hull ≤ S ideal × E hull > S ideal × Distribution of GS/LS-NASICONs ≤ N Na ≤
4 60/1199 5.0% 1< N Na <
3 127/1576 8.1% N Na ≤
1 454/1106 41.0% Pure phosphate/silicate/sulfate 225/683 32.9% Mixed polyanions 416/3198 13.0%
To obtain more statistical information about the distribution of E hull among the calculated compositions, all the NASICONs are separated into three categories: thermodynamic ground states (GS) ( 𝐸 ,-)) ≤ 0 ), likely synthesizable (LS) ( ,-)) ≤ 𝑆 × 1000𝐾 ) and unlikely synthesizable (US) ( 𝐸 ,-)) > 𝑆 × 1000𝐾 ). The fractions of the three NASICON groups are shown in Table 1 , with the GS- and LS-NASICONs only making up 6.3% and 10.2% of the computed compositions, respectively. The distribution of GS/LS-NASICONs across Na content and polyanion compositions is also shown in Table 1 . In total, 41.1% of the NASICONs with 𝑁 !" ≤ 1 are in the GS/LS group, which is more than eight times the GS/LS percentage in the !" ≤ 4 region. In addition, the percentage of GS/LS-NASICONs with a single type of polyanion is three times higher than that of mixed polyanions NASICONs. Experimental validation of the predicted Na-rich mixed-polyanion NASICONs
To identify which NASICON compositions may have a higher probability of being interesting fast-ion conductors, we extracted the experimentally measured conductivity vs. Na content x for more than two hundred NASICONs from the literature and plot the result in Fig. 3A ( ). In general, the conductivity increases with Na content until x ≈ 3. For this reason, we focus our synthesis attempts on novel NASICONs with Na content. The top panel of
Fig. 3B shows how many predicted GS/LS-Na -NASICONs contain a specific M-M’ pair . A comparison with the bottom panel, which shows the same for experimentally reported NASICONs, indicates that there appears to be many promising candidates beyond the well-known Zr-containing NASICONs ( ). Therefore, in the subspace of Na -NASICONs, there is a large compositional space that remains undiscovered, within which 60 GS/LS compounds are predicted. To further narrow down targets for experimental validation, we focused on the Na -NASICONs with mixed SiO and PO because such NASICON compositions are 1) more likely to be superionic conductors ( ) and 2) more challenging to synthesize as they have both high Na content and mix polyanions, ( Fig. 2B , Table 1 ) which makes them good candidates to test the accuracy of our prediction. This selection criterion led us to attempt to synthesize the six NASICONs highlighted with a red rectangle in the bottom panel of
Fig. 3B . Detailed information on these phases, such as E hull and competing phases are listed in Table S1 . Fig. 3. Solid-state synthesis exploration of predicted NASICONs within Na M y M’ Si z P O subspace. (A) experimentally measured conductivity in log scale vs. Na content of NASICONs at room temperature and 300 ℃ as extracted from the experimental literature. The upper curve shows the relative distribution of experimentally reported NASICONs as a function of Na content. (B) Distribution of DFT predicted GS/LS-NASICON (top panel) and experimentally reported NASICON compositions (bottom panel) with Na stoichiometry. The six experimentally attempted compounds in this work are highlighted with a red box. The number of NASICONs in each block is indicated by the color bar and specified by the number. (C) X-ray diffraction patterns of the six experimentally attempted NASICON compounds, five of which are pure NASICON phase or with a trace amount of impurity phases such as (Hf/Zr)O (triangle), NaInP O (circle). Fig. 3C shows the diffraction pattern for the six NASICON compounds for which we attempted synthesis. Five of the six target NASICONs were obtained phase-pure or with only a trace amount of impurity phase. Na HfSn(SiO ) (PO ) was synthesized with a large fraction of SnO impurity, indicating that the solubility of Sn may be limited under the synthetic conditions attempted. However, the presence of a NASICON phase in the diffraction pattern suggests that this composition may be synthesizable under different conditions (i.e. higher temperature, longer calcination time) or by an alternative synthesis route (
60, 61 ). Physical trends of stability in NASICON
The stability trend of NASICONs across composition space can be rationalized by analyzing the bonding compatibility and site miscibility within the two bond complexes in the structure, as schematically illustrated in
Fig. 4A . In the M -O- A bond complex, where M is the metal and A is the cation in the polyanion, both M and A ions compete for covalence with oxygen. In the M -O-Na bond complex, the electropositive nature of Na causes electrons to be mainly localized within M -O bonds. Because NASICONs can exist with mixed species on both the metal site ( M ) and polyanion central site ( A ), site miscibility (e.g, bond distortion due to size effect) is an additional factor influencing stability. Three electronegativities are considered to capture the bonding compatibility: the average electronegativity in the metal site 𝑋 + , the average electronegativity of the central cation in the polyanion 𝑋 , and the sodium electronegativity 𝑋 !" . The distribution of 𝑋 + − 𝑋 !" and 𝑋 − 𝑋 + for all calculated NASICONs are plotted in Fig. 4B . The orange triangles mark the computed NASICONs, while the blue lines define contours along which the fraction of GS/LS-NASICONs is constant, with higher fraction of stable NASICONs towards the center of the contour lines. Two groups of NASICONs are highlighted in Fig. 4B to illustrate the strong correlation between electronegativity and phase stability. The first group contains three NASICONs with stoichiometry Na , i.e., A (NaHfMg(SO ) (PO )), A1 (NaHfCa(SO ) (PO )) and A2 (NaGeMg(SO ) (PO )); while the second group contains three NASICONs with stoichiometry Na , i.e., B (Na HfZr(SiO ) (PO )), B1 (Na HfGe(SiO ) (PO )) and B2 (Na Ge (SiO ) (PO )). Compositions A and B are representatives of NASICONs with optimal metal electronegativity for stoichiometry Na and stoichiometry Na , respectively, and appear at relatively central positions of the contour lines in Fig. 4B . A1/A2 and B1/B2 are obtained by varying the electronegativities of the metal sites. In the first comparison group, the metal site is made more electropositive as Mg (A) is replaced by Ca (A1) (A à A1), while it becomes more electronegative with a substitution of Hf (A) to Ge (A2) (A à A2). Both substitutions increase the E hull value. More specifically, A à A1 increases E hull from 1.0 to 48.7 meV/atom, while A-A2 increases E hull from 1.0 to 54.8 meV/atom. In the second group, substitution of Zr by Ge (B à B1), or both Hf and Zr by Ge (B à B2) makes the metal site more electronegative. As a result, E hull rises from 4.0 meV/atom in B, to 56.0 and 90.3 meV/atom in B1 and B2. These results show that deviation from the optimal electronegativity in the metal site leads to decreased NASICON stability. Fig. 4. Physical factors that influence NASICON stability , (A) Schematic of two main factors responsible for NASICON stability, i.e., bond compatibility and site miscibility. (B)
Distribution of the computed NASICON compositions (orange triangles) in the space defined by two relevant electronegativity differences. x-axis: the difference between the average electronegativity of metal ( 𝑋 ) and Na ( 𝑋 ); y-axis: the difference between the average electronegativity of metal ( 𝑋 ) and polyanion central cation ( 𝑋 ). The probability distribution of GS/LS-NASICONs are plotted with contour lines. The contour at the center is the higher probability one. Two groups of representative NASICONs are also highlighted, i.e., A (NaHfMg(SO ) (PO )) ( E hull =1.0 meV/atom), A1 (NaHfCa(SO ) (PO )) ( E hull = 48.7 meV/atom), A2 (NaGeMg(SO ) (PO )) ( E hull =54.8 meV/atom). The other group contains three Na rich NASICONs, i.e., B (Na HfZr(SiO ) (PO )) ( E hull =4.0 meV/atom), B1 (Na HfGe(SiO ) (PO )) ( E hull =56.0 meV/atom), B2 (Na Ge (SiO ) (PO )) ( E hull =90.3 meV/atom); (C) The corresponding Integrated Crystal Orbital Hamilton Population (-ICOHP) per metal for M-O bonds of two groups in (B) . (D) The distribution of NASICONs as a function of metal-oxygen bond deviation ( 𝜀 ) and central cation-oxygen bond deviation ( 𝜀 ) in the polyanions; The color code applied in all NASICONs and GS/LS-NASICONs are the same as (B) . ( E) The distribution of bond deviations 𝜀 and 𝜀 for four groups of NASICONs, i.e., NASICONs with no mixing, with only cation mixing, with only anion mixing and with both mixing of cations and anions. To understand the electronic origin of the destabilization that occurs when moving away from the ideal electronegativity differences, the -ICOHP (Integrated Crystal Orbital Hamilton Population) ( ) values for M -O bonds were calculated for both groups of NASICONs. As shown in Fig. 4C , with cation substitution A à A1, A à A2, or B à B1 à B2, all the -ICOHP values drop, indicating weakened M -O bonds and thus higher E hull values as mentioned above. It is worth noting that the calculated -ICOHP values of A-O bonds in the polyanion group do not change much upon cation substitution ( Supplementary Note 3 ), which can be understood by the strong covalency within the tetrahedron polyanion group. The ICOHP data strengthens the idea that the NASICONs are stabilized by bond compatibility, as illustrated by the cartoon in
Fig. 4A : The species in the M site should bond to O with reasonable covalency (setting a lower limit to the electronegativity) to stabilize the polarized M -O-Na complex; however, M-O cannot be too covalent (setting an upper limit for electronegativity) otherwise it is not compatible with the A -O bonds in the polyanion group. Size differences of the cations that are mixed will also influence phase stability. In general, such a size effect is destabilizing when it becomes too large ( ). In order to understand the influence of site miscibility on the NASICON stability, the bond length deviations extracted from the DFT computation are analyzed in
Fig. 4D . The variations of the bond lengths d , ( 𝑑 +67"8 −𝑑 +6719 )/𝑑 +67"8 , are denoted as 𝜀 +6:";1"%1<9 and 𝜀 for M-O and
A-O bonds respectively. The fraction of GS/LS-NASICONs are plotted with contour lines. It can be observed in
Fig. 4D that GS/LS-NASICONs are mostly found in the region with less bond deviations, consistent with the general idea that the size difference of ions mixed on a site contributes a positive amount (favoring unmixing) to the energy of formation ( ). The mixing effects introduced by cation and polyanion are further distinguished and plotted in Fig. 4E . Clearly, the maximum lattice distortion is observed when both cation and polyanion mixing are present; while with only polyanion mixing, both 𝜀 +6:";1"%1<9 and 𝜀 are generally larger than when only M cations are mixed. This explains the distinct destabilization observed in Fig. 2 which is much more pronounced for mixed polyanion NASICONs than for mixed metal NASICONs. It is possible that the A -O bond is less tolerant to length variation and therefore has to pass size effects on as distortions in the M -O bonds. The ability of in particular d metals to be accommodated in distorted octahedra has recently been explained ( ). Machine-learned tolerance factor for NASICONs
Having established that bond compatibility and site miscibility are two major physical factors responsible for the stability of NASICONs, a useful phenomenological model that is able to capture those physical factors and quantify them without additional DFT calculations can further facilitate the search for new NASICON compositions. Such model should involve all the important physical information, while at the same time be simple enough to interpret and analyze. Therefore, we combine sure independent screening (SIS) ( ) and machine-learned ranking (MLR) (
38, 39 ) to search for the best mathematical representation of NASICON stability (
37, 69 ). We limit the model to a 2D representation so that it can easily be visualized and interpreted. Starting with 24 basic physical properties (primary features) that include elemental, compositional and structural information of specific NASICONs, we iteratively apply 17 mathematical operators to generate increasing complex features until we have a set of 1,895,020 (details in details in the
Methods and
Supplementary Note 3 ) candidate features. SIS is performed to select the top 2000 features that correlate the most with E hull values. The 2D descriptors pool is then formed by the combinatorial space of these features (1,999,000 in all). We perform MLR to identify the 2D descriptors that has the best capability to rank the stability ( E hull values) of NASICONs. Such SIS+MLR process can be regarded as a variation of the SISSO method developed by Ouyang et. al.( ) except that in our case MLR would optimize the ranking of relative stability, which is more important than the absolute 𝐸 ,-)) values for evaluating synthesizability of materials. The SIS+MLR-learned 2D descriptors for NASICON stability are plotted in Fig. 5A . The first dimension of the feature is formulated as 𝑡1 = %𝑁 !" ! + 𝑋 , where 𝑁 !" is the Na content and 𝑋 is the standard deviation of electronegativity of the central cation in polyanion. This term thus evolves with respect to Na content and anion mixing. The second dimension is formulated as 𝑡2 = 𝐸 ’(")&* ∙ 𝑋 +!" ∙ 𝑅 +$%& , where 𝐸 ’(")& is the Ewald energy of the lattice model normalized by that in Na Zr (SiO ) (PO ) (details in Supplementary Note 4 ), 𝑋 +!" is the metal site electronegativity taking the Na electronegativity as a reference, and 𝑅 +$%& is the standard deviation of the cation radii normalized by the ionic radius of Na + . This term thus captures the cation mixing effect. Here we used non-depersonalized 𝐸 ’(")& and 𝑅 +$%& to keep the MLR and later classification process physically meaningful. These two features extracted by the ML algorithm incorporate the physical factors identified in Fig. 4 . As shown in
Fig. 5A , most of the US-NASICONs can be well separated from the GS/LS-NASICONs in the identified feature space. To provide practical guidance for experimental synthesis, a linear decision boundary to distinguish the GS/LS-NASICONs from US-NASICONs can be drawn with the assistance of support vector classification (SVC), as illustrated in
Fig. 5A (details in
Methods ). More specifically, a binary classification model with this boundary indicates that with the machine-learned t1 and t2, a synthetically accessible NASICON at 1000K should satisfy the condition . Applied to all 3881 compositions, such stability model has an overall classification accuracy of 81.8%, as well as a recall value of 96.6% and F1 score of 88.2% for capturing GS/LS-NASICONs. To test whether the aforementioned stability model can capture the stability trend with diverse composition variations, we evaluated three measures of its performance: 1) its ability to capture experimentally reported NASICONs ( Fig. 5B ); 2) its robustness in identifying GS/LS-NASICONs in both Na rich and Na poor regions (
Fig. 5C ); 3) its robustness in identifying GS/LS-NASICONs with different metal chemistries (
Fig. 5D ). To test the capability of capturing experimentally synthesized NASICONs, we applied the stability model to classify 27 previously reported and 5 newly synthesized NASICONs (shown in
Fig. 5B ). The stability model correctly captured 20 of the 32 experimentally synthesized NASICONs as synthetically accessible (accessible probability >50%). In addition, the compounds classified as false negatives are still given a relatively large accessible probability by the model( ). To be more specific, seven out of the twelve false negative predictions have ≥
40% probability to be accessible while the most far-off prediction still has 18% probability to be accessible (
Supplementary Note 5 ). Those observations reveal that the stability of all experimentally reported NASICONs is reasonably captured by such a stability boundary. As the fraction of GS/LS-NASICONs are quite different for Na-poor and Na-rich compositions (
Fig. 2B and
Table 1 ), we also tested the ability of the model to identify GS/LS-NASICONs for both Na and Na compositions separately. The distribution of GS/LS Na and Na compounds and corresponding accessible probabilities (P accessible , indicated by the colormap) predicted by the above-mentioned SVC model are demonstrated in Fig. 5C (The distributions of all Na and Na NASICONs are shown in
Supplementary Note 6 ). The P accessible is derived from the Platt scaling probability calculated from trained SVC model. It is clear that most of the GS/LS-NASICONs in both Na and Na group are well-captured by the demarcation boundary. Moreover, the false negative classifications still remain likely to be synthetically accessible, as the average accessible probability of the false negative predictions is 30.9% for Na NASICONs and 39.0% for Na NASICONs.
Fig. 5. Machine learned tolerance factor of NASICONs and distribution of different groups of NASICONs in this feature space.
The white and grey shaded regions in each plot are spaces predicted to be GS/LS-NASICONs and US-NASICONs respectively from SVC. The symbols in x and y represent the cube root of Na contents ( ,𝑁 ! ), standard deviation of anion electronegativity ( 𝑋 ), square of Ewald summation ( 𝐸 =>/+-? ), electronegativity difference between average cation electronegativity and Na electronegativity ( 𝑋 ) and standard deviation of cation radius 𝑅 ). (A) Distribution of all NASICONs in the machine-learned 2D space. The computed E hull -S max T is represented by the color of the data point. (B) Reported NASICONs and newly synthesized NASICONs plotted in the machine-learned 2D space; (C)
Na rich (Na /f.u.) and Na poor (Na /f.u.) GS/LS-NASICONs plotted in the machine- learned 2D feature space, with the accessible probability predicted by SVC indicated by the color of the data point; (D) Ca/Ge or Zr/Hf containing GS/LS-NASICONs in the machine-learned 2D space, with the accessible probability predicted by SVC indicated by the color of the data point.
Similarly, we have plotted the GS/LS-NASICONs for Ca/Ge or Zr/Hf containing compositions in
Fig. 5D (the distribution of all NASICONs containing these two metals is also plotted in
Supplementary Note 6 ). The Hf/Zr group is mostly distributed in the space where NASICONs are predicted to be stable, thereby serving as a good contrast to the Ca/Ge group which has much lower GS/LS frequency (
Fig. 2 ). As shown by
Fig. 5D , the GS/LS-NASICONs of both groups are well captured, even the misclassified ones are not far off given that the average probability of being accessible is 36.6% for Hf/Zr group and 29.5% for Ca/Ge group. Based on the performance as shown in
Fig. 5 , it can therefore be concluded that the stability model developed should be able to offer reasonably stability prediction for NASICONs regardless of Na contents and metal chemistries.
Discussion
The efficient discovery of inorganic solids with complexity in both composition and crystal structure is a challenging task. While simple modifications to a basic compound are prevalent in the literature, as is the case for NASICONs (
32, 59, 71, 72 ), there is a lack of thorough understanding of the stability rules to help one navigate through the unexplored chemical space. In this work, we conducted high-throughput phase diagram calculations to evaluate 3881 possible NASICON compositions. We not only successfully synthesized five out of six Na-rich silicon phosphates predicted to be synthesizable, but also rationalized the variation in thermodynamic stability as a function of Na content, metal species and polyanion mixing. Our analys together with a machine-learned predictive phenomenological model points out several key factors for NASICON stability:
1) The electronegativity of a metal is bounded on both sides for it to be stable in a NASICON due to its covalency competition with the cation in the polyanion and the necessity to maintain Na ionicity. A metal that is too electronegative favors competing phases where those metals are present in a more covalent bonding environment. To give one example, Ge-containing NASICONs, such as Na HfGe(SiO ) (PO ) (B1) and Na Ge (SiO ) (PO ) (B2), tend to decompose into other phases with the Ge accommodated in Na Ge O . The calculated -ICOHP per Ge-O bond is 1.65 and 1.38 in B1 and B2 compared with 3.38 in Na Ge O (details in Supplementary Note 7 ), which supports the strengthening of Ge-O by decomposition. On the other hand, a metal that is too electropositive metal is also unstable in NASICON, as demonstrated by the Ca-containing NASICONs such as NaHfCa(SO ) (PO ) (A1). The small -ICOHP per Ca-O bond, 0.22, indicates a low covalency of the Ca-O bond, and its decomposition into CaSO removes the Na-Ca competition for ionization and enhances the strength of Ca-O bonds (-ICOHP = 0.64 per bond) (details in Supplementary Note 8 ). This concept of bond competition can be extended to other polyanionic ceramic materials which contain multiple metal sites. One important group of such materials is the alkali (earth) containing polyanion materials, which are widely studied as both electrodes and solid electrolytes in the energy storage field ( ). Those materials generally have two metal sites, one of which is the alkali (earth) metal site, while the other is typically a transition metal site whose polyhedron corner-shares with polyanions. In such materials, the alkali (earth) metal site would set the lower bound of electronegativity, while the upper bound is set by the central cation of polyanion in order to maintain the strong covalency within the polyanion group. 2) The size difference of the metal and polyanion cation also influences stability and the bond length variation analysis (
Fig. 4D ) demonstrates that the polyanion site is generally less tolerant to lattice distortion. Therefore, most of the misfit strain due to mixing polyanions has to be passed on to the cation site. If that is not sufficiently possible, the compound will form other phases to release the strain energy. Our analysis on NASICON materials thus explains why it is in general more difficult to create mixed polyanion systems. Our work also demonstrates a computational paradigm to evaluate the stability of complex ceramic materials with phenomenological models derived from machine learning. The complicated compositional variety and bond topology of such materials makes the development of a phenomenological stability model similar to Goldschmidt’s tolerance factor remarkably difficult. However, a relatively simple and quantitatively predictive model was obtained in this work by applying the sure independence screening combined with machine learning ranking algorithms to search through billions of possible mathematical representations. The successful derivation of a “tolerance factor” for NASICON provides an example of how effective such a toolset can be for establishing stability rules for complex oxides starting from a certain amount of DFT data. Once a simple machine-learned model is established, no additional DFT calculations are required for exploring the broader chemical space. Such a framework would be remarkably useful for multi-component ceramic materials with complex structures (
75, 76 ). Methods: Experiments
For the solid-state synthesis of NASICON compounds, typical metal oxides or hydroxides (HfO (Aldrich, 99.8%), MgO (Aldrich, ≥ O (Sigma, 99.9%), Zr(OH) (Aldrich, 97%), SnO (Alfa Aesar, 99.9%), CaO (Sigma-Aldrich, 99.9%)) were used as precursors to introduce metal cations. SiO (Sigma-Aldrich, nanopowder) and NaH PO (Sigma, >99%) were used as silicate and phosphate sources. Na CO (Sigma-Aldrich, >99%) was used as an extra sodium source. In addition, 10% excess NaH PO was introduced to compensate for the possible sodium and phosphate loss during the high-temperature treatment. The stoichiometric powder mixtures were wet ball-milled for 12h using a Planetary Ball Mill PM200 (Retsch) to achieve a thorough mixing before being pressed into pellets. The pelletized samples were first annealed at 900 ℃ under Ar flow, then grounded with mortar and pestle, wet ball-milled, pelletized and re-annealed at 1100 ℃ . The crystal structures of the obtained materials were analyzed using X-ray diffraction (Rigaku Miniflex 600 and Bruker D8 Diffractometer) with Cu Kα radiation. First-principles calculations
First-principles total energies calculations were performed with the Vienna ab initio simulation package (VASP) with a plane-wave basis set( ). Projector augmented-wave potentials( ) with a kinetic energy cutoff of 520 eV and the exchange-correlation form in the Perdew–Burke–Ernzerhof generalized gradient approximation (GGA-PBE)( ) were employed in all the structural optimizations and total-energy calculations. For all the calculations, a reciprocal space discretization of 25 k-points per Å −1 was applied, and the convergence criteria were set as 10 −6 eV for electronic iterations and 0.02 eV/Å for ionic iterations. A rhombohedral conventional cell was used for each of the NASICON structures, and the Na-vacancy ordering, cation ordering, and anion ordering were set as the one with lowest electrostatic energy. In addition to the 3881 NASICON structures, all the competing phases in the relevant chemical spaces that are given in The Materials Project( ) were also calculated to construct the phase diagrams. To estimate the synthetic accessibility of the NASICON at finite temperature, we assumed that a synthetic accessible NASICON could be made due to a configurational entropy stabilization, i.e., 𝐸 ,-)) − 𝑇𝑆 ≤ 0 . Therefore, the ideal entropy was calculated by assuming a fully disordered distribution of Na, cation and anion sites ( Supplementary Note 2 ). Machine learning setup
To consider as many physical factors as possible, 24 basic features covering the basic element properties of chemical species are constructed (a list of features is presented in
Supplementary Note 4 ). Combination of these features with 17 mathematical operators listed in
Supplementary Note 4 with a complexity of two leads to 1,895,020 possible features. SIS was then performed on 80% of the data, using the SISSO package developed by Ouyang et.al.( ) to identify the top 2000 features with highest Pearson correlation to E hull values. 3104 out of 3881 data points are generated from stratified sampling due to unbalanced sample amount between the GS/LS-NASICONs and US-NASICONs. With the selected 2000 features, the possible 2D models will be 𝐶 *===* = E hull values, as demonstrated in Fig. 5 . The objective function of the SVM is set as the accuracy of predicting the pairwise ranking matrix for model optimization. During the MLR process, we used the stratified sampling training set applied at SIS step for optimizing the hyper parameters. The validation accuracy was used as the target function to optimize through hyper parameter searching. With the SIS+MLR learned 2D descriptor space, a linear decision boundary of the synthetic accessible NASICONs was drawn with the assistance of support vector classification (SVC) with balanced weight given that the populations of GS/LS samples are not equivalent. References
1. A. Jain et al. , Commentary: The Materials Project: A materials genome approach to accelerating materials innovation.
Apl Materials , 011002 (2013). 2. C. C. Fischer, K. J. Tibbetts, D. Morgan, G. Ceder, Predicting crystal structure by merging data mining with quantum mechanics. Nature materials , 641-646 (2006). 3. G. Hautier et al. , Phosphates as lithium-ion battery cathodes: an evaluation based on high-throughput ab initio calculations. Chemistry of Materials , 3495-3508 (2011). 4. S. Curtarolo et al. , The high-throughput highway to computational materials design. Nature materials , 191-201 (2013). 5. W. Sun et al. , A map of the inorganic ternary metal nitrides. Nature materials , 732 (2019). 6. A. Jain, Y. Shin, K. A. Persson, Computational predictions of energy materials using density functional theory. Nature Reviews Materials , 1-13 (2016). 7. D. W. Davies et al. , Computational screening of all stoichiometric inorganic materials. Chem , 617-627 (2016). 8. C. J. Bartel, A. W. Weimer, S. Lany, C. B. Musgrave, A. M. Holder, The role of decomposition reactions in assessing first-principles predictions of solid stability. npj Computational Materials , 4 (2019). 9. S. P. Ong, L. Wang, B. Kang, G. Ceder, Li−Fe−P−O2 Phase Diagram from First Principles Calculations. Chemistry of Materials , 1798-1807 (2008). 10. C. J. Bartel et al. , A critical examination of compound stability predictions from machine-learned formation energies. npj Computational Materials , 97 (2020). 11. V. M. Goldschmidt, Die gesetze der krystallochemie. Naturwissenschaften , 477-485 (1926). 12. Z. Song, Q. Liu, Tolerance Factor and Phase Stability of the Normal Spinel Structure. Crystal Growth & Design , 2014-2018 (2020). 13. V. Stevanovic, M. d’Avezac, A. Zunger, Universal electrostatic origin of cation ordering in A2BO4 spinel oxides. Journal of the American Chemical Society , 11649-11654 (2011). 14. S. Müller, L.-W. Wang, A. Zunger, C. Wolverton, Coherent phase stability in Al-Zn and Al-Cu fcc alloys: The role of the instability of fcc Zn.
Physical Review B , 16448 (1999). 15. V. Stevanović, M. d’Avezac, A. Zunger, Simple point-ion electrostatic model explains the cation distribution in spinel oxides. Physical review letters , 075501 (2010). 16. L. M. Ghiringhelli, J. Vybiral, S. V. Levchenko, C. Draxl, M. Scheffler, Big data of materials science: critical role of the descriptor.
Physical review letters , 105503 (2015). 17. C. J. Bartel et al. , New tolerance factor to predict the stability of perovskite oxides and halides.
Science advances , eaav0693 (2019). 18. L. Pauling, The principles determining the structure of complex ionic crystals. Journal of the american chemical society , 1010-1026 (1929). 19. Y. Lu et al. , A High - Performance Monolithic Solid - State Sodium Battery with Ca2+ Doped Na3Zr2Si2PO12 Electrolyte.
Advanced Energy Materials , 1901205 (2019). 20. Q. Ma et al. , Scandium-substituted Na3Zr2 (SiO4) 2 (PO4) prepared by a solution-assisted solid-state reaction method as sodium-ion conductors. Chemistry of Materials , 4821-4828 (2016). 26 21. Y. Yue, W. Pang, Hydrothermal synthesis and characterization of NaSn2 (PO4) 3. Journal of materials science letters , 148-149 (1992). 22. P. Slater, C. Greaves, Synthesis and Conductivities of Sulfate/Selenate Phases Related to Nasicon: NaxM ′ (II) xM ″ (III) 2-x (SO4) 3-y (SeO4) y. Journal of Solid State Chemistry , 12-18 (1993). 23. J.-M. Winand, A. Rulmont, P. Tarte, Nouvelles solutions solides LI (MIV) 2− x (NIV) x (PO4) 3 (L= Li, Na M, N= Ge, Sn, Ti, Zr, Hf) synthèse et étude par diffraction x et conductivité ionique.
Journal of Solid State Chemistry , 341-349 (1991). 24. R. Cava, E. M. Vogel, D. W. Johnson Jr, Effect of homovalent framework cation substitutions on the sodium ion conductivity in Na3Zr2Si2PO12. Journal of the American Ceramic Society , c157-c159 (1982). 25. M. Zahir, R. Olazcuaga, P. Hagenmuller, Crystal chemistry and ionic conductivity in Nasicon-type phases Na1+ xZr2− xYbx (AsO4) 3 with 0 ⩽ x ⩽ Materials Letters , 234-236 (1984). 26. S. Fujitsu, M. Nagai, T. Kanazawa, I. Yasui, Conduction paths in sintered ionic conductive material Na1+ xYxZr2− x (PO4) 3. Materials Research Bulletin , 1299-1309 (1981). 27. J. B. Goodenough, H.-P. Hong, J. Kafalas, Fast Na+-ion transport in skeleton structures. Materials Research Bulletin , 203-220 (1976). 28. F. Cherkaoui, J. Viala, C. Delmas, P. Hagenmuller, Crystal chemistry and ionic conductivity of a new Nasicon-related solid solution Na1+ xZr2− x2Mgx2 (PO4) 3. Solid State Ionics , 333-337 (1986). 29. Y. Saito, K. Ado, T. Asai, H. Kageyama, O. Nakamura, Ionic conductivity of NASICON-type conductors Na1. 5M0. 5Zr1. 5 (PO4) 3 (M: Al3+, Ga3+, Cr3+, Sc3+, Fe3+, In3+, Yb3+, Y3+). Solid State Ionics , 327-331 (1992). 30. T. Takahashi, K. Kuwabara, M. Shibata, Solid-state ionics-conductivities of Na+ ion conductors based on NASICON. Solid State Ionics , 163-175 (1980). 31. Y. Shimizu, N. Yamashita, Solid electrolyte CO2 sensor using NASICON and perovskite-type oxide electrode. Sensors and Actuators B: Chemical , 102-106 (2000). 32. N. Anantharamulu et al. , A wide-ranging review on Nasicon type materials. Journal of materials science , 2821-2837 (2011). 33. Z. Jian, Y. S. Hu, X. Ji, W. Chen, Nasicon - structured materials for energy storage. Advanced Materials , 1601925 (2017). 34. P. Fabry, J. Gros, J. Million-Brodaz, M. Kleitz, NASICON, an ionic conductor for solid-state Na+-selective electrode. Sensors and Actuators , 33-49 (1988). 35. J. Fan, R. Samworth, Y. Wu, Ultrahigh dimensional feature selection: beyond the linear model. The Journal of Machine Learning Research , 2013-2038 (2009). 36. J. Fan, J. Lv, Sure independence screening for ultrahigh dimensional feature space. Journal of the Royal Statistical Society: Series B (Statistical Methodology) , 849-911 (2008). 37. R. Ouyang, S. Curtarolo, E. Ahmetcik, M. Scheffler, L. M. Ghiringhelli, SISSO: A compressed-sensing method for identifying the best low-dimensional descriptor in an immensity of offered candidates. Physical Review Materials , 083802 (2018). 38. T.-Y. Liu, Learning to rank for information retrieval. Foundations and Trends in Information Retrieval , (2009). 27 39. T. Joachims, in Proceedings of the eighth ACM SIGKDD international conference on Knowledge discovery and data mining . (2002), pp. 133-142. 40. Y. Deng et al. , Crystal structures, local atomic environments, and ion diffusion mechanisms of scandium-substituted sodium superionic conductor (NASICON) solid electrolytes.
Chemistry of Materials , 2618-2630 (2018). 41. E. R. Losilla et al. , Understanding Na mobility in NASICON materials: a Rietveld, 23Na and 31P MAS NMR, and impedance study. Chemistry of materials , 665-673 (1998). 42. Z. Khakpour, Influence of M: Ce4+, Gd3+ and Yb3+ substituted Na3+ xZr2-xMxSi2PO12 solid NASICON electrolytes on sintering, microstructure and conductivity. Electrochimica Acta , 337-347 (2016). 43. P. Yadav, M. Bhatnagar, Preparation, structure and conductivity of Sn modified NASICON material.
Journal of Electroceramics , 145-151 (2013). 44. S. Sigaryov, A. Vasiliev, Disorder of PO4 tetrahedra and ionic conductivity of Na3In2 (PO4) 3. Journal of Physics and Chemistry of Solids , 467-470 (1991). 45. Y. Yue, F. Zhou, W. Pang, Hydrothermal crystallization and structural investigation of Na1+ 2xZr2− xMgx (PO4) 3 systems (0< x< 1). Materials Chemistry and Physics , 28-30 (1993). 46. S. Song, H. M. Duong, A. M. Korsunsky, N. Hu, L. Lu, A Na+ superionic conductor for room-temperature sodium batteries. Scientific reports , 32330 (2016). 47. W. Baur, J. Dygas, D. Whitmore, J. Faber, Neutron powder diffraction study and ionic conductivity of Na2Zr2SiP2O12 and Na3Zr2Si2PO12. Solid State Ionics , 935-943 (1986). 48. A. Feltz, S. Barth, Preparation and conductivity behaviour of Na3MIIZr (PO4) 3,(MII: Mn, Mg, Zn). Solid State Ionics , 817-821 (1983). 49. P. Maldonado-Manso, M. A. Aranda, S. Bruque, J. Sanz, E. R. Losilla, Nominal vs. actual stoichiometries in Al-doped NASICONs: A study of the Na1. 4Al0. 4M1. 6 (PO4) 3 (M= Ge, Sn, Ti, Hf, Zr) family. Solid State Ionics , 1613-1625 (2005). 50. S. Naqash, F. Tietz, O. Guillon, Synthesis and characterization of equimolar Al/Y-substituted NASICON solid solution Na1+ 2x+ yAlxYxZr2− 2xSiyP3− yO12.
Solid State Ionics , 13-21 (2018). 51. V. Nicholas, P. Johnson, A. Kingon, Conductivity measurements in the Nasicon system.
Solid state ionics , 351-357 (1985). 52. T. Okura, M. Saimaru, H. Monma, K. Yamasita, Ionic conductivities of Nasicon-type glass-ceramic superionic conductors in the system Na2O–Y2O3–XO2–SiO2 (X= Ti, Ge, Te). Solid State Ionics , 537-540 (2009). 53. P. J. Squattrito et al. , Sodium and oxygen disorder in a scandium-substituted nasicon: A time of flight neutron powder diffraction study of Na2. 5Zr1. 8Sc0. 2Si1. 3P1. 7O12.
Solid state ionics , 31-40 (1988). 54. S. Susman, C. Delbecq, T. Brun, E. Prince, Fast ion transport in the NASICON analog Na3Sc2 (PO4) 3: Structure and conductivity. Solid State Ionics , 839-844 (1983). 55. C. Veríssimo et al. , Ionic conductivity and structural characterization of Na1. 5Nb0. 3Zr1. 5 (PO4) 3 with NASICON-type structure. Solid State Ionics , 127-134 (1997). 56. W. Wang, D. Li, J. Zhao, Solid phase synthesis and characterization of Na3Zr2− yNb0. 8ySi2PO12 system.
Solid state ionics , 97-100 (1992). 28 57. W. Wang, X. Liu, Study on sodium fast ion conductors in the Na1+ 3xAlxHf2− xSi2xP3− 2xO12 system. Solid state ionics , 165-168 (1996). 58. Y. Zhu, L. Li, C. Li, L. Zhou, Y. Wu, Na1+ xAlxGe2− xP3O12 (x= 0.5) glass–ceramic as a solid ionic conductor for sodium ion. Solid State Ionics , 113-117 (2016). 59. M. Guin, F. Tietz, Survey of the transport properties of sodium superionic conductor materials for use in sodium batteries.
Journal of power sources , 1056-1064 (2015). 60. M. Bianchini et al. , The interplay between thermodynamics and kinetics in the solid-state synthesis of layered oxides.
Nature Materials , 1-8 (2020). 61. A. Miura et al. , Selective metathesis synthesis of MgCr 2 S 4 by control of thermodynamic driving forces.
Materials Horizons , 1310-1316 (2020). 62. R. Dronskowski, P. E. Blöchl, Crystal orbital Hamilton populations (COHP): energy-resolved visualization of chemical bonding in solids based on density-functional calculations. The Journal of Physical Chemistry , 8617-8624 (1993). 63. V. L. Deringer, A. L. Tchougréeff, R. Dronskowski, Crystal orbital Hamilton population (COHP) analysis as projected from plane-wave basis sets. The journal of physical chemistry A , 5461-5466 (2011). 64. S. Maintz, V. L. Deringer, A. L. Tchougr é eff, R. Dronskowski, Analytic projection from plane - wave and PAW wavefunctions and application to chemical - bonding analysis in solids. Journal of computational chemistry , 2557-2567 (2013). 65. A. Kohan, G. Ceder, Tight-binding calculation of formation energies in multicomponent oxides: Application to the MgO-CaO phase diagram. Physical Review B , 805 (1996). 66. G. Ceder, P. Tepesch, A. Kohan, A. Van der Ven, A Model to Predict Ionic Disorder and Phase Diagrams: Application to CaO-MgO, Gd 2 O 3-ZrO 2, and to Sodium β′′ -alumina. Journal of Electroceramics , 15-26 (1997). 67. P. D. Tepesch et al. , A model to compute phase diagrams in oxides with empirical or first-principles energy methods and application to the solubility limits in the CaO-MgO system. J. Am. Ceram. Soc. , 2033-2040 (1996). 68. A. Urban, A. Abdellahi, S. Dacek, N. Artrith, G. Ceder, Electronic-structure origin of cation disorder in transition-metal oxides. Physical review letters , 176402 (2017). 69. S. Kotz, N. Balakrishnan, N. L. Johnson,
Continuous multivariate distributions, Volume . (John Wiley & Sons, 2004), vol. 1. 70. J. Platt, Probabilistic outputs for support vector machines and comparisons to regularized likelihood methods.
Advances in large margin classifiers , 61-74 (1999). 71. C. Zhou, S. Bag, V. Thangadurai, Engineering materials for progressive all-solid-state Na batteries. ACS Energy Letters , 2181-2198 (2018). 72. J. Wang et al. , A High - Energy NASICON - Type Cathode Material for Na - Ion Batteries.
Advanced Energy Materials , 1903968 (2020). 73. Y. Tian et al. , Promises and Challenges of Next-Generation “Beyond Li-ion” Batteries for Electric Vehicles and Grid Decarbonization. Chemical Reviews , (2020). 74. P. Canepa et al. , Odyssey of multivalent cathode materials: open questions and future challenges.
Chemical reviews , 4287-4341 (2017). 75. C. Zhao, F. Ding, Y. Lu, L. Chen, Y. S. Hu, High - Entropy Layered Oxide Cathodes for Sodium - Ion Batteries.
Angewandte Chemie International Edition , 264-269 (2020). 29 76. Z. Lun et al. , Cation-disordered rocksalt-type high-entropy cathodes for Li-ion batteries. Nature Materials , 1-8 (2020). 77. G. Kresse, J. Furthmüller, Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set.
Physical Review B , 11169-11186 (1996). 78. J. P. Perdew, K. Burke, M. Ernzerhof, Generalized Gradient Approximation Made Simple. Physical Review Letters , 3865-3868 (1996). Acknowledgments
We also thank Zijian Cai and Leeann Sun at University of California, Berkeley for their assistances during the NASICON synthesis.
Funding
This work was supported by the Samsung Advanced Institute of Technology. The computational analysis was performed using computational resources sponsored by the Department of Energy’s Office of Energy Efficiency and Renewable Energy at the National Renewable Energy Laboratory. Computational resources were also provided by the Extreme Science and Engineering Discovery Environment (XSEDE), which is supported by the National Science Foundation grant number ACI1053575 and the National Energy Research Scientific Computing Center (NERSC), a DOE Office of Science User Facility supported by the Office of Science and the U.S. Department of Energy under contract no. DE-AC02- 05CH11231.
Author contributions
B. O., J.W. and G.C. initiated and designed the project. G.C. supervised all aspects of the research. B.O. performed the high-throughput DFT calculations. B.O. and J.W. analyzed the high-throughput data and conceived the phase stability trend. J. W. synthesized all proposed compounds. B.O. conducted the ICOHP calculations. B.O. performed the machine learning work with the help of T. H., H. H., C. J. B. and J. W. Y. W., V. L. and H. K. contributed valuable discussion and insights. B. O., J. W. and G. C. wrote the manuscript. The manuscript was revised by all authors. Competing interests
The authors declare no competing interests. Synthetic accessibility and stability rules of NASICONs SUPPLEMENTARY MATERIALS
Bin Ouyang † , Jingyang Wang † , Tanjin He , Christopher J. Bartel , Haoyan Huo , Yan Wang , Valentina Lacivita , Haegyeom Kim and Gerbrand Ceder * Materials Sciences Division, Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA Department of Materials Science and Engineering, University of California, Berkeley, CA 94720, USA Advanced Material Lab, Samsung Research America, 10 Wilson Rd., Cambridge, MA 02138, USA † These authors contributed equally to this work * corresponding author: Prof. Gerbrand Ceder (Email: [email protected])
This PDF file includes:
Supplementary note 1 to 8 Figs. S1 to S9 Tables S1 to S3 Supplementary note 1. E hull and calculated competing phases of six synthesis attempts of NASICONs Table S1.
Six newly predicted and experimentally synthesized Na-rich silicon phosphate NASICON compounds hull (meV/atom) Competing Phases
1 Na HfZr(SiO ) (PO ) 1.9 Na P O + HfSiO + Na Zr (SiO ) (PO ) + Na ZrSi O + HfO
2 Na Hf Mg (SiO )(PO ) NaMgPO + NaHf (PO ) + Na P O + HfSiO + HfO
3 Na HfSn(SiO ) (PO ) 12.6 HfSiO + Na P O + Na Si O + SnO
4 Na HfSc(SiO )(PO ) HfSiO + Na P O + Na Sc (PO ) + NaSc(SiO ) + HfO
5 Na Ca Hf (SiO )(PO ) HfSiO + NaCaPO + Na P O + NaHf (PO ) + HfO
6 Na Mg Zr (SiO )(PO ) Na Mg P O + Na P O + NaZr (PO ) + ZrSiO + ZrO Supplementary note 2. Evaluation of ideal entropy of mixing ( 𝑺 𝑰𝒅𝒆𝒂𝒍 ) The partially occupied sites of NASICON are considered for ideal entropy calculations. For a NASICON with the chemical formula Na x M y M’ (AO ) z (BO ) , the ideal entropy is calculated as below with the unit of eV/atom: 𝑆 &’()* = −𝑘 +
1𝑥 + 17 *14 ,𝑥4 ln ,𝑥4/ + (1 − 𝑥4) ln ,1 − 𝑥4// + 12 ,𝑦2 log ,𝑦2/ + (1 − 𝑦2) log ,1 − 𝑦2//+ 13 ,𝑧3 log ,𝑧3/ + (1 − 𝑧3) log ,1 − 𝑧3//8 Supplementary note 3. Calculated COHP and kernel density for analyzing bond compatibility and site miscibility.
Figure S1:
Calculated -COHP normalized by number of metal-O bonds for four comparison groups, i.e., A/A1, A/A2, B/B1 and B/B2 from left to right.
Figure S2:
Kernel distribution of NASICONs as a function of cation/anion mixing factor. The mixing factor is defined as (𝑦/2,1 − 𝑦/2) for cation mixing and 𝑚𝑖𝑛(𝑧/3,1 − 𝑧/3) for anion mixing. Figure S3:
Calculated -COHP normalized by number of S-O/P-O bonds for A/A1 and A/A2 comparisons.
Figure S4:
Calculated -COHP normalized by number of Si-O/P-O bonds for A/A1 and A/A2 comparisons.
Table S2:
Calculated ICOHP normalized by number of ions for A/A1/A2 and B/B1/B2 comparisons.
ICOHP M-O A-O
A (NaHfMg(SO ) (PO )) A1 (NaHfCa(SO ) (PO )) A2 (NaGeMg(SO ) (PO )) B (Na HfZr(SiO ) (PO )) B1 (Na HfGe(SiO ) (PO )) B2 (Na Ge (SiO ) (PO )) Supplementary note 4. Basic Machine Learning Features.
The properties related to stability of NASICON would be electronegativity (X), ionic radius (R), ionic charge states (Q), Na content (N Na ), formation energy of competing phases (E Form ), Ewald summation energy (E
Ewald ). These basic properties can be obtained from different sites, and with different mathematical form. For electronegativity (X), ionic radius (R) and ionic charge states (Q), both metal sites and anion sites are investigated, leading to X
M/A , R
M/A and Q
M/A , respectively. In addition, the average across all sites, standard deviation, maximum deviation were considered, which are noted as 𝑋 ,-./ , 𝑋 ,01’ , 𝑋 ,01’ , 𝑋 --./ , 𝑋 -01’ , 𝑋 -01’ , 𝑅 ,-./ , 𝑅 ,01’ , 𝑅 ,01’ , 𝑅 --./ , 𝑅 -01’ , 𝑅 -01’ , 𝑄 ,-./ , 𝑄 ,01’ , 𝑄 ,01’ , 𝑄 --./ , 𝑄 -01’ , 𝑄 -01’ . For ionic radius and electronegativity, three additional features were considered, i.e., 𝑅 ,-./ /𝑅 , 𝑋 ,-./ − 𝑋 and 𝑋 ,-./ − 𝑋 - . For the anion site, the property of center ion i.e., Si , P and S , was also considered Therefore, using the aforementioned E Form, N Na and E Ewald , a total of 24 basic features are constructed. To enlarge the physically meaningful searching space for mathematical representations, all ionic radii are normalized with the radius of Na + , while all charge states are normalized by the charge of an electron. It is worth mentioning that the real 𝐸 energy is related to not only Na contents and average cation/anion charge states, but also cation/anion ordering and lattice parameter. Here we did not want to involve the impacts of cation/anion ordering and lattice parameter, as those properties cannot be estimated without first principal calculations. The 𝐸 we calculated was on base of a structural model with the lattice parameter of Na Zr (SiO ) (PO ) and electrostatic ground state ionic ordering. Similarly, the 𝐸 value is normalized by the Ewald energy of Na Zr (SiO ) (PO ). Therefore, this parameter can be interpreted as a dimensionaless physical term that measures the electrostatic ground-state electrostatic interaction as a function of Na content and average cation/anion charge states. It can then be used as a descriptor that couples the influenced Na contents and cation/anion charge states. Supplementary note 5. Predicted probability for 32 experimentally synthesized NAISCONs
The calculated probabilities of the 32 experimentally synthesized NASICONs being GS/LS according to Platt scaling are indicated Table S3.
Table S3. Calculated Platt scaling probability for 32 experimentally synthesized NASICONs.
Composition P accessible
Classified Composition P accessible
Classified Na Zr (SiO ) (PO ) 0.449 No NaZnIn(SO ) Zr (SiO ) ) (PO ) HfTi(SiO ) (PO ) 0.328 No NaSn (PO ) (PO ) (PO ) ) MgZr(PO ) ) Zr (SiO )(PO ) ) ZrTi(SiO ) (PO ) 0.314 No Na ZrAl(PO ) Ti (SiO ) (PO ) 0.449 No Na In (PO ) Sc (PO ) (PO ) ZrZn(PO ) Al (PO ) Hf (SiO ) (PO ) 0.449 No Na HfZr(SiO ) (PO ) 0.440 No Na Y (PO ) Hf Mg (SiO )(PO ) HfAl(SiO ) (PO ) 0.177 No Na HfSc(SiO )(PO ) ) Ca Hf (SiO )(PO ) ) Mg Zr (SiO )(PO ) Supplementary note 6. Predicted probability map for Na Poor/Na Rich NASICONs and Ca/Ge/Hf/Zr NASICONs
The predicted probability distributions of Na Poor/Na Rich NASICONs and Ca/Ge/Hf/Zr NASICONs are demonstrated in
Figure S5 : Figure S5: (a) The Platt scaling probability to be accessible (P accessible ) distribution of Na rich and Na poor NASICONs; (b) The Platt scaling accessible probability (P accessible ) distribution of Ca/Ge/Zr/Hf contained NASICONs.
Likely Stable RegionLikely Unstable Region (a) (b) P !" P !" Supplementary note 7. Comparisons between -ICOHP of M-O bonds in B, B1, B2 and the -ICOHP of the M-O bonds in the competing phases of B1 and B2
Figure S6.
The calculated -ICOHP for different M-O bonds in NASICONs composition B, B1 and B2. As B1 and B2 are predicted to be unlikely synthesizable ( E hull = 56.0 meV/atom and 90.3 meV/atom respectively), the calculated -ICOHP per cation of the M-O bonds in competing phase is also plotted. Figure S7.
The calculated -COHP comparisons (normalized by M-O bond number) for (a) Ge-O in B1 and Na Ge O ; (b) Ge-O in B2 and Na Ge O . (a) (b) Supplementary note 8. Comparisons between -ICOHP in A, A1, and the Ca contained competing phase of A1
Figure S8.
The calculated -ICOHP for different M-O bonds in NASICONs composition A and A1. As A1 is predicted to be unlikely synthesizable by DFT ( E hull = 48.7 meV/atom), the calculated -ICOHP per cation of the M-O bonds in competing phase is also plotted. Figure S9.
The calculated -COHP comparisons (normalized by M-O bond number) for Ca-O in A1 and CaSO4