Open Access Article

This Open Access Article is licensed under a Creative Commons Attribution-Non Commercial 3.0 Unported Licence

Joshua
Ojih
^{a},
Mohammed
Al-Fahdi
^{a},
Yagang
Yao
^{b},
Jianjun
Hu
^{c} and
Ming
Hu
*^{a}
^{a}Department of Mechanical Engineering, University of South Carolina, SC 29208, USA. E-mail: hu@sc.edu
^{b}National Laboratory of Solid State Microstructures, College of Engineering and Applied Sciences, Jiangsu Key Laboratory of Artificial Functional Materials, and Collaborative Innovation Center of Advanced Microstructures, Nanjing University, Nanjing 210093, China
^{c}Department of Computer Science and Engineering, University of South Carolina, Columbia, SC 29208, USA

Received
11th October 2023
, Accepted 23rd February 2024

First published on 6th March 2024

Prediction of crystal structures with desirable material properties is a grand challenge in materials research, due to the enormous search space of possible combinations of elements and their countless arrangements in 3D space. Despite the recent progress of a few crystal structure prediction algorithms, most of those methods only target a few specific material families or are restricted to simple systems with limited element diversity. Moreover, these algorithms are usually coupled with first principles calculations and thus are computationally expensive and very time consuming. Therefore, establishing a workflow that can generate a large number of hypothetical structures with diverse elements and quickly optimize and screen out stable structures is urgently needed for the crystal structure prediction field. In this study, we take 17277 compositions involving 63 elements across the periodic table from the open quantum materials database (OQMD) and use a graph theory assisted universal structure searcher (MAGUS) to generate more than 3.4 million hypothetical structures. We employ a pre-trained universal interatomic potential named Crystal Hamiltonian Graph neural Network (CHGNet) to rapidly optimize this large number of hypothetical structures. Subsequently, we validate these optimized structures using density functional theory (DFT) and find 4145 structures are successfully optimized and 2368 structures have energy lower than those from the original OQMD. 647 structures are further identified to be dynamically stable using CHGNet. Moreover, the stability of 123 out of 200 randomly chosen structures are validated by DFT, corresponding to a high success rate of 61.5%. We further use 4706 DFT data points to train 3 graph neural network models to predict lattice thermal conductivity (LTC) and heat capacity. Numerous structures with ultralow LTC and high heat capacity, which are promising for advanced energy conversion and energy storage, have been identified. The success of our workflow demonstrates that combining graph theory with pre-trained universal interatomic potential is highly expected to accelerate the search for new both thermodynamically and dynamically stable structures with target material properties.

In the past decade, several algorithms have been developed to search for new hypothetical structures. For example, CALYPSO^{7} uses particle swarm optimization (PSO) to predict the high pressure phase of many different combinations of compositions. USPEX^{8} uses an evolutionary algorithm (EA) to find many stable structures as well. Some recent CSP methods or algorithms in this line include crystal diffusion variational autoencoder^{9} and machine learning (ML) enabled or accelerated general inverse design of inorganic crystals^{10} and porous materials.^{11} These newer methods aim to use ML algorithms to capture the physical inductive bias of material stability and hidden rules of periodic arrangements of atoms among all possible combinations of elements. Despite the huge success, there is at least one major issue associated with these approaches: the number of local minima in the vast material space grows exponentially with the number of atoms, and the cost of density functional theory (DFT) calculations used for local structure optimization grows rapidly with the number of atoms and huge amount of hypothetical structures. Hence, the first motivation of this work is to think about an alternative approach that can bypass the tedious optimization of hypothetical structures. In recent years, machine learning techniques have gained rapid adoption as a precise and efficient tool for atomistic simulations, enabling the resolution of a wide array of fundamental problems.^{12} This trend has given rise to a novel class of interatomic potentials known as machine learning potentials (MLPs). Notable examples include Gaussian approximation potential (GAP),^{13–15} higher order equivariant message passing neural networks (MACE),^{16,17} high-dimensional neural network potential (HDNNP),^{18} moment tensor potential (MTP),^{19} spectral neighbor analysis potential (SNAP),^{20} spatial density neural network force fields (SDNNFFs)^{21} and deep potential for molecular dynamics (DeepMD).^{22} MLPs have found applications across a diverse spectrum of materials systems, consistently delivering high accuracy compared to the conventional DFT approaches. Based on these studies, we expect that the state-of-the-art MLPs would be straightforwardly applied to solve the issue in CSP by speeding up energy and force evaluation necessary for structure optimization.

