Unleashing the power of artificial intelligence in phonon thermal transport: Current challenges and prospects

Phonons, the quanta of atomic vibrations, dominate the thermal transport properties in semiconductors and insulators and contribute non-negligibly to energy transport in metals and metallic systems. Because of their ubiquity, fast and accurately predicting phonon properties and phonon-mediated heat transfer process in ordered structures is enormously important for the development of innovative energy oriented technologies (Fig. 1), including but not limited to energy conversion,^1–3 thermal management,^4–7 superconductivity,^8–10 photovoltaics,^11–13 quantum computing,^14–16 etc. Discovery of advanced thermal materials with exceptional phonon properties drives technological advancements and impacts innovations. Thermal conductivity is a basic property of materials dictating the rate of heat transfer. For crystalline solids, the lattice thermal conductivity (LTC) is subdivided into two primary contributions, namely, by phonons and electrons. In metals, both electron and phonon thermal conductivity coexist due to the free conduction of electrons, whereas in covalently bonded structures such as ceramics, phonons dominate the heat transport. For example, thermoelectric devices are crucial for reversibly converting waste heat into electricity. These devices use a special material that, when under a temperature gradient, generates an electric current. Thermoelectric materials scale in efficiency with the dimensionless figure of merit, , whereby high $ZT$ is needed for competitive electricity generation. The thermal conductivity $κ$ is further split into phonon $κ p h$ and electron $κ e l$ contributions, whereby the sum of these two yields the net thermal conductance. To increase $ZT, κ e l$ must be low enough despite owning strong coupling with the electrical conductivity $σ$⁠. In addition, the phonon contribution $κ p h$ is required to remain negligibly low which is strongly correlated with the lattice anharmonicity, i.e., the deviation from vibrational harmonics providing finite phonon lifetimes. Provided with this understanding, a scientific basis for phonon transport and its contributions to the LTC of ordered materials is necessary for the discovery of new thermal materials for applications such as high performance thermoelectrics.

Fundamentally, phonons exist as a function of the thermodynamic state (i.e., temperature and pressure) and are dependent on the entropy of the crystal. The third law of thermodynamics states that the entropy of a perfect crystal approaches zero as the temperature approaches zero, thereby suppressing the vibrations and phonon transport. As the temperature increases, the resulting vibrations drive the entropy and dictate both the stability (free energy) of the crystal and its capacity to transfer heat. Like a spring-mass system, the displacements of oscillating atoms due to lattice vibrations are correlated with the total potential energy and the atomic forces. An approximation of the potential energy is made by considering a Taylor expansion about the atomic displacements up to a designated order. In a real solid at finite temperatures, the order of the expansion is infinitely large, but in application the interaction may be truncated for computational considerations. Relating the potential energy and the atomic displacements are the interatomic force constants (IFCs) which analogously act as the spring constants in a crystalline system. IFCs for second-order interactions, or equivalently the second-order IFCs, correspond with the quadratic nature of harmonic oscillations, whereby the phonons are non-interacting and have infinite lifetimes. Third- and higher-order IFCs are the source of anharmonicity in the crystal and provide phonon scattering mechanisms as well as finite properties such as LTC and thermal expansion.

Currently, density functional theory (DFT) is the primary method for atomic force calculation and corresponding phonon property prediction as a realistic chemical modeling tool.^17–22 The general workflow of DFT coupled with phonon Boltzmann transport equation (BTE) for predicting comprehensive phonon properties including LTC of crystalline structures will be detailed in Sec. II. Despite the high accuracy of such DFT + BTE method, given the costly nature of DFT, the calculation of atomic forces for the IFCs is an extremely time-consuming process as it requires evaluation of several hundred supercells containing displaced atoms by high resolution DFT calculations. For instance, based on the widely employed finite displacement method (FDM), the required configurations for DFT evaluation for third-order IFCs scale proportionally to $N 1× N 2× N 3× n 2$⁠, where $N 1$⁠, $N 2$⁠, and $N 3$ is the number of duplicated unit cells in each of the crystallographic directions forming the supercell and n is the number of atoms per unit cell. For instance, a $3×3×3$ primitive cell for boron arsenide requires ∼360 supercell evaluations from Phono3py code, which may take several days or even weeks depending on the computing system. Even with a certain cutoff radius, in which supercells with pairs of interacting atoms beyond a designated distance are neglected, the demand is still significant in the context of hundreds to tens of thousands of materials for LTC prediction. In general, the limitations imposed by computationally expensive DFT paired with the demanding requirements of commonplace BTE solvers currently limits high-throughput phonon properties including but not limited to LTC for large-scale material screening or database development. Therefore, there is an urgent need for circumventing the time-consuming nature of DFT + BTE calculations of phonons in the context of large potential structures.

