Disentangling the Discrepancy Between Theoretical and Experimental Curie Temperatures in Ferroelectric PbTiO$_3$

TL;DR

The paper investigates why theoretical predictions of the Curie temperature $T_c$ for PbTiO$_3$ systematically fall short of experimental values. By conducting the largest constant-pressure AIMD simulations and benchmarking against MLFF-based MD (DP, NEP, and qNEP), the authors show that the primary source of error is the exchange–correlation functional (PBEsol), not the MLFF fitting. They reveal a nuanced interplay between finite-size effects and long-range electrostatics: short-range ML descriptors can artificially stiffen the lattice and overestimate $T_c$, while explicit long-range corrections bring the predicted $T_c$ down toward a converged limit around $600$ K. The work underscores the necessity of long-range electrostatics and improved exchange-correlation functionals for reliable finite-temperature predictions in ferroelectrics and provides a high-quality dataset and methodological guidance for developing robust ML force fields in strongly anharmonic systems.

Abstract

Accurately predicting the Curie temperature ($T_c$) of ferroelectrics from first principles remains a major challenge, as theoretical estimates often fall significantly below experimental values. In this work, we investigate the origin of these discrepancies in the prototypical ferroelectric PbTiO$_3$ by performing extensive constant-pressure ab initio molecular dynamics (AIMD) simulations and benchmarking them against classical molecular dynamics (MD) using machine learning force fields (MLFFs) derived from first-principles data. Our results show that the underestimation of $T_c$ primarily stems from the limitations of the exchange-correlation functional, rather than inaccuracies in the MLFF fitting. We uncover a critical interplay between finite-size effects and the range of interatomic interactions: although short-range MLFFs appear to yield better agreement with experimental $T_c$, this improvement results from a fortuitous cancellation of errors. Incorporating explicit long-range interactions improves accuracy for larger supercells but ultimately leads to lower predicted $T_c$ values. These findings highlight that accurate finite-temperature predictions require not only high-quality training data and sufficiently large simulation cells, but also the explicit treatment of long-range interactions and improved exchange-correlation functionals.

Figures (6)

• Figure 1: Time evolution of lattice constants and polarization fluctuations near the Curie temperature ($T_c$), obtained from $NPT$ ab initio molecular dynamics simulations using a $4 \times 4 \times 4$ (320-atom) supercell. Panels (a)–(c) present the lattice constants at 450, 500, and 550 K, respectively, while panels (d)–(f) show the corresponding polarization components along the Cartesian axes. A switching of the polar axis is observed at 500 K, indicating a transition temperature ($T_c$) of approximately 500 K.
• Figure 2: Temperature dependence of (a) lattice constants and (b) polarization components along the Cartesian axes. Error bars represent the standard deviation, with significantly larger fluctuations observed above $T_c\approx 500$ K. Panel (c) shows the distribution of the polarization along the polar axis ($P_z$) of the supercell at various temperatures. At $T_c \approx 500$ K, a substantial population emerges near $P_z = 0$. The analysis at each temperature is based on configurations sampled from the final 20 ps of the equilibrium AIMD trajectories.
• Figure 3: Comparison of DPMD and AIMD results. Parity plots comparing (a) energies and (b) atomic forces predicted by the DP model against reference values from AIMD configurations not included in the training dataset. The fitting errors are comparable to those from the training set, indicating that the model is not overfitting and demonstrates good generalizability. Panels (c) and (d) compare the temperature dependence of lattice constants, tetragonality, and polarization component $P_z$ as predicted by DPMD and AIMD. At each temperature, ensemble-averaged properties are computed using the same number of supercell configurations sampled over an equivalent duration of equilibrium trajectories.
• Figure 4: Comparison of configuration diversity between datasets from DPGEN and AIMD simulations. (a) Distributions of energies (left) and atomic forces (right). The DPGEN dataset (red) spans a significantly broader range compared to the AIMD dataset (blue). (b) Principal component analysis (PCA) of structural descriptors. The AIMD configurations (blue), obtained from a 310 ps trajectory, are confined to a narrow region near the tetragonal ($P4mm$) ground state. In contrast, the DPGEN dataset (red) explores a much broader region of configuration space. (c) Comparison of the temperature dependence of the polarization component $P_z$ (top panel) and the $c/a$ ratio (bottom panel), showing results from a DP model trained only on AIMD configurations (DP$_{\text{AIMD}}$) versus reference AIMD results. (d) Comparison of radial distribution functions $g(r)$ obtained using different methods.
• Figure 5: Supercell size effect on the Curie temperature. (a) Temperature dependence of the $c/a$ ratio for supercells of size $L \times L \times L$. The Curie temperature ($T_c$) converges at $L = 8$. (b) Temperature evolution of the distribution of unit-cell-resolved local polarization components $p_z$, obtained from DPMD simulations using a $10 \times 10 \times 10$ supercell.