Another issue in high-throughput CSP with desired material properties is to quickly screen or accurately predict the physical properties of hypothetical structures after successfully optimizing structures. Some physical properties are easy to train with ML models while some are not. For instance, most mechanical properties such as Young's modulus, bulk modulus, etc. can be well trained and predicted with high accuracy by both traditional ML models^{23} and newly developed models such as graph neural networks.^{24} Unfortunately, the lattice thermal conductivity (LTC) we are interested in in this work belongs to the latter case. Screening and predicting the LTC of materials is crucial in the design of thermal management, energy conversion, and energy storage systems such as thermoelectrics, electronic cooling, solid-state phase change, etc. A recent study^{25} established a strong correlation between the thermal conductivity of amorphous GaO_{x} and mass density, O to Ga composition ratio, and structural similarity factor (SSF) which inherits the sensitivity of the smooth overlap of atomic positions (SOAP) descriptor to density and composition. These results demonstrate that combining SOAP and SSF could be very promising for identifying more hidden relationships or correlations between thermal transport and structural features. Despite the ubiquitous application of thermal materials, it is hard to train a good ML model to quantitatively predict LTC across diverse material types and families or sometimes the developed ML model involves material descriptors that are hard to compute in a high-throughput manner due to the nonlinear relationship between atomic structure and thermal transport properties.^{26,27} Furthermore, another obstacle before the struggling LTC training and prediction is to judge whether a crystal is dynamically stable, i.e., there should be no imaginary phonon frequencies in the full Brillouin zone. Without ensuring the dynamical stability of a crystalline structure, predicting its thermal transport properties would be physically meaningless. To this end, the second motivation of this work is to use ML models to screen dynamically stable materials with high speed and high accuracy so that the prediction on the LTC of the final promising structures will have a high success rate.

In this paper, we established a workflow consisting of hypothetical structure generation, fast structure optimization, and quick screening of dynamical stability for searching for ultralow LTC for thermoelectric energy conversion and high heat capacity for thermal energy storage. We first use an ML algorithm called “machine learning and graph theory assisted universal structure searcher (MAGUS)^{4}” for the generation of large-scale hypothetical structures, followed by fast structure optimization and quick screening of dynamical stability using a pre-trained MLP named “Crystal Hamiltonian Graph neural Network (CHGNet)”.^{28} The LTC and heat capacity of the filtered structures in the previous steps are predicted by our separately trained graph neural network (GNN) models. Selected promising candidates were validated by full DFT calculations. In all, after generating 200 hypothetical structures each for the 17277 unstable structures from the OQMD, we found 4145 of those structures to be successfully optimized by DFT. Furthermore, 647 structures were found to be dynamically stable by our CHGNet model, and 123 out of 200 randomly selected structures were verified to be indeed dynamically stable by full DFT calculations, showing a high success rate of screening dynamically stable structures by using pre-trained MLPs.

Here, we would like to discuss the accuracy of the CHGNet MLP which plays a significant role in obtaining more structures compared to DFT calculations. However, it is extremely hard to train a highly accurate universal MLP that can cover the vast space of different combinations of elements across the periodic table and various atomic arrangements in space to the largest extent. To the best of our knowledge, so far there are only two such types of pre-trained universal MLPs that are publicly accessible, namely CHGNet (the one used herein) and M3GNEt.^{41} These pre-trained universal MLPs are expected to accelerate discovery of novel materials through fast evaluation and screening of large previously unexplored materials space. Both MLPs were trained on the raw data from the Materials Project database gathered in the past 10 years or even longer. The biggest advantage of these pre-trained universal MLPs is that they can cover different combinations of elements across the periodic table and various atomic arrangements in space, such that one can immediately use them to quickly evaluate and screen a very large number of hypothetical systems (3.4 million in this work) in an acceptable timeframe. We would also like to point out that the development of these MLPs is still in their infancy and their accuracy needs more comprehensive testing on different structures in the vast materials space. Fine tuning or even re-training a universal MLP would require a large amount of new DFT calculations, which is not an easy task. The main barrier is to identify the unseen and to-be-poorly-predicted configurations by the pre-trained CHGNet model out of the large amount of hypothetical structures and subsequently perform DFT calculations to serve as new training data. Very recently, we also noticed that the Google DeepMind team^{42} claimed that their newly developed universal MLP outperforms CHGNet. This might be true considering two facts: (1) huge amounts of high quality DFT data might be generated in that work, and thus it is not surprising to train an ML model better than previous ones. (2) The statistics or distribution of materials for the 384938 materials^{42} is fundamentally different from the Materials Project database. For example, our quick screening found that at least 78% of the materials reported by the DeepMind team have lanthanide or actinide elements, while this percentage in the Materials Project database is only 25.5% (data downloaded on Jan. 3rd, 2022). Therefore, it is not surprising that the CHGNet model that was trained on different materials datasets has poor performance on predicting the lanthanide/actinide dominated material pool. Comprehensive evaluation of the accuracy of CHGNet would require more careful work and a fair materials dataset for comparison.

