Large-scale molecular dynamics simulations of negative thermal expansion in Bi1-x Lax CoO3 using a pretrained machine learning force field

Large-scale molecular dynamics simulations of negative thermal expansion in Bi1-x Lax CoO3 using a pretrained machine learning force field | Research Square window.SnipcartSettings = { analytics: { enabled: false } }; (function() { var accessVector = localStorage.getItem('access_vector') || ''; window.dataLayer = window.dataLayer || []; if (accessVector) { window.dataLayer.push({ user: { profile: { profileInfo: { snid: accessVector } } } }); } })(); (function(w,d,s,l,i){w[l]=w[l]||[];w[l].push({'gtm.start':new Date().getTime(),event:'gtm.js'});var f=d.getElementsByTagName(s)[0],j=d.createElement(s),dl=l!='dataLayer'?'&l='+l:'';j.async=true;j.src='https://www.googletagmanager.com/gtm.js?id='+i+dl;f.parentNode.insertBefore(j,f);})(window,document,'script','dataLayer','GTM-K279D39R'); Browse Preprints In Review Journals COVID-19 Preprints AJE Video Bytes Research Tools Research Promotion AJE Professional Editing AJE Rubriq About Preprint Platform In Review Editorial Policies Our Team Advisory Board Help Center Sign In Submit a Preprint Cite Share Download PDF Article Large-scale molecular dynamics simulations of negative thermal expansion in Bi1-x Lax CoO3 using a pretrained machine learning force field Shogo Wakazaki, Yuki Sakai, Kazuki Takahashi, Hena Das, Masaki Azuma This is a preprint; it has not been peer reviewed by a journal. https://doi.org/ 10.21203/rs.3.rs-8425406/v1 This work is licensed under a CC BY 4.0 License Status: Published Journal Publication published 10 Apr, 2026 Read the published version in Scientific Reports → Version 1 posted 14 You are reading this latest preprint version Abstract Negative Thermal Expansion (NTE) is a key property for achieving dimensional stability in advanced technologies, yet understanding its underlying mechanisms, particularly the dynamic processes of phase transitions, remains a formidable challenge. In this study, we applied large-scale molecular dynamics simulations using the Crystal Hamilton Graph Network (CHGNet), a pretrained universal machine learning force field (MLFF), to elucidate the atomic-level mechanisms of the NTE phenomenon in the perovskite oxide system Bi 1 − x La x CoO 3 . CHGNet achieved a computational speed approximately 20,000 times faster than first-principles calculations and qualitatively reproduced the chemical trend of the NTE transition temperature systematically decreasing as La substitution increased. Leveraging this computational efficiency, we successfully visualized the dynamic process of the NTE phase transition at the atomic level, which had previously been unobservable. The results revealed a nucleation and propagation mechanism wherein the phase transition initiates heterogeneously in La-rich regions and propagates to Bi-rich regions. These findings demonstrate that universal MLFFs are powerful tools for more quickly elucidating complex phase transition phenomena and open up new possibilities for computational science-driven exploration of new materials. Physical sciences/Chemistry Physical sciences/Materials science Physical sciences/Mathematics and computing Physical sciences/Physics Figures Figure 1 Figure 2 Figure 3 Figure 4 Figure 5 Figure 6 Figure 7 Figure 8 Figure 9 Figure 10 Figure 11 Introduction While most materials expand upon heating, a few materials exhibit the peculiar property of contracting, known as "Negative Thermal Expansion" (NTE). This property is highly sought after for applications in advanced technologies requiring dimensional stability under extreme temperature changes, such as precision optical instruments and semiconductor manufacturing equipment, as it enables composite materials to be created that experience near-zero overall thermal expansion when combined with conventional positive thermal expansion materials. The mechanisms driving NTE are diverse, but a type driven by a first-order phase transition from a large-volume, low-temperature stable phase to a small-volume, high-temperature phase is of particular interest due to the colossal volume changes often involved. 1 Research on high-pressure synthesized perovskite oxides provides representative examples of this phase-transition-driven NTE. For instance, in the PbVO₃ derivatives with PbTiO 3 type structure with a giant c / a ratio, a colossal NTE exceeding 8% has been reported, induced by a phase transition from a polar ferroelectric structure to a nonpolar paraelectric structure. 2 , 3 In the BiNiO 3 derivatives, NTE occurs through a phase transition from a divalent Ni phase with a charge-disproportionated state on the Bi site (Bi 3+ /Bi 5+ , Ni 2+ ) to a uniform trivalent Ni phase (Bi 3+ , Ni 3+ ) involving charge transfer between Bi and Ni. By optimizing elemental substitution in this system, the phase transition temperature has been controlled to achieve giant NTE near room temperature. 4 Thus, phase-transition-driven NTE is a complex phenomenon where multiple degrees of freedom—lattice, charge, and spin—are intricately coupled, making the atomic-level elucidation of its mechanism a critical challenge from both fundamental and applied perspectives. The perovskite oxide BiCoO 3 , the focus of this study, belongs to this family of phase-transition-type NTE materials. Under ambient pressure, the stereochemically active 6 s 2 lone pair of Bi 3+ and the Jahn-Teller effect of high-spin Co 3+ lead to a large polar displacement of Co within the CoO 6 octahedra as schematically illustrated as the inset of Fig. 1 (b), resulting in a pyramid-like coordination and a ferroelectric tetragonal structure (space group P 4 mm ) with a significantly elongated c -axis as c / a = 1.27. 5 When a pressure of approximately 3 GPa is applied, it undergoes a spin-crossover phenomenon where Co 3+ ions transform from high-spin to low-spin states. 6 This spin-crossover is accompanied by a transition to a nonpolar paraelectric orthorhombic phase with tilted, isotropic CoO 6 octahedra, resulting in volume contracting colossally by about 13%. This small-volume phase, stable under high pressure, can be stabilized at ambient pressure through elemental substitution. Although irreversible, temperature-induced NTE has been reported in systems where Bi and Co are substituted for Ba and Ti, respectively. 7 Similarly, in the Bi 1 − x La x CoO 3 system studied here, substituting La 3+ for Bi 3+ is expected to stabilize the high-pressure phase, allowing thermal energy to drive the transition and induce NTE. First-principles calculations, particularly those based on Density Functional Theory (DFT), are a powerful tool for elucidating the mechanisms of such complex phase transitions at the atomic level. 8 However, the computational cost of DFT scales as O ( N 3 ) with the number of atoms N , where " O " denotes the order of complexity in Big O notation, primarily due to the cubic scaling of core algorithms like matrix diagonalization required to be solved for the electronic wavefunctions. This makes simulations extremely difficult for phenomena occurring on time and length scales greater than a few 100 atoms and nanoseconds. To resolve this trade-off between accuracy and cost, Machine Learning Force Fields (MLFFs) have recently emerged as a new paradigm. 9 MLFFs are trained on vast datasets of energies and forces obtained from DFT calculations to rapidly reproduce the potential energy surface from atomic configurations. The development of "universal force fields" pretrained on large-scale databases is further accelerating this trend. M3GNet 10 , MACE 11 , and CHGNet 12 (Crystal Hamilton Graph Network, which is used in this study) are prominent examples aimed at applicability to a wide range of unknown materials across the periodic table without needing material-specific fine-tuning. The effectiveness of these universal force fields is already being demonstrated in phase transition simulations of other perovskite oxides, representing a major step toward high-throughput screening of computational materials. 13 Recently, the application of MLFFs to NTE materials has gained attraction, offering atomic-scale insights that are difficult to access experimentally. For instance, MLFFs based on the Deep Potential scheme have successfully elucidated the role of polyhedral vibrations in the pressure-induced amorphization of ZrW 2 O 8 14 , the "lateral and tilt" vibrational modes in Zn(CN) 2 15 , and the impact of octahedral tilting on the NTE magnitude in fluorides such as ScF 3 and CaZrF 6 16 . Furthermore, modular development strategies for MLFFs have extended these capabilities to solid solutions like Pb x Sr 1− x TiO 3 , reproducing complex topological textures 17 . However, these pioneering studies typically rely on constructing system-specific models, which necessitate the generation of extensive DFT datasets and dedicated training processes for each target material system. This "train-from-scratch" paradigm remains a bottleneck for high-throughput exploration, particularly when searching for new compositions or investigating substitution effects across a wide range of elements (e.g., exploring various dopants in BiCoO 3 ). This is where "universal" MLFFs pretrained on massive databases, such as CHGNet, offer a distinct advantage by eliminating the need for ad hoc DFT calculations. Yet, the applicability of these universal models to transition metal oxides exhibiting colossal volume changes (> 5%) and discontinuous first-order phase transitions, phenomena deeply coupled with electronic and spin states, has not been sufficiently verified. Unlike the continuous NTE mechanisms seen in framework structures like ZrW 2 O 8 , the phase-transition-driven NTE in Bi-based perovskites involves drastic structural rearrangements and electronic structure changes that challenge the transferability of pretrained models. Therefore, demonstrating that a universal MLFF can reproduce these complex behaviors in a solid solution system without additional training represents a critical step toward the democratization of computational materials design. In this study, we focus on the Bi 1 − x La x CoO 3 system and perform Molecular Dynamics (MD) simulations using both CHGNet and first-principles calculations to reproduce the finite-temperature phase transition and elucidate its underlying principles at the atomic level. Results Figure 1 shows the pressure dependence of the lattice parameters and the resulting unit cell volume of BiCoO₃, as calculated by the universal MLFF, CHGNet. In order to reproduce the structural transition from P 4 mm to Pnma phases and the random distribution of substituted elements, the symmetry is lowered to P1 as indicated in the cif file found in the supplementary information. To clarify the origin of the pressure-induced phase transition, we first analyzed the changes in lattice parameters in detail. In the low-pressure phase, the structure, which is experimentally known to be tetragonal ( P 4 mm , a = b ), showed a slight deviation from this symmetry in our simulation, resulting in distinct a - and b -axis lattice parameters (see Fig. 1 a). To clarify the origin of the symmetry breaking ( a ≠ b ) observed in the MLFF simulation at 0 GPa, we performed static energy calculations for pure BiCoO 3 using both CHGNet and DFT (VASP). We compared two structures: the experimentally reported tetragonal phase ( a = b ) and the slightly distorted monoclinic-like phase ( a ≠ b ) observed in our MLFF simulation. The DFT energies were corrected using the MaterialsProject2020 Compatibility scheme implemented in pymatgen to ensure consistency with the training data. 18 The calculations revealed a subtle discrepancy in the ground state stability. While DFT (VASP) correctly predicts the tetragonal phase to be more stable by 13.31 meV/f.u., CHGNet predicts the distorted phase to be more stable by 18.44 meV/f.u. compared to the tetragonal one. This small energy inversion (≒18 meV) in the MLFF's potential energy surface explains why the simulation settles into a lower-symmetry structure with different a - and b -lattice parameters. However, since the energy difference is minute, the global trend of the pressure-induced phase transition is preserved. Despite this slight symmetry breaking, a significantly elongated c -axis ( c / a = 1.27 at 0 GPa) was clearly observed, reflecting the effect of the 6 s 2 lone pair in Bi 3+ . A very sharp first-order phase transition is observed at 3.7 GPa. At this transition pressure, the c -axis length dramatically contracts from about 4.5 to 4.1 Å, while the a - and b -axes lengths slightly expand. This anisotropic lattice deformation is qualitatively consistent with the experimentally observed structural phase transition to a GdFeO 3 -type structure involving tilting of the CoO 6 octahedra. 6 This selective compression along the c -axis is the primary cause of the colossal volume contraction. As a result of these lattice parameter changes, the volume discontinuously decreases from approximately 62 Å 3 to about 56 Å 3 at the transition pressure. The total volume change, Δ V / V 0 , reaches − 9.7%, In the initial structure optimization, a slight deviation from the experimentally observed tetragonal symmetry ( a = b ) was noted, resulting in a minor discrepancy between the a - and b -axis lattice parameters. This is likely due to the inherent limitations of the universal MLFF. While CHGNet is remarkably powerful, the potential energy surface it predicts for a complex system like BiCoO 3 , which is strongly influenced by Bi 3+ lone pairs and the Jahn-Teller effect of Co 3+ , may contain shallow local minima. This can lead to “systematic softening,” which will be discussed later. Compared to experimental results, the transition pressure obtained in our simulation (3.7 GPa) overestimates the experimental value (approx. 3.0 GPa), while the volume change (-9.7%) underestimates the experimental value (approx. -13%). 6 These quantitative discrepancies may be attributed to inherent errors in the universal MLFF and the simulation conditions, as discussed later. However, the crucial point is that the universal MLFF, without any additional training for this complex system, was able to qualitatively reproduce the essential physics of a pressure-induced first-order phase transition in a solid accompanied by volume greatly shrinking. This result strongly suggests that CHGNet can describe the structural energy landscape of BiCoO 3 -derivatives reasonably accurately, providing a solid foundation for investigating the more complex temperature-induced NTE phenomenon. Next, to verify the occurrence of a phase transition upon La substitution, we performed Ab Initio Molecular Dynamics (AIMD) simulations for Bi 0.5 La 0.5 CoO 3 . Table 1 shows the average volume and the c / a ratio, an indicator of tetragonality, after a 2-ps simulation at each temperature to ensure structural equilibration. In the temperature range from 300 to 450 K, the volume slightly increased from 71.21 to 71.34 Å 3 , exhibiting normal positive thermal expansion. During this interval, the c / a ratio remained large, increasing from 1.265 to 1.284, indicating the stability of the low-temperature tetragonal phase. However, a dramatic change is observed between 450 and 500 K. The volume sharply decreases by approximately 7.0%, from 71.34 to 66.34 Å 3 , confirming a clear NTE. Simultaneously, the c / a ratio approaches 1, decreasing from 1.284 to 1.069, which suggests a transformation to a more isotropic crystal structure. Table 1 Volume and c/a ratio of at finite temperatures from AIMD calculations. Temperature (K) Volume (Å 3 ) c/a ratio 300 71.21 1.265 400 71.23 1.269 450 71.34 1.284 500 66.34 1.069 This abrupt change in lattice parameters is due to a local atomic coordination change, as shown in Fig. 2 . In the low-temperature phase at 450 K, the Jahn-Teller distortion of Co 3+ persists, resulting in a polar, near-tetragonal structure where Co ions adopt a five-coordinated pyramidal configuration. In contrast, in the high-temperature phase at 500 K, the thermal energy overcomes the stereochemical activity of the 6s 2 lone pairs of Bi 3+ and the Hund coupling stabilizing the high-spin Co 3+ , and the structure transforms into a smaller-volume, higher-symmetry phase where Co ions are centered in six-coordinated octahedral environments. These AIMD calculations confirm at a first-principles level that Bi 0.5 La 0.5 CoO 3 undergoes a first-order phase transition accompanied by distinct local structural changes in the 450–500 K temperature range, resulting in NTE. Next, using the NTE phenomenon confirmed by AIMD as a benchmark, we extended our investigation to a systematic compositional search of Bi 1 − x La x CoO 3 using MLFF-MD. Figure 3 shows the temperature dependence of the lattice volume from MLFF-MD simulations with the La substitution systematically varied from 0 to 100%. This comprehensive simulation clearly revealed three important trends: Occurrence of NTE : For compositions with La substitution ranging from 25% to 75%, a first-order phase transition, i.e., NTE, characterized by volume discontinuously and sharply decreasing upon heating, is clearly reproduced. Systematic change in transition temperature : The most significant finding is that the NTE transition temperature systematically shifts to a lower side as the La substitution increases. This suggests that La substitution dilutes and weakens the effect of the Bi 3+ lone pair, enabling the transition to the smaller-volume high-temperature phase to be induced with lower thermal energy. Disappearance of NTE : When the La substitution was 12.5% or less, the volume no longer sharply contracts upon heating; only monotonic positive thermal expansion is seen. Furthermore, these compounds were found to decompose at elevated temperatures instead of undergoing the NTE transition, highlighting the critical role of La substitution in increased stability. To identify the specific symmetry and octahedral tilting pattern of the high-temperature phase, we performed a detailed structural analysis using the MLFF-MD trajectories. Direct observation of snapshots at high temperatures is often hindered by significant thermal fluctuations. Therefore, we calculated the time-averaged structure over a window of 100 fs (100 frames) to extract the equilibrium atomic positions. Figure 4 presents the projected structures of Bi 0.75 La 0.25 CoO 3 before and after the transition. The low-temperature phase (Fig. 4 (a)) clearly maintains the pyramidal coordination of Co with no observable tilting of the polyhedra ( a 0 a 0 c 0 ), preserving the tetragonal P 4 mm -like local symmetry. In contrast, the high-temperature phase (Fig. 4 (b)) exhibits a cooperative rotation of the CoO 6 octahedra. By analyzing the rotation patterns along each axis, we identified an antiphase rotation along the a - and c -axes ( a − , c − ) and an in-phase rotation along the b -axis (b + ). This a − b + c − (or simplified as a − b + a − ) tilting system is the distinct signature of the GdFeO 3 -type orthorhombic perovskite structure (space group Pnma ). This result provides robust computational evidence that the thermally induced phase transition in La-substituted BiCoO 3 involves a symmetry change from P 4 mm to Pnma , consistent with the high-pressure phase of the parent BiCoO 3 . Furthermore, a detailed analysis of the temperature and composition dependence of the c -lattice parameter (Fig. 5 ) reveals that, similar to the pressure-induced phase transition, the temperature-induced NTE is also primarily driven by the selective contraction of the c -axis. Notably, a clear relationship is observed between the La substitution amount and the lattice parameters of the low-temperature phase. In compositions exhibiting NTE, the a - and b -axis length did not significantly change after La substitution (see Figures S1 and S2 in the Supplementary Information), whereas the c -axis length markedly tended to contract even below the transition as La content increased. This suggests that substitution with La 3+ relaxes the structural anisotropy caused by the 6 s 2 lone pair of Bi 3+ . To clarify the electronic origin of the observed NTE, we focused on the magnetic moments of the Co ions, a secondary output from the simulations. The CHGNet model used in this study can predict not only the forces and energy but also the magnetic moment at each atomic site. Figure 6 shows the temperature and composition dependence of the average magnetic moment on the Co sites obtained from the MLFF-MD simulations. The results clearly show a strong correlation between the NTE and the spin state of the Co ions. In compositions where NTE is clearly observed (0.25 ≤ x ≤ 0.75), the magnetic moment decreases discontinuously from approximately 3.0 µ B to about 2.6–2.7 µ B , perfectly coinciding with the transition temperatures where the volume sharply contracts. This behavior corresponds to the spin-crossover phenomenon, where the high-spin state Co 3+ ions, stable at low temperatures, transition to a low-spin or intermediate state upon heating. The difference in the magnetic moment below and above the transition temperature decreases as the La substitution decreases, most probably because the thermally activated low- (high-) spin state is mixed to the ground state high- (low-) spin state. The magnetic moment of approximately 2.6–2.7 µ B in the post-transition phase is too large to be attributed to a low-spin (LS, S = 0) state. Instead, this value strongly suggests a spin-crossover to a phase dominated by the intermediate-spin (IS, S = 1) configuration of Co³⁺. This interpretation is consistent with detailed studies on the related perovskite LaCoO₃. Through GGA + U calculations, Pandey et al. concluded that the electronic structure of LaCoO₃ at room temperature is best described by a dominant IS state configuration. 19 They also highlighted the strong hybridization between Co 3 d and O 2 p orbitals within the perovskite octahedra, which differentiates the system from a purely ionic picture. Therefore, we interpret the observed magnetic moment not as a covalent LS state, but as clear evidence of a transition to a high-temperature phase with a dominant IS character. This change in spin state directly impacts the lattice structure. High-spin Co 3+ ion is Jahn-Teller active, which distorts the octahedron and stabilizes the pyramidal coordination. In the intermediate- or low-spin state, however, the Jahn-Teller effect vanishes, allowing for a more symmetric octahedral coordination. This change in local structure is the driving force behind the lattice volume contraction, i.e., NTE. On the other hand, in compositions where NTE was not observed ( x ≤ 0.125 and x = 1.0), the absence of such a clear spin-crossover transition supports this mechanism. From the analysis above, our MLFF-MD simulations strongly suggest that the NTE in the Bi 1 − x La x CoO 3 system is a phenomenon in which the lattice and spin degrees of freedom are strongly coupled, accompanied by a temperature-induced spin-crossover of Co 3+ ions. This demonstrates that a universal MLFF is not merely a tool for structural prediction but also an extremely effective method for capturing the essential physics, such as changes in electronic states, behind complex phase transitions. To investigate the reversibility of the phase transition and the hysteresis behavior, we performed cooling simulations for the x = 0.25 composition following the heating run. Figure 7 displays the evolution of total energy and cell volume during the thermal cycle. Interestingly, the simulation revealed a "one-way" irreversible phase transition. As shown in Fig. 7 (a) , once the system transforms into the small-volume phase (GdFeO 3 -type) around 600 K during heating, it does not revert to the original large-volume phase upon cooling. Instead, it maintains the compact, octahedral-coordinated structure down to low temperatures. The energy profile (Fig. 7 (b) ) provides a clear thermodynamic explanation for this behavior: the potential energy of the post-transition phase is lower than that of the initial pre-transition phase at the same temperature. This indicates that, within the potential energy surface described by CHGNet, the experimentally observed low-temperature PbTiO 3 -type phase (tetragonal) corresponds to a metastable local minimum, while the GdFeO 3 -type phase (orthorhombic) represents the global minimum for this composition. This irreversibility may be partly exaggerated by the finite simulation time and supercell size, which can hinder the nucleation of the reverse transition. As will be discussed in the next section, experimental measurements on La-substituted samples ( x = 0.2) exhibit a reversible NTE, while the fraction of the orthorhombic phase changes only 25 to 57% on heating and to 32% on cooling, not from 0 to 100%. To verify the "one-way" transition behavior suggested by the simulation, we compared these results with experimental synchrotron X-ray diffraction data for a compositionally similar sample, Bi 0.8 La 0.2 CoO 3 ( x = 0.20). Figure 8 shows the temperature dependence of the unit cell volumes and the phase fraction. The experimental system exists as a two-phase mixture of the large-volume tetragonal phase (PbTiO 3 -type) and the small-volume orthorhombic phase (GdFeO 3 -type) even at 100 K. Upon heating, the fraction of the orthorhombic phase (black line) increases, from 25 to 57%, leading to the NTE behavior observed in the average volume (green line). Crucially, during the subsequent cooling process, the phase fraction exhibits significant hysteresis. The system does not fully revert to its initial state; instead, the fraction of the orthorhombic phase (32%) remains higher than it was before heating. This experimental observation qualitatively aligns with the irreversibility predicted by our MLFF-MD simulation. The simulation indicated that the orthorhombic phase is energetically more stable than the tetragonal phase over a wide temperature range for La-substituted compositions. In the idealized simulation box with periodic boundary conditions, this led to a complete lock-in to the stable orthorhombic phase. In the real polycrystalline sample, the transformation is spatially constrained by domain walls and inter-granular strain, preventing a complete transformation and resulting in the observed hysteresis and mixed-phase state. Nevertheless, both simulation and experiment point to the same conclusion: La substitution fundamentally stabilizes the collapsed GdFeO 3 -type structure, making the NTE transition partially or wholly irreversible depending on the microstructural constraints. Figure 9 shows the dependence of the NTE transition temperature on the La substitution amount, as obtained from MLFF-MD. Experimental results shown in the Fig. 8 indicate that for a Bi 0.8 La 0.2 CoO 3 , the NTE transition is observed over a wide temperature range from approximately 480 K to 700 K. The closest composition in our simulation, x = 0.25, predicted a transition temperature of 725 K, which is considerably higher than the onset temperature observed experimentally. For the x = 0.5 composition, the MLFF prediction (approx. 350 K) was 100–150 K lower than the AIMD prediction (450–500 K). These results indicate that the universal MLFF has quantitative challenges in predicting the absolute values of transition temperatures. Nevertheless, it can be said to possess sufficient capability to correctly capture the chemical trends, such as whether NTE will occur with La substitution and how the transition temperature changes as the substitution amount increases. The true value of the MLFF-MD simulation used in this study extends beyond the prediction of thermodynamic properties. It functions as a "computational microscope" that enables the dynamic process of how a phase transition proceeds to be directly observed, with the movement of each individual atom tracked. Figure 10 shows representative snapshots of atomic configurations within ± 20 fs of the moment the NTE phase transition occurs in the heating simulation of Bi 0.5 La 0.5 CoO 3 . Before the transition : In the low-temperature phase, the entire system adopts a near-tetragonal structure elongated along the c -axis, with Co 3+ ions in a 5-coordinated pyramidal configuration. During the transition : Most notably, the snapshot during the transition clearly shows that the transformation to the small-volume high-temperature phase (with Co 3+ in 6-coordinated octahedral structures) does not occur uniformly throughout the system. It clearly captures the transition propagating from the La 3+ -rich regions (which lack lone pairs) to the Bi 3+ -rich regions. After the transition : Ultimately, the entire system completely transforms into the high-temperature phase with octahedral coordination. This visualization provides a new perspective for understanding complex phase transition mechanisms. Experimental techniques, typified by in-situ transmission electron microscopy (TEM) observation, are powerful methods able to capture structural changes under heating with atomic resolution. However, the time resolution of conventional in-situ TEM observations, which is typically limited by the video frame rate to the order of milliseconds, is insufficient to resolve the femtosecond-to-picosecond timescale on which the elementary processes of a phase transition, such as atomic vibrations, occur. Therefore, the information obtained from such techniques is essentially a time-averaged structure over these ultrafast events. MLFF-MD simulation can track the dynamic behavior of individual atoms in this ultrafast time domain, which is difficult to observe experimentally. The heterogeneous propagation behavior, nucleating from La-rich regions, as revealed here, is a transient phenomenon that cannot be captured from the averaged information. Furthermore, the ability to perform such systematic dynamic analyses for numerous compositions at a realistic computational cost is a significant advantage over experiments. Therefore, MLFF-MD not only enables high-speed screening but also plays a complementary role to experimental observations, demonstrating its potential to become a new analytical method for elucidating the fundamental elementary processes of complex solid-state reactions. Discussion The quantitative deviation in the transition temperature observed in this simulation is thought to stem from the inherent limitations of universal MLFFs. Several recent studies have pointed out that many universal MLFFs, including CHGNet, exhibit a tendency known as "systematic softening," where they underestimate the interatomic forces in high-energy regions. 20 The primary cause of this phenomenon is the training dataset for the universal model (in CHGNet's case, the Materials Project DFT relaxation trajectories 21 ) biased towards near-equilibrium structures around energy minima, with insufficient sampling of high-energy regions such as transition states and highly strained structures. The phase transition temperature is determined by the point where the free energies of the two phases become equal, and the contribution of vibrational entropy is significant to the free energy. If the potential energy surface (PES) is soft (i.e., has a small curvature), phonon frequencies are underestimated, which in turn affects the evaluation of vibrational entropy and can shift the predicted transition temperature. The underestimation of the transition temperatures seen in this study can be regarded as a direct manifestation of this systematic softening. To further validate the reliability of the MLFF-MD simulation and assess the impact of this systematic softening, we performed single-point DFT calculations on 400 snapshots extracted from the MD trajectory of Bi 0.5 La 0.5 CoO 3 during the heating process. Figure 11 compares the potential energies predicted by CHGNet with those calculated by DFT (VASP). Crucially, both methods exhibit a clear and simultaneous drop in total energy around snapshot 4,000, which corresponds to the experimentally and computationally observed phase transition temperature. This agreement confirms that the phase transition reproduced by CHGNet is not a mere artifact of the machine learning model but is physically driven by the stabilization of the high-temperature phase, consistent with first-principles energetics. However, as shown in the high-temperature region (around and after snapshot 4,000), the energy difference between CHGNet and DFT gradually increases. This trend provides that MLFF tends to under- or overestimate the energy of structures with large thermal displacements. Nevertheless, the fact that the distinct energy jump associated with the first-order transition is qualitatively captured supports the validity of using this universal MLFF to elucidate the transition mechanism, even if quantitative discrepancies in the transition temperature persist. On the other hand, a likely reason for CHGNet's ability to capture the qualitatively correct physics is that its architecture explicitly incorporates magnetic moment (spin) information. 12 Since the phase transition in the BiCoO 3 system is closely related to the change in the spin state of the Co ions, this feature likely enabled the model to somewhat learn the coupling between the electronic state and the structure, leading to its qualitative success. Furthermore, the computational efficiency of the MLFF-MD method employed in this study is remarkable. An AIMD calculation using a 40-atom 2 × 2 × 2 supercell required approximately 20 hours to run a 2-ps simulation at a single fixed temperature. In contrast, an MLFF-MD calculation using a 320-atom 4 × 4 × 4 supercell (8 times larger than AIMD) and simulating a continuous temperature ramp from 1 K to 1000 K (1000 points) took only 8 hours. A simple throughput comparison shows that the MLFF achieved a computational speed over 20,000 times faster than AIMD. This overwhelming speed-up makes comprehensive simulations (screenings) of material behavior over wide compositional and temperature ranges, which were previously impractical due to computational cost, a tangible reality. However, the current universal training data is insufficient to quantitatively and accurately describe the subtle energy changes associated with spin-crossover, which is considered another source of quantitative error. For future high-precision predictions, an approach that uses the universal model as a base and then fine-tunes it with a small amount of high-accuracy DFT data specific to the target system (such as the AIMD trajectories generated in this study) is suggested to be effective. 22 Conclusion In this study, we applied large-scale molecular dynamics (MD) simulations using the Crystal Hamilton Graph Network (CHGNet), a pretrained universal machine learning force field (MLFF), to elucidate the negative thermal expansion (NTE) phenomenon in the perovskite oxide Bi 1 − x La x CoO 3 . The main achievements are as follows: Methodology Validation and NTE Reproduction : The universal MLFF, CHGNet, qualitatively reproduced the giant volume contraction induced by pressure in BiCoO 3 . Furthermore, it successfully reproduced the experimentally observed pressure-induced phase transition in BiCoO 3 and NTE in the Bi 1 − x La x CoO 3 , and indicated the chemical trend, where the NTE transition temperature decreases as La substitution increases. Elucidation of Atomic-Level Dynamic Mechanism : By leveraging the overwhelming computational efficiency of MLFF-MD (over 20,000 times faster than ab initio molecular dynamics (AIMD)), we visualized the dynamic process of the NTE phase transition at the atomic level. As a result, we uncovered the nucleation and propagation behavior, where the phase transition initiates heterogeneously in La-rich regions and propagates to Bi-rich regions. These achievements concretely demonstrate that universal MLFF-MD simulation can serve as a powerful "computational microscope" for predicting the behavior of functional inorganic materials with complex phase transitions and for elucidating atomic-level mechanisms that are unobservable in experiments. This research not only embodies a new paradigm for materials science research driven by computational and data sciences but also provides important feedback for developing next-generation, higher-precision foundational models. Methods Machine Learning Force Field Molecular Dynamics Simulations In this study, the Crystal Hamilton Graph Network (CHGNet) (ver. 0.4.0) was employed as the machine learning force field (MLFF). 12 CHGNet is a universal MLFF pretrained on the extensive first-principles calculation database of the Materials Project. 21 It was used in this research without any additional training (fine-tuning). Molecular dynamics (MD) simulations were performed using the open-source physics simulation library Atomic Simulation Environment (ASE) (ver. 3.26.0) 23 with PyTorch (ver. 2.7.1) as its backend. All MLFF-MD calculations were executed on a single NVIDIA RTX 4090 GPU. A 4 × 4 × 4 supercell (320 atoms) created by expanding the unit cell of BiCoO 3 and Bi 1 − x La x CoO 3 solid solutions was used as the simulation cell. Prior to the MD simulations, the initial structure for each composition was optimized by using the BFGS algorithm, a quasi-Newton method, until the residual force on each atom was less than 0.001 eV/Å. For the subsequent evaluation of thermal expansion behavior, the NpT ensemble (isothermal-isobaric ensemble) was employed to maintain constant pressure and temperature. Temperature was controlled using a Nosé-Hoover thermostat (time constant = 25 fs) 24 , and pressure was controlled with a Nosé-Hoover barostat. A time step of 1.0 fs was used. In the pressure-induced phase transition simulations, the temperature was kept at 10 K while the pressure was incrementally increased from ambient pressure (1 bar) to 5 GPa in steps of 0.1 GPa. The structure was held for 2 ps at each pressure point to ensure sufficient relaxation. For the temperature-induced phase transition simulations, the pressure was maintained at ambient pressure (1 bar), and the temperature was ramped up and down from 0 K to 1000 K in 1 K intervals. At each temperature point, an MD calculation of 2000 steps (2 ps) was performed to equilibrate the system, and the lattice parameters and volume were recorded. In the heating simulations, the final structure from the previous temperature step was used as the initial structure for the next step to ensure continuity. First-Principles Molecular Dynamics Simulations To compare and validate the MLFF-MD results, ab initio molecular dynamics (AIMD) calculations were also performed using the VASP (Vienna Ab initio Simulation Package). 25 A 2 × 2 × 2 supercell of Bi 0.5 La 0.5 CoO 3 (40 atoms: Bi 4 La 4 Co 8 O 24 ) was used for the simulation cell. The Perdew-Burke-Ernzerhof (PBE) functional was used to describe the exchange-correlation energy of electrons. To account for the strong electron correlation among the 3d electrons of Co atoms, the GGA + U method (Dudarev's approach) 26 , 27 was applied, with an effective U parameter ( U eff ) of 3.32 eV set for the Co d-orbitals. The cutoff energy for the plane-wave basis set was 520 eV, and the first Brillouin zone was sampled using a 3×3×2 k-point mesh centered at the Gamma point. The simulations were conducted in the NpT ensemble, with temperature and pressure controlled by using Langevin dynamics. Calculations were performed at temperatures of 300, 400, 450, and 500 K for 2 ps (2000 steps) with a time step of 1.0 fs. Structural information for the equilibrium state was obtained from the trajectory of the final 1 ps. Declarations Ethics declarations All authors declare no financial or non-financial competing interests. Funding This work was partially supported by JSPS KAKENHI (grant no. JP24H00374), JST-CREST (JPMJCR22O1), the Kanagawa Institute of Industrial Science and Technology and Design and Engineering by Joint Inverse Innovation for Materials Architecture project of the Ministry of Education, Culture, Sports, Science and Technology. Author Contribution S.W. performed the MD simulations and analyzed the data. Y.S, K.T. and M.A performed the high-pressure synthesis and synchrotron X-ray diffraction experiments. H.D. and M.A. conceived the project and supervised the work. All authors contributed to writing the manuscript. Acknowledgement This work was partially supported by JSPS KAKENHI (grant no. JP24H00374), JST-CREST (JPMJCR22O1), the Kanagawa Institute of Industrial Science and Technology and Design and Engineering by Joint Inverse Innovation for Materials Architecture project of the Ministry of Education, Culture, Sports, Science and Technology. Data Availability All data generated or analyzed during this study are included in this published article and its supplementary information files Code Availability The underlying code for this study is not publicly available but may be made available to qualified researchers on reasonable requests from the corresponding author. References Azuma, M. Phase transition type giant negative thermal expansion materials. J. Ceram. Soc. Jpn . 133 , 450–454 (2025). Yamamoto, H. et al. Colossal negative thermal expansion in electron-doped PbVO₃. Angew Chem. Int. Ed. 57 , 8170–8173 (2018). Nishikubo, T. et al. Polar–Nonpolar Transition-Type Negative Thermal Expansion with 11.1% Volume Shrinkage by Design. Chem. Mater. 35 , 870–878 (2023). Azuma, M. et al. Giant negative thermal expansion in BiNiO₃. Nat. Commun. 2 , 347 (2011). Belik, A. A. et al. Neutron Powder Diffraction Study on the Crystal and Magnetic Structures of BiCoO 3 . Chem. Mater. 18 , 798–803 (2006). Oka, K. et al. Pressure-induced spin-state transition in BiCoO₃. J. Am. Chem. Soc. 132 , 9438–9443 (2010). Pan, Z. et al. Pronounced negative thermal expansion in lead-free BiCoO₃-based ferroelectrics triggered by the stabilized perovskite structure. Chem. Mater. 31 , 6187–6192 (2019). Becke, A. D. & Perspective Fifty years of density-functional theory in chemical physics. J. Chem. Phys. 140 , 18A301 (2014). Unke, O. T. et al. Machine learning force fields. Chem. Rev. 121 , 10142–10186 (2021). Chen, C. & Ong, S. P. A universal graph deep learning interatomic potential for the periodic table. Nat. Comput. Sci. 2 , 718–728 (2022). Batatia, I. et al. Higher order equivariant message passing neural networks for fast and accurate force fields. Adv. Neural Inf. Process. Syst. 35 , 11487–11501 (2022). Deng, B. et al. CHGNet as a pretrained universal neural network potential for charge-informed atomistic modelling. Nat. Mach. Intell. 5 , 1031–1041 (2023). Li, D., Yang, J., Chen, X., Yu, L. & Liu, S. To Use or Not to Use a Universal Force Field. Preprint at (2025). https://doi.org/10.48550/ARXIV.2503.08207 He, R., Wu, H., Lu, Y. & Zhong, Z. Origin of negative thermal expansion and pressure-induced amorphization in zirconium tungstate from a machine-learning potential. Phys. Rev. B . 106 , 174101 (2022). Shi, G., Zhong, Z., Su, Y., He, R. & Zhang, C. Investigating the atomic-level mechanisms of negative thermal expansion in Zn(CN)2 using machine learning-enabled molecular dynamics. Phys. Scr. 100 , 066002 (2025). Zhao, K. et al. Octahedral tilt distortion in negative thermal expansion in the fluorides CaZrF 6 and ScF 3 . Phys. Rev. B . 110 , 064322 (2024). Wu, J., Yang, J., Ma, L., Zhang, L. & Liu, S. Modular development of deep potential for complex solid solutions. Phys. Rev. B . 107 , 144102 (2023). Jain, A. et al. Formation enthalpies by mixing GGA and GGA calculations. Phys. Rev. B . 84 , 045115 (2011). Pandey, S. K. et al. Investigation of the spin state of Co in LaCoO 3 at room temperature: Ab initio calculations and high-resolution photoemission spectroscopy of single crystals. Phys Rev. B 77 , (2008). Deng, B. et al. Systematic softening in universal machine learning interatomic potentials. npj Comput. Mater 11 , (2025). Jain, A. et al. The materials project: A materials genome approach to accelerating materials innovation. APL Mater. 1 , 011002 (2013). Chen, P. Y., Shibata, K. & Mizoguchi, T. High precision machine learning force field development for BaTiO3 phase transitions, amorphous, and liquid structures. APL Mach. Learn. 3 , 036115 (2025). Larsen, A. H. et al. The atomic simulation environment—a Python library for working with atoms. J. Phys. : Condens. Matter . 29 , 273002 (2017). Nosé, S. A molecular dynamics method for simulations in the canonical ensemble. J. Chem. Phys. 81, 511–519 ; Hoover, W. G. Canonical dynamics: Equilibrium phase-space distributions. Phys. Rev. A 31, 1695–1697 (1985). (1984). Kresse, G. & Furthmüller, J. Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set. Phys. Rev. B . 54 , 11169–11186 (1996). Perdew, J. P., Burke, K. & Ernzerhof, M. Generalized gradient approximation made simple. Phys. Rev. Lett. 77 , 3865–3868 (1996). Dudarev, S. L., Botton, G. A., Savrasov, S. Y., Humphreys, C. J. & Sutton, A. P. Electron-energy-loss spectra and the structural stability of nickel oxide: An LSDA + U study. Phys. Rev. B . 57 , 1505–1509 (1998). Additional Declarations No competing interests reported. Supplementary Files SupplementaryInformations.docx BiCoO3tetra444op.cif BiCoO3tetra444op3.125.cif BiCoO3tetra444op6.25.cif BiCoO3tetra444op12.5.cif BiCoO3tetra444op37.5.cif BiCoO3tetra444op50.cif BiCoO3tetra444op100.cif BiCoO3tetra444op75.cif BiCoO3tetra444op62.5.cif Cite Share Download PDF Status: Published Journal Publication published 10 Apr, 2026 Read the published version in Scientific Reports → Version 1 posted Editorial decision: Revision requested 19 Jan, 2026 Reviews received at journal 16 Jan, 2026 Reviews received at journal 13 Jan, 2026 Reviewers agreed at journal 12 Jan, 2026 Reviewers agreed at journal 12 Jan, 2026 Reviewers agreed at journal 12 Jan, 2026 Reviews received at journal 04 Jan, 2026 Reviewers agreed at journal 24 Dec, 2025 Reviewers agreed at journal 24 Dec, 2025 Reviewers invited by journal 24 Dec, 2025 Editor invited by journal 24 Dec, 2025 Editor assigned by journal 22 Dec, 2025 Submission checks completed at journal 22 Dec, 2025 First submitted to journal 22 Dec, 2025 You are reading this latest preprint version Research Square lets you share your work early, gain feedback from the community, and start making changes to your manuscript prior to peer review in a journal. As a division of Research Square Company, we’re committed to making research communication faster, fairer, and more useful. We do this by developing innovative software and high quality services for the global research community. Our growing team is made up of researchers and industry professionals working together to solve the most critical problems facing scientific publishing. Also discoverable on Platform About Our Team In Review Editorial Policies Advisory Board Help Center Resources Author Services Accessibility API Access RSS feed Manage Cookie Preferences © Research Square 2026 | ISSN 2693-5015 (online) Privacy Policy Terms of Service Do Not Sell My Personal Information {"props":{"pageProps":{"initialData":{"identity":"rs-8425406","acceptedTermsAndConditions":true,"allowDirectSubmit":false,"archivedVersions":[],"articleType":"Article","associatedPublications":[],"authors":[{"id":565449839,"identity":"360b5705-3a79-4699-ae2e-cc4a6b5fdaad","order_by":0,"name":"Shogo Wakazaki","email":"","orcid":"","institution":"Institute of Science Tokyo","correspondingAuthor":false,"prefix":"","firstName":"Shogo","middleName":"","lastName":"Wakazaki","suffix":""},{"id":565449840,"identity":"4e976e36-f7c2-4ed2-97d7-8c4ddbe55844","order_by":1,"name":"Yuki Sakai","email":"","orcid":"","institution":"Institute of Science Tokyo","correspondingAuthor":false,"prefix":"","firstName":"Yuki","middleName":"","lastName":"Sakai","suffix":""},{"id":565449841,"identity":"81ad340f-c05f-4b24-8f00-5a1cf9ea1554","order_by":2,"name":"Kazuki Takahashi","email":"","orcid":"","institution":"Institute of Science Tokyo","correspondingAuthor":false,"prefix":"","firstName":"Kazuki","middleName":"","lastName":"Takahashi","suffix":""},{"id":565449842,"identity":"a39f7d17-5529-4914-9ed4-375e3f2d7167","order_by":3,"name":"Hena Das","email":"","orcid":"","institution":"Kanagawa Institute of Industrial Science and Technology","correspondingAuthor":false,"prefix":"","firstName":"Hena","middleName":"","lastName":"Das","suffix":""},{"id":565449843,"identity":"152f7684-26f7-4e06-aa55-a76c403b6325","order_by":4,"name":"Masaki Azuma","email":"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAZAAAAAyAQMAAABI0h/eAAAABlBMVEX///8AAABVwtN+AAAACXBIWXMAAA7EAAAOxAGVKw4bAAAA9klEQVRIiWNgGAWjYPACGwYGdiCVAMQGYAEeglrSGBiYGRgbEhKI13IYooUBrgUP0G0//kzqRs35aPlmBvYHD3/YyJszMD/8wCBzB6cWszMJadI5x27nbjgMdlia4c4GNmMJBp5nuLUcSDgmncMG1ALxy+EEgwMMZkC/HMat5fzDNumcf+dy5zeDtfwHamH/hl/LjWQ26dy2A7kNEIcdAGrhIWDLjWfM1rl9yUC/MDbOSEhLNtxwmKdYIgGfX86nP7yd880ud35784GPP2zs5A2Ot2/88LEHd4ghAWC0gAEzECf2HCBGCwr4QbqWUTAKRsEoGLYAAAxMVzkEqbyBAAAAAElFTkSuQmCC","orcid":"","institution":"Institute of Science Tokyo","correspondingAuthor":true,"prefix":"","firstName":"Masaki","middleName":"","lastName":"Azuma","suffix":""}],"badges":[],"createdAt":"2025-12-22 14:07:17","currentVersionCode":1,"declarations":"","doi":"10.21203/rs.3.rs-8425406/v1","doiUrl":"https://doi.org/10.21203/rs.3.rs-8425406/v1","draftVersion":[],"editorialEvents":[{"content":"https://doi.org/10.1038/s41598-026-47983-9","type":"published","date":"2026-04-10T15:57:49+00:00"}],"editorialNote":"","failedWorkflow":false,"files":[{"id":99314377,"identity":"64a9f216-e073-4482-a72b-90b933ed445b","added_by":"auto","created_at":"2025-12-31 16:21:18","extension":"docx","order_by":0,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":6825187,"visible":true,"origin":"","legend":"","description":"","filename":"Bi1xLaxCoO3CHGNetmanuscriptrevised.docx","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/7dce3bb0540ed3d900b14f96.docx"},{"id":99057864,"identity":"358b761e-d7c2-4981-9928-1556706f2792","added_by":"auto","created_at":"2025-12-26 19:42:27","extension":"json","order_by":1,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":6779,"visible":true,"origin":"","legend":"","description":"","filename":"10fcf57c72cf4c29b0ce8a977ac7982b.json","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/eb85cd32aa0de68ada61067c.json"},{"id":99314044,"identity":"6163f256-c467-4504-916a-67f65b9d9456","added_by":"auto","created_at":"2025-12-31 16:20:47","extension":"cif","order_by":2,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":29090,"visible":true,"origin":"","legend":"","description":"","filename":"BiCoO3tetra444op.cif","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/eb93cf4aa071cea0f841dc27.cif"},{"id":99314534,"identity":"a1af6209-61e1-4586-9ac8-a0e67515560e","added_by":"auto","created_at":"2025-12-31 16:21:46","extension":"cif","order_by":3,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":29090,"visible":true,"origin":"","legend":"","description":"","filename":"BiCoO3tetra444op100.cif","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/d9fd2d77925ed8a26918ce49.cif"},{"id":99313999,"identity":"70e78754-d974-4ef4-8486-ec453b2ef6cf","added_by":"auto","created_at":"2025-12-31 16:20:43","extension":"cif","order_by":4,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":29090,"visible":true,"origin":"","legend":"","description":"","filename":"BiCoO3tetra444op12.5.cif","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/cb59c36f57e8484e7e855196.cif"},{"id":99314478,"identity":"540862b3-01f4-48ee-8d05-5cf4f395b87e","added_by":"auto","created_at":"2025-12-31 16:21:34","extension":"cif","order_by":5,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":29090,"visible":true,"origin":"","legend":"","description":"","filename":"BiCoO3tetra444op3.125.cif","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/c617b35e8cb0d45c68f98c2a.cif"},{"id":99313667,"identity":"b1c39b9a-80c5-4056-8805-18e9c48e310c","added_by":"auto","created_at":"2025-12-31 16:20:24","extension":"cif","order_by":6,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":29090,"visible":true,"origin":"","legend":"","description":"","filename":"BiCoO3tetra444op37.5.cif","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/2eb2bd9419ea1457f52bd1ad.cif"},{"id":99057889,"identity":"5880c225-4d80-4c52-a510-41a8ff2217a5","added_by":"auto","created_at":"2025-12-26 19:42:27","extension":"cif","order_by":7,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":29090,"visible":true,"origin":"","legend":"","description":"","filename":"BiCoO3tetra444op50.cif","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/6ffe9eefe83dbd5b3ed264a6.cif"},{"id":99314212,"identity":"0f4c4a27-a662-40e5-a2f6-6c2e6628ec9a","added_by":"auto","created_at":"2025-12-31 16:21:00","extension":"cif","order_by":8,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":29090,"visible":true,"origin":"","legend":"","description":"","filename":"BiCoO3tetra444op6.25.cif","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/6a9ca574f82238219ddc9111.cif"},{"id":99314570,"identity":"b8d73973-09fd-435e-9f94-c40d297455a0","added_by":"auto","created_at":"2025-12-31 16:21:54","extension":"cif","order_by":9,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":29090,"visible":true,"origin":"","legend":"","description":"","filename":"BiCoO3tetra444op62.5.cif","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/ceaabd27e2e0149ff1295311.cif"},{"id":99057888,"identity":"bfbbde65-dc27-4f17-af90-f4dee38db876","added_by":"auto","created_at":"2025-12-26 19:42:27","extension":"cif","order_by":10,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":29090,"visible":true,"origin":"","legend":"","description":"","filename":"BiCoO3tetra444op75.cif","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/bf2a841da90cd54dc4b2da7c.cif"},{"id":99314148,"identity":"b682d04c-0b3a-4993-8f09-2866ed9322dc","added_by":"auto","created_at":"2025-12-31 16:20:54","extension":"docx","order_by":11,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":2999738,"visible":true,"origin":"","legend":"","description":"","filename":"SupplementaryInformations.docx","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/4ae3a94bb64ba3296cea803f.docx"},{"id":99314286,"identity":"260caee5-aef4-48ef-934b-cc613c37e838","added_by":"auto","created_at":"2025-12-31 16:21:07","extension":"xml","order_by":12,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":104026,"visible":true,"origin":"","legend":"","description":"","filename":"10fcf57c72cf4c29b0ce8a977ac7982b1enriched.xml","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/2a1035b49445dfee4939586e.xml"},{"id":99314552,"identity":"eb9a8a65-9e21-4b90-8c6b-2068e24cac43","added_by":"auto","created_at":"2025-12-31 16:21:51","extension":"png","order_by":13,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":156166,"visible":true,"origin":"","legend":"","description":"","filename":"floatimage1.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/0cc8833251815d089717f291.png"},{"id":99314497,"identity":"b81c214e-aa97-4b2f-b967-bf657a1e4f6a","added_by":"auto","created_at":"2025-12-31 16:21:36","extension":"png","order_by":14,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":728466,"visible":true,"origin":"","legend":"","description":"","filename":"floatimage10.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/e611f6c7e91b27c2f21bff49.png"},{"id":99314892,"identity":"fdfd2b15-ff41-4d5f-b40b-714022099e15","added_by":"auto","created_at":"2025-12-31 16:24:21","extension":"png","order_by":15,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":172987,"visible":true,"origin":"","legend":"","description":"","filename":"floatimage11.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/35a705aa9b48419cd49e68c3.png"},{"id":99314607,"identity":"615695d1-d027-4a8e-b822-2e91fa34fa63","added_by":"auto","created_at":"2025-12-31 16:22:00","extension":"png","order_by":16,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":593175,"visible":true,"origin":"","legend":"","description":"","filename":"floatimage2.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/ada54a8ae80ff00b1c246d84.png"},{"id":99314609,"identity":"d786abe3-6196-4fb3-b77a-fa7406a6fc05","added_by":"auto","created_at":"2025-12-31 16:22:00","extension":"png","order_by":17,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":583592,"visible":true,"origin":"","legend":"","description":"","filename":"floatimage3.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/c9a841af80423f82a2d753d9.png"},{"id":99057892,"identity":"79658375-cbcc-4cb3-8ec7-b721756cb7c8","added_by":"auto","created_at":"2025-12-26 19:42:28","extension":"png","order_by":18,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":1948123,"visible":true,"origin":"","legend":"","description":"","filename":"floatimage4.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/e8320c6653bf350e6c0f3b5c.png"},{"id":99057895,"identity":"7765506a-7cd2-4300-88da-f8c78ab9459f","added_by":"auto","created_at":"2025-12-26 19:42:28","extension":"png","order_by":19,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":496605,"visible":true,"origin":"","legend":"","description":"","filename":"floatimage5.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/1db1d4155eea6b4eb1c5a2f8.png"},{"id":99057890,"identity":"e4304fce-f6e8-4769-ada0-e91588332995","added_by":"auto","created_at":"2025-12-26 19:42:28","extension":"png","order_by":20,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":887563,"visible":true,"origin":"","legend":"","description":"","filename":"floatimage6.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/e7c2640b7a1718375b56f9a5.png"},{"id":99057893,"identity":"32740c2a-0a9d-4320-9968-de5d0598a16d","added_by":"auto","created_at":"2025-12-26 19:42:28","extension":"png","order_by":21,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":329899,"visible":true,"origin":"","legend":"","description":"","filename":"floatimage7.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/e995371859bb676e01a885c3.png"},{"id":99057913,"identity":"cd240b87-e767-4e97-bca8-92de6417ddfc","added_by":"auto","created_at":"2025-12-26 19:42:28","extension":"png","order_by":22,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":393393,"visible":true,"origin":"","legend":"","description":"","filename":"floatimage8.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/40c471b421c805dadb797112.png"},{"id":99057900,"identity":"7754b769-9082-4b16-8447-a9d33349d173","added_by":"auto","created_at":"2025-12-26 19:42:28","extension":"png","order_by":23,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":113257,"visible":true,"origin":"","legend":"","description":"","filename":"floatimage9.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/d027c77698d9379e918b51e6.png"},{"id":99057906,"identity":"dccea94b-4374-44df-8fa8-8e8a4e38d7c5","added_by":"auto","created_at":"2025-12-26 19:42:28","extension":"png","order_by":24,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":27830,"visible":true,"origin":"","legend":"","description":"","filename":"Onlinefloatimage1.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/fd746299b861e4ff758e9598.png"},{"id":99057911,"identity":"4a412568-d0d7-4c01-8cd3-8675dfc9b8b7","added_by":"auto","created_at":"2025-12-26 19:42:28","extension":"png","order_by":25,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":114291,"visible":true,"origin":"","legend":"","description":"","filename":"Onlinefloatimage10.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/8a43655e5b4789a99d649183.png"},{"id":99313859,"identity":"535d789c-4ce0-4cf9-a676-af3e2843f035","added_by":"auto","created_at":"2025-12-31 16:20:33","extension":"png","order_by":26,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":54878,"visible":true,"origin":"","legend":"","description":"","filename":"Onlinefloatimage11.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/687993c364471bafa3846872.png"},{"id":99057918,"identity":"48682bb0-59c3-4598-814e-d93c33f1b507","added_by":"auto","created_at":"2025-12-26 19:42:28","extension":"png","order_by":27,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":87696,"visible":true,"origin":"","legend":"","description":"","filename":"Onlinefloatimage2.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/3eee2044e02c7ef0d748d879.png"},{"id":99315108,"identity":"9ca8c9cd-d0ea-49cc-9e4f-db0f57a991b2","added_by":"auto","created_at":"2025-12-31 16:26:22","extension":"png","order_by":28,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":148144,"visible":true,"origin":"","legend":"","description":"","filename":"Onlinefloatimage3.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/e0af290a523d5fdc4fc4db47.png"},{"id":99057916,"identity":"cccd6347-f820-4f08-8dd5-15a4f1705a9f","added_by":"auto","created_at":"2025-12-26 19:42:28","extension":"png","order_by":29,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":307383,"visible":true,"origin":"","legend":"","description":"","filename":"Onlinefloatimage4.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/6980e2667c861ea69c50d515.png"},{"id":99057897,"identity":"50b85441-2285-441f-8d68-3398e21d5130","added_by":"auto","created_at":"2025-12-26 19:42:28","extension":"png","order_by":30,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":138057,"visible":true,"origin":"","legend":"","description":"","filename":"Onlinefloatimage5.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/4056dfc786d10f672bedae00.png"},{"id":99314679,"identity":"e8d8c11e-ab13-4849-a4b7-d42a6781b5ce","added_by":"auto","created_at":"2025-12-31 16:22:21","extension":"png","order_by":31,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":200914,"visible":true,"origin":"","legend":"","description":"","filename":"Onlinefloatimage6.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/fb58d3a0f0dbc76e13106ead.png"},{"id":99313796,"identity":"fc6e2855-df92-4a9a-b0d4-e798a33faf62","added_by":"auto","created_at":"2025-12-31 16:20:31","extension":"png","order_by":32,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":86678,"visible":true,"origin":"","legend":"","description":"","filename":"Onlinefloatimage7.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/900bbc08be872e69cb8bc85e.png"},{"id":99057891,"identity":"b031c54d-2743-44f7-bc5d-3ff5e96b0370","added_by":"auto","created_at":"2025-12-26 19:42:28","extension":"png","order_by":33,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":80750,"visible":true,"origin":"","legend":"","description":"","filename":"Onlinefloatimage8.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/2ceceaca7ac055a7e1f04739.png"},{"id":99314559,"identity":"0ac7315d-93a9-4601-ac05-73820f8c3ddf","added_by":"auto","created_at":"2025-12-31 16:21:52","extension":"png","order_by":34,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":43807,"visible":true,"origin":"","legend":"","description":"","filename":"Onlinefloatimage9.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/83bc72f73a65481cb82ec485.png"},{"id":99057919,"identity":"0e30cb18-ebf9-4e17-b26d-c81afc77b4cd","added_by":"auto","created_at":"2025-12-26 19:42:29","extension":"xml","order_by":35,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":103367,"visible":true,"origin":"","legend":"","description":"","filename":"10fcf57c72cf4c29b0ce8a977ac7982b1structuring.xml","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/c790505ea28502591f41d111.xml"},{"id":99314521,"identity":"88e1462d-386f-4183-9e71-bc4d73d5fec3","added_by":"auto","created_at":"2025-12-31 16:21:44","extension":"html","order_by":36,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":116545,"visible":true,"origin":"","legend":"","description":"","filename":"earlyproof.html","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/e357f121527d42c8703d024a.html"},{"id":99057861,"identity":"cf33f3a2-0d21-4fd9-a43b-78c5680f87c2","added_by":"auto","created_at":"2025-12-26 19:42:27","extension":"png","order_by":1,"title":"Figure 1","display":"","copyAsset":false,"role":"figure","size":156166,"visible":true,"origin":"","legend":"\u003cp\u003ePressure dependence of the lattice parameters (a) and volume (b) and of BiCoO\u003csub\u003e3\u003c/sub\u003e calculated by machine learning force field MD reproducing the experimentally observed polar to non-polar transition accompanied by a volume collapse.\u003c/p\u003e","description":"","filename":"floatimage1.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/d7e57cbe397f9ba49f5eaa76.png"},{"id":99314975,"identity":"b6959ca3-6dae-4238-aac3-0c8e43da7cee","added_by":"auto","created_at":"2025-12-31 16:25:16","extension":"png","order_by":2,"title":"Figure 2","display":"","copyAsset":false,"role":"figure","size":593175,"visible":true,"origin":"","legend":"\u003cp\u003eCrystal structures of Bi\u003csub\u003e0.5\u003c/sub\u003eLa\u003csub\u003e0.5\u003c/sub\u003eCoO\u003csub\u003e3\u003c/sub\u003e from AIMD simulations (at 450 and 500 K).\u003c/p\u003e","description":"","filename":"floatimage2.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/90003dc353f37c801e372a0e.png"},{"id":99057867,"identity":"042c7f78-c4cd-4887-9d20-84b544c86dc9","added_by":"auto","created_at":"2025-12-26 19:42:27","extension":"png","order_by":3,"title":"Figure 3","display":"","copyAsset":false,"role":"figure","size":583592,"visible":true,"origin":"","legend":"\u003cp\u003eTemperature dependence of the lattice volume of Bi\u003csub\u003e1-x\u003c/sub\u003eLa\u003csub\u003ex\u003c/sub\u003eCoO\u003csub\u003e3\u003c/sub\u003e calculated by MLFF-MD.\u003c/p\u003e","description":"","filename":"floatimage3.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/196d6e1448da83a2adb4d4ed.png"},{"id":99057869,"identity":"8eb6f209-153f-46a9-8eb7-1739a52b7d65","added_by":"auto","created_at":"2025-12-26 19:42:27","extension":"png","order_by":4,"title":"Figure 4","display":"","copyAsset":false,"role":"figure","size":1948123,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eVisualization of octahedral tilting and structural phase transition in Bi\u003c/strong\u003e\u003csub\u003e\u003cstrong\u003e0.75\u003c/strong\u003e\u003c/sub\u003e\u003cstrong\u003eLa\u003c/strong\u003e\u003csub\u003e\u003cstrong\u003e0.25\u003c/strong\u003e\u003c/sub\u003e\u003cstrong\u003eCoO\u003c/strong\u003e\u003csub\u003e\u003cstrong\u003e3\u003c/strong\u003e\u003c/sub\u003e\u003cstrong\u003e analyzed by time-averaged MLFF-MD trajectories.\u003c/strong\u003e (a) The low-temperature phase exhibits a polar structure with 5-coordinated pyramidal Co ions and no octahedral tilting (Glazer notation \u003cem\u003ea\u003c/em\u003e\u003csup\u003e0\u003c/sup\u003e\u003cem\u003ea\u003c/em\u003e\u003csup\u003e0\u003c/sup\u003e\u003cem\u003ec\u003c/em\u003e\u003csup\u003e0\u003c/sup\u003e), consistent with the tetragonal \u003cem\u003eP\u003c/em\u003e4\u003cem\u003emm\u003c/em\u003e symmetry. (b) The high-temperature phase transforms into a non-polar structure with 6-coordinated octahedra with tilting pattern of \u003cem\u003ea\u003c/em\u003e\u003csup\u003e-\u003c/sup\u003e\u003cem\u003eb\u003c/em\u003e\u003csup\u003e+\u003c/sup\u003e\u003cem\u003ec\u003c/em\u003e\u003csup\u003e-\u003c/sup\u003e (or \u003cem\u003ea\u003c/em\u003e\u003csup\u003e-\u003c/sup\u003e\u003cem\u003eb\u003c/em\u003e\u003csup\u003e+\u003c/sup\u003e\u003cem\u003ea\u003c/em\u003e\u003csup\u003e-\u003c/sup\u003e assuming \u003cem\u003ea\u003c/em\u003e ≒ \u003cem\u003ec\u003c/em\u003e), corresponds to the \u003cem\u003ePnma\u003c/em\u003e space group (GdFeO\u003csub\u003e3\u003c/sub\u003e-type)\u003c/p\u003e","description":"","filename":"floatimage4.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/0d878d8bdd78502724c3dc49.png"},{"id":99314031,"identity":"0abbf5b3-5de5-49d8-8aa7-3fa61627df35","added_by":"auto","created_at":"2025-12-31 16:20:46","extension":"png","order_by":5,"title":"Figure 5","display":"","copyAsset":false,"role":"figure","size":496605,"visible":true,"origin":"","legend":"\u003cp\u003eTemperature dependence of the \u003cem\u003ec\u003c/em\u003e-axis length of Bi\u003csub\u003e1-x\u003c/sub\u003eLa\u003csub\u003ex\u003c/sub\u003eCoO\u003csub\u003e3 \u003c/sub\u003ecalculated by MLFF-MD.\u003c/p\u003e","description":"","filename":"floatimage5.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/7e6713a5434b04bc19afbf7b.png"},{"id":99057886,"identity":"0b955b0d-12dd-49f1-90b6-be631229c872","added_by":"auto","created_at":"2025-12-26 19:42:27","extension":"png","order_by":6,"title":"Figure 6","display":"","copyAsset":false,"role":"figure","size":887563,"visible":true,"origin":"","legend":"\u003cp\u003eTemperature and composition dependence of the average magnetic moment on Co sites obtained from MLFF-MD simulations.\u003c/p\u003e","description":"","filename":"floatimage6.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/70b8d58a2338541ae651b150.png"},{"id":99313818,"identity":"0d81b4ad-55af-4ada-8542-643d543a50a3","added_by":"auto","created_at":"2025-12-31 16:20:31","extension":"png","order_by":7,"title":"Figure 7","display":"","copyAsset":false,"role":"figure","size":329899,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eIrreversibility of the phase transition in B\u003c/strong\u003e\u003csub\u003e\u003cstrong\u003e0.75\u003c/strong\u003e\u003c/sub\u003e\u003cstrong\u003eLa\u003c/strong\u003e\u003csub\u003e\u003cstrong\u003e0.25\u003c/strong\u003e\u003c/sub\u003e\u003cstrong\u003eCoO\u003c/strong\u003e\u003csub\u003e\u003cstrong\u003e3\u003c/strong\u003e\u003c/sub\u003e\u003cstrong\u003e revealed by heating and cooling MLFF-MD simulations. \u003c/strong\u003e(a) unit cell volume and (b) total potential energy for the \u003cem\u003ex\u003c/em\u003e = 0.25 composition during heating (orange solid line) and cooling (blue dashed line) cycles. While a sharp first-order transition from the large-volume tetragonal phase to the small-volume orthorhombic phase is observed around 600 K upon heating, the reverse transition does not occur upon cooling. The system remains in the small-volume phase down to low temperatures. The energy profile (b) indicates that the post-transition phase is energetically more stable (lower energy) than the initial phase at 0 K, suggesting the irreversible nature of this transition in the simulation.\u003c/p\u003e","description":"","filename":"floatimage7.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/de7dc2d54fec573487e31222.png"},{"id":99314458,"identity":"a2cdc103-9e82-44a3-be65-6f1877173255","added_by":"auto","created_at":"2025-12-31 16:21:31","extension":"png","order_by":8,"title":"Figure 8","display":"","copyAsset":false,"role":"figure","size":393393,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eExperimental thermal evolution of lattice parameters and phase fraction in Bi\u003c/strong\u003e\u003csub\u003e\u003cstrong\u003e0.8\u003c/strong\u003e\u003c/sub\u003e\u003cstrong\u003eLa\u003c/strong\u003e\u003csub\u003e\u003cstrong\u003e0.2\u003c/strong\u003e\u003c/sub\u003e\u003cstrong\u003eCoO\u003c/strong\u003e\u003csub\u003e\u003cstrong\u003e3\u003c/strong\u003e\u003c/sub\u003e\u003cstrong\u003e. \u003c/strong\u003eTemperature dependence of the unit cell volume and phase fraction for Bi\u003csub\u003e0.8\u003c/sub\u003eLa\u003csub\u003e0.2\u003c/sub\u003eCoO\u003csub\u003e3\u003c/sub\u003e obtained from synchrotron X-ray diffraction experiments. The blue squares and red triangles represent the volumes of the tetragonal (PbTiO\u003csub\u003e3\u003c/sub\u003e-type) and orthorhombic (GdFeO\u003csub\u003e3\u003c/sub\u003e-type) phases, respectively. The black circles indicate the phase fraction of the orthorhombic phase. \u0026nbsp;The green circles show the weighted average volume (\u003cem\u003eV\u003c/em\u003e\u003csub\u003ePC\u003c/sub\u003e). Arrows indicate the heating and cooling processes.\u003c/p\u003e","description":"","filename":"floatimage8.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/7aa20d87b0e752d381f34930.png"},{"id":99314619,"identity":"7e624470-fc84-40fa-9a1a-5a265d3dada2","added_by":"auto","created_at":"2025-12-31 16:22:02","extension":"png","order_by":9,"title":"Figure 9","display":"","copyAsset":false,"role":"figure","size":113257,"visible":true,"origin":"","legend":"\u003cp\u003eNTE transition temperatures calculated by MLFF-MD.\u003c/p\u003e","description":"","filename":"floatimage9.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/801c1fbdebb9d21a3d09eb87.png"},{"id":99057884,"identity":"f3877deb-ed75-43dc-9ce0-407b9e978aed","added_by":"auto","created_at":"2025-12-26 19:42:27","extension":"png","order_by":10,"title":"Figure 10","display":"","copyAsset":false,"role":"figure","size":728466,"visible":true,"origin":"","legend":"\u003cp\u003eVisualization of the phase transition dynamics in the MLFF-MD simulation on heating. Green, purple, blue, and red atoms show La, Bi, Co, and O, respectively. (a) Before, (b) during, and (c) after the phase transition.\u003c/p\u003e","description":"","filename":"floatimage10.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/298f0e8570adcdaa067d7ea4.png"},{"id":99057899,"identity":"4af9f8fe-96d1-4c4a-932b-67fe143e6fe4","added_by":"auto","created_at":"2025-12-26 19:42:28","extension":"png","order_by":11,"title":"Figure 11","display":"","copyAsset":false,"role":"figure","size":172987,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eValidation of the potential energy surface along the MD trajectory.\u003c/strong\u003e Comparison of the total energies predicted by CHGNet (blue circles) and calculated by DFT (red squares) for 400 snapshots extracted from the AIMD trajectory of Bi\u003csub\u003e0.5\u003c/sub\u003eLa\u003csub\u003e0.5\u003c/sub\u003eCoO\u003csub\u003e3\u003c/sub\u003e during the heating process. The DFT energies were corrected using the Pymatgen MaterialsProject2020Compatibility module..\u003c/p\u003e","description":"","filename":"floatimage11.png","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/2437aa0546a47d0119d7f7bb.png"},{"id":106809005,"identity":"43a95560-6867-4a6f-9f88-1bf42513582c","added_by":"auto","created_at":"2026-04-13 16:05:42","extension":"pdf","order_by":0,"title":"","display":"","copyAsset":false,"role":"manuscript-pdf","size":7229786,"visible":true,"origin":"","legend":"","description":"","filename":"manuscript.pdf","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/03508b91-8d1a-455f-95cd-3433ebbf3700.pdf"},{"id":99313747,"identity":"a96885c6-1a90-4445-8baa-19ae543896ea","added_by":"auto","created_at":"2025-12-31 16:20:28","extension":"docx","order_by":0,"title":"","display":"","copyAsset":false,"role":"supplement","size":2999738,"visible":true,"origin":"","legend":"","description":"","filename":"SupplementaryInformations.docx","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/1c53079e84ca1630413af0d0.docx"},{"id":99313835,"identity":"10eb326e-cda7-4363-8748-63a0993ee50e","added_by":"auto","created_at":"2025-12-31 16:20:32","extension":"cif","order_by":1,"title":"","display":"","copyAsset":false,"role":"supplement","size":29090,"visible":true,"origin":"","legend":"","description":"","filename":"BiCoO3tetra444op.cif","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/eae5a4c258fa4bba365412f6.cif"},{"id":99314477,"identity":"cd1ebaaf-0b88-4d59-af0e-a9ebbeb01615","added_by":"auto","created_at":"2025-12-31 16:21:33","extension":"cif","order_by":2,"title":"","display":"","copyAsset":false,"role":"supplement","size":29090,"visible":true,"origin":"","legend":"","description":"","filename":"BiCoO3tetra444op3.125.cif","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/83a2cf87a652a0307330bcf9.cif"},{"id":99315034,"identity":"e1c0ef18-3ebb-47e7-89e3-001314b7971f","added_by":"auto","created_at":"2025-12-31 16:26:07","extension":"cif","order_by":3,"title":"","display":"","copyAsset":false,"role":"supplement","size":29090,"visible":true,"origin":"","legend":"","description":"","filename":"BiCoO3tetra444op6.25.cif","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/f9572205aa0feafc039d504b.cif"},{"id":99314630,"identity":"253da64e-aab9-4f75-a000-8e20b3858dcc","added_by":"auto","created_at":"2025-12-31 16:22:04","extension":"cif","order_by":4,"title":"","display":"","copyAsset":false,"role":"supplement","size":29090,"visible":true,"origin":"","legend":"","description":"","filename":"BiCoO3tetra444op12.5.cif","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/b792f18de9db0d1e7fdc7386.cif"},{"id":99057875,"identity":"a9406039-ac60-407f-860a-a5f3e4b885eb","added_by":"auto","created_at":"2025-12-26 19:42:27","extension":"cif","order_by":5,"title":"","display":"","copyAsset":false,"role":"supplement","size":29090,"visible":true,"origin":"","legend":"","description":"","filename":"BiCoO3tetra444op37.5.cif","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/03b4df9cb38fd7e76128778b.cif"},{"id":99057885,"identity":"19fbd0a8-f316-460e-bb28-ba97ca629bb1","added_by":"auto","created_at":"2025-12-26 19:42:27","extension":"cif","order_by":6,"title":"","display":"","copyAsset":false,"role":"supplement","size":29090,"visible":true,"origin":"","legend":"","description":"","filename":"BiCoO3tetra444op50.cif","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/d00c8f8b1ce7104661861d02.cif"},{"id":99057877,"identity":"9932ef19-4d58-4c04-ae1a-b883c6b6bd7d","added_by":"auto","created_at":"2025-12-26 19:42:27","extension":"cif","order_by":7,"title":"","display":"","copyAsset":false,"role":"supplement","size":29090,"visible":true,"origin":"","legend":"","description":"","filename":"BiCoO3tetra444op100.cif","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/eaf1062b788bac4ab5180fa9.cif"},{"id":99057878,"identity":"a2eaebed-c27a-4f79-a1aa-251505c5db17","added_by":"auto","created_at":"2025-12-26 19:42:27","extension":"cif","order_by":8,"title":"","display":"","copyAsset":false,"role":"supplement","size":29090,"visible":true,"origin":"","legend":"","description":"","filename":"BiCoO3tetra444op75.cif","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/fd934ca3f5b5ad6ec76a2fb3.cif"},{"id":99057883,"identity":"5cda92de-9801-4214-b9e8-11346b424567","added_by":"auto","created_at":"2025-12-26 19:42:27","extension":"cif","order_by":9,"title":"","display":"","copyAsset":false,"role":"supplement","size":29090,"visible":true,"origin":"","legend":"","description":"","filename":"BiCoO3tetra444op62.5.cif","url":"https://assets-eu.researchsquare.com/files/rs-8425406/v1/f99ece0abf16d7e995914ece.cif"}],"financialInterests":"No competing interests reported.","formattedTitle":"Large-scale molecular dynamics simulations of negative thermal expansion in Bi1-x Lax CoO3 using a pretrained machine learning force field","fulltext":[{"header":"Introduction","content":"\u003cp\u003eWhile most materials expand upon heating, a few materials exhibit the peculiar property of contracting, known as \"Negative Thermal Expansion\" (NTE). This property is highly sought after for applications in advanced technologies requiring dimensional stability under extreme temperature changes, such as precision optical instruments and semiconductor manufacturing equipment, as it enables composite materials to be created that experience near-zero overall thermal expansion when combined with conventional positive thermal expansion materials.\u003c/p\u003e \u003cp\u003eThe mechanisms driving NTE are diverse, but a type driven by a first-order phase transition from a large-volume, low-temperature stable phase to a small-volume, high-temperature phase is of particular interest due to the colossal volume changes often involved.\u003csup\u003e\u003cspan citationid=\"CR1\" class=\"CitationRef\"\u003e1\u003c/span\u003e\u003c/sup\u003e Research on high-pressure synthesized perovskite oxides provides representative examples of this phase-transition-driven NTE. For instance, in the PbVO₃ derivatives with PbTiO\u003csub\u003e3\u003c/sub\u003e type structure with a giant \u003cem\u003ec\u003c/em\u003e/\u003cem\u003ea\u003c/em\u003e ratio, a colossal NTE exceeding 8% has been reported, induced by a phase transition from a polar ferroelectric structure to a nonpolar paraelectric structure.\u003csup\u003e\u003cspan citationid=\"CR2\" class=\"CitationRef\"\u003e2\u003c/span\u003e,\u003cspan citationid=\"CR3\" class=\"CitationRef\"\u003e3\u003c/span\u003e\u003c/sup\u003e In the BiNiO\u003csub\u003e3\u003c/sub\u003e derivatives, NTE occurs through a phase transition from a divalent Ni phase with a charge-disproportionated state on the Bi site (Bi\u003csup\u003e3+\u003c/sup\u003e/Bi\u003csup\u003e5+\u003c/sup\u003e, Ni\u003csup\u003e2+\u003c/sup\u003e) to a uniform trivalent Ni phase (Bi\u003csup\u003e3+\u003c/sup\u003e, Ni\u003csup\u003e3+\u003c/sup\u003e) involving charge transfer between Bi and Ni. By optimizing elemental substitution in this system, the phase transition temperature has been controlled to achieve giant NTE near room temperature.\u003csup\u003e\u003cspan citationid=\"CR4\" class=\"CitationRef\"\u003e4\u003c/span\u003e\u003c/sup\u003e Thus, phase-transition-driven NTE is a complex phenomenon where multiple degrees of freedom\u0026mdash;lattice, charge, and spin\u0026mdash;are intricately coupled, making the atomic-level elucidation of its mechanism a critical challenge from both fundamental and applied perspectives.\u003c/p\u003e \u003cp\u003eThe perovskite oxide BiCoO\u003csub\u003e3\u003c/sub\u003e, the focus of this study, belongs to this family of phase-transition-type NTE materials. Under ambient pressure, the stereochemically active 6\u003cem\u003es\u003c/em\u003e\u003csup\u003e2\u003c/sup\u003e lone pair of Bi\u003csup\u003e3+\u003c/sup\u003e and the Jahn-Teller effect of high-spin Co\u003csup\u003e3+\u003c/sup\u003e lead to a large polar displacement of Co within the CoO\u003csub\u003e6\u003c/sub\u003e octahedra as schematically illustrated as the inset of Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003e (b), resulting in a pyramid-like coordination and a ferroelectric tetragonal structure (space group \u003cem\u003eP\u003c/em\u003e4\u003cem\u003emm\u003c/em\u003e) with a significantly elongated \u003cem\u003ec\u003c/em\u003e-axis as \u003cem\u003ec\u003c/em\u003e/\u003cem\u003ea\u003c/em\u003e\u0026thinsp;=\u0026thinsp;1.27.\u003csup\u003e5\u003c/sup\u003e When a pressure of approximately 3 GPa is applied, it undergoes a spin-crossover phenomenon where Co\u003csup\u003e3+\u003c/sup\u003e ions transform from high-spin to low-spin states.\u003csup\u003e\u003cspan citationid=\"CR6\" class=\"CitationRef\"\u003e6\u003c/span\u003e\u003c/sup\u003e This spin-crossover is accompanied by a transition to a nonpolar paraelectric orthorhombic phase with tilted, isotropic CoO\u003csub\u003e6\u003c/sub\u003e octahedra, resulting in volume contracting colossally by about 13%. This small-volume phase, stable under high pressure, can be stabilized at ambient pressure through elemental substitution. Although irreversible, temperature-induced NTE has been reported in systems where Bi and Co are substituted for Ba and Ti, respectively.\u003csup\u003e\u003cspan citationid=\"CR7\" class=\"CitationRef\"\u003e7\u003c/span\u003e\u003c/sup\u003e Similarly, in the Bi\u003csub\u003e1\u0026thinsp;\u0026minus;\u0026thinsp;\u003cem\u003ex\u003c/em\u003e\u003c/sub\u003eLa\u003csub\u003e\u003cem\u003ex\u003c/em\u003e\u003c/sub\u003eCoO\u003csub\u003e3\u003c/sub\u003e system studied here, substituting La\u003csup\u003e3+\u003c/sup\u003e for Bi\u003csup\u003e3+\u003c/sup\u003e is expected to stabilize the high-pressure phase, allowing thermal energy to drive the transition and induce NTE.\u003c/p\u003e \u003cp\u003eFirst-principles calculations, particularly those based on Density Functional Theory (DFT), are a powerful tool for elucidating the mechanisms of such complex phase transitions at the atomic level.\u003csup\u003e\u003cspan citationid=\"CR8\" class=\"CitationRef\"\u003e8\u003c/span\u003e\u003c/sup\u003e However, the computational cost of DFT scales as \u003cem\u003eO\u003c/em\u003e(\u003cem\u003eN\u003c/em\u003e\u003csup\u003e3\u003c/sup\u003e) with the number of atoms \u003cem\u003eN\u003c/em\u003e, where \"\u003cem\u003eO\u003c/em\u003e\" denotes the order of complexity in Big \u003cem\u003eO\u003c/em\u003e notation, primarily due to the cubic scaling of core algorithms like matrix diagonalization required to be solved for the electronic wavefunctions. This makes simulations extremely difficult for phenomena occurring on time and length scales greater than a few 100 atoms and nanoseconds. To resolve this trade-off between accuracy and cost, Machine Learning Force Fields (MLFFs) have recently emerged as a new paradigm.\u003csup\u003e\u003cspan citationid=\"CR9\" class=\"CitationRef\"\u003e9\u003c/span\u003e\u003c/sup\u003e MLFFs are trained on vast datasets of energies and forces obtained from DFT calculations to rapidly reproduce the potential energy surface from atomic configurations. The development of \"universal force fields\" pretrained on large-scale databases is further accelerating this trend. M3GNet\u003csup\u003e\u003cspan citationid=\"CR10\" class=\"CitationRef\"\u003e10\u003c/span\u003e\u003c/sup\u003e, MACE\u003csup\u003e\u003cspan citationid=\"CR11\" class=\"CitationRef\"\u003e11\u003c/span\u003e\u003c/sup\u003e, and CHGNet\u003csup\u003e\u003cspan citationid=\"CR12\" class=\"CitationRef\"\u003e12\u003c/span\u003e\u003c/sup\u003e (Crystal Hamilton Graph Network, which is used in this study) are prominent examples aimed at applicability to a wide range of unknown materials across the periodic table without needing material-specific fine-tuning. The effectiveness of these universal force fields is already being demonstrated in phase transition simulations of other perovskite oxides, representing a major step toward high-throughput screening of computational materials.\u003csup\u003e\u003cspan citationid=\"CR13\" class=\"CitationRef\"\u003e13\u003c/span\u003e\u003c/sup\u003e\u003c/p\u003e \u003cp\u003eRecently, the application of MLFFs to NTE materials has gained attraction, offering atomic-scale insights that are difficult to access experimentally. For instance, MLFFs based on the Deep Potential scheme have successfully elucidated the role of polyhedral vibrations in the pressure-induced amorphization of ZrW\u003csub\u003e2\u003c/sub\u003eO\u003csub\u003e8\u003c/sub\u003e\u003csup\u003e14\u003c/sup\u003e, the \"lateral and tilt\" vibrational modes in Zn(CN)\u003csub\u003e2\u003c/sub\u003e\u003csup\u003e15\u003c/sup\u003e, and the impact of octahedral tilting on the NTE magnitude in fluorides such as ScF\u003csub\u003e3\u003c/sub\u003e and CaZrF\u003csub\u003e6\u003c/sub\u003e\u003csup\u003e16\u003c/sup\u003e. Furthermore, modular development strategies for MLFFs have extended these capabilities to solid solutions like Pb\u003csub\u003e\u003cem\u003ex\u003c/em\u003e\u003c/sub\u003eSr\u003csub\u003e1\u0026minus;\u003cem\u003ex\u003c/em\u003e\u003c/sub\u003eTiO\u003csub\u003e3\u003c/sub\u003e, reproducing complex topological textures\u003csup\u003e\u003cspan citationid=\"CR17\" class=\"CitationRef\"\u003e17\u003c/span\u003e\u003c/sup\u003e. However, these pioneering studies typically rely on constructing system-specific models, which necessitate the generation of extensive DFT datasets and dedicated training processes for each target material system. This \"train-from-scratch\" paradigm remains a bottleneck for high-throughput exploration, particularly when searching for new compositions or investigating substitution effects across a wide range of elements (e.g., exploring various dopants in BiCoO\u003csub\u003e3\u003c/sub\u003e). This is where \"universal\" MLFFs pretrained on massive databases, such as CHGNet, offer a distinct advantage by eliminating the need for ad hoc DFT calculations. Yet, the applicability of these universal models to transition metal oxides exhibiting colossal volume changes (\u0026gt;\u0026thinsp;5%) and discontinuous first-order phase transitions, phenomena deeply coupled with electronic and spin states, has not been sufficiently verified. Unlike the continuous NTE mechanisms seen in framework structures like ZrW\u003csub\u003e2\u003c/sub\u003eO\u003csub\u003e8\u003c/sub\u003e, the phase-transition-driven NTE in Bi-based perovskites involves drastic structural rearrangements and electronic structure changes that challenge the transferability of pretrained models. Therefore, demonstrating that a universal MLFF can reproduce these complex behaviors in a solid solution system without additional training represents a critical step toward the democratization of computational materials design.\u003c/p\u003e \u003cp\u003eIn this study, we focus on the Bi\u003csub\u003e1\u0026thinsp;\u0026minus;\u0026thinsp;\u003cem\u003ex\u003c/em\u003e\u003c/sub\u003eLa\u003csub\u003e\u003cem\u003ex\u003c/em\u003e\u003c/sub\u003eCoO\u003csub\u003e3\u003c/sub\u003e system and perform Molecular Dynamics (MD) simulations using both CHGNet and first-principles calculations to reproduce the finite-temperature phase transition and elucidate its underlying principles at the atomic level.\u003c/p\u003e"},{"header":"Results","content":"\u003cp\u003eFigure \u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003e shows the pressure dependence of the lattice parameters and the resulting unit cell volume of BiCoO₃, as calculated by the universal MLFF, CHGNet. In order to reproduce the structural transition from \u003cem\u003eP\u003c/em\u003e4\u003cem\u003emm\u003c/em\u003e to \u003cem\u003ePnma\u003c/em\u003e phases and the random distribution of substituted elements, the symmetry is lowered to P1 as indicated in the cif file found in the supplementary information. To clarify the origin of the pressure-induced phase transition, we first analyzed the changes in lattice parameters in detail. In the low-pressure phase, the structure, which is experimentally known to be tetragonal (\u003cem\u003eP\u003c/em\u003e4\u003cem\u003emm\u003c/em\u003e, \u003cem\u003ea\u003c/em\u003e\u0026thinsp;=\u0026thinsp;\u003cem\u003eb\u003c/em\u003e), showed a slight deviation from this symmetry in our simulation, resulting in distinct \u003cem\u003ea\u003c/em\u003e- and \u003cem\u003eb\u003c/em\u003e-axis lattice parameters (see Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003ea).\u003c/p\u003e \u003cp\u003eTo clarify the origin of the symmetry breaking (\u003cem\u003ea\u003c/em\u003e\u0026thinsp;\u0026ne;\u0026thinsp;\u003cem\u003eb\u003c/em\u003e) observed in the MLFF simulation at 0 GPa, we performed static energy calculations for pure BiCoO\u003csub\u003e3\u003c/sub\u003e using both CHGNet and DFT (VASP). We compared two structures: the experimentally reported tetragonal phase (\u003cem\u003ea\u003c/em\u003e\u0026thinsp;=\u0026thinsp;\u003cem\u003eb\u003c/em\u003e) and the slightly distorted monoclinic-like phase (\u003cem\u003ea\u003c/em\u003e\u0026thinsp;\u0026ne;\u0026thinsp;\u003cem\u003eb\u003c/em\u003e) observed in our MLFF simulation. The DFT energies were corrected using the MaterialsProject2020 Compatibility scheme implemented in pymatgen to ensure consistency with the training data.\u003csup\u003e\u003cspan citationid=\"CR18\" class=\"CitationRef\"\u003e18\u003c/span\u003e\u003c/sup\u003e\u003c/p\u003e \u003cp\u003eThe calculations revealed a subtle discrepancy in the ground state stability. While DFT (VASP) correctly predicts the tetragonal phase to be more stable by 13.31 meV/f.u., CHGNet predicts the distorted phase to be more stable by 18.44 meV/f.u. compared to the tetragonal one. This small energy inversion (≒18 meV) in the MLFF's potential energy surface explains why the simulation settles into a lower-symmetry structure with different \u003cem\u003ea\u003c/em\u003e- and \u003cem\u003eb\u003c/em\u003e-lattice parameters. However, since the energy difference is minute, the global trend of the pressure-induced phase transition is preserved. Despite this slight symmetry breaking, a significantly elongated \u003cem\u003ec\u003c/em\u003e-axis (\u003cem\u003ec\u003c/em\u003e/\u003cem\u003ea\u003c/em\u003e\u0026thinsp;=\u0026thinsp;1.27 at 0 GPa) was clearly observed, reflecting the effect of the 6\u003cem\u003es\u003c/em\u003e\u003csup\u003e2\u003c/sup\u003e lone pair in Bi\u003csup\u003e3+\u003c/sup\u003e. A very sharp first-order phase transition is observed at 3.7 GPa. At this transition pressure, the \u003cem\u003ec\u003c/em\u003e-axis length dramatically contracts from about 4.5 to 4.1 \u0026Aring;, while the \u003cem\u003ea\u003c/em\u003e- and \u003cem\u003eb\u003c/em\u003e-axes lengths slightly expand. This anisotropic lattice deformation is qualitatively consistent with the experimentally observed structural phase transition to a GdFeO\u003csub\u003e3\u003c/sub\u003e-type structure involving tilting of the CoO\u003csub\u003e6\u003c/sub\u003e octahedra.\u003csup\u003e\u003cspan citationid=\"CR6\" class=\"CitationRef\"\u003e6\u003c/span\u003e\u003c/sup\u003e\u003c/p\u003e \u003cp\u003eThis selective compression along the \u003cem\u003ec\u003c/em\u003e-axis is the primary cause of the colossal volume contraction. As a result of these lattice parameter changes, the volume discontinuously decreases from approximately 62 \u0026Aring;\u003csup\u003e3\u003c/sup\u003e to about 56 \u0026Aring;\u003csup\u003e3\u003c/sup\u003e at the transition pressure. The total volume change, Δ\u003cem\u003eV\u003c/em\u003e/\u003cem\u003eV\u003c/em\u003e\u003csub\u003e0\u003c/sub\u003e, reaches \u0026minus;\u0026thinsp;9.7%,\u003c/p\u003e \u003cp\u003eIn the initial structure optimization, a slight deviation from the experimentally observed tetragonal symmetry (\u003cem\u003ea\u003c/em\u003e\u0026thinsp;=\u0026thinsp;\u003cem\u003eb\u003c/em\u003e) was noted, resulting in a minor discrepancy between the \u003cem\u003ea\u003c/em\u003e- and \u003cem\u003eb\u003c/em\u003e-axis lattice parameters. This is likely due to the inherent limitations of the universal MLFF. While CHGNet is remarkably powerful, the potential energy surface it predicts for a complex system like BiCoO\u003csub\u003e3\u003c/sub\u003e, which is strongly influenced by Bi\u003csup\u003e3+\u003c/sup\u003e lone pairs and the Jahn-Teller effect of Co\u003csup\u003e3+\u003c/sup\u003e, may contain shallow local minima. This can lead to \u0026ldquo;systematic softening,\u0026rdquo; which will be discussed later.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eCompared to experimental results, the transition pressure obtained in our simulation (3.7 GPa) overestimates the experimental value (approx. 3.0 GPa), while the volume change (-9.7%) underestimates the experimental value (approx. -13%).\u003csup\u003e6\u003c/sup\u003e These quantitative discrepancies may be attributed to inherent errors in the universal MLFF and the simulation conditions, as discussed later. However, the crucial point is that the universal MLFF, without any additional training for this complex system, was able to qualitatively reproduce the essential physics of a pressure-induced first-order phase transition in a solid accompanied by volume greatly shrinking. This result strongly suggests that CHGNet can describe the structural energy landscape of BiCoO\u003csub\u003e3\u003c/sub\u003e-derivatives reasonably accurately, providing a solid foundation for investigating the more complex temperature-induced NTE phenomenon.\u003c/p\u003e \u003cp\u003eNext, to verify the occurrence of a phase transition upon La substitution, we performed Ab Initio Molecular Dynamics (AIMD) simulations for Bi\u003csub\u003e0.5\u003c/sub\u003eLa\u003csub\u003e0.5\u003c/sub\u003eCoO\u003csub\u003e3\u003c/sub\u003e.\u003c/p\u003e \u003cp\u003eTable\u0026nbsp;\u003cspan refid=\"Tab1\" class=\"InternalRef\"\u003e1\u003c/span\u003e shows the average volume and the \u003cem\u003ec\u003c/em\u003e/\u003cem\u003ea\u003c/em\u003e ratio, an indicator of tetragonality, after a 2-ps simulation at each temperature to ensure structural equilibration. In the temperature range from 300 to 450 K, the volume slightly increased from 71.21 to 71.34 \u0026Aring;\u003csup\u003e3\u003c/sup\u003e, exhibiting normal positive thermal expansion. During this interval, the \u003cem\u003ec\u003c/em\u003e/\u003cem\u003ea\u003c/em\u003e ratio remained large, increasing from 1.265 to 1.284, indicating the stability of the low-temperature tetragonal phase. However, a dramatic change is observed between 450 and 500 K. The volume sharply decreases by approximately 7.0%, from 71.34 to 66.34 \u0026Aring;\u003csup\u003e3\u003c/sup\u003e, confirming a clear NTE. Simultaneously, the \u003cem\u003ec\u003c/em\u003e/\u003cem\u003ea\u003c/em\u003e ratio approaches 1, decreasing from 1.284 to 1.069, which suggests a transformation to a more isotropic crystal structure.\u003c/p\u003e \u003cp\u003e \u003cdiv class=\"gridtable\"\u003e\u003ctable float=\"Yes\" id=\"Tab1\" border=\"1\"\u003e \u003ccaption language=\"En\"\u003e \u003cdiv class=\"CaptionNumber\"\u003eTable 1\u003c/div\u003e \u003cdiv class=\"CaptionContent\"\u003e \u003cp\u003eVolume and c/a ratio of at finite temperatures from AIMD calculations.\u003c/p\u003e \u003c/div\u003e \u003c/caption\u003e \u003ccolgroup cols=\"3\"\u003e \u003cdiv align=\"left\" class=\"colspec\" colname=\"c1\" colnum=\"1\"\u003e\u003c/div\u003e \u003cdiv align=\"char\" char=\".\" class=\"colspec\" colname=\"c2\" colnum=\"2\"\u003e\u003c/div\u003e \u003cdiv align=\"char\" char=\".\" class=\"colspec\" colname=\"c3\" colnum=\"3\"\u003e\u003c/div\u003e \u003cthead\u003e \u003ctr\u003e \u003cth align=\"left\" colname=\"c1\"\u003e \u003cp\u003eTemperature (K)\u003c/p\u003e \u003c/th\u003e \u003cth align=\"left\" colname=\"c2\"\u003e \u003cp\u003eVolume (\u0026Aring;\u003csup\u003e3\u003c/sup\u003e)\u003c/p\u003e \u003c/th\u003e \u003cth align=\"left\" colname=\"c3\"\u003e \u003cp\u003ec/a ratio\u003c/p\u003e \u003c/th\u003e \u003c/tr\u003e \u003c/thead\u003e \u003ctbody\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cp\u003e300\u003c/p\u003e \u003c/td\u003e \u003ctd align=\"char\" char=\".\" colname=\"c2\"\u003e \u003cp\u003e71.21\u003c/p\u003e \u003c/td\u003e \u003ctd align=\"char\" char=\".\" colname=\"c3\"\u003e \u003cp\u003e1.265\u003c/p\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cp\u003e400\u003c/p\u003e \u003c/td\u003e \u003ctd align=\"char\" char=\".\" colname=\"c2\"\u003e \u003cp\u003e71.23\u003c/p\u003e \u003c/td\u003e \u003ctd align=\"char\" char=\".\" colname=\"c3\"\u003e \u003cp\u003e1.269\u003c/p\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cp\u003e450\u003c/p\u003e \u003c/td\u003e \u003ctd align=\"char\" char=\".\" colname=\"c2\"\u003e \u003cp\u003e71.34\u003c/p\u003e \u003c/td\u003e \u003ctd align=\"char\" char=\".\" colname=\"c3\"\u003e \u003cp\u003e1.284\u003c/p\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cp\u003e500\u003c/p\u003e \u003c/td\u003e \u003ctd align=\"char\" char=\".\" colname=\"c2\"\u003e \u003cp\u003e66.34\u003c/p\u003e \u003c/td\u003e \u003ctd align=\"char\" char=\".\" colname=\"c3\"\u003e \u003cp\u003e1.069\u003c/p\u003e \u003c/td\u003e \u003c/tr\u003e \u003c/tbody\u003e \u003c/colgroup\u003e \u003c/table\u003e\u003c/div\u003e \u003c/p\u003e \u003cp\u003eThis abrupt change in lattice parameters is due to a local atomic coordination change, as shown in Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003e. In the low-temperature phase at 450 K, the Jahn-Teller distortion of Co\u003csup\u003e3+\u003c/sup\u003epersists, resulting in a polar, near-tetragonal structure where Co ions adopt a five-coordinated pyramidal configuration. In contrast, in the high-temperature phase at 500 K, the thermal energy overcomes the stereochemical activity of the 6s\u003csup\u003e2\u003c/sup\u003e lone pairs of Bi\u003csup\u003e3+\u003c/sup\u003e and the Hund coupling stabilizing the high-spin Co\u003csup\u003e3+\u003c/sup\u003e, and the structure transforms into a smaller-volume, higher-symmetry phase where Co ions are centered in six-coordinated octahedral environments.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eThese AIMD calculations confirm at a first-principles level that Bi\u003csub\u003e0.5\u003c/sub\u003eLa\u003csub\u003e0.5\u003c/sub\u003eCoO\u003csub\u003e3\u003c/sub\u003e undergoes a first-order phase transition accompanied by distinct local structural changes in the 450\u0026ndash;500 K temperature range, resulting in NTE.\u003c/p\u003e \u003cp\u003eNext, using the NTE phenomenon confirmed by AIMD as a benchmark, we extended our investigation to a systematic compositional search of Bi\u003csub\u003e1\u0026thinsp;\u0026minus;\u0026thinsp;x\u003c/sub\u003eLa\u003csub\u003ex\u003c/sub\u003eCoO\u003csub\u003e3\u003c/sub\u003e using MLFF-MD.\u003c/p\u003e \u003cp\u003eFigure \u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003e shows the temperature dependence of the lattice volume from MLFF-MD simulations with the La substitution systematically varied from 0 to 100%. This comprehensive simulation clearly revealed three important trends:\u003c/p\u003e \u003cp\u003e \u003cul\u003e \u003cli\u003e \u003cp\u003e \u003cb\u003eOccurrence of NTE\u003c/b\u003e: For compositions with La substitution ranging from 25% to 75%, a first-order phase transition, i.e., NTE, characterized by volume discontinuously and sharply decreasing upon heating, is clearly reproduced.\u003c/p\u003e \u003c/li\u003e \u003cli\u003e \u003cp\u003e \u003cb\u003eSystematic change in transition temperature\u003c/b\u003e: The most significant finding is that the NTE transition temperature systematically shifts to a lower side as the La substitution increases. This suggests that La substitution dilutes and weakens the effect of the Bi\u003csup\u003e3+\u003c/sup\u003e lone pair, enabling the transition to the smaller-volume high-temperature phase to be induced with lower thermal energy.\u003c/p\u003e \u003c/li\u003e \u003cli\u003e \u003cp\u003e \u003cb\u003eDisappearance of NTE\u003c/b\u003e: When the La substitution was 12.5% or less, the volume no longer sharply contracts upon heating; only monotonic positive thermal expansion is seen. Furthermore, these compounds were found to decompose at elevated temperatures instead of undergoing the NTE transition, highlighting the critical role of La substitution in increased stability.\u003c/p\u003e \u003c/li\u003e \u003c/ul\u003e \u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eTo identify the specific symmetry and octahedral tilting pattern of the high-temperature phase, we performed a detailed structural analysis using the MLFF-MD trajectories. Direct observation of snapshots at high temperatures is often hindered by significant thermal fluctuations. Therefore, we calculated the time-averaged structure over a window of 100 fs (100 frames) to extract the equilibrium atomic positions.\u003c/p\u003e \u003cp\u003eFigure \u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003e presents the projected structures of Bi\u003csub\u003e0.75\u003c/sub\u003eLa\u003csub\u003e0.25\u003c/sub\u003eCoO\u003csub\u003e3\u003c/sub\u003e before and after the transition. The low-temperature phase (Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003e(a)) clearly maintains the pyramidal coordination of Co with no observable tilting of the polyhedra (\u003cem\u003ea\u003c/em\u003e\u003csup\u003e0\u003c/sup\u003e\u003cem\u003ea\u003c/em\u003e\u003csup\u003e0\u003c/sup\u003e\u003cem\u003ec\u003c/em\u003e\u003csup\u003e0\u003c/sup\u003e), preserving the tetragonal \u003cem\u003eP\u003c/em\u003e4\u003cem\u003emm\u003c/em\u003e-like local symmetry. In contrast, the high-temperature phase (Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003e(b)) exhibits a cooperative rotation of the CoO\u003csub\u003e6\u003c/sub\u003e octahedra. By analyzing the rotation patterns along each axis, we identified an antiphase rotation along the \u003cem\u003ea\u003c/em\u003e- and \u003cem\u003ec\u003c/em\u003e-axes (\u003cem\u003ea\u003c/em\u003e\u003csup\u003e\u0026minus;\u003c/sup\u003e, \u003cem\u003ec\u003c/em\u003e\u003csup\u003e\u0026minus;\u003c/sup\u003e) and an in-phase rotation along the \u003cem\u003eb\u003c/em\u003e-axis (b\u003csup\u003e+\u003c/sup\u003e). This \u003cem\u003ea\u003c/em\u003e\u003csup\u003e\u0026minus;\u003c/sup\u003e\u003cem\u003eb\u003c/em\u003e\u003csup\u003e+\u003c/sup\u003e\u003cem\u003ec\u003c/em\u003e\u003csup\u003e\u0026minus;\u003c/sup\u003e (or simplified as \u003cem\u003ea\u003c/em\u003e\u003csup\u003e\u0026minus;\u003c/sup\u003e\u003cem\u003eb\u003c/em\u003e\u003csup\u003e+\u003c/sup\u003e\u003cem\u003ea\u003c/em\u003e\u003csup\u003e\u0026minus;\u003c/sup\u003e) tilting system is the distinct signature of the GdFeO\u003csub\u003e3\u003c/sub\u003e-type orthorhombic perovskite structure (space group \u003cem\u003ePnma\u003c/em\u003e). This result provides robust computational evidence that the thermally induced phase transition in La-substituted BiCoO\u003csub\u003e3\u003c/sub\u003e involves a symmetry change from \u003cem\u003eP\u003c/em\u003e4\u003cem\u003emm\u003c/em\u003e to \u003cem\u003ePnma\u003c/em\u003e, consistent with the high-pressure phase of the parent BiCoO\u003csub\u003e3\u003c/sub\u003e.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eFurthermore, a detailed analysis of the temperature and composition dependence of the \u003cem\u003ec\u003c/em\u003e-lattice parameter (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003e) reveals that, similar to the pressure-induced phase transition, the temperature-induced NTE is also primarily driven by the selective contraction of the \u003cem\u003ec\u003c/em\u003e-axis. Notably, a clear relationship is observed between the La substitution amount and the lattice parameters of the low-temperature phase. In compositions exhibiting NTE, the \u003cem\u003ea\u003c/em\u003e- and \u003cem\u003eb\u003c/em\u003e-axis length did not significantly change after La substitution (see Figures \u003cspan refid=\"MOESM1\" class=\"InternalRef\"\u003eS1\u003c/span\u003e and S2 in the Supplementary Information), whereas the \u003cem\u003ec\u003c/em\u003e-axis length markedly tended to contract even below the transition as La content increased. This suggests that substitution with La\u003csup\u003e3+\u003c/sup\u003e relaxes the structural anisotropy caused by the 6\u003cem\u003es\u003c/em\u003e\u003csup\u003e2\u003c/sup\u003e lone pair of Bi\u003csup\u003e3+\u003c/sup\u003e.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eTo clarify the electronic origin of the observed NTE, we focused on the magnetic moments of the Co ions, a secondary output from the simulations. The CHGNet model used in this study can predict not only the forces and energy but also the magnetic moment at each atomic site. Figure\u0026nbsp;\u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003e shows the temperature and composition dependence of the average magnetic moment on the Co sites obtained from the MLFF-MD simulations.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eThe results clearly show a strong correlation between the NTE and the spin state of the Co ions. In compositions where NTE is clearly observed (0.25\u0026thinsp;\u0026le;\u0026thinsp;\u003cem\u003ex\u003c/em\u003e\u0026thinsp;\u0026le;\u0026thinsp;0.75), the magnetic moment decreases discontinuously from approximately 3.0 \u0026micro;\u003csub\u003eB\u003c/sub\u003e to about 2.6\u0026ndash;2.7 \u0026micro;\u003csub\u003eB\u003c/sub\u003e, perfectly coinciding with the transition temperatures where the volume sharply contracts. This behavior corresponds to the spin-crossover phenomenon, where the high-spin state Co\u003csup\u003e3+\u003c/sup\u003e ions, stable at low temperatures, transition to a low-spin or intermediate state upon heating. The difference in the magnetic moment below and above the transition temperature decreases as the La substitution decreases, most probably because the thermally activated low- (high-) spin state is mixed to the ground state high- (low-) spin state.\u003c/p\u003e \u003cp\u003eThe magnetic moment of approximately 2.6\u0026ndash;2.7 \u003cem\u003e\u0026micro;\u003c/em\u003e\u003csub\u003eB\u003c/sub\u003e in the post-transition phase is too large to be attributed to a low-spin (LS, \u003cem\u003eS\u003c/em\u003e\u0026thinsp;=\u0026thinsp;0) state. Instead, this value strongly suggests a spin-crossover to a phase dominated by the intermediate-spin (IS, \u003cem\u003eS\u003c/em\u003e\u0026thinsp;=\u0026thinsp;1) configuration of Co\u0026sup3;⁺. This interpretation is consistent with detailed studies on the related perovskite LaCoO₃. Through GGA\u0026thinsp;+\u0026thinsp;\u003cem\u003eU\u003c/em\u003e calculations, Pandey et al. concluded that the electronic structure of LaCoO₃ at room temperature is best described by a dominant IS state configuration.\u003csup\u003e\u003cspan citationid=\"CR19\" class=\"CitationRef\"\u003e19\u003c/span\u003e\u003c/sup\u003e They also highlighted the strong hybridization between Co 3\u003cem\u003ed\u003c/em\u003e and O 2\u003cem\u003ep\u003c/em\u003e orbitals within the perovskite octahedra, which differentiates the system from a purely ionic picture. Therefore, we interpret the observed magnetic moment not as a covalent LS state, but as clear evidence of a transition to a high-temperature phase with a dominant IS character.\u003c/p\u003e \u003cp\u003eThis change in spin state directly impacts the lattice structure. High-spin Co\u003csup\u003e3+\u003c/sup\u003e ion is Jahn-Teller active, which distorts the octahedron and stabilizes the pyramidal coordination. In the intermediate- or low-spin state, however, the Jahn-Teller effect vanishes, allowing for a more symmetric octahedral coordination. This change in local structure is the driving force behind the lattice volume contraction, i.e., NTE. On the other hand, in compositions where NTE was not observed (\u003cem\u003ex\u003c/em\u003e\u0026thinsp;\u0026le;\u0026thinsp;0.125 and \u003cem\u003ex\u003c/em\u003e\u0026thinsp;=\u0026thinsp;1.0), the absence of such a clear spin-crossover transition supports this mechanism.\u003c/p\u003e \u003cp\u003eFrom the analysis above, our MLFF-MD simulations strongly suggest that the NTE in the Bi\u003csub\u003e1\u0026thinsp;\u0026minus;\u0026thinsp;x\u003c/sub\u003eLa\u003csub\u003ex\u003c/sub\u003eCoO\u003csub\u003e3\u003c/sub\u003e system is a phenomenon in which the lattice and spin degrees of freedom are strongly coupled, accompanied by a temperature-induced spin-crossover of Co\u003csup\u003e3+\u003c/sup\u003e ions. This demonstrates that a universal MLFF is not merely a tool for structural prediction but also an extremely effective method for capturing the essential physics, such as changes in electronic states, behind complex phase transitions.\u003c/p\u003e \u003cp\u003eTo investigate the reversibility of the phase transition and the hysteresis behavior, we performed cooling simulations for the \u003cem\u003ex\u003c/em\u003e\u0026thinsp;=\u0026thinsp;0.25 composition following the heating run. Figure\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003e displays the evolution of total energy and cell volume during the thermal cycle. Interestingly, the simulation revealed a \"one-way\" irreversible phase transition. As shown in Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003e\u003cb\u003e(a)\u003c/b\u003e, once the system transforms into the small-volume phase (GdFeO\u003csub\u003e3\u003c/sub\u003e-type) around 600 K during heating, it does not revert to the original large-volume phase upon cooling. Instead, it maintains the compact, octahedral-coordinated structure down to low temperatures. The energy profile (Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003e\u003cb\u003e(b)\u003c/b\u003e) provides a clear thermodynamic explanation for this behavior: the potential energy of the post-transition phase is lower than that of the initial pre-transition phase at the same temperature. This indicates that, within the potential energy surface described by CHGNet, the experimentally observed low-temperature PbTiO\u003csub\u003e3\u003c/sub\u003e-type phase (tetragonal) corresponds to a metastable local minimum, while the GdFeO\u003csub\u003e3\u003c/sub\u003e-type phase (orthorhombic) represents the global minimum for this composition. This irreversibility may be partly exaggerated by the finite simulation time and supercell size, which can hinder the nucleation of the reverse transition. As will be discussed in the next section, experimental measurements on La-substituted samples (\u003cem\u003ex\u003c/em\u003e\u0026thinsp;=\u0026thinsp;0.2) exhibit a reversible NTE, while the fraction of the orthorhombic phase changes only 25 to 57% on heating and to 32% on cooling, not from 0 to 100%.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eTo verify the \"one-way\" transition behavior suggested by the simulation, we compared these results with experimental synchrotron X-ray diffraction data for a compositionally similar sample, Bi\u003csub\u003e0.8\u003c/sub\u003eLa\u003csub\u003e0.2\u003c/sub\u003eCoO\u003csub\u003e3\u003c/sub\u003e (\u003cem\u003ex\u003c/em\u003e\u0026thinsp;=\u0026thinsp;0.20). Figure\u0026nbsp;\u003cspan refid=\"Fig8\" class=\"InternalRef\"\u003e8\u003c/span\u003e shows the temperature dependence of the unit cell volumes and the phase fraction.\u003c/p\u003e \u003cp\u003eThe experimental system exists as a two-phase mixture of the large-volume tetragonal phase (PbTiO\u003csub\u003e3\u003c/sub\u003e-type) and the small-volume orthorhombic phase (GdFeO\u003csub\u003e3\u003c/sub\u003e-type) even at 100 K. Upon heating, the fraction of the orthorhombic phase (black line) increases, from 25 to 57%, leading to the NTE behavior observed in the average volume (green line). Crucially, during the subsequent cooling process, the phase fraction exhibits significant hysteresis. The system does not fully revert to its initial state; instead, the fraction of the orthorhombic phase (32%) remains higher than it was before heating.\u003c/p\u003e \u003cp\u003eThis experimental observation qualitatively aligns with the irreversibility predicted by our MLFF-MD simulation. The simulation indicated that the orthorhombic phase is energetically more stable than the tetragonal phase over a wide temperature range for La-substituted compositions. In the idealized simulation box with periodic boundary conditions, this led to a complete lock-in to the stable orthorhombic phase. In the real polycrystalline sample, the transformation is spatially constrained by domain walls and inter-granular strain, preventing a complete transformation and resulting in the observed hysteresis and mixed-phase state. Nevertheless, both simulation and experiment point to the same conclusion: La substitution fundamentally stabilizes the collapsed GdFeO\u003csub\u003e3\u003c/sub\u003e-type structure, making the NTE transition partially or wholly irreversible depending on the microstructural constraints.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eFigure \u003cspan refid=\"Fig9\" class=\"InternalRef\"\u003e9\u003c/span\u003e shows the dependence of the NTE transition temperature on the La substitution amount, as obtained from MLFF-MD. Experimental results shown in the Fig.\u0026nbsp;\u003cspan refid=\"Fig8\" class=\"InternalRef\"\u003e8\u003c/span\u003e indicate that for a Bi\u003csub\u003e0.8\u003c/sub\u003eLa\u003csub\u003e0.2\u003c/sub\u003eCoO\u003csub\u003e3\u003c/sub\u003e, the NTE transition is observed over a wide temperature range from approximately 480 K to 700 K. The closest composition in our simulation, \u003cem\u003ex\u003c/em\u003e\u0026thinsp;=\u0026thinsp;0.25, predicted a transition temperature of 725 K, which is considerably higher than the onset temperature observed experimentally. For the \u003cem\u003ex\u003c/em\u003e\u0026thinsp;=\u0026thinsp;0.5 composition, the MLFF prediction (approx. 350 K) was 100\u0026ndash;150 K lower than the AIMD prediction (450\u0026ndash;500 K). These results indicate that the universal MLFF has quantitative challenges in predicting the absolute values of transition temperatures. Nevertheless, it can be said to possess sufficient capability to correctly capture the chemical trends, such as whether NTE will occur with La substitution and how the transition temperature changes as the substitution amount increases.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eThe true value of the MLFF-MD simulation used in this study extends beyond the prediction of thermodynamic properties. It functions as a \"computational microscope\" that enables the dynamic process of how a phase transition proceeds to be directly observed, with the movement of each individual atom tracked. Figure\u0026nbsp;\u003cspan refid=\"Fig10\" class=\"InternalRef\"\u003e10\u003c/span\u003e shows representative snapshots of atomic configurations within \u0026plusmn;\u0026thinsp;20 fs of the moment the NTE phase transition occurs in the heating simulation of Bi\u003csub\u003e0.5\u003c/sub\u003eLa\u003csub\u003e0.5\u003c/sub\u003eCoO\u003csub\u003e3\u003c/sub\u003e.\u003c/p\u003e \u003cp\u003e \u003cul\u003e \u003cli\u003e \u003cp\u003e \u003cb\u003eBefore the transition\u003c/b\u003e: In the low-temperature phase, the entire system adopts a near-tetragonal structure elongated along the \u003cem\u003ec\u003c/em\u003e-axis, with Co\u003csup\u003e3+\u003c/sup\u003e ions in a 5-coordinated pyramidal configuration.\u003c/p\u003e \u003c/li\u003e \u003cli\u003e \u003cp\u003e \u003cb\u003eDuring the transition\u003c/b\u003e: Most notably, the snapshot during the transition clearly shows that the transformation to the small-volume high-temperature phase (with Co\u003csup\u003e3+\u003c/sup\u003e in 6-coordinated octahedral structures) does not occur uniformly throughout the system. It clearly captures the transition propagating from the La\u003csup\u003e3+\u003c/sup\u003e-rich regions (which lack lone pairs) to the Bi\u003csup\u003e3+\u003c/sup\u003e-rich regions.\u003c/p\u003e \u003c/li\u003e \u003cli\u003e \u003cp\u003e \u003cb\u003eAfter the transition\u003c/b\u003e: Ultimately, the entire system completely transforms into the high-temperature phase with octahedral coordination.\u003c/p\u003e \u003c/li\u003e \u003c/ul\u003e \u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eThis visualization provides a new perspective for understanding complex phase transition mechanisms. Experimental techniques, typified by in-situ transmission electron microscopy (TEM) observation, are powerful methods able to capture structural changes under heating with atomic resolution. However, the time resolution of conventional in-situ TEM observations, which is typically limited by the video frame rate to the order of milliseconds, is insufficient to resolve the femtosecond-to-picosecond timescale on which the elementary processes of a phase transition, such as atomic vibrations, occur. Therefore, the information obtained from such techniques is essentially a time-averaged structure over these ultrafast events. MLFF-MD simulation can track the dynamic behavior of individual atoms in this ultrafast time domain, which is difficult to observe experimentally. The heterogeneous propagation behavior, nucleating from La-rich regions, as revealed here, is a transient phenomenon that cannot be captured from the averaged information. Furthermore, the ability to perform such systematic dynamic analyses for numerous compositions at a realistic computational cost is a significant advantage over experiments. Therefore, MLFF-MD not only enables high-speed screening but also plays a complementary role to experimental observations, demonstrating its potential to become a new analytical method for elucidating the fundamental elementary processes of complex solid-state reactions.\u003c/p\u003e"},{"header":"Discussion","content":"\u003cp\u003eThe quantitative deviation in the transition temperature observed in this simulation is thought to stem from the inherent limitations of universal MLFFs. Several recent studies have pointed out that many universal MLFFs, including CHGNet, exhibit a tendency known as \"systematic softening,\" where they underestimate the interatomic forces in high-energy regions. \u003csup\u003e20\u003c/sup\u003e The primary cause of this phenomenon is the training dataset for the universal model (in CHGNet's case, the Materials Project DFT relaxation trajectories\u003csup\u003e\u003cspan citationid=\"CR21\" class=\"CitationRef\"\u003e21\u003c/span\u003e\u003c/sup\u003e) biased towards near-equilibrium structures around energy minima, with insufficient sampling of high-energy regions such as transition states and highly strained structures. The phase transition temperature is determined by the point where the free energies of the two phases become equal, and the contribution of vibrational entropy is significant to the free energy. If the potential energy surface (PES) is soft (i.e., has a small curvature), phonon frequencies are underestimated, which in turn affects the evaluation of vibrational entropy and can shift the predicted transition temperature. The underestimation of the transition temperatures seen in this study can be regarded as a direct manifestation of this systematic softening.\u003c/p\u003e \u003cp\u003eTo further validate the reliability of the MLFF-MD simulation and assess the impact of this systematic softening, we performed single-point DFT calculations on 400 snapshots extracted from the MD trajectory of Bi\u003csub\u003e0.5\u003c/sub\u003eLa\u003csub\u003e0.5\u003c/sub\u003eCoO\u003csub\u003e3\u003c/sub\u003e during the heating process. Figure\u0026nbsp;\u003cspan refid=\"Fig11\" class=\"InternalRef\"\u003e11\u003c/span\u003e compares the potential energies predicted by CHGNet with those calculated by DFT (VASP). Crucially, both methods exhibit a clear and simultaneous drop in total energy around snapshot 4,000, which corresponds to the experimentally and computationally observed phase transition temperature. This agreement confirms that the phase transition reproduced by CHGNet is not a mere artifact of the machine learning model but is physically driven by the stabilization of the high-temperature phase, consistent with first-principles energetics. However, as shown in the high-temperature region (around and after snapshot 4,000), the energy difference between CHGNet and DFT gradually increases. This trend provides that MLFF tends to under- or overestimate the energy of structures with large thermal displacements. Nevertheless, the fact that the distinct energy jump associated with the first-order transition is qualitatively captured supports the validity of using this universal MLFF to elucidate the transition mechanism, even if quantitative discrepancies in the transition temperature persist.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eOn the other hand, a likely reason for CHGNet's ability to capture the qualitatively correct physics is that its architecture explicitly incorporates magnetic moment (spin) information.\u003csup\u003e\u003cspan citationid=\"CR12\" class=\"CitationRef\"\u003e12\u003c/span\u003e\u003c/sup\u003e Since the phase transition in the BiCoO\u003csub\u003e3\u003c/sub\u003e system is closely related to the change in the spin state of the Co ions, this feature likely enabled the model to somewhat learn the coupling between the electronic state and the structure, leading to its qualitative success. Furthermore, the computational efficiency of the MLFF-MD method employed in this study is remarkable. An AIMD calculation using a 40-atom 2 \u0026times; 2 \u0026times; 2 supercell required approximately 20 hours to run a 2-ps simulation at a single fixed temperature. In contrast, an MLFF-MD calculation using a 320-atom 4 \u0026times; 4 \u0026times; 4 supercell (8 times larger than AIMD) and simulating a continuous temperature ramp from 1 K to 1000 K (1000 points) took only 8 hours. A simple throughput comparison shows that the MLFF achieved a computational speed over 20,000 times faster than AIMD. This overwhelming speed-up makes comprehensive simulations (screenings) of material behavior over wide compositional and temperature ranges, which were previously impractical due to computational cost, a tangible reality.\u003c/p\u003e \u003cp\u003eHowever, the current universal training data is insufficient to quantitatively and accurately describe the subtle energy changes associated with spin-crossover, which is considered another source of quantitative error. For future high-precision predictions, an approach that uses the universal model as a base and then fine-tunes it with a small amount of high-accuracy DFT data specific to the target system (such as the AIMD trajectories generated in this study) is suggested to be effective.\u003csup\u003e\u003cspan citationid=\"CR22\" class=\"CitationRef\"\u003e22\u003c/span\u003e\u003c/sup\u003e\u003c/p\u003e"},{"header":"Conclusion","content":"\u003cp\u003eIn this study, we applied large-scale molecular dynamics (MD) simulations using the Crystal Hamilton Graph Network (CHGNet), a pretrained universal machine learning force field (MLFF), to elucidate the negative thermal expansion (NTE) phenomenon in the perovskite oxide Bi\u003csub\u003e1\u0026thinsp;\u0026minus;\u0026thinsp;x\u003c/sub\u003eLa\u003csub\u003ex\u003c/sub\u003eCoO\u003csub\u003e3\u003c/sub\u003e.\u003c/p\u003e \u003cp\u003eThe main achievements are as follows:\u003c/p\u003e \u003cp\u003e \u003col\u003e \u003cspan\u003e \u003cli\u003e \u003cp\u003e \u003cb\u003eMethodology Validation and NTE Reproduction\u003c/b\u003e: The universal MLFF, CHGNet, qualitatively reproduced the giant volume contraction induced by pressure in BiCoO\u003csub\u003e3\u003c/sub\u003e. Furthermore, it successfully reproduced the experimentally observed pressure-induced phase transition in BiCoO\u003csub\u003e3\u003c/sub\u003e and NTE in the Bi\u003csub\u003e1\u0026thinsp;\u0026minus;\u0026thinsp;\u003cem\u003ex\u003c/em\u003e\u003c/sub\u003eLa\u003csub\u003e\u003cem\u003ex\u003c/em\u003e\u003c/sub\u003eCoO\u003csub\u003e3\u003c/sub\u003e, and indicated the chemical trend, where the NTE transition temperature decreases as La substitution increases.\u003c/p\u003e \u003c/li\u003e \u003c/span\u003e \u003cspan\u003e \u003cli\u003e \u003cp\u003e \u003cb\u003eElucidation of Atomic-Level Dynamic Mechanism\u003c/b\u003e: By leveraging the overwhelming computational efficiency of MLFF-MD (over 20,000 times faster than ab initio molecular dynamics (AIMD)), we visualized the dynamic process of the NTE phase transition at the atomic level. As a result, we uncovered the nucleation and propagation behavior, where the phase transition initiates heterogeneously in La-rich regions and propagates to Bi-rich regions.\u003c/p\u003e \u003c/li\u003e \u003c/span\u003e \u003c/ol\u003e \u003c/p\u003e \u003cp\u003eThese achievements concretely demonstrate that universal MLFF-MD simulation can serve as a powerful \"computational microscope\" for predicting the behavior of functional inorganic materials with complex phase transitions and for elucidating atomic-level mechanisms that are unobservable in experiments. This research not only embodies a new paradigm for materials science research driven by computational and data sciences but also provides important feedback for developing next-generation, higher-precision foundational models.\u003c/p\u003e"},{"header":"Methods","content":"\u003cp\u003eMachine Learning Force Field Molecular Dynamics Simulations\u003c/p\u003e \u003cp\u003eIn this study, the Crystal Hamilton Graph Network (CHGNet) (ver. 0.4.0) was employed as the machine learning force field (MLFF). \u003csup\u003e12\u003c/sup\u003e CHGNet is a universal MLFF pretrained on the extensive first-principles calculation database of the Materials Project.\u003csup\u003e\u003cspan citationid=\"CR21\" class=\"CitationRef\"\u003e21\u003c/span\u003e\u003c/sup\u003e It was used in this research without any additional training (fine-tuning).\u003c/p\u003e \u003cp\u003eMolecular dynamics (MD) simulations were performed using the open-source physics simulation library Atomic Simulation Environment (ASE) (ver. 3.26.0)\u003csup\u003e23\u003c/sup\u003e with PyTorch (ver. 2.7.1) as its backend. All MLFF-MD calculations were executed on a single NVIDIA RTX 4090 GPU. A 4 \u0026times; 4 \u0026times; 4 supercell (320 atoms) created by expanding the unit cell of BiCoO\u003csub\u003e3\u003c/sub\u003e and Bi\u003csub\u003e1\u0026thinsp;\u0026minus;\u0026thinsp;x\u003c/sub\u003eLa\u003csub\u003ex\u003c/sub\u003eCoO\u003csub\u003e3\u003c/sub\u003e solid solutions was used as the simulation cell. Prior to the MD simulations, the initial structure for each composition was optimized by using the BFGS algorithm, a quasi-Newton method, until the residual force on each atom was less than 0.001 eV/\u0026Aring;.\u003c/p\u003e \u003cp\u003eFor the subsequent evaluation of thermal expansion behavior, the NpT ensemble (isothermal-isobaric ensemble) was employed to maintain constant pressure and temperature. Temperature was controlled using a Nos\u0026eacute;-Hoover thermostat (time constant\u0026thinsp;=\u0026thinsp;25 fs)\u003csup\u003e\u003cspan citationid=\"CR24\" class=\"CitationRef\"\u003e24\u003c/span\u003e\u003c/sup\u003e, and pressure was controlled with a Nos\u0026eacute;-Hoover barostat. A time step of 1.0 fs was used.\u003c/p\u003e \u003cp\u003eIn the pressure-induced phase transition simulations, the temperature was kept at 10 K while the pressure was incrementally increased from ambient pressure (1 bar) to 5 GPa in steps of 0.1 GPa. The structure was held for 2 ps at each pressure point to ensure sufficient relaxation. For the temperature-induced phase transition simulations, the pressure was maintained at ambient pressure (1 bar), and the temperature was ramped up and down from 0 K to 1000 K in 1 K intervals. At each temperature point, an MD calculation of 2000 steps (2 ps) was performed to equilibrate the system, and the lattice parameters and volume were recorded. In the heating simulations, the final structure from the previous temperature step was used as the initial structure for the next step to ensure continuity.\u003c/p\u003e\n\u003ch3\u003eFirst-Principles Molecular Dynamics Simulations\u003c/h3\u003e\n\u003cp\u003eTo compare and validate the MLFF-MD results, ab initio molecular dynamics (AIMD) calculations were also performed using the VASP (Vienna Ab initio Simulation Package).\u003csup\u003e\u003cspan citationid=\"CR25\" class=\"CitationRef\"\u003e25\u003c/span\u003e\u003c/sup\u003e\u003c/p\u003e \u003cp\u003eA 2 \u0026times; 2 \u0026times; 2 supercell of Bi\u003csub\u003e0.5\u003c/sub\u003eLa\u003csub\u003e0.5\u003c/sub\u003eCoO\u003csub\u003e3\u003c/sub\u003e (40 atoms: Bi\u003csub\u003e4\u003c/sub\u003eLa\u003csub\u003e4\u003c/sub\u003eCo\u003csub\u003e8\u003c/sub\u003eO\u003csub\u003e24\u003c/sub\u003e) was used for the simulation cell. The Perdew-Burke-Ernzerhof (PBE) functional was used to describe the exchange-correlation energy of electrons. To account for the strong electron correlation among the 3d electrons of Co atoms, the GGA\u0026thinsp;+\u0026thinsp;U method (Dudarev's approach)\u003csup\u003e\u003cspan citationid=\"CR26\" class=\"CitationRef\"\u003e26\u003c/span\u003e,\u003cspan citationid=\"CR27\" class=\"CitationRef\"\u003e27\u003c/span\u003e\u003c/sup\u003e was applied, with an effective U parameter (\u003cem\u003eU\u003c/em\u003e\u003csub\u003eeff\u003c/sub\u003e) of 3.32 eV set for the Co d-orbitals. The cutoff energy for the plane-wave basis set was 520 eV, and the first Brillouin zone was sampled using a 3\u0026times;3\u0026times;2 k-point mesh centered at the Gamma point. The simulations were conducted in the NpT ensemble, with temperature and pressure controlled by using Langevin dynamics. Calculations were performed at temperatures of 300, 400, 450, and 500 K for 2 ps (2000 steps) with a time step of 1.0 fs. Structural information for the equilibrium state was obtained from the trajectory of the final 1 ps.\u003c/p\u003e"},{"header":"Declarations","content":"\u003cp\u003e \u003ch2\u003eEthics declarations\u003c/h2\u003e \u003cp\u003eAll authors declare no financial or non-financial competing interests.\u003c/p\u003e \u003c/p\u003e\u003ch2\u003eFunding\u003c/h2\u003e \u003cp\u003eThis work was partially supported by JSPS KAKENHI (grant no. JP24H00374), JST-CREST (JPMJCR22O1), the Kanagawa Institute of Industrial Science and Technology and Design and Engineering by Joint Inverse Innovation for Materials Architecture project of the Ministry of Education, Culture, Sports, Science and Technology.\u003c/p\u003e\u003ch2\u003eAuthor Contribution\u003c/h2\u003e\u003cp\u003eS.W. performed the MD simulations and analyzed the data. Y.S, K.T. and M.A performed the high-pressure synthesis and synchrotron X-ray diffraction experiments. H.D. and M.A. conceived the project and supervised the work. All authors contributed to writing the manuscript.\u003c/p\u003e\u003ch2\u003eAcknowledgement\u003c/h2\u003e\u003cp\u003eThis work was partially supported by JSPS KAKENHI (grant no. JP24H00374), JST-CREST (JPMJCR22O1), the Kanagawa Institute of Industrial Science and Technology and Design and Engineering by Joint Inverse Innovation for Materials Architecture project of the Ministry of Education, Culture, Sports, Science and Technology.\u003c/p\u003e\u003ch2\u003eData Availability\u003c/h2\u003e\u003cp\u003eAll data generated or analyzed during this study are included in this published article and its supplementary information files\u003c/p\u003e\u003cp\u003e\u003cstrong\u003eCode Availability\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eThe underlying code for this study is not publicly available but may be made available to qualified researchers on reasonable requests from the corresponding author.\u003c/p\u003e"},{"header":"References","content":"\u003col\u003e\u003cli\u003e\u003cspan\u003eAzuma, M. Phase transition type giant negative thermal expansion materials. \u003cem\u003eJ. Ceram. Soc. Jpn\u003c/em\u003e. \u003cb\u003e133\u003c/b\u003e, 450\u0026ndash;454 (2025).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eYamamoto, H. et al. Colossal negative thermal expansion in electron-doped PbVO₃. \u003cem\u003eAngew Chem. Int. Ed.\u003c/em\u003e \u003cb\u003e57\u003c/b\u003e, 8170\u0026ndash;8173 (2018).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eNishikubo, T. et al. Polar\u0026ndash;Nonpolar Transition-Type Negative Thermal Expansion with 11.1% Volume Shrinkage by Design. \u003cem\u003eChem. Mater.\u003c/em\u003e \u003cb\u003e35\u003c/b\u003e, 870\u0026ndash;878 (2023).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eAzuma, M. et al. Giant negative thermal expansion in BiNiO₃. \u003cem\u003eNat. Commun.\u003c/em\u003e \u003cb\u003e2\u003c/b\u003e, 347 (2011).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBelik, A. A. et al. Neutron Powder Diffraction Study on the Crystal and Magnetic Structures of BiCoO\u003csub\u003e3\u003c/sub\u003e. \u003cem\u003eChem. Mater.\u003c/em\u003e \u003cb\u003e18\u003c/b\u003e, 798\u0026ndash;803 (2006).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eOka, K. et al. Pressure-induced spin-state transition in BiCoO₃. \u003cem\u003eJ. Am. Chem. Soc.\u003c/em\u003e \u003cb\u003e132\u003c/b\u003e, 9438\u0026ndash;9443 (2010).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003ePan, Z. et al. Pronounced negative thermal expansion in lead-free BiCoO₃-based ferroelectrics triggered by the stabilized perovskite structure. \u003cem\u003eChem. Mater.\u003c/em\u003e \u003cb\u003e31\u003c/b\u003e, 6187\u0026ndash;6192 (2019).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBecke, A. D. \u0026amp; Perspective Fifty years of density-functional theory in chemical physics. \u003cem\u003eJ. Chem. Phys.\u003c/em\u003e \u003cb\u003e140\u003c/b\u003e, 18A301 (2014).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eUnke, O. T. et al. Machine learning force fields. \u003cem\u003eChem. Rev.\u003c/em\u003e \u003cb\u003e121\u003c/b\u003e, 10142\u0026ndash;10186 (2021).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eChen, C. \u0026amp; Ong, S. P. A universal graph deep learning interatomic potential for the periodic table. \u003cem\u003eNat. Comput. Sci.\u003c/em\u003e \u003cb\u003e2\u003c/b\u003e, 718\u0026ndash;728 (2022).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBatatia, I. et al. Higher order equivariant message passing neural networks for fast and accurate force fields. \u003cem\u003eAdv. Neural Inf. Process. Syst.\u003c/em\u003e \u003cb\u003e35\u003c/b\u003e, 11487\u0026ndash;11501 (2022).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eDeng, B. et al. CHGNet as a pretrained universal neural network potential for charge-informed atomistic modelling. \u003cem\u003eNat. Mach. Intell.\u003c/em\u003e \u003cb\u003e5\u003c/b\u003e, 1031\u0026ndash;1041 (2023).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLi, D., Yang, J., Chen, X., Yu, L. \u0026amp; Liu, S. To Use or Not to Use a Universal Force Field. Preprint at (2025). \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://doi.org/10.48550/ARXIV.2503.08207\u003c/span\u003e\u003cspan address=\"10.48550/ARXIV.2503.08207\" targettype=\"DOI\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHe, R., Wu, H., Lu, Y. \u0026amp; Zhong, Z. Origin of negative thermal expansion and pressure-induced amorphization in zirconium tungstate from a machine-learning potential. \u003cem\u003ePhys. Rev. B\u003c/em\u003e. \u003cb\u003e106\u003c/b\u003e, 174101 (2022).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eShi, G., Zhong, Z., Su, Y., He, R. \u0026amp; Zhang, C. Investigating the atomic-level mechanisms of negative thermal expansion in Zn(CN)2 using machine learning-enabled molecular dynamics. \u003cem\u003ePhys. Scr.\u003c/em\u003e \u003cb\u003e100\u003c/b\u003e, 066002 (2025).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eZhao, K. et al. Octahedral tilt distortion in negative thermal expansion in the fluorides CaZrF\u003csub\u003e6\u003c/sub\u003e and ScF\u003csub\u003e3\u003c/sub\u003e. \u003cem\u003ePhys. Rev. B\u003c/em\u003e. \u003cb\u003e110\u003c/b\u003e, 064322 (2024).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eWu, J., Yang, J., Ma, L., Zhang, L. \u0026amp; Liu, S. Modular development of deep potential for complex solid solutions. \u003cem\u003ePhys. Rev. B\u003c/em\u003e. \u003cb\u003e107\u003c/b\u003e, 144102 (2023).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eJain, A. et al. Formation enthalpies by mixing GGA and GGA calculations. \u003cem\u003ePhys. Rev. B\u003c/em\u003e. \u003cb\u003e84\u003c/b\u003e, 045115 (2011).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003ePandey, S. K. et al. Investigation of the spin state of Co in LaCoO\u003csub\u003e3\u003c/sub\u003e at room temperature: Ab initio calculations and high-resolution photoemission spectroscopy of single crystals. \u003cem\u003ePhys Rev. B\u003c/em\u003e \u003cb\u003e77\u003c/b\u003e, (2008).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eDeng, B. et al. Systematic softening in universal machine learning interatomic potentials. \u003cem\u003enpj Comput. Mater\u003c/em\u003e \u003cb\u003e11\u003c/b\u003e, (2025).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eJain, A. et al. The materials project: A materials genome approach to accelerating materials innovation. \u003cem\u003eAPL Mater.\u003c/em\u003e \u003cb\u003e1\u003c/b\u003e, 011002 (2013).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eChen, P. Y., Shibata, K. \u0026amp; Mizoguchi, T. High precision machine learning force field development for BaTiO3 phase transitions, amorphous, and liquid structures. \u003cem\u003eAPL Mach. Learn.\u003c/em\u003e \u003cb\u003e3\u003c/b\u003e, 036115 (2025).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLarsen, A. H. et al. The atomic simulation environment\u0026mdash;a Python library for working with atoms. \u003cem\u003eJ. Phys. : Condens. Matter\u003c/em\u003e. \u003cb\u003e29\u003c/b\u003e, 273002 (2017).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eNos\u0026eacute;, S. A molecular dynamics method for simulations in the canonical ensemble. J. Chem. Phys. 81, 511\u0026ndash;519 ; Hoover, W. G. Canonical dynamics: Equilibrium phase-space distributions. \u003cem\u003ePhys. Rev. A\u003c/em\u003e 31, 1695\u0026ndash;1697 (1985). (1984).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eKresse, G. \u0026amp; Furthm\u0026uuml;ller, J. Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set. \u003cem\u003ePhys. Rev. B\u003c/em\u003e. \u003cb\u003e54\u003c/b\u003e, 11169\u0026ndash;11186 (1996).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003ePerdew, J. P., Burke, K. \u0026amp; Ernzerhof, M. Generalized gradient approximation made simple. \u003cem\u003ePhys. Rev. Lett.\u003c/em\u003e \u003cb\u003e77\u003c/b\u003e, 3865\u0026ndash;3868 (1996).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eDudarev, S. L., Botton, G. A., Savrasov, S. Y., Humphreys, C. J. \u0026amp; Sutton, A. P. Electron-energy-loss spectra and the structural stability of nickel oxide: An LSDA\u0026thinsp;+\u0026thinsp;U study. \u003cem\u003ePhys. Rev. B\u003c/em\u003e. \u003cb\u003e57\u003c/b\u003e, 1505\u0026ndash;1509 (1998).\u003c/span\u003e\u003c/li\u003e\u003c/ol\u003e"}],"fulltextSource":"","fullText":"","funders":[],"hasAdminPriorityOnWorkflow":false,"hasManuscriptDocX":true,"hasOptedInToPreprint":true,"hasPassedJournalQc":"","hasAnyPriority":false,"hideJournal":false,"highlight":"","institution":"","isAcceptedByJournal":true,"isAuthorSuppliedPdf":false,"isDeskRejected":"","isHiddenFromSearch":false,"isInQc":false,"isInWorkflow":false,"isPdf":false,"isPdfUpToDate":true,"isWithdrawnOrRetracted":false,"journal":{"display":true,"email":"[email protected]","identity":"scientific-reports","isNatureJournal":false,"hasQc":true,"allowDirectSubmit":false,"externalIdentity":"scirep","sideBox":"Learn more about [Scientific Reports](http://www.nature.com/srep/)","snPcode":"","submissionUrl":"","title":"Scientific Reports","twitterHandle":"","acdcEnabled":true,"dfaEnabled":true,"editorialSystem":"stoa","reportingPortfolio":"Scientific Reports","inReviewEnabled":true,"inReviewRevisionsEnabled":true},"keywords":"","lastPublishedDoi":"10.21203/rs.3.rs-8425406/v1","lastPublishedDoiUrl":"https://doi.org/10.21203/rs.3.rs-8425406/v1","license":{"name":"CC BY 4.0","url":"https://creativecommons.org/licenses/by/4.0/"},"manuscriptAbstract":"\u003cp\u003eNegative Thermal Expansion (NTE) is a key property for achieving dimensional stability in advanced technologies, yet understanding its underlying mechanisms, particularly the dynamic processes of phase transitions, remains a formidable challenge. In this study, we applied large-scale molecular dynamics simulations using the Crystal Hamilton Graph Network (CHGNet), a pretrained universal machine learning force field (MLFF), to elucidate the atomic-level mechanisms of the NTE phenomenon in the perovskite oxide system Bi\u003csub\u003e1\u0026thinsp;\u0026minus;\u0026thinsp;\u003cem\u003ex\u003c/em\u003e\u003c/sub\u003eLa\u003csub\u003e\u003cem\u003ex\u003c/em\u003e\u003c/sub\u003eCoO\u003csub\u003e3\u003c/sub\u003e. CHGNet achieved a computational speed approximately 20,000 times faster than first-principles calculations and qualitatively reproduced the chemical trend of the NTE transition temperature systematically decreasing as La substitution increased. Leveraging this computational efficiency, we successfully visualized the dynamic process of the NTE phase transition at the atomic level, which had previously been unobservable. The results revealed a nucleation and propagation mechanism wherein the phase transition initiates heterogeneously in La-rich regions and propagates to Bi-rich regions. These findings demonstrate that universal MLFFs are powerful tools for more quickly elucidating complex phase transition phenomena and open up new possibilities for computational science-driven exploration of new materials.\u003c/p\u003e","manuscriptTitle":"Large-scale molecular dynamics simulations of negative thermal expansion in Bi1-x Lax CoO3 using a pretrained machine learning force field","msid":"","msnumber":"","nonDraftVersions":[{"code":1,"date":"2025-12-26 19:42:22","doi":"10.21203/rs.3.rs-8425406/v1","editorialEvents":[{"type":"communityComments","content":0},{"type":"decision","content":"Revision requested","date":"2026-01-19T08:44:22+00:00","index":"","fulltext":""},{"type":"editorInvitedReview","content":"","date":"2026-01-16T14:19:09+00:00","index":"hide","fulltext":""},{"type":"editorInvitedReview","content":"","date":"2026-01-13T16:20:53+00:00","index":"hide","fulltext":""},{"type":"reviewerAgreed","content":"331539483552123604971899834200561880440","date":"2026-01-13T01:47:04+00:00","index":"hide","fulltext":""},{"type":"reviewerAgreed","content":"91075228569335716281505361242254166482","date":"2026-01-12T20:06:52+00:00","index":"hide","fulltext":""},{"type":"reviewerAgreed","content":"7296019319972653084968396017081727435","date":"2026-01-12T07:02:43+00:00","index":"hide","fulltext":""},{"type":"editorInvitedReview","content":"","date":"2026-01-04T08:40:47+00:00","index":"hide","fulltext":""},{"type":"reviewerAgreed","content":"244301348387655879105190633491269249084","date":"2025-12-25T04:16:41+00:00","index":"hide","fulltext":""},{"type":"reviewerAgreed","content":"140793826969848144414753115097489725673","date":"2025-12-24T19:45:06+00:00","index":"hide","fulltext":""},{"type":"reviewersInvited","content":"","date":"2025-12-24T13:35:33+00:00","index":"","fulltext":""},{"type":"editorInvited","content":"","date":"2025-12-24T12:44:27+00:00","index":"","fulltext":""},{"type":"editorAssigned","content":"","date":"2025-12-23T03:06:03+00:00","index":"","fulltext":""},{"type":"checksComplete","content":"","date":"2025-12-23T03:04:42+00:00","index":"","fulltext":""},{"type":"submitted","content":"Scientific Reports","date":"2025-12-22T13:06:08+00:00","index":"","fulltext":""}],"status":"published","journal":{"display":true,"email":"[email protected]","identity":"scientific-reports","isNatureJournal":false,"hasQc":true,"allowDirectSubmit":false,"externalIdentity":"scirep","sideBox":"Learn more about [Scientific Reports](http://www.nature.com/srep/)","snPcode":"","submissionUrl":"","title":"Scientific Reports","twitterHandle":"","acdcEnabled":true,"dfaEnabled":true,"editorialSystem":"stoa","reportingPortfolio":"Scientific Reports","inReviewEnabled":true,"inReviewRevisionsEnabled":true}}],"origin":"","ownerIdentity":"d591f92c-7b0d-4586-9f1c-16f44c99d4a1","owner":[],"postedDate":"December 26th, 2025","published":true,"recentEditorialEvents":[],"rejectedJournal":[],"revision":"","amendment":"","status":"published-in-journal","subjectAreas":[{"id":60201996,"name":"Physical sciences/Chemistry"},{"id":60201997,"name":"Physical sciences/Materials science"},{"id":60201998,"name":"Physical sciences/Mathematics and computing"},{"id":60201999,"name":"Physical sciences/Physics"}],"tags":[],"updatedAt":"2026-04-13T16:02:25+00:00","versionOfRecord":{"articleIdentity":"rs-8425406","link":"https://doi.org/10.1038/s41598-026-47983-9","journal":{"identity":"scientific-reports","isVorOnly":false,"title":"Scientific Reports"},"publishedOn":"2026-04-10 15:57:49","publishedOnDateReadable":"April 10th, 2026"},"versionCreatedAt":"2025-12-26 19:42:22","video":"","vorDoi":"10.1038/s41598-026-47983-9","vorDoiUrl":"https://doi.org/10.1038/s41598-026-47983-9","workflowStages":[]},"version":"v1","identity":"rs-8425406","journalConfig":"researchsquare"},"__N_SSP":true},"page":"/article/[identity]/[[...version]]","query":{"redirect":"/article/rs-8425406","identity":"rs-8425406","version":["v1"]},"buildId":"8U1c8b4HqxoKbykW_rLl7","isFallback":false,"isExperimentalCompile":false,"dynamicIds":[84888],"gssp":true,"scriptLoader":[]}