With the development of electronic structure theory, a new class of materials—quantum ones—has been recognized by the community. Traditionally, it has been believed that the properties of such compounds cannot be described within the framework of modern density functional theory, and indeed, more advanced post-mean-field theory methods are needed. Motivated by this, herein, we develop a fundamental understanding of such complex materials using the example of paramagnetic YNiO3, which is experimentally known to exhibit metal-to-insulator phase transition. We show that this material has a temperature-dependent distribution of local motifs. Thus, while at low temperatures, YNiO3 has distinct structural disproportionation with the formation of large and small octahedra, as the temperature increases, this disproportionation is suppressed. We also explain the paramagnetic monoclinic to paramagnetic orthorhombic phase transition within the double-well to single-well energy profile, predicting the variation in the corresponding energy profile as a function of octahedral size distribution. In this way, we demonstrate a fundamental understanding of structural phase transitions in quantum materials, giving insights into how they can be used for different applications and what minimum level of theory is needed to describe such types of complex materials at finite temperatures.
I. INTRODUCTION
Atomic identities, composition, and structure define the electronic properties of solid-state materials1 and have been used to develop open-access databases (e.g., Materials Project,2 AFLOW,3 and OQMD4) containing electronic structures for thousands of materials. In recent years, it has also become clear that in addition to “classical” solid-state compounds (e.g., Si, ZnO), many quantum materials (especially paramagnetic phases)5–14 exhibit a range of spin and/or structural symmetry breaking. For instance, in the case of magnetic materials, systems are considered to be magnetically long-range ordered until a certain critical temperature—the Néel temperature, above which the spin long-range order breaks, leading to the formation of a spin polymorphous system. Similarly, many quantum materials exhibit significant changes in local structural motifs with temperature, i.e., one associated with a monoclinic to orthorhombic or orthorhombic to cubic crystal structure phase transition.5 The common tendency observed experimentally is that the crystal structure symmetry increases with temperature,15 leading to a transition from multiple local environments to a single local environment.
For a long time, the theoretical community has largely ignored the distribution of local motifs in quantum materials. This is primarily attributed to the deeply entrenched beliefs that the physics of these compounds can be explained through strong correlation,16–19 where the properties of materials are often defined by the ratio of electron–electron repulsion (U) to bandwidth (W) without considering the distribution of local motifs. This perspective presumes that the interplay of electrons within a material—their correlations—dominates the behaviors we observe, such as conductivity, magnetism, and superconductivity. Recently, however, it has become evident that the role of strong correlation (which is different than the conventional Coulomb correlation20) has been overestimated in many materials, and indeed, local motifs play a key/important role in determining material properties.5–14,21,22 This also, in large part, challenges a correlation picture of bandwidth-driven insulator-to-metal phase transition (i.e., insulator-to-metal phase transition caused by the reduction in U/W ratio) and highlights the possible role of the distribution of local motifs in materials properties. For instance, in strongly correlated materials, dynamical mean field theory has recently provided a means to account for structural energy lowering21–25 and even phonon calculations,26 effectively addressing temperature effects in some degree of accuracy.21,27 These studies provide a thorough explanation of the metal–insulator transition in paramagnetic quantum materials, emphasizing charge disproportionation, suppression of orbital fluctuation, and bond disproportionation. However, these methods are not yet routinely used for all studies, and moreover, for paramagnets, even such seminal works are still usually based on density functional theory (DFT) for a global average structure where local and total magnetic moments are set to 0 (which is a naïve approximation of the paramagnets5–7) as the basis for more advanced theories. Moreover, the understanding of the temperature-dependent distribution of local motifs and its impact on material properties remains largely unexplored, especially in the strongly correlated compounds.
Herein, we use the example of YNiO3, which exists in three known experimentally observed phases having different distributions of local motifs:
The low-temperature α phase has a monoclinic crystal structure and an antiferromagnetic spin arrangement.28 It exhibits structural disproportionation in the form of the NiO6 octahedral breathing mode (i.e., the defining feature of the monoclinic distortion) and magnetic disproportionation in the form of low and high spin sites.5,7
Above the Néel temperature of 145 K,28 YNiO3 exists in the β phase, which retains the crystal structure and structural disproportionation of the α phase, but the magnetic long-range order is reduced to a paramagnetic arrangement with a net-zero magnetic moment. Importantly, it shows a distribution of local spins, which is crucial to maintain its insulating properties.5–7
At temperatures exceeding 582 K, the system undergoes a transition to an orthorhombic lattice system in the γ phase, concomitant with an insulator-to-metal transition,29 while both of the aforementioned disproportionation factors are washed out due to thermal motion.5
While the polymorphous theory of some representative quantum materials has been developed5–14 and some fundamental physics of such compounds has been understood, the temperature-dependent physics of local motifs is rarely explored. Motivated by this, herein, we apply a spin polymorphous description of the paramagnetic phase within static and ab initio molecular dynamics (AIMD) to study the temperature-dependent distribution of local motifs in the β and γ phases. This approach allows us to observe the gradual softening of the structural disproportionation tendency, revealing the physics of phase transition and helping us to understand the temperature-induced changes in the distribution of local motifs and the origin of β and γ phase transition. Using these results, we further develop a correlation of the distribution of local motifs to thermodynamic potentials, providing insights into the possible detection of local motifs in quantum materials. We note that such exploration has fundamental applications, as recent experimental observations for different quantum materials demonstrate the apparent existence of piezoelectric signals and second harmonic generation signals in nominally cubic phases,30 the magnetic response of paramagnetic quantum materials above the Néel temperature,31 the experimentally detected presence of spin short-range order in paramagnets,32,33 and symmetry breaking in high symmetry phases (even at high temperatures, which are expected to be symmetrical without the memory of local motifs present at low temperature)34,35 that all cannot be explained using the global average symmetric structure.
II. METHODS
All first-principles calculations were performed utilizing the Vienna Ab Initio Simulation Package (VASP) software36–39 with projector augmented-wave pseudopotentials40 and SCAN functional.41 SCAN has been shown to improve the prediction of properties over other exchange-correlation functionals [e.g., Perdew–Burke–Ernzerhof (PBE)42 and even DFT + U methods43,44] in the handling of symmetry breaking,45 internal structures,46,47 and energies46–48 for a diverse set of molecules and solids. Importantly, the utilization of this functionality allows for the reproduction of the insulating nature of a wide range of 3d oxides7,11,45 without utilizing any additional corrections. The cutoff energy was set to 500 and 550 eV for final relaxation and volume optimization/molecular dynamics calculations, respectively. Volume optimization began with a 1000 k-point grid per reciprocal atom and was increased to 10 000 for the final internal structure relaxation. The AIMD calculations utilized a Г point-only k-grid. The relaxation calculations were iteratively performed until the forces acting on each atom in the system were less than 0.01 eV/Å. The spin polymorphous model13 of the PM phase with a 160-atom cell was used to describe the phase, where spin distribution corresponds to a special quasi-random structure (SQS)49 in contrast to the starting point in many other works that consider paramagnet as a global average structure with both local and total magnetic moments of 0. To investigate the impact of temperature on structural characteristics, AIMD simulation using the NPT ensemble has been performed. In this ensemble, the number of atoms (N), pressure (P), and temperature (T) remain constant. The Langevin thermostat with Parrinello–Rahman barostat50,51 was employed to maintain both temperature and pressure at consistent levels. The time step was set to 1 fs. The Langevin thermostat was utilized to establish a friction coefficient of 3 ps−1 for the atomic degree of freedom of all atoms, while the lattice degree of freedom was set to 10 ps−1. These calculations were performed with a cutoff energy of 550 eV. It should be noted that the method discussed above does not allow accounting for the temperature-dependent distribution of spin short-range order and, hence, cannot be used, for instance, to describe accurately the phase transition from magnetically ordered AFM monoclinic to the paramagnetic orthorhombic phase of NdNiO3. However, in principle, with different degrees of accuracy, such phase transitions may be described using the superposition methods,52 structural dynamics,10 or utilizing dynamic mean field theory22,25,27 (which, for instance, is capable of describing electronic entropy with a higher degree of accuracy than DFT). The data analysis and visualization were done using pymatgen53 and Vesta software.54
III. RESULTS AND DISCUSSION
A. Distribution of local motifs in paramagnetic YNiO3 phases
The paramagnetic β phase exhibits the structural and spin symmetry breaking, with the formation of two different octahedra, resulting in the insulating state. This is well captured within the spin polymorphous theory,5–7 where the paramagnetic compound is considered within the spin SQS approximation.49 This model aligns with the high-temperature limit of the paramagnetic phase,10,52 where the total magnetization is effectively zero, and each magnetic species exhibits an uncorrelated random spin distribution. Utilizing this spin SQS framework and the large supercell size needed to converge materials properties (see discussion in Ref. 55), our SCAN calculations reveal a distribution of local NiO6 octahedral volumes, with the resulting system having an insulating nature and a bandgap energy of 0.62 eV. Specifically, we observe the coexistence of large and small NiO6 octahedra, measuring 10.61 ± 0.05 and 9.16 ± 0.05 Å3, respectively. Here, large and small octahedra are identified by the sublattice corresponding to the high (≈1.4 μB) and low (≈0 μB) spin states, respectively. We note that accounting for such distribution of local motifs is essential for describing the insulating nature of the β phase without accounting for any dynamic correlation effects, as discussed in our previous work.7 We recognize the feasibility of considering intermediate-temperature spin configurations. Yet, the spin SQS model, which is commonly adopted in physics, seems adept at effectively representing the core physics of such systems.10 It is also worth noting that in a spin polymorphous system, the orthorhombic YNiO3 phase tends to disproportion unless the temperature is factored in.5,56 This, thus, implies the properties of orthorhombic YNiO3 cannot be properly analyzed within the spin polymorphous theory that only relies on the minimization of internal energy. Indeed, this behavior is not surprising as at finite temperature, the behavior of a system is not defined anymore only by the internal energy minimization. Hence, the understanding of the β and γ phase transition does require post-static calculations.
To examine the impact of temperature on the distribution of local structural motifs, we performed AIMD simulations using the NPT ensemble. It is important to mention that an increasing number of studies have performed similar simulations for different types of materials,10,57–64 but critical input parameters are not explicitly defined in the literature (for instance, the distribution of local motifs is often discussed for extremely short AIMD simulations with a total simulation time of a couple of picoseconds or even fewer). Moreover, YNiO3 represents a unique set of compounds where both distributions of structural and spin motifs coexist simultaneously. For clarity, herein, our primary focus is to analyze structural motifs, which can be robustly studied experimentally using techniques like synchrotron x-ray diffraction (s-XRD)29 or local probe analysis of the pair distribution function.65,66 While local probe measurements for spin motif distribution are feasible,32,33 they are less common and necessitate an advanced theoretical framework to account for various magnetic sublattices. Hence, the distribution of local spin motifs will not be discussed in detail here.
When given atomic identities, composition, and structure are introduced to AIMD, the atomic velocities are usually randomly assigned in accordance with the corresponding temperature. This, thus, initializes a system that is not explicitly adjusted for applied external conditions (e.g., temperature, pressure) and interatomic or intermolecular interactions. Because of this, there is a finite time for the system to reach equilibrium. Traditionally, such equilibration is verified by monitoring the time dependence of internal energy, temperature, or lattice parameters (Figs. 1–3). However, for quantum materials, the distribution of local motifs becomes pivotal, and, hence, achieving its convergence is critical. We should note that the time needed for convergence is influenced by the system proximity to the transition point and, indeed, may require sufficiently long simulation time and the use of post-AIMD methods (e.g., active learning MD67–70). Because of this, the exact temperature for the phase transition estimated from the AIMD simulation should be taken with caution; the most significant part that one can learn here is the temperature dependence of the distribution of local motifs and its impact on material properties.
First-principles molecular dynamics of the PM phase of YNiO3 at 450 K using a spin polymorphous model. Time-dependent variation in (a) lattice constants, (b) total free energy, and (c) temperature during NPT molecular dynamics simulation at 450 K. (d) Distribution of local octahedral volumes at different moments of time.
First-principles molecular dynamics of the PM phase of YNiO3 at 450 K using a spin polymorphous model. Time-dependent variation in (a) lattice constants, (b) total free energy, and (c) temperature during NPT molecular dynamics simulation at 450 K. (d) Distribution of local octahedral volumes at different moments of time.
First-principles molecular dynamics of the PM phase of YNiO3 at 500 K using a spin polymorphous model. Time-dependent variation in (a) lattice constants, (b) total free energy, and (c) temperature during NPT molecular dynamics simulation at 500 K. (d) Distribution of local octahedral volumes at different moments of time.
First-principles molecular dynamics of the PM phase of YNiO3 at 500 K using a spin polymorphous model. Time-dependent variation in (a) lattice constants, (b) total free energy, and (c) temperature during NPT molecular dynamics simulation at 500 K. (d) Distribution of local octahedral volumes at different moments of time.
First-principles molecular dynamics of the PM phase of YNiO3 at 550 K using a spin polymorphous model. Time-dependent variation in (a) lattice constants, (b) total free energy, and (c) temperature during NPT molecular dynamics simulation at 550 K. (d) Distribution of local octahedral volumes at different moments of time.
First-principles molecular dynamics of the PM phase of YNiO3 at 550 K using a spin polymorphous model. Time-dependent variation in (a) lattice constants, (b) total free energy, and (c) temperature during NPT molecular dynamics simulation at 550 K. (d) Distribution of local octahedral volumes at different moments of time.
B. Temperature-dependent distribution of local motifs in paramagnetic YNiO3 phases
While at T = 0 K calculations, the distribution of local motifs is not very significant. When the temperature is introduced, the intrinsic tendency of the system (defined by the minimization of internal energy) can be significantly modified by the temperature contribution (e.g., driven by entropy increase). This is, indeed, reflected in the distribution of local structural motifs, for simplicity calculated here as volumes of NiO6 octahedra. Thus, at low temperatures, one can see the bimodal Gaussian distribution around two distinct types of octahedral volumes (i.e., small and large). For instance, at 200 K, the average volumes of large and small octahedra are 10.67 and 9.19 Å3, respectively. These results align with the static calculations for the β phase and indicate that at lower temperatures, the structure oscillates around the T = 0 K spin polymorphous configuration. This observation underscores the effectiveness of the spin SQS model we developed. It also confirms that the T = 0 K spin polymorphous configuration serves as a kernel for low-temperature AIMD simulations. As the temperature increases (Fig. 4), the two distinct octahedral peaks begin to shift slightly toward each other, and the distribution of local motifs becomes broader, with the two distributions starting to overlap. Eventually, at high temperatures (550 K), a broad-shifted Gaussian-like distribution of octahedral volumes emerges, and the average size of the octahedra is 10.03 Å3. This shift in distribution, especially in higher temperature regimes, provides a deeper understanding of the phase transition. Indeed, the transition from double-peak volume distribution to single peak volume distribution is responsible for the temperature-induced fall of the insulating state.5 Taking into account that the monoclinic phase has a distinct disproportionation, it also becomes clear that the transition from the bimodal to single Gaussian distribution is closely related to the monoclinic to orthorhombic structure phase transition, which is experimentally observed at approximately 582 K.29 We emphasize that all these results are obtained with mean-field DFT without accounting for any additional dynamic correlation effect or using U correction, demonstrating the temperature-dependent suppression of the structural disproportionation. It is also important to note that while different Hamiltonian models usually consider the existence of only a single local environment, the current results clearly demonstrate that at all temperatures, there is a wide distribution of local structural environments, which, however, noticeably changes with temperature.
Correlation of distribution of local motifs to the Gibbs free energy profile. (a) Distribution of local octahedra volumes for different temperatures extracted over 20 ps after the initial equilibration of 5 ps. (b) Relative Gibbs free energy profile computed based on the distribution of octahedral volumes. The results are presented based on the thermal equilibration of the monoclinic (β) PM phase (as starting configuration) modeled as a 160-atom supercell with a polymorphous spin distribution.
Correlation of distribution of local motifs to the Gibbs free energy profile. (a) Distribution of local octahedra volumes for different temperatures extracted over 20 ps after the initial equilibration of 5 ps. (b) Relative Gibbs free energy profile computed based on the distribution of octahedral volumes. The results are presented based on the thermal equilibration of the monoclinic (β) PM phase (as starting configuration) modeled as a 160-atom supercell with a polymorphous spin distribution.
C. Potential surface profile extracted from first-principles theory at finite temperature
Knowing the distribution of local structural motifs, one may wonder whether extracting the effective coordinate-dependent potential surface energy is possible. This is critical as it can potentially demonstrate the origin of the potential phase transition and, more generally, explain the role of temperature in the distribution of local motifs. We note that describing the presence of a double local environment leads to a double-well free energy potential.56 The relationship between statistical physics and thermodynamics allows for the computation of macroscopic properties like the Gibbs free energy from the microscopic behavior of a system. AIMD simulations provide a series of “snapshots” representing different microstates, offering a statistical framework for understanding these thermodynamic quantities. Since different snapshots can be considered non-interactive microstates of the systems, it is possible to apply statistical physics to extract energy profiles; specifically, the Gibbs free energy can be expressed as , where is a probability distribution of microstates and corresponds to the highest probability among all microstates. Assuming that each microstate of the system can be characterized by the distribution of octahedra volumes extracted over time (we note, however, that such an assumption requires the converged distribution of local motifs with time [Fig. 4(a)]), we calculated the relative Gibbs free energy profile shown in Fig. 4(b). Based on the computed results, one can see the following details. At low temperatures, there is a distinct double-well picture showing that the system tends to structurally disproportionation (i.e., form small and large NiO6 octahedra), which is consistent with that described by Mercy et al.56 It should be noted here that imperfection in the potential surface energy profile is mainly caused by the rather short sampling time (about 15 ps for each temperature after the initial equilibration) and may completely converge with increasing time. As temperature increases to about 500–550 K, there is the suppression of the double-well distribution, resulting in a single-well Gibbs free energy profile. These results, thus, demonstrate how the temperature can significantly modify the potential energy profile going from double well to single well, and, indeed, it is the temperature responsible for the suppression of spontaneous energy lowering disproportionation. The important part of these results is that the depth of the double-well energy profile slowly reduces with increasing temperature. This picture, thus, demonstrates that the phase transition in YNiO3 is, indeed, similar to that observed in some other materials.57 Given the temperature-dependent distribution of local motifs and comprehension of the Gibbs free energy profile, it is imperative to understand the novel insights derived from these data. First, this study illustrates the power of the spin polymorphous model in accurately representing the transition from the monoclinic to orthorhombic phase in paramagnetic YNiO3. Despite observing anticipated patterns at reduced temperatures, a crucial observation is that certain defining features dissipate earlier [e.g., already at low temperatures, experimentally observed s-XRD patterns29 show only the barely distinguishable peak splits (2θ ∼ 22.6°) used to identify the presence of monoclinic phase] than the suppression of the disproportionation temperature. Thus, concentrating on only one aspect might not provide an accurate representation of phase transition temperatures, especially if they are defined based on the suppression of structural disproportionation. This observation aligns with the experimental results, reinforcing the effectiveness of the spin polymorphous model in considering temperature-dependent factors. Second, we also observe that while heating, as the system approaches the transition temperature (we note, however, that the exact transition temperature may, indeed, depend on some computational factors such as the time of molecular dynamic equilibration, the density of k-points, etc.), some interesting physics can be observed. For instance, in Fig. 2, the time-resolved distribution of octahedral motifs fluctuates between single and double-peak octahedral volume distributions. This opens up interesting avenues in the physics of quantum materials close to the transition temperature. These observations demonstrate the complexity of phase transitions near the phase transition temperature. We note that experimental investigation of such complex symmetry breaking requires additional local probe analysis as has been done for NaMnO2,71 GeTe,34 BaTiO3,35 and Sr1−xNaxFe2As2,72 all of which, otherwise, cannot be easily detected with standard x-ray diffraction techniques. This further emphasizes the importance of using methods that can distinguish local symmetry breaking, such as pair distribution analysis.
IV. CONCLUSIONS
In summary, using a combination of density functional theory and ab initio molecular dynamics within the polymorphous framework, we provide a fundamental understanding of the distribution of local spin and structural motifs in paramagnetic YNiO3 phases, showing the existence of two types of octahedra (large and small) that play a key role in the insulating nature of this compound. Specifically, we focus on the temperature effect on the distribution of local motifs and provide comprehensive insights into the thermally induced transformation of these motifs and its correlation to the potential surface profile. For instance, we observe that the energy potential profile changes from a double-well at low temperatures to a single-well profile at around 500–550 K. These results, thus, show how temperature modifies the potential energy profile, reflecting the temperature-induced phase transition. Taking into account that the properties of the given system should be calculated as some superposition of properties of local motifs (and not as properties of global average structure), we, thus, highlight the need to account for the temperature-induced distribution of local motifs in describing properties of correlated materials. Moreover, we also reaffirmed the limitations of relying solely on x-ray diffraction peaks to determine phase transitions. Our findings hint that considering the degree of structural disproportionation could offer a more precise indication of a phase transition. However, despite the promising results, the study recognizes that the simulations and models utilized have limitations and may not precisely capture the intricate behavior of these systems.
ACKNOWLEDGMENTS
The authors acknowledge the “ENSEMBLE3—Centre of Excellence for nanophotonics, advanced materials and novel crystal growth-based technologies” project (GA No. MAB/2020/14) carried out within the International Research Agendas program of the Foundation for Polish Science co-financed by the European Union under the European Regional Development Fund and the European Union's Horizon 2020 research and innovation program Teaming for Excellence (GA No. 857543) for the support of this work. We gratefully acknowledge Poland’s high-performance computing infrastructure PLGrid (HPC Centers: ACK Cyfronet AGH) for providing computer facilities and support within computational Grant No. PLG/2023/016228. The authors acknowledge fruitful discussions with Professor Andriy Gusak.
AUTHOR DECLARATIONS
Conflict of Interest
The authors have no conflicts to disclose.
Author Contributions
Himanshu Joshi: Data curation (lead); Formal analysis (equal); Investigation (lead); Methodology (equal); Validation (equal); Writing – original draft (supporting). Mateusz Wlazło: Data curation (supporting); Formal analysis (equal); Investigation (supporting); Methodology (equal); Validation (supporting); Writing – original draft (supporting). Harshan Reddy Gopidi: Data curation (supporting); Formal analysis (supporting); Investigation (supporting); Methodology (equal); Validation (supporting); Writing – original draft (supporting). Oleksandr I. Malyi: Data curation (supporting); Formal analysis (equal); Investigation (supporting); Project administration (lead); Supervision (lead); Writing – original draft (lead).
DATA AVAILABILITY
The data that support the findings of this study are available from the corresponding author upon reasonable request.