Regarding the difference in the number of structures between CHGNet screening and DFT confirmation, one of the main reasons is there could be tons of new local atomic environments that have not been seen in the training data of the CHGNet MLP and thus will be poorly predicted by the ML model. Such phenomenon is very common when dealing with huge numbers of hypothetical structures. A very standard procedure is to use the ML model to do screening first, due to its fast speed. In this step, accuracy is not a big concern since the finally obtained “good” structures will still be large. Nevertheless, despite some discrepancies as compared with full DFT calculations, CHGNet is a great success in predicting the energy of a given cell and atomic forces therein with fast speed. In fact, our success rate of 28.3% is way higher than in a recent study using transfer learning to predict stable materials.^{43} The relatively higher success rate is due to the different screening/filtering procedure and the different ML model used.

We then combine the CHGNet MLP and PHONOPY package to obtain the second order IFCs for all 2368 MAGUS structures. After checking phonon frequencies in the full Brillouin zone, we identify 647 structures to be dynamically stable (a positive dispersion percentage of 27.3%), i.e., no negative phonon frequencies are found. We would rate the percentage of positive frequency as normal, although the total number of dynamically stable structures is not impressively large. This judgement is based on our recent observation of generally positive dispersion percentages between 23% and 28% on some other material datasets and families. We would like to emphasize that it is much stricter for a crystal structure to be dynamically stable than simply being optimized. Negative formation energy and even energy above hull close to zero does not necessarily guarantee dynamical stability. In other words, the low percentage of predicted dynamically stable structures is more likely due to the small amount of truly stable structures found. Nevertheless, we are still able to use the workflow developed herein to discover some new stable structures. A higher success rate of finding dynamically stable structures would require a more sophisticated ML model to understand the intrinsic nature and physical/chemical feature of dynamic stability, and an equally important thing is to generate structures in the previously unexplored vast space of compositions and arrangements of atoms.

To validate the success rate of dynamic stability screening by CHGNet, we randomly selected 200 structures and performed full DFT calculations for their phonon dispersions. Due to the expensive computational cost by DFT calculations of full interatomic force constants required for evaluating phonon frequencies, in particular for complex structures with low material symmetry, we only randomly selected 200 structures for DFT validation. We found 123 out of 200 structures to be dynamically stable, corresponding to a success rate of 61.5%. Despite one-time sampling, the success rate is comparable to our separate validation on another material family with prototype ABC_{2}D_{6} and noncubic symmetry (results have not been published yet). Therefore, we believe that the CHGNet MLP has decent accuracy in predicting dynamic stability. The 61.5% success rate of positive dispersions as validated by DFT calculations is good enough for high-throughput screening of a large number of previously unexplored structures if the initial screening pool is very large (say >10^{4} predicted positive dispersions). The DFT validation result confirms the high success rate of using the pre-trained CHGNet MLP for quickly screening dynamic stability, which would have taken much longer simulation time and resources if running by full DFT calculations. Fig. 3a–f show the phonon dispersion comparison between CHGNet and DFT for 6 selected structures. All 6 structures show no imaginary frequency along the high symmetry paths. It should be noted that there are indeed some discrepancies between the phonon dispersions by CHGNet and DFT. For instance, it seems that there is a systematic underestimation of frequencies by CHGNet as compared with DFT. The discrepancy is understandable considering that the DFT data used for training CHGNet is quite different from the DFT calculations performed by us, in terms of both structure and element diversity. This leaves some room for further improvement of the CHGNet model. It is well-known that the size and quality of the training data determine the accuracy of an MLP. When more DFT data is added into training in the future, the CHGNet MLP will be more accurate. Nevertheless, generally speaking, the frequency near the Γ-point is well predicted by CHGNet, while the high frequency optical phonons are under-predicted as compared with DFT results, and the pre-trained CHGNet is efficient in predicting a wide range of structures with elements across the entire periodic table and material symmetry across all 230 space groups.