Artificial intelligence (AI) coupled with promising machine learning (ML) techniques well known from computer science is transforming science and engineering and poised to transform discovery and innovation. This era offers a unique opportunity to establish a new paradigm for the discovery of advanced materials by leveraging databases, simulations, and accumulated knowledge, venturing into experimental frontiers, and incorporating cutting-edge AI technologies. On the one hand, a primary route that recently has taken researchers by storm is using machine learning potentials (MLPs) to replace the computationally costly DFT calculations in the phonon transport workflow. Several representative MLP algorithms and approaches will be given out in Sec. III. In general, these MLP algorithms share the features of predicting total energy and atomic forces of the systems with computational speed usually orders of magnitude faster than full DFT calculations and in the meantime with accuracy comparable or reasonably comparable to DFT. Such alternative approach for evaluating the required atomic forces and the subsequent IFCs would significantly speed up the entire workflow of phonon transport calculations. However, care must be taken for the accuracy and extrapolation ability of the trained MLPs, in particular, when processing a large number of diverse structures with a broad spectrum of constituent elements, material symmetries, etc. On the other hand, predicting phonon properties of crystalline structures can also be understood in the context of the general structure–property relationship. In fact, understanding the intricate relationship between composition, structure, and material properties has been a long standing and knotty problem for material physicists and chemists for decades, i.e., the well-known processing–structure–property–performance (PSPP) relationship. Exploring innovative materials involves navigating vast design spaces and considering chemical and structural factors on multiple scales and modalities. AI/ML algorithms and approaches provide unique strength in mapping the highly nonlinear or hidden relationships between atomic structures and material properties. In the past decade, tremendous ML techniques have been developed to analyze high-throughput data with a view to obtaining useful physical/chemical insights, categorizing, predicting, and making evidence-based decisions in novel ways, which will promote the growth of novel applications and fuel the sustainable booming of AI. While traditional materials science studies depend heavily on the knowledge of individual experts (called domain knowledge in the language of computer science), the ML models are increasingly used in materials science field because of their exceptional accuracy and efficiency. Deployment of AI/ML approaches in thermal science, in particular phonon, thermal transport field just took off recently. There are also some comprehensive review/perspective papers in this field.^23–25 In Sec. III A, we will discuss some early studies of using AI/ML methods for predicting phonon properties of crystals and their limitations, followed by some perspective on current challenges and future trends in Sec. IV. In short, there is plenty of room to unleash the power of AI in phonon thermal transport, which is the main topic of this perspective.

The phonon scattering rate and LTC can be obtained by solving the semi-classical phonon Boltzmann transport equation such as ShengBTE package,²⁸ which requires the inputs of second-order harmonic and third order anharmonic interatomic force constants (IFCs). The second order harmonic IFCs can be used to determine the phonon frequencies and eigenvectors using the Phonopy package.²⁹ For the third-order IFCs, a script of thirdorder.py can be applied to generate different supercell configurations with the consideration of both point-group symmetry and translational invariance, and the interaction cutoff up to a certain nearest-neighbor distance is usually implemented.²⁶ There are also some other efficient and accurate approaches to fit the force constants, for example, compressive sensing lattice dynamics (CSLD)³⁰ and the temperature-dependent effective potential (TDEP).^31,32

In 2018, three independent experiments^33–35 confirmed the significance of four-phonon scattering in BAs. Subsequent studies pursuing four-phonon effect have proved its importance in a broader range of systems, including insulators, semiconductors, and semimetals with topics covering thermal conductivity predictions, radiative properties, and phonon linewidths. FourPhonon module was also built within ShengBTE package and its execution is fully compatible with ShengBTE.³⁶ It uses the second- and third-order IFCs which are obtained from the methods described above. CONTROL file is also needed as an input, which specifies settings and parameters, including crystal structural information, temperature, q-mesh, broadening factor, etc. On the basis of this workflow, FourPhonon requires an additional input file with fourth-order IFC. Four-phonon scattering calculation is generally computationally expensive, since it involves huge amount of possible four-phonon process. A simple metric for strong intrinsic fourth-order phonon anharmonicity was then proposed,³⁷ which is promising to be combined with AI/ML methods to predict large-scale materials in a high-throughput manner. To obtain more accurate phonon frequencies and LTC, self-consistent phonon (SCPH) renormalization are introduced into the calculations for analyzing the effects of temperature and higher order lattice anharmonicity on phonon thermal transport properties.^38–40 For instance, as implemented in the HIPHIVE package,⁴¹ the structures for the calculations of third and fourth-order IFCs are generated using the random displacement method. The second order IFCs considering the temperature effects based on the SCPH theory are obtained by training the data using the HIPHIVE package.^41,42 Tadano et al. proposed the method of incorporating frequency renormalization effects by the bubble self-energy within the quasiparticle approximation, and then applied the developed methodology to the strongly anharmonic α-CsPbBr₃.⁴³ Tadano et al. also analyzed phonon transport including a coherent interbranch component by using first-principles based self-consistent phonon theory and solving off diagonal components of group velocity operators. Temperature-dependent harmonic interatomic force constants were then obtained for quadruple-well potential of guest atoms in type-I Ba₈Ga₁₆Sn₃₀.⁴⁴ The relevant methodology has been implemented into the Alamode package.^45,46

Note that the phonon lifetime $τ λ$ is replaced with the scattering rate $Γ q s=1/ τ λ$⁠. In addition, the group velocity is generalized to the off diagonal form $V q s , s ′$⁠. Since then, there are several studies that have explicitly quantified the contributions from off diagonal terms in the heat-flux operator^44,50–56 and even combined with a complete treatment of quartic anharmonicity for both phonon frequencies and lifetimes at finite temperatures.⁵⁷ Despite significant attention has been given to the off diagonal contributions to the overall thermal transport, currently a systematic understanding and analysis across diverse material families, complexity, symmetry, elements, bonding nature, etc. are lacking. AI/ML approaches can exactly fit in this niche as we will see below.

