Subterahertz dielectric relaxation in lead-free Ba(Zr,Ti)O3 relaxor ferroelectrics

Relaxors are complex materials with unusual properties that have been puzzling the scientific community since their discovery. The main characteristic of relaxors, that is, their dielectric relaxation, remains unclear and is still under debate. The difficulty to conduct measurements at frequencies ranging from ≃1 GHz to ≃1 THz and the challenge of developing models to capture their complex dynamical responses are among the reasons for such a situation. Here, we report first-principles-based molecular dynamic simulations of lead-free Ba(Zr0.5Ti0.5)O3, which allows us to obtain its subterahertz dynamics. This approach reproduces the striking characteristics of relaxors including the dielectric relaxation, the constant-loss behaviour, the diffuse maximum in the temperature dependence of susceptibility, the substantial widening of dielectric spectrum on cooling and the resulting Vogel–Fulcher law. The simulations further relate such features to the decomposed dielectric responses, each associated with its own polarization mechanism, therefore, enhancing the current understanding of relaxor behaviour.

(HN) relaxation function, which determine the shape of relaxation spectrum at frequencies below and above ν HN , respectively. Both of these coefficients are found to deviate from unity for temperatures below~250 K, which indicates that this relaxation differs from a pure Debye mechanism.
One can also see that these two parameters exhibit some fluctuations with temperature, which are most probably due to the fact that the fitting model described by Eq. (1)  as well as the coupling constant, δ, inherent to the COR mechanism. We have found that, for T < 70 K, the use of a heavily damped harmonic oscillator term instead of the COR term in Eq.
(1) of the manuscript leads to better fitting. The change from COR to DHO model results in jumps in the temperature dependences of the related parameters at 70 K as observed in Supplementary   Figure 1 (note that the DHO model is equivalent to the COR model with δ = 0).
Furthermore, Supplementary Figure 1(d) shows the resonance frequencies of these DHOs and COR. Similar temperature behaviors were previously reported in Ref. [3] for the resonant frequencies of the Zr-and Ti-related DHOs (ν 1 , ν 2 , respectively), while the present addition of the COR mechanism (via its ν CO and ν CR frequencies) is found necessary here to further improve the fittings of the MD data. 4 Finally, Supplementary Figure 1(e) shows separate contributions of the HN relaxation, χ R0 , and of the COR mechanism, χ COR to the total static susceptibility. It indicates that χ R0 vanishes abovẽ 350 K, and narrowly peaks around 100 K. On the other hand, the temperature dependent peak of χ COR is rather broad, centering around~300 K and χ CO is still significant for temperatures as high as 1000K. Note, however, and as cautioned in the manuscript, that the HN and COR contributions significantly overlap, especially at comparatively high temperatures. Consequently, some of their parameters can only be estimated rather than uniquely determined, which is also consistent with the fluctuations of (1 − α) γ and ν CR in Supplementary Figure 1.  [4] which predicts that the T VF parameter of the imaginary part is smaller than that of the real part. Note also that T VF can be negative in this theory as it is not a temperature that has physical significance; instead, it is simply a parameter that is determined by the temperature dependence of the dielectric strength and the parameters of the relaxation spectrum.
Let us now discuss in more detail the relation of our results to the published experimental data on Ba(Zr x Ti 1−x )O 3 system. While no high-frequency measurements of the samples of exactly the same composition (x = 0.5) we theoretically investigate have been reported, the comparison with the available data of similar compositions (x = 0.4 and x = 0.6) is meaningful, as it is known that the characteristic relaxation frequencies and their temperature behaviour are similar for all the compositions with 0.2 ≤ x ≤ 0.8 (see Figure 2 of Ref. [7] which shows the experimental results obtained for BZT ceramics). Significant differences only concern the magnitude and temperature variation of the static susceptibility [5,6]. In Supplementary Figure 3, the available experimental spectra of BZT ceramics with x = 0.4 are compared to our MD predictions.
Similar to the MD results, two relaxation-related excitations were revealed in the THz and sub-THz ranges experimentally. The first one is an excitation with the characteristic frequency of ∼ 3 THz (∼ 100 cm −1 ) which is practically temperature-independent (upper curve in Figure 2 of Ref. [7]). It was successfully simulated with the coupled oscillator-relaxator (COR) at T > 70 K and the overdamped harmonic oscillator at lower temperatures (their frequencies ν CR and ν CO are shown in Supplementary Figure 1(d) ). Another relaxation with a frequency of~100 GHz, which we found to dominate in the subterahertz range for T < T * , has also been observed experimentally.
The characteristic frequency of this relaxation has not been determined accurately because of insufficiency of available experimental data points, but it has been found to remain within the range of 10 -100 GHz at all studied temperatures (see Figure 2 of Ref. [7]). This experimental result is in agreement with our prediction that reveals the Arrhenius law at T > 70 K and the saturation of ν m at lower temperatures (see Figure 3 (c) in the main text). The predicted dielectric strength (contribution to the static susceptibility) of this mode, χ R0 , passes through a maximum upon cooling and practically vanishes at temperatures close to zero (see Supplementary Figure   1(e)). Similar behaviour has also been observed in experiments: as shown in Supplementary Figure   3(e), almost the same susceptibility is measured at 1 THz and at 1 MHz (i.e., well below and well above the experimentally determined characteristic frequency, respectively), which means that the relaxation strength at 10 K is much smaller than the relaxation strength at the temperature for which the real part of the static dielectric response is maximum.

Supplementary Note 3: Detailed comparison of MD simulation results with experiment
As shown in Supplementary Figure 3, in the THz and subterahertz parts of the spectrum, the results obtained from MD simulations agree well with experiments, not only for the characteristic frequencies, but also for the magnitude of the susceptibility. The most significant discrepancy is found for the magnitude at 10 K (85 and 60 from experiments and simulations, respectively).
However, at ν ∼ 1 GHz (close to the lowest possible frequency available via MD simulations), the situation is more complex. While at 150 K the agreement is fairly good (cf. Supplementary   Figure 3 (c)), at 300 K (cf. Supplementary Figure 3 (a)) the measured susceptibility is twice as large as the simulated one. The large increase of the measured susceptibility observed below ∼ 10 GHz, which is not seen in simulations, is evidently due to the third relaxation mode having a characteristic frequency of ∼ 1 GHz at 200-300 K (cf. Figure 2 of Ref. [7]). Upon cooling this mode moves to lower frequencies and disappears from the frequency window available in our simulation. We conclude, therefore, that this relaxation mode, which dominates the BZT's spectrum for frequencies lower than 1 GHz according to experiments, cannot be reproduced in our 6