We further show the relationship between the material properties P_{3} parameter and mean square displacement (MSD) with LTC in Fig. 6a and b, respectively. The P_{3} parameter corresponds to the availability of three-phonon scattering opportunities throughout the entire Brillouin zone. A high P_{3} parameter value indicates the presence of numerous pathways for phonon–phonon interaction (scattering) within the crystal lattice, which in turn usually results in a low LTC.^{45,46} The MSD offers information about how atomic positions in crystals change with temperature, shedding light on how atoms deviate from their equilibrium positions in a harmonic phonon lattice. The large MSD value leads to the identification of rattling atoms within the crystal.^{47,48} Both the P_{3} parameter and MSD serve as valuable material descriptors, allowing for a rapid screening of crystalline materials with exceptionally ultralow LTC. These plots depict the relationship between LTC and the P_{3} parameter and MSD for the 123 DFT verified structures. It is evident that this observation is consistent with the established knowledge that the LTC generally decreases as the MSD or P_{3} parameter increases, showcasing an inverse proportionality.

To gain deep insight into the electronic level understanding of ultralow LTC materials, we conducted a more in-depth analysis of our predicted LTC using Crystal Orbital Hamilton Population (COHP)^{49} to measure the contribution of bonding and antibonding states. The single values of bonding and antibonding for each structure are obtained by performing integration over COHP curves for each atomic pair as evaluated by the LOBSTER package and then taking the average for all pairs in the primitive cell. This analysis helps us to understand how chemical bonding influences thermal transport properties such as LTC. In Fig. 7a, we illustrate bonding and antibonding states, highlighting the log value of predicted LTC. From the kinetic theory of phonon transport, the lattice thermal conductivity of a crystal can be expressed as , where C_{v}, ν, and τ are the volumetric heat capacity, average phonon group velocity, and average phonon relaxation time, respectively. In principle, the group velocity and relaxation time are determined by the harmonic and anharmonic characteristics of the lattice, respectively, which can be roughly correlated with the bonding and antibonding COHP, respectively. Generally speaking, low bonding COHP and high antibonding COHP would more likely lead to low LTC, as reflected by the lower-right corner of Fig. 7a. When bonding values are less than 10 and antibonding greater than 0.01, the structures exhibit LTC less than 1 W m^{−1} K^{−1}, i.e., the LTC is ultralow. In this region, the interatomic bonding strengths are weak while phonon anharmonicity is strong, resulting in exceptionally low LTC. This observation is consistent with our recent big data analysis of bonding–antibonding states on 13718 stable structures.^{29} In contrast, high antibonding COHP does not necessarily lead to low LTC. This is because the phonon group velocity (the harmonic term in the kinetic theory) competes with the phonon relaxation time (the anharmonic term). Therefore, medium to high LTC materials can happen in a broad range of antibonding. Also note that the absolute value of high LTC materials shown in Fig. 7a only ranges from 16 W m^{−1} K^{−1} to 32 W m^{−1} K^{−1}, which is not very high at all compared with the really high values of many other LTC materials. Thus, in Fig. 7a it is easy to find some materials with harmonic terms overwhelming the anharmonic terms and then showing relatively high LTC. Such trend can be clearly seen from more data (13718 cubic structures in total) in our recent work.^{29}

Fig. 7 (a) Insight into the relationship between bonding and antibonding with color coded by the log-value of predicted LTC. (b) Pearson correlation between material descriptors and LTC. |