In the past decade, plenty of physical properties have been well trained and predicted by AI/ML models, provided that there is significant amount of training data that can be obtained from either high-throughput calculations or experiments. These successful cases of material properties span from simple thermodynamics (e.g., formation energy and energy above hull^58–60 to basic and easy-to-calculate mechanical properties (e.g., elastic constants, Young's modulus, hardness),^61–65 to complex transport properties (e.g., thermoelectrics,^66–70 superconductivity).^71–76 For properties that are hard to compute, a common strategy used in many previous studies is that, in order to get large enough number of training data, complex physics usually needs to be simplified by some empirical models, or estimates by some rough theories, or refining or improving previous semi-empirical models.^77,78 For instance, in thermoelectrics, electronic conductivity is determined by the relaxation time of electrons, which, in principle, requires complicated first-principles calculations to explicitly consider electron–electron and electron–phonon scattering, but such calculations are extremely computationally demanding. To bypass this obstacle, an empirical or arbitrary (more importantly adjustable) parameter is usually adopted, i.e., the well-known constant relaxation time of electrons (in practice, this value is chosen as about 10 fs). Despite the success of such assumptions in many cases, the rude treatment of single constant relaxation time of electrons could lead to unacceptable error if the number of to-be-predicted structures is significant, and more importantly nobody knows which materials will have how much error.

Following the procedure of using AI/ML algorithms to predict other physical/chemical properties, an immediate approach the researchers have recently pursued in thermal science is data-driven approach through extraction of vital information from already existing DFT data and/or experimental measurements to explore previously unseen structures. Arguably, the most common ML strategy used in data-driven materials science is supervised learning, whereby a vector of input variables $x ⇀$ are mapped to the desired output(s) $y$ via a tunable, continuous function $f$ returning the model $y= f( x ⇀)$⁠. As such, several papers have been published with supervised learning of the LTC with material descriptors as the input which are significantly cheaper to generate than conventional DFT-LTC calculations. Here, a critical step of training a high quality and explainable ML model is to carefully define and select the inputs for the model, the so-called feature or descriptor engineering. Indeed, material descriptors are the distinguishing factor in recent works with high-throughput LTC due to their importance in physically representing crystalline systems, dictating the prediction efficiency and accuracy. There are a few basic standards for defining material descriptors:

1. The dimensions of the descriptor should be as low as possible. The lower dimension of the descriptors, the smaller amount of data will be required for training.

2. The descriptor uniquely characterizes the material as well as the property-relevant elementary process. The ideal case is to ensure one-to-one correspondence between input (descriptor) and output (material property) to the largest extent. In such cases, it will be much easier to train a ML model.

3. Materials that are very different (similar) should be characterized and selected by very different (similar) descriptor values. This property will help ML model itself adjust the internal parameters to recognize the hidden relationship or pattern.

4. The determination of the descriptor must not involve calculations as intensive as those needed for the evaluation of the property to be predicted (simple calculation). This is extremely important for phonon property prediction, since the traditional DFT + BTE approach is very time and resource consuming.

The first big category that is straightforward to execute is to train traditional ML models with elemental descriptors and/or simple features inspired by existing phonon transport physics and empirical models, and then deploy trained ML models to predict phonon properties of new structures. This method is also called forward prediction (see the orange arrow in Fig. 2). Conventionally, the ML models are descriptor-based, where the key descriptors representing the system must first be defined prior to fitting a suitable ML model for prediction. By analyzing the impact of different choice of descriptors and four different ML algorithms including support vector regression (SVR), fully connected neural networks (FC), kernel ridge regression (KRR), and eXtreme Gradient Boosting (XGBoost), Wang et al.⁸¹ found that the XGBoost algorithm based on the descriptors of crystal structural and compositional information can accurately predict the LTC with an average mean absolute error (MAE) of 2.13 W/m K. The training dataset of LTC contains 5486 materials with the Automatic GIBBS Library (AGL) method using a quasi-harmonic Debye model and was taken from AFLOW repositories. Therefore, the accuracy of LTC training data itself is a concern. Nevertheless, feature importance analysis yielded top five important features to impact the prediction of LTC, being the compositional weighted atomization enthalpy $( avg(Δ H atomic))$⁠, the average period [avg(period)], the volume per atom (⁠ $V p a$⁠), the mass density (⁠ $ρ$⁠), the volume per primitive cell (⁠ $V p c$⁠), and the average atomic number [avg(N)]. The LTC distribution was then projected to the $avg(Δ H atomic)− V p c$ and $V p a− V p c$ dual-descriptor space, and low LTC materials can be found in the bottom-left corner of the plots, i.e., a large volume of the primitive cell and atoms. Such features learned from LTC data analysis are consistent with already known empirical physical models for LTC, such as Slack model⁸² and Debye–Callaway model.^83–86 For example, both models emphasize the importance of atomic mass, cell volume, etc.

Using a benchmark data set of experimentally measured LTC of about 100 inorganic materials, Chen et al.⁸⁷ built a Gaussian process regression (GPR) based ML model to predict LTC. They compared the model performance between the original 63-dimensional features (using the Matminer package) and the reduced 29-dimensional features by recursive feature elimination (RFE) using the linear support vector regression algorithm with fivefold cross-validation. The accuracy of the developed ML models was found to be comparable to the past semi-empirical Slack and Debye–Callaway models. Since very limited number of data was used in this study (the accuracy of the data is high though, as they come from experiments) and there is no validation on more diverse structures, the realistic accuracy and transferability of the trained ML model to broader material families, plus the potential overfitting problem (the dimension of input descriptors is on the same order as the data size), are unknown. In addition, it is worth noting that an important feature used in the model training and prediction is not easy to compute, i.e., bulk modulus, which will be a drawback for fast material screening. Similarly, using experimentally measured LTC for 350 different materials, Qin et al.⁸⁸ constructed and trained 15 traditional ML models for LTC prediction. During the training process, eight basic properties of the materials were used as descriptors (inputs) and the experimentally measured LTC values were used as targets (output). The trained deep learning models showed the highest performance of accurately predicting LTC spanning four orders of magnitude, which have a great advantage over the widely used empirical Slack model. They continued to deploy the trained four deep learning models, combined with semi-supervised learning strategy, to predict LTC of 3716 materials, and the results were validated by the optimized Slack model.⁷⁸ 

Juneja et al.⁸⁹ carried high-throughput ab initio calculations on a dataset containing 195 compounds composed of 60 binary, 85 ternary, and 50 quaternaries. Out of these 195 structures, 120 compounds were found to be dynamically stable with LTC spanning over three orders of magnitude. Pearson correlation analysis on this dataset reveals a strong dependence of LTC on some descriptors, namely, maximum phonon frequency, integrated Grüneisen parameter up to 3 THz, average atomic mass, and volume of the unit cell. Using these descriptors, they then trained a GPR-based ML model for LTC and found better performance than the Slack model. Again, they did not do further validation or study on more unseen structures by deploying their trained ML model. Same as the previous study, their dataset size is small and the two important features (namely, maximum phonon frequency, integrated Grüneisen parameter up to 3 THz) are hard to obtain from high-throughput point of view.

Similar study was conducted by Jaafreh et al.,⁹⁰ who calculated LTCs of 119 compounds at various temperatures (100−1000 K) based on DFT and then built a predictive model using various ML algorithms including decision tree, random forest (RF), gradient boosting regression tree, and extreme gradient boosting. Unlike previous studies, the accuracy of their models was validated using new cases of four compounds, which was not seen for the model before. They continued to use their model to screen 32 116 compounds in the Inorganic Crystal Structure Database (ICSD). Cs₂SnI₆ and SrS were selected to validate the ML prediction with fairly consistent results with DFT calculations. This is a comprehensive study that includes all necessary components of data-driven approach for phonon transport, including data generation, model training, prediction, and validation.

The second big category is to develop ML models without explicitly defining descriptors. As ML algorithms have been further developed in recent years, graph neural networks (GNNs) have received intense interest as a rapidly expanding class of ML models remarkably well-suited for materials applications. It is well known that finding effective descriptors in conventional ML models as described above could be very challenging for problems with a large amount of compositionally and structurally diverse materials. Instead, GNNs could potentially overcome the limitations of static descriptors by learning the representations on flexible, graph-based inputs. To date, a large number of successful GNNs have been proposed and demonstrated for systems ranging from crystal stability to electronic transport property prediction and even to catalyst chemistry. This brings a new opportunity for phonon property prediction. Zhu et al.⁹¹ used their recently prepared high-throughput LTC data and trained crystal graph convolutional neural network (CGCNN) to map with the LTC and achieved MAE of 0.14 on log-scale of LTC and R² score of 0.85 on the same log–log scale when compared with DFT-LTC. The GNN performance is comparable to that of the RF ML model with a 154-dimensional descriptor as the input (MAE of 0.12 and R² score of 0.87). The trained model was then deployed to predict the LTC of 92 919 materials taken from the ICSD. Further feature analysis shows three most important features for LTC prediction, being the average volume per atom in the ground state, average bond length, and volume per atom, which again reproduce the previous knowledge from phonon physics and empirical models. These relevant physical mechanisms can be seen from Slack's model for lattice thermal conductivity,⁸⁸ namely, $κ=A M ¯ δ n 1 / 3 Θ 3 γ 2 T$⁠, where $M ¯$ is the averaged atomic mass, $δ$ is the cubic root of the average volume per atom, n is the number of atoms in the primitive cell that determines the number of phonon branches, $Θ$ is the acoustic Debye temperature, $γ$ is the Grüneisen parameter, and T is the absolute temperature. For instance, the $δ$ parameter in Slack's model exactly means the average volume per atom, while the average bond length can be linked to the effect of $δ$ and n. They also extended CGCNN with transfer learning (TL-CGCNN) to predict experimental LTC values, with the purpose of finding new knowledge content. This attempt is interesting as it will lead researchers to search and design new materials in previously unexplored regions if new knowledge is successfully extracted.

My group recently trained various GNN models on a few thousand high-quality LTC data from our own full DFT + BTE calculations. The second- and third-order IFCs required for phonon band structure and LTC calculations were calculated using the compressive sensing lattice dynamics (CSLD) method, which extracts the IFCs from the Taylor-expanded interatomic forces in terms of atomic displacements via the advanced compressive sensing technique. Generally, for each supercell, 16–50 randomly displaced configurations were used for obtaining IFCs, depending on the symmetry of the materials. In a recent work by Ojih et al.⁹² various ML models including the state-of-the-art CGCNN, and global attention graph neural network (deeperGATGNN) were trained on 3377 high-quality DFT data, and then potential materials with high heat capacity were efficiently searched. The deeperGATGNN model exhibits high prediction accuracy and is used for predicting the heat capacity of 32 026 structures screened from the open quantum material database. Deep insight into the correlation between heat capacity and structure descriptors was gained, such as space group, prototype, lattice volume, atomic weight, etc. We also identified one structure, namely, MnIn₂Se₄, with space group No. 227, that exhibits heat capacity even higher than that of the Dulong−Petit limit at room temperature.

Many mechanical properties have been found to have strong correlations with LTC. From Eq. (1) based on the kinetic theory of phonon transport⁹³ and using the single mode relaxation time approximation of the Boltzmann equation, we know that, if the phonon group velocity is significantly reduced, the thermal conductivity is anticipated to be very low. Therefore, materials with low group velocity should have low LTC. From the physics law, we also know that the phonon group velocity can be roughly estimated as $v∝ E / ρ$⁠, where E is Young's modulus characterizing the bonding strength and $ρ$ is the mass density of the material. Therefore, Young's modulus or bulk modulus could be a good strategy for screening low LTC materials. Inspired by this idea, using 10 158 elastic constants as training data, Ojih et al.⁹⁴ first trained deeperGATGNN model on five mechanical properties, namely, bulk modulus, shear modulus, Young's modulus, Poisson's ratio, and hardness. With the trained model, they then predicted 775 947 data in search of materials with low bulk modulus and potentially low LTC. 338 structures were finally verified with first principles. The results demonstrate that one can find materials with extreme mechanical properties recommended by high fidelity GNN models and low LTC material from bulk modulus prediction with minimal first-principles calculations of the structures (only 0.04% in the case study) in the large-scale materials pool.

A critical point that has been neglected in most of existing ML models for phonon property prediction is the dynamical stability, i.e., there should be no imaginary phonon frequencies in the full Brillouin zone. It is worth pointing out the rigor of using imaginary phonon frequency as a criterion for determining dynamical stability, which could be a good recipe for high throughput of a large number of structures. However, care must be taken in some special cases. For example, some recent studies have used the self-consistent phonon theory or temperature-dependent effective potentials (TEDPs) and obtained positive dispersions at finite or elevated temperatures even if imaginary frequencies occur at 0 K phonon calculations. Moreover, even with fully positive phonon frequencies, structures may lack dynamical stability under certain elevated temperatures, especially for materials undergoing phase transitions or melting at low temperatures. Nevertheless, without confirming the dynamical stability of a crystalline structure, predicting its thermal transport property would be physically meaningless. To this end, more efforts should be dedicated to screening dynamically stable materials with high speed and high accuracy so that the prediction on LTC of final promising structures will have a high success rate. This task can also be done by ML models. Ojih et al.⁹⁵ reported an efficient workflow combining high-throughput DFT computing and two different types of ML models for fast and accurately screening ultralow LTC from large-scale inorganic crystals. First, seven classification ML models on 8077 data obtained from high-throughput full DFT calculations were trained to classify 50 574 structures into positive and negative dispersions, among which 22 899 structures are predicted to be dynamically stable. Second, with 4041 high-quality DFT-LTC data, three GNN models were trained and used to predict LTC, with 359 randomly selected structures verified by full DFT calculations. The result showed the ML model successfully predicted 90% of 359 structures to possess ultralow LTC (lower than 1 W/m K). An additional 3218 structures with ultralow LTC are also predicted and provided. This workflow integrating dual ML models offers a new route to accelerate the discovery of novel dynamically stable materials with a high success rate for predicting effective LTC. Moreover, the correlation analysis reproduced the phonon transport physics for ultralow LTC from two aspects: (1) the large P₃ parameter represents a large number of three-phonon scattering channels; (2) the large thermal mean squared displacement (MSD) reflects the soft phonon modes in the lattice usually resulting in strong phonon anharmonicity.

Despite the powerful strength of ML models in phonon transport prediction, either traditional ML or state-of-the-art GNN model, direct establishing relationship between atomic structures and target material property (e.g., LTC, heat capacity, etc.) may lose the comprehensive information for phonons as would have been obtained in the intermediate steps if traditional DFT + BTE workflow is applied (see the red arrows in Fig. 2). As a matter of fact, the comprehensive properties of phonons are beneficial for many relevant emergent applications. For example, materials with target ultrahigh cutoff phonon frequencies would be potential candidates for high temperature phonon-mediated superconductivity.^96–101 Large negative Grüneisen parameter in the acoustic phonon modes, which are accessible after knowing the anharmonic IFCs, would help design negative thermal expansion materials near room temperature^102–104 and more interestingly zero thermal expansion materials over a wide temperature range.^105–107 Quantitative characterization of vibrational entropy due to the contribution of different phonon modes will facilitate design of advanced high entropy alloys (HEAs), in particular, vibrational entropy stabilized HEAs.^108–111 This motivates us to re-think the whole design of ML models for phonon properties. In particular, one step backward to the atomic force level would be a good try.

Regardless of the ML models, again keep in mind that our central task is to circumvent the time-consuming nature of DFT + BTE approach for phonon-based computation. One approach, in particular, without using AI/ML is to minimize the number of displaced supercells and corresponding DFT calculations to compute the IFCs, as this step is most costly and time-consuming in the pipeline. A few methodologies have been proposed in recent years. For example, the Hiphive package in which a cluster representation yields a reduction of the number of required configurations for LTC of monolayer MoS₂ from 571 to only 20 to 25 configurations while maintaining the accuracy of IFCs.¹¹² Similar approaches include aforementioned CSLD method and temperature-dependent effective potential (TDEP).^31,32 These methods can reduce the total number of displaced supercells that need to be evaluated by DFT by from several folders to up to one order of magnitude. While they are efficient for predicting phonon properties of single material or a limited number of materials or material families, the required DFT computational demand is still unbearable when dealing with large-scale materials (say, a few thousand or more). Nevertheless, these methodologies are very promising to accelerate phonon calculations when combined with MLPs, as we will see shortly.

To get the desired LTC values significantly faster than traditional DFT + BTE approach using AI/ML algorithms, a very early study was conducted by Carrete et al.¹¹³ who combined ML algorithms, physical insights, and automatic ab initio calculations and scanned approximately 79 000 half-Heusler entries in the AFLOWLIB database. Specifically relevant to ML, part of their LTC values were obtained by random forest regression (denoted as $κ forest$⁠) by leveraging the LTC of 32 Heusler materials fully calculated by DFT as a training set. They then employed the fitted model to predict the remaining LTCs. The random forest algorithm works well in this study partially because the crystal structures are the same half-Heusler type and the only changes among different structures are the constituent elements occupying the lattice sites. This method can be extended to similar problems such as predicting thermal transport properties of high entropy alloys or alloying structures (e.g., Al_xGa_1−xN), where the main lattice symmetry and the relative position of atoms and their neighbors keeps the same while only the atom type and mass on the lattice sites changes. Certainly, this method cannot be used to predict structure pool with diverse material families and symmetries and mixed atomic environment, etc. A fundamentally new approach should be developed in those cases.

While training aforementioned ML models with material descriptors offers physical insights toward feature importance for LTC prediction, limitations are present when facing high throughput. Mainly, these models are still required to generate reference LTC data to serve as the target during training. ML models such as artificial neural networks (ANNs) depend on data diversity due to their interpolative nature, i.e., they cannot perform well when provided data outside the training set. As such, the data generation for a sufficiently robust model is expensive and may limit the predictions to a small subset of materials. Additionally, because these models are usually trained on one temperature designed to output a single value of LTC, they are unable to provide the plethora of information that comes with DFT + BTE LTC calculations. This means that the plethora of information such as phonon dispersion, mode-dependent scattering rates, temperature-dependent LTC, phonon frequency-dependent accumulative LTC, and even off diagonal contribution to overall thermal conductivity are inaccessible, all of which however are standard outputs in DFT + BTE LTC calculations.

To circumvent these issues, the AI/ML model for phonon thermal transport may be re-designed from a lower level, more specifically through the atomic forces. As mentioned previously, the IFCs are critical for modeling harmonic and anharmonic properties and their contributions to the LTC and are derived from the atomic displacements in supercells and their corresponding forces. Consequentially, approaching BTE solvers with already computed atomic forces from ML maintains the rich output of information that comes standard from phonon calculators. More importantly, training for the atomic forces has the potential to reduce the costly demand for training set generation. As the primary cost for model development is the DFT calculation needed for training, the goal should be to maximize the ratio of datapoints to DFT calculations while maintaining accuracy. For models with material descriptors, many DFT calculations are required for one single LTC, serving as a single datapoint for training. Meanwhile, training on atomic descriptors is advantageous in terms of data abundance per DFT run because each simulation provides $( 3 N + 1)$ data corresponding to N atoms worth of force vectors and one total energy. Models trained on such data are referred to as machine learning potentials (MLPs). Namely, MLPs capture the explicit electronic-level features from DFT by implicit representation of the potential energy surfaces as functions of the atomic nuclei positions. Due to the purely mathematical nature of ML models, to capture the appropriate physics the accuracy of MLPs strongly depends on the description of the atomic environment surrounding central atoms. It is well known that cartesian coordinates are not well-suited for the description of atomic environments, since (1) they are not rotationally invariant and cannot physically represent the potential energy, (2) the order in which the cartesian coordinates are fed into MLPs is ill-defined and could effectively “confuse” models, and (3) atoms may leave/enter the designated cutoff range and is incompatible with MLPs due to the allocation of weights at the beginning of training. As such, several MLPs have been proposed over the past decade by using various improved descriptors. Many studies have shown excellent representation of DFT-level energetics and realistic property prediction with MLPs. Some representative MLPs in this line include the high dimensional neural network potential (HDNNP),¹¹⁴ deep potential molecular dynamics (DeePMD),¹¹⁵ and GAP with SOAP descriptors,¹¹⁶ higher order equivariant message passing neural networks (MACE).¹¹⁷ Those MLPs have been applied for a wide array of materials, including but not limited to molten salt,¹¹⁸ metals,¹¹⁹ semiconductors, organic molecules,¹²⁰ superionic conductors,^121,122 and even amorphous structures.¹²³ In the context of phonon property prediction, several studies have been published using MLPs as the force calculator. For instance, Marques et al. reproduce the phonon dispersions for cubic Si and Ge in comparison to those found from DFT,¹²⁴ and Minamitani et al. capture both Si and GaN dispersions and LTC within 5.4%.¹²⁵ Typically, all of these studies own a root mean square error (RMSE) of the force predictions from the MLP within a few tens or above $100 meV/Å$ which is at least one or two orders of magnitude less accurate than DFT pseudopotentials but with approximately 1000× faster evaluation time.

Undoubtedly, the robustness of the recently developed MLPs has the potential to mitigate the current speed-related bottlenecks in the DFT + BTE LTC workflow. However, to date a majority of studies using MLPs also share a common denominator in that the models are limited to a material-to-material basis or a fractional material family.^126–128 This is primarily due to the exponential scaling of model parameters with the number of atomic species or elements (denoted as $N elem$⁠) contained in the training data. For instance, the HDNNP requires $N elem$ element-specific networks each containing approximately $N elem$ radial and $N elem( N elem + 1)$ angular symmetry functions. When faced with data containing elements spanning the periodic table, the training efficiency and evaluation time is reduced significantly due to the $∼ N elem 2$ scaling of the input descriptor. Additionally, training of each element-specific network requires central atoms dedicated only to said element, e.g., copper networks are only trained on copper data, meaning that little to no knowledge of atomic environments from other central atom species are shared. For instance, an atom sharing the same column on the periodic table may provide additional electronic-level information for another species due to their similar valences. In general, recent MLPs represent atomic positions numerically and the elements with sub-models and/or specific input descriptors, in turn diminishing prediction quality with ten or more elements. Furthermore, descriptors should be designed to scale independently of the elements without needing to reconstruct and train the model from scratch. This opens the opportunity to develop MLPs without the need to re-fit ab initio data when new chemical element is introduced, allowing for MLPs to grow rather than require retraining from the ground-up. Overall, independent elemental scaling and centralized ML training are two major factors necessary for evaluation of theoretical materials databases containing a plethora of structures and atomic species that would otherwise be too difficult to handle with modern MLPs. Computing forces across many atomic environments is especially challenging for high-throughput LTC considering the notoriously strict force accuracy requirements for the IFCs and the resulting LTC.

In my group's recent studies,^79,80,129 we focus on the development of accurate force calculators including high scope of transferability between atomic structures and elements. We have pipelined the first-principles calculations for atomic structure optimization, high precision DFT calculations, and postprocessing DFT results to generate training data, and automated training of deep neural network potential for obtaining phonon properties (Fig. 3). The unique feature of our recently developed machine learning algorithm, dubbed “spatial density neural network force field (SDNNFF)” (also see the green box in Fig. 2),⁸⁰ is the usage of a three-dimensional mesh of density functions, which together act to map the atomic environment and provide a physical representation of the forces acting on the central atom. The high efficiency and accuracy of our SDNNFF method benefits from several notable advantages, including: (1) avoiding the chain rule of derivatives of the total energy (a big source of error in obtaining atomic forces from machine learning models); (2) using a single atom and its neighboring environment to train the model, which gives us ever increasing $n×m$ scaling of the training data, where n is the number of atoms in a supercell and m is the number of first-principles-evaluated structures; (3) significantly reduced number of parameters and human effort needed to successfully train a per-atom property-converged neural network model. We first tested our SDNNFF method on thermal transport in simple crystals, namely, bulk diamond, Si, SiC, and BAs.¹²⁹ Our new SDNNFF method yields unprecedently low root-mean-square-error (RMSE) of atomic forces of 1.5 meV/Å for bulk Si with phonon dispersions and temperature-dependent LTC excellently matching the DFT and experimental results, which outperform existing MLPs reported so far. We also did a benchmark study on amorphous Si and the results show that our SDNNFF method has much faster convergence in training process than other neural network potentials with several folders to even one order of magnitude smaller RMSE in predicting atomic forces, confirming its high efficiency and accuracy.

Following our initial SDNNFF model, the newest solution for descriptor development including both atomic positions and elements requires only a single neural network with independent scaling of input size with respect to the available species in the training set. The basic idea of the new model, dubbed the Elemental-SDNNFF,^79,80 is to simultaneously capture the previously accurate spatial mapping of neighbors and how atomic elements influence the signals measured at allocated 3D grid points surrounding the central atom. Two advantages arise from the Elemental-SDNNFF descriptor: (1) the summation of density functions multiplied by the neighboring atomic weights eliminates the need for designated slots in the descriptor vector for each element and removes the scaling of input size with respect to number of elements; (2) by providing the central atom atomic weight in the input descriptor vector, the network can distinguish central atoms whereby individual element-specific SDNNFFs are not required. The result is a neural network force field capable of modeling atomic systems spanning the periodic table without sacrificing efficiency and power of deep neural network training. To demonstrate the effectiveness and fast speed of our new Elemental-SDNNFF algorithm, we apply it to train on a large dataset containing a mixture of ∼80 000 cubic materials, totaling 63 unique elements across the periodic table. With individual atoms in each supercell as separate input data entry and additional data augmentation technique, the total amount of training data size reaches ∼4 × 10⁷, which is 3–4 orders of magnitude larger than traditional MLPs that only treat the entire supercell system as single input data. Benefiting from such an extremely large training dataset, the Elemental-SDNNFF can capture both atomic structure and species with significant improvement to training speed and accuracy. The RMSE of atomic forces for the latest training round reaches a level as low as 29.3 meV/Å, which is a few folders to one order of magnitude lower than most of existing MLPs. The RMSE of predicted phonon frequencies for >90% of tested structures is within only 3%, indicating that our Elemental-SDNNFF algorithm has predicted the interatomic force constants very accurately (Fig. 4). We further used the Elemental-SDNNFF to predict LTC of these structures, and good agreement has been achieved as compared to the full first-principles calculations. These codes and workflow are expected to be used for more complex and larger amount of new local atomic environment learning and thus precisely predict phonon properties of sub-million-scale materials.

With predicted second- and third-order IFCs and then further predicted phonon properties for 80 198 cubic structures,⁷⁹ we found an outstanding performance with the maximum mean square displacement (MSD) of atoms and incorporate this as a descriptor for lattice thermal conductivity (⁠ $κ L$⁠) of crystals. This is most likely due to the all-encompassing nature of thermal MSD in describing other vibrational properties. Additionally, MSD may be computed as a function of temperature and is more useful to observe temperature-dependent trends of almost all vibrational dominated or related properties. A generally inverse-linear relationship is observed between the log $κ L$ and the log(max MSD) (Fig. 5), which provides evidence of linearity through predicted $κ L$ and the max MSD. Structures with extremely high maximum MSD are indicative of rattling atoms, where strong phonon–phonon scattering and ultralow LTC is prevalent. These predictions remain highly beneficial for quickly marking structures with strong phonon anharmonicity, which is closely related to mechanical and transport properties. Specifically, we found 9306 total structures with norm. log(max MSD) higher than 0.464, among which 8873 (95.4%) structures possess LTC below 1 W/m K. Thus, the max MSD is a reliable approach for indicating highly unique structures with out-of-trend values of atomic displacement and corresponding strong phonon anharmonicity, which provides a route to narrow down potential candidates in searching ultralow LTC. The results on 80 198 cubic structures show great promise of our Elemental-SDNNFF approach for accurate, high-throughput, and comprehensive phonon property prediction including phonon dispersion, scattering rates, and LTC at a fraction of the expensive computational cost of traditional DFT + BTE workflow for phonon property calculation.

Despite the data-driven AI/ML methods have shown great potential in reducing computational cost of well executed DFT + BTE approach and accurate phonon properties prediction, there are still a lot of challenges and issues that limit their further widespread applications:

1. Quality and quantity of phonon data: Many phonon properties are complex properties that are hard to either be calculated by DFT or train by AI/ML models, and the structurephonon property relationship is highly nonlinear and implicit. To let the AI/ML models cover as many material families and symmetries as possible, a large amount of training data are required, in terms of both quality and quantity. The DFT + BTE approach will still be the major way to obtain high accuracy phonon data in the near future. To address this challenge, further work can be pursued through database infrastructure construction. Creating a public phonon database (see the blue arrow in Fig. 2) containing comprehensive phonon properties will not only benefit the thermal science field, but also become an excellent complement to the currently widely used material databases, such as Materials Project,¹³⁰ OQMD.¹³¹ In such phonon database, a web server for non-profit researchers to effectively retrieve, quickly navigate, visualize, and compare large pools of phonon band structures to pinpoint the materials of interest is expected. More importantly, the new phonon database will provide multi-channels for communications between users and developers. For example, users around the world can calculate their own DFT data with standard input and output files compatible with the new phonon database, and then upload or contribute to the phonon database. Such user-interactive toolsets are expected to accelerate the accumulation speed of high-quality phonon data and will benefit broad communities of material physics, chemistry, and engineering with new collaborative opportunities for novel materials discovery in many societally important areas.

2. Machine learning of multi-resolution data: Currently, the majority of accurate phonon calculations are all based on first-principles calculations, which is too computational demanding for high-throughput screening. That is the main reason why most existing ML studies only used a few hundred LTC data. On the other hand, there are some theories that have captured the correct phonon physics and can be used to generate large amounts of data in a quick way, but the predictions are less accurate. It is then intuitive to develop deep learning models that can take advantage of heterogeneous multi-resolution phonon data and scarce data for training efficient and accurate phonon prediction models. This approach allows one to combine the strength of data-driven modeling to learn physically consistent phonon transport prediction models, which may have higher generalization capability as well to extend to previously unexplored material region. The latter is always a challenge and issue for the current AI/ML algorithms, i.e., the so called extrapolation problem (see details below). There are a few ML algorithms in computer science that researchers from thermal transport field can borrow.

3. High fidelity universal MLPs: MLPs have been proved to accelerate the phonon property predictions by order of magnitude faster. However, current MLPs are usually well designed and trained on small amounts of data, and thus those MLPs can be only applied to a limited number of materials, or a material-to-material basis, or a fractional material family. It is still a grand challenge to train ultra-large MLPs that can cover broad material families, diverse symmetries, various atom species, and constituent elements. Initial attempt of developing such universal or pre-trained MLPs for elements across the periodic table have been made in recent years, such as CHGNet¹³² and M3GNet.¹³³ These universal MLPs are successful in some degree primarily benefited from the millions of high resolution DFT calculations accumulated over the development of the Materials Project database.¹³⁰ While it is impossible to enumerate all possible atomic configurations to do DFT and then train universal MLPs, an effective approach could be using ML algorithms to detect the “holes” in the current potential “net” so that small amount of DFT calculations can be targeted to conduct and then the fidelity and robustness of the MLPs can be further improved. The speed of MLPs is only slightly slower than empirical theories, so universal MLPs hold great promise to quickly screen ever-large-scale unknown structures.

4. Generalization ability or extrapolation of trained ML models and inverse design: Current ML models are good at intrinsic interpolation, i.e., they are capable of make accurate predictions to the areas that are in the vicinity of original training data, while it is challenging to extrapolate or predict material properties beyond the training scope. In fact, materials with exceptional properties are more interesting to researchers but such outliers are usually “alone” in the vast material space and thus it is hard to find them. Therefore, it is desirable that future ML models should be designed as generative rather than interpolative and also should be able to generate new structures with target phonon properties, i.e., the so-called inverse design (see the yellow arrow in Fig. 2). This requires the ML model to understand hidden laws from the training data rather than relying on prior intuition, and then be able to generate hypothetical structures or “first materials” with specific properties. Ongoing efforts have been made to explore and advance in these directions. Encouraging news were presented by a new AI tool, GNoME,¹³⁴ developed by the Google DeepMind team, where the authors claimed that > 380 000 stable materials have been predicted using a strategy of active learning combined with GNN. Despite a big question mark on whether those predicted materials are truly stable, in particular, dynamically stable (absence of negative frequencies) from phonons point of view, the predictions were claimed to span orders of magnitude beyond human knowledge.

5. Extend ML models beyond phonon properties of single crystals and even beyond phonons: Thermal transport in heterogeneous materials other than single crystals is equally important for many applications, such as thermoelectric energy conversion. ML models that are trained on phonon properties of single crystals can be extended to lots of systems where structural heterogeneity occurs or even dominates. For example, predicting thermal transport across interfaces requires the knowledge of phonon properties of not only the two bulk materials in contact, but also the interface they form. The application of AL/ML in predicting the phonon thermal transport in the presence of interfaces is currently an important open research topic. An interface loses its material symmetry and thus cannot be treated as single (periodic) crystal anymore, which disables the continuing applications of the well deployed recipe of ML algorithms. However, the aforementioned universal MLPs can be well-suited to study interfacial thermal transport, e.g., by nonequilibrium molecular dynamics (NEMD) simulations.^135,136 To save computational cost of NEMD simulations (possible interfaces or material pairs scale up as $N 2$⁠, where N is the total number of single crystals), one can think of predicting phonon density of states (DOS) of two contacting materials first and then calculating their phonon DOS overlap as a rough estimate to narrow down the pool of promising material pairs. Quantitative predictions of interfacial thermal conductance (ITC) can be also made by using DMM model as implemented in the almaBTE package¹³⁷ to help filter out unnecessary candidates further. All these calculations and screening process can be accelerated by incorporating ML models, such as training phonon DOS and ITC, designing new substrates with both high ITC and high LTC for cooling specific heat source material. In the presence of multiple interfaces, one of the main scientific questions is the role of disorder in interfacial phonon transport. Several studies have been devoted to optimizing aperiodic superlattices to minimize phonon thermal conduction with AL/ML.^138–144 ML-enabled research have showed the possibility of minimization of phonon thermal conduction for moderate disorder in heterogeneous superlattices against previous physics intuition for width-modulated homogeneous superlattices. These recent developments on the power of AI in detecting underlying physics rules governing phonon thermal transport will inspire more and deeper studies in this field.

With further development of improving AI/ML algorithms for phonon thermal transport, ultimately, the AI/ML-based approaches are poised to revolutionize our understanding of new phonon thermal transport phenomena in vast material space and enable rational design of novel thermal materials with desired or controllable thermal transport properties, ushering in a new era of accelerated innovation and advancement in thermal science.

This project was supported in part by the NSF (Award Nos. 2030128, 2110033, 2311202, an 2320292), SC EPSCoR Program under Award No. 23-GC01, and an ASPIRE grant from the Office of the Vice President for Research at the University of South Carolina (Project No. 80005046).

The author has no conflicts to disclose.

Ming Hu: Conceptualization (lead); Funding acquisition (lead); Project administration (lead); Supervision (lead); Writing – original draft (lead); Writing – review & editing (lead).

Data sharing is not applicable to this article as no new data were created or analyzed in this study.