Furthermore, we explore the correlation between our predicted LTC values and material descriptors which are simple characteristics that can help quickly screen structures with low LTC. This serves as a bridge between DFT accuracy and classical level efficiency. In Fig. 7b, we present the Pearson correlation between some representative material descriptors with LTC. The Pearson correlation coefficient is the most common way of measuring a linear correlation between two variables. The Pearson correlation number is between −1 and +1 that measures the strength and direction of the relationship between two variables, with −1 meaning strong negative correlation and +1 meaning strong positive correlation. Volume of the primitive cell, bond length, total weight and average number of electrons exhibit negative correlations with LTC, with Pearson correlation values of −0.49, −0.45, −0.21 and −0.18, respectively. This means higher values of these descriptors correspond to structures with low LTC. Physically speaking, this trend is understandable since a larger cell volume or longer bond length usually corresponds to loose bonding between atoms, and thus leads to low phonon group velocities. On the other hand, atom number density, the number of unpaired electrons, and mass density are positively correlated with LTC, with Pearson correlation values of 0.54, 0.41 and 0.16, respectively. Higher atom number density means denser packing of atoms in space and thus usually corresponds to strong interatomic bonding and thus high phonon group velocities, which facilitates phonon energy transport in the lattice. Meanwhile, descriptors like the total number of atoms, maximum principal quantum number, and Pauling electronegativity show almost no correlations or weak correlations with LTC. The positive and negative correlated descriptors identified herein are expected to accelerate the screening process of even larger-scale hypothetical structures in the future.

Fig. 9 DFT validation of deeperGATGNN prediction for the 123 structures randomly selected from the final dynamically stable structure pool (647) as screened by the workflow in this work. |

MAGUS ID | Formula | Symmetry | Space group number | Formation energy (eV per atom) | LTC (W m^{−1} K^{−1}) |
Heat capacity (J mol^{−1} K^{−1}) |
Ratio of heat capacity to Dulong–Petit limit |
---|---|---|---|---|---|---|---|

313161_167 | YTl_{3} |
Pmm | 221 | −0.23374 | 1.780 | 98.7352 | 0.9896 |

546287_12 | CaTlBi_{2} |
P4/mmm | 123 | −0.47185 | 1.078 | 98.7321 | 0.9896 |

310222_80 | CaBi_{3} |
Pmm | 221 | −0.50218 | 1.471 | 98.7256 | 0.9896 |

313873_43 | CaPb_{3} |
Pmm | 221 | −0.32828 | 4.031 | 98.5031 | 0.9873 |

371322_34 | YTl_{2}In |
P4/mmm | 123 | −0.28082 | 2.228 | 98.4791 | 0.9871 |

453343_147 | YIn_{2}Bi |
P4/mmm | 123 | −0.41288 | 1.512 | 98.4355 | 0.9866 |

499515_125 | LiBiPb_{2} |
P4/mmm | 123 | −0.19219 | 0.742 | 98.3930 | 0.9862 |

468937_47 | BiPd_{2}Pt |
P4/mmm | 123 | −0.20358 | 1.003 | 98.1140 | 0.9834 |

454575_137 | ScIn_{2}Pb |
P4/mmm | 123 | −0.21603 | 1.474 | 98.0507 | 0.9828 |

1548987_90 | K2YTlBr_{6} |
P | 147 | −1.67237 | 0.096 | 245.111 | 0.9827 |

314264_194 | Ag_{3}Ge |
Amm2 | 38 | 0.113334 | 0.029 | 97.9828 | 0.9821 |

416611_76 | Hf_{2}CdAu |
P | 2 | −0.19819 | 0.177 | 97.9752 | 0.9820 |

308064_41 | TaAu_{3} |
P4/mmm | 123 | −0.05549 | 0.628 | 97.8864 | 0.9811 |

312610_45 | CaIn_{3} |
Pmm | 221 | −0.31326 | 1.804 | 97.7639 | 0.9799 |

314157_46 | RhAu_{3} |
Pm1 | 164 | 0.138769 | 0.150 | 97.7346 | 0.9796 |

515810_43 | CdPdAu_{2} |
P4/mmm | 123 | −0.26516 | 2.601 | 97.7262 | 0.9795 |

308890_145 | TaPt_{3} |
Pmm | 221 | −0.54038 | 2.168 | 97.7192 | 0.9795 |

463270_167 | Pd_{2}PtPb |
P4/mmm | 123 | −0.24499 | 3.143 | 97.5625 | 0.9779 |

311765_44 | ScIn_{3} |
Pmm | 221 | −0.28857 | 3.058 | 97.4223 | 0.9765 |

308939_100 | Pd_{3}Pb |
Pmm | 221 | −0.30703 | 3.652 | 97.3603 | 0.9759 |

- C. C. Pantelides, C. S. Adjiman and K. Andrei V, Prediction and Calculation of Crystal Structures, ed. S. A.-E. A. Aspuru-Guzik, Springer Cham, 2014, pp. 25–58 Search PubMed .
- H. Callen, Thermodynamics and an Introduction to Thermostatistics, John Wiley & Sons, 1985 Search PubMed .
- B. C. Revard, W. W. Tipton and R. G. Hennig, Structure and Stability Prediction of Compound with Evolutionary Algorithm, 2014 Search PubMed .
- J. Wang, H. Gao, Y. Han, C. Ding, S. Pan, Y. Wang, Q. Jia, H.-T. Wang, D. Xing and J. Sun, Natl. Sci. Rev., 2023, 10(7), nwad128 CrossRef CAS PubMed .
- R. O. Jones and O. Gunnarsson, Rev. Mod. Phys., 1989, 61, 689–746 CrossRef CAS .
- J. Ojih, M. Al-fahdi, A. D. Rodriguez and K. Choudhary, npj Comput. Mater., 2022, 1–12 Search PubMed .
- Y. Wang, J. Lv, L. Zhu and Y. Ma, Comput. Phys. Commun., 2012, 183, 2063–2070 CrossRef CAS .
- C. W. Glass, A. R. Oganov and N. Hansen, Comput. Phys. Commun., 2006, 175, 713–720 CrossRef CAS .
- T. Xie, X. Fu, O.-E. Ganea, R. Barzilay and T. Jaakkola, Crystal Diffusion Variational Autoencoder for Periodic Material Generation, 2021, pp. 1–20 Search PubMed .
- Z. Ren, S. I. P. Tian, J. Noh, F. Oviedo, G. Xing, J. Li, Q. Liang, R. Zhu, A. G. Aberle, S. Sun, X. Wang, Y. Liu, Q. Li, S. Jayavelu, K. Hippalgaonkar, Y. Jung and T. Buonassisi, Matter, 2022, 5, 314–335 CrossRef CAS .
- B. Kim, S. Lee and J. Kim, Sci. Adv., 2020, 6, 1–8 Search PubMed .
- Q. Tong, P. Gao, H. Liu, Y. Xie, J. Lv, Y. Wang and J. Zhao, J. Phys. Chem. Lett., 2020, 11, 8710–8720 CrossRef CAS PubMed .
- A. P. Bartõk and G. Csányi, Int. J. Quantum Chem., 2015, 115, 1051–1057 CrossRef .
- A. P. Bartók, M. C. Payne, R. Kondor and G. Csányi, Phys. Rev. Lett., 2010, 104, 1–4 CrossRef PubMed .
- A. P. Bartók, J. Kermode, N. Bernstein and G. Csányi, Phys. Rev. X, 2018, 8, 41048 Search PubMed .
- I. Batatia, D. P. Kovács, G. N. C. Simm, C. Ortner and G. Csányi, Adv. Neural Inf. Process. Syst., 2022, 35, 11423–11436 Search PubMed .
- D. P. Kovács, I. Batatia, E. S. Arany and G. Csányi, J. Chem. Phys., 2023, 159(4), 044118 CrossRef PubMed .
- J. Behler and M. Parrinello, Phys. Rev. Lett., 2007, 98, 1–4 CrossRef PubMed .
- A. V. Shapeev, Multiscale Model. Simul., 2016, 14, 1153–1173 CrossRef .
- M. A. Wood and A. P. Thompson, J. Chem. Phys., 2018, 148, 241721 CrossRef PubMed .
- A. Rodriguez, Y. Liu and M. Hu, Phys. Rev. B, 2020, 102, 35203 CrossRef CAS .
- H. Wang, L. Zhang, J. Han and W. E, Comput. Phys. Commun., 2018, 228, 178–184 CrossRef CAS .
- J. Ojih, U. Onyekpe, A. Rodriguez, J. Hu, C. Peng and M. Hu, ACS Appl. Mater. Interfaces, 2022, 14, 43277–43289 CrossRef CAS PubMed .
- J. Ojih, A. Rodriguez, J. Hu and M. Hu, Energy AI, 2023, 14, 100286 CrossRef .
- Y. Liu, H. Liang, L. Yang, G. Yang, H. Yang, S. Song, Z. Mei, G. Csányi and B. Cao, Adv. Mater., 2023, 35, 1–13 Search PubMed .
- G. Qin, Y. Wei, L. Yu, J. Xu, J. Ojih, A. D. Rodriguez, H. Wang, Z. Qin and M. Hu, J. Mater. Chem. A, 2023, 5801–5810 RSC .
- C. Loftis, K. Yuan, Y. Zhao, M. Hu and J. Hu, J. Phys. Chem. A, 2021, 125, 435–450 CrossRef CAS PubMed .
- B. Deng, P. Zhong, K. Jun, J. Riebesell, K. Han, C. J. Bartel and G. Ceder, Nat. Mach. Intell., 2023, 5, 1031–1041 CrossRef .
- A. Rodriguez, C. Lin, C. Shen, K. Yuan, M. Al-Fahdi, X. Zhang, H. Zhang and M. Hu, Commun. Mater., 2023, 4, 61 CrossRef CAS .
- A. Rodriguez, C. Lin, H. Yang, M. Al-Fahdi, C. Shen, K. Choudhary, Y. Zhao, J. Hu, B. Cao, H. Zhang and M. Hu, npj Comput. Mater., 2023, 9, 20 CrossRef CAS .
- A. Togo and I. Tanaka, Scr. Mater., 2015, 108, 1–5 CrossRef CAS .
- G. Kresse and J. Hafner, Phys. Rev. B: Condens. Matter Mater. Phys., 1993, 47, 558–561 CrossRef CAS PubMed .
- D. Joubert, Phys. Rev. B: Condens. Matter Mater. Phys., 1999, 59, 1758–1775 CrossRef .
- R. A. Vargas-Hernández, J. Phys. Chem. A, 2020, 124, 4053–4061 CrossRef PubMed .
- J. P. Perdew, K. Burke and M. Ernzerhof, Phys. Rev. Lett., 1996, 77, 3865–3868 CrossRef CAS PubMed .
- K. Choudhary and B. DeCost, npj Comput. Mater., 2021, 7, 1–8 CrossRef .
- M. Karamad, R. Magar, Y. Shi, S. Siahrostami, I. D. Gates and A. Barati Farimani, Phys. Rev. Mater., 2020, 4, 1–6 Search PubMed .
- S. S. Omee, S.-Y. Louis, N. Fu, L. Wei, S. Dey, R. Dong, Q. Li and J. Hu, Patterns, 2021, 3, 100491 CrossRef PubMed .
- F. Zhou, W. Nielson, Y. Xia and V. Ozoliņš, Phys. Rev. B, 2019, 100, 1–15 Search PubMed .
- W. Li, J. Carrete, N. A. Katcho and N. Mingo, Comput. Phys. Commun., 2014, 185, 1747–1758 CrossRef CAS .
- C. Chen and S. P. Ong, Nat. Comput. Sci., 2022, 2, 718–728 CrossRef PubMed .
- A. Merchant, S. Batzner, S. S. Schoenholz, M. Aykol, G. Cheon and E. D. Cubuk, Nature, 2023, 624, 80–85 CrossRef CAS PubMed .
- J. Schmidt, N. Hoffmann, H. C. Wang, P. Borlido, P. J. M. A. Carriço, T. F. T. Cerqueira, S. Botti and M. A. L. Marques, Adv. Mater., 2023, 35, 2210788 CrossRef CAS PubMed .
- T. Zhu, R. He, S. Gong, T. Xie, P. Gorai, K. Nielsch and J. C. Grossman, Energy Environ. Sci., 2021, 14, 3559–3566 RSC .
- L. Elalfy, D. Music and M. Hu, Phys. Rev. B, 2021, 103, 75203 CrossRef CAS .
- T. Ouyang and M. Hu, Phys. Rev. B: Condens. Matter Mater. Phys., 2015, 92, 1–13 CrossRef .
- S. Ju, R. Yoshida, C. Liu, S. Wu, K. Hongo, T. Tadano and J. Shiomi, Phys. Rev. Mater., 2021, 5, 1–2 Search PubMed .
- A. Jain, H. P. Veeravenkata, S. Godse and Y. Srivastava, arXiv, 2022, arXiv:2204.03628, DOI:10.48550/arXiv.2204.03628.
- R. Dronskowski and P. E. Blochl, J. Phys. Chem., 1993, 97, 8617–8624 CrossRef CAS .

## Footnote |

† Electronic supplementary information (ESI) available. See DOI: https://doi.org/10.1039/d3ta06190f |

This journal is © The Royal Society of Chemistry 2024 |