History-dependent and frequency-dependent dielectric nonlinearities induced by polar nanoregions in a thin film of 0.5 ( Ba 0.7 Ca 0.3 TiO 3 ) – 0.5 ( BaZr 0.2 Ti 0.8 O 3 )

Kevin Nadaud, Guillaume Nataf, Nazir Jaber, Edgar Chaslin, Béatrice
Negulescu, Jérôme Wolfman

► To cite this version:

Kevin Nadaud, Guillaume Nataf, Nazir Jaber, Edgar Chaslin, Béatrice Negulescu, et al.. History-dependent and frequency-dependent dielectric nonlinearities induced by polar nanoregions in a thin film of 0.5 ( Ba 0.7 Ca 0.3 TiO 3 ) – 0.5 ( BaZr 0.2 Ti 0.8 O 3 ). Physical Review B, 2025, 112 (1), pp.014202. 10.1103/fhdr-wn22. hal-05149785

HAL Id: hal-05149785

https://hal.science/hal-05149785v1

Submitted on 10 Jul 2025

**HAL** is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
History- and Frequency-Dependent Dielectric Nonlinearities
Induced by Polar Nanoregions in
0.5 (Ba<sub>0.7</sub>Ca<sub>0.3</sub>TiO<sub>3</sub>)-0.5 (BaZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub>) Thin Films

Kevin Nadaud,\* Guillaume F. Nataf, Nazir Jaber, Edgar

Chaslin, Béatrice Negulescu, and Jérôme Wolfman

GREMAN UMR 7347, Université de Tours, CNRS,
INSA-CVL, 16 rue Pierre et Marie Curie, 37071 Tours, France
Abstract

In this article, dielectric non-linearities in 0.5 (Ba<sub>0.7</sub>Ca<sub>0.3</sub>TiO<sub>3</sub>)–0.5 (BaZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub>) thin film
are studied using impedance spectroscopy and harmonic measurements, as a function of the AC
measuring field, at different frequencies and upon cycling. The measurements reveal that the
pinching of the hysteresis loop, characterized by a phase angle of the third harmonic close to –270°,
is stronger for low frequencies. This confirms that the pinching is induced by the presence of polar
nanoregions (PNRs), whose responses are also strongly dependent on frequency. When repeating
the measurement, the PNR contribution changes since the phase angle of the third harmonic
response evolves from pinched to conventional relaxor. This shows that strong changes in the PNR
configuration can be induced, even for low AC fields. First Order Reversal Curves (FORC) confirm
the presence of a pinched hysteresis loop. When repeating the FORC measurement a second time,
the distribution drastically changes and corresponds to a soft ferroelectric. The asymmetry of the
Preisach plane measured using FORC is confirmed by a proposed novel measurement strategy:
unipolar impedance measurements.

# I. INTRODUCTION

The high dielectric permittivity and piezoelectric coefficient of relaxor ferroelectrics make them promising materials for energy storage and actuator applications [1–4]. Relaxors exhibit polar nanoregions (PNRs), clusters in which the polarization is randomly aligned in absence of external electric field [5–8]. For relaxors and conventional ferroelectrics, the polar cluster boundaries, PNRs and ferroelectric domain walls (DWs) respectively, contribute significantly to the macroscopic dielectric and piezoelectric responses [9–12]. Those PNRs can subsist several kelvins above the Curie temperature (in the centrosymmetric phase) [13–15], which results in a dielectric permittivity peak much broader in temperature than that of conventional ferroelectric [6] and frequency dependent [16–18].

Due to their non-linear response to the applied electric field, ferroelectrics exhibit specific
dielectric responses, even below the coercive field. The measurement of the dielectric per-
mittivity as function of the AC field can be used to determine the reversible and irreversible
domain wall motion contributions. According to the Rayleigh law, the relative permittivity

<footnote>* Author to whom correspondence should be addressed: kevin.nadaud@univ-tours.fr</footnote>
linearly increases: [19, 20]

$$ \varepsilon_r = \varepsilon_{r-l} + \alpha_r E_{AC} \quad (1) $$

$\varepsilon_{r-l}$ correspond to the low field permittivity and $\alpha_r$ to the slope and represents the irreversible domain wall motion contribution. The Rayleigh law assumes that domain walls interact with a homogeneous distribution of the pinning centers.

However, in real ferroelectric materials, the distribution of pinning centers is not homogeneous and the linear increase is only visible after a given threshold [21–23]. For low field, domain walls can only vibrate around an equilibrium position and the contribution is called reversible. In that case, the relative permittivity evolution with the measuring field can be described using the hyperbolic law: [24–26]

$$ \varepsilon_r = \varepsilon_{r-l} + \sqrt{\varepsilon_{r-rev}^2 + (\alpha_r E_{AC})^2} \quad (2) $$

with $\varepsilon_{r-rev}$ the reversible domain wall motion contribution, proportional to the domain wall density [27–29].

In addition to the conventional impedance spectroscopy, which considers only the first harmonic of the polarization, it is possible to study higher order harmonics. When the applied electric field is $E(t) = E_0 \sin(\omega t)$, the polarization can be expressed as : [30]

$$ P(t) = P^{(0)} + \sum_{n=1}^{\infty} P^{(n)} \sin(n\omega t + \delta_n) \quad (3) $$

$$ = P^{(0)} + \sum_{n=1}^{\infty} P'^{(n)} \sin(n\omega t) + P''^{(n)} \cos(n\omega t) \quad (4) $$

$P^{(0)}$ is an eventual DC contribution, $P^{(n)}$ is the magnitude of the polarization of the n-th order and $\delta_n$ is its phase. $P'^{(n)}$ and $P''^{(n)}$ represent the in-phase and quadrature component of the polarization of the n-th order. Depending on the analysis, either magnitude/phase or in-phase/quadrature representations are used. With that decomposition, we can define the n-th order relative real and imaginary part of the permittivity:

$$ \varepsilon_r'^{(n)} = \frac{P'^{(n)}}{E_0} \quad \text{and} \quad \varepsilon_r''^{(n)} = \frac{P''^{(n)}}{E_0} \quad (5) $$

For ferroelectrics, the phase angle of the third harmonics is particularly studied [24, 31, 32]. For an ideal ferroelectric, having domain walls that interact with a homogeneous distribution of pinning centers, a phase angle of the third harmonic $\delta_3 = -90^{\circ}$ is found as stated by the Rayleigh law [33]. For many ferroelectric materials, the phase angle of the
third harmonic differs, at least in some range of AC field, from the ideal value of -90°. A deviation from this angle gives information regarding the deformation of the polarization loop.

(i) $\delta_3 = 0°$, the third harmonic is of opposite sign with respect to the first harmonic, when the field is maximum. This corresponds to saturation-like response (Fig. 1a). A saturation like response is found at large electric field for ferroelectrics and relaxors.

(ii) $\delta_3 = -90°$, the third harmonic opens the center of the loop. This corresponds to the ideal Rayleigh behavior (Fig. 1b), which is found for $BaTiO_3$ [33] and weakly doped (Pb, Zr)TiO$_3$ [34].

(iii) $\delta_3 = -180°$, the third harmonic is of same sign with respect to the first harmonic, when the field is maximum. This corresponds to a divergent-like response (Fig. 1c). This response is found at low field for almost all ferroelectrics and relaxors.

(iv) $\delta_3 = -270°$, the third harmonic constricts the center of the loop. This corresponds to the pinching of the loop (Fig. 1d). Pinching of the loop can occur for fresh state for hard ferroelectrics [31] or relaxor above global phase transition temperature [24].

For relaxors, the non linearities are also studied but using the susceptibility decomposition of the polarization:

$$P = \varepsilon_0 (\chi_1 E + \chi_2 E^2 + \chi_3 E^3 + \dots) \qquad (6)$$

where $\chi_i$ is the i-th order susceptibility.

The third order susceptibility is widely studied but often it is only its real part [16, 35-38], since for many models dealing with relaxors, the imaginary part is null [39]. The sign of this contribution can be used to determine the framework explaining the PNR contribution to the non linearities [16, 40-42].

$\chi_i$ can be obtained from the measurement, by applying a sinusoidal electric field, and are a linear combination of $\varepsilon_r^{(n)}$, with $n \ge i$ [36, 39]. For example, we have for the third order susceptibility: [36, 39]

$$\chi_3 = -\frac{4}{E_0^2} (\varepsilon_r^{(3)} + 5\varepsilon_r^{(5)} + 14\varepsilon_r^{(7)}) \qquad (7)$$

Even if the consideration of the higher order harmonics allows gaining a higher precision on the extraction of the $\chi_3$ value [36], it can lead to a significant increase of the noise [16].
<figure><img src="image_2.png" /><figcaption>Figure 1. Hysteresis deformation due to third-harmonic response shown for four limiting cases: $\delta_3 = 0°$ (a), $\delta_3 = -90°$ (b), $\delta_3 = -180°$ (c) and $\delta_3 = -270°$ (d). Dashed curves correspond to the contribution of the first harmonic and solid black curve to the contribution of the third harmonic. Gray arrows indicate the rotational sense of the third harmonic and colored arrows the deformation of the hysteresis. Adapted from [30, 32].</figcaption></figure>

It is thus possible to consider as a first approximation (especially for weak non-linearity at low driving field $E_0$):

$$ \chi_3 = -\frac{4}{E_0^2} \varepsilon_r^{(3)} \qquad (8) $$

The dynamics of the PNRs can be thus obtained by studying the sign of $\varepsilon_r^{(3)}$, more specifically if $\delta_3$ is close to 0° or -180°. In relaxor ferroelectrics, the total dielectric response is the sum of different contributions (PNRs, ferroelectric domain walls, polar boundaries, etc.) and depending on the measurement conditions, frequency or magnitude of the AC field, temperature, their respective contribution may change [39]. The decorrelation of the respective contributions is thus made by measuring the dielectric properties at different magnitude and frequency of the measuring AC field.

In addition to impedance spectroscopy, First Order Reversal Curve (FORC) measurements, which allow the extraction of the Preisach distribution [43–45] can be used. Preisach model has been initially developed for the description of the ferromagnetic hysteresis loop [46, 47] and has been extended to ferroelectric materials [43–45, 48]. In this model, the $P(E)$ hysteresis loop can be seen as the sum of a large number of elementary hysteresis loop called
hysterons. Each hysteron has given values of of up-switching field *E* and down switching *E*<sub>r</sub>. The knowledge of the hysteron decomposition of the *P*(*E*) cycle gives information on the fatigue and polarization of the material [43, 49–51]. It also allows the decomposition of a material in different phases [52–55]. Depending on the shape of the hysteresis loop, peculiar FORC distribution can be found [56].

In this article we study the effect of the frequency of the AC measuring signal and the number of measurement cycles on the dielectric response in BCTZ relaxor thin film. The objective is to show that the pinching response of the studied material is strongly correlated with the presence of PNRs. In addition, we show that the pinching, and thus the presence of PNR, is strongly influenced by the sub-coercive applied electric field used for the measurement. FORC measurements are also used to study the pinching response and a novel unipolar impedance measurement method is proposed to probe specific portions of the Preisach plane.

## II. EXPERIMENTS

Polycristalline 380 nm thick BCTZ film has been grown on Pt/TiO<sub>2</sub>/MgO substrate by pulsed laser deposition. Details on growth conditions can be found elsewhere [24, 57]. Top circular Au/Ti electrodes (150 µm radius) were deposited through a shadow mask. The metal-insulator-metal topology has been chosen for the simple extraction of the dielectric properties using the parallel plate capacitor formula. This composition is close to the morphotropic phase boundary, and a conventional ferroelectric behavior is found when the material is in ceramic form [58, 59]. In the present case, a relaxor behavior has been observed (see Fig. S1 in Supplementary material [60], see also references [7, 10, 61–70] therein).

The dielectric characterizations presented in this article have been acquired using a lock-in amplifier (MFLI with MD option, Zurich Instrument) connected to a probe station (Summit 12000, Cascade Microtech). For the impedance and harmonics analysis, the AC measuring signal has been generated using the embedded lock-in amplifier generator. Its amplitude has been swept from 10 mV<sub>rms</sub> to 1 V<sub>rms</sub> at a frequency from 31 Hz to 100 kHz. The applied voltage and current through the capacitor are measured by the lock-in amplifier and are demodulated simultaneously. The capacitance and the harmonics measurements are described elsewhere [24, 55, 71, 72]. Dielectric permittivity is obtained using parallel plate
formula. Since measurement of the third-harmonic phase-angle can be tricky, careful setting of the lock-in amplifier has been made. Filters with an order of 8 have been selected. The bandwidth of the low-pass filter has been set to have a minimum of 80 dB of rejection at the demodulated frequency, corresponding to a factor 10 between the excitation frequency and the –3 dB cut-off frequency of the filter. To further increase the sensitivity of the measurement a minimum of $15\tau/0.1$ s for the averaging is considered.

FORC measurement have been acquired using the same lock-in amplifier, but with a de-modulation frequency of 0 Hz, which allows having the time evolution of the applied voltage and the resulting current. The applied voltage has been generated using an Agilent 33250A Waveform Generator. Numerical integration has been done to obtain the polarization. To avoid poling effects inevitable with conventional FORC measurement, starting from a satu-ration state [43–45], we use a bipolar waveform [55, 73, 74], allowing to access of the FORC distribution of the fresh state. This method is then called symmetric FORC. The detailed procedure of the FORC distribution extraction is given in refs. [55, 73]. For the FORC distributions, $E$ is the switching field for increasing field and $E_r$ is the switching field for the decreasing one. Unless otherwise specified, all the characterization have been conducted at 300 K.

### III. RESULTS AND DISCUSSION

The first part is dedicated to measurements as function of frequency and allows us to identify the frequency dependency of the pinching behavior. The second part contains the study as function of the number of measurements, and thus its effect on the non-linear response. In the third part, FORC measurements are used to gain a deeper insight into the pinching response. Finally, unipolar impedance measurement is introduced as a novel method to probe only a desired part of the FORC distribution at low field.

#### A. Influence of measurement frequency

Fig. 2 shows the first and third harmonics contribution to the permittivity, as a function of the applied AC electric field, for different frequencies. Measurement is slightly noisier for low frequencies, which is common [16, 24, 38]. This comes from the capacitive behavior of the
<figure><img src="image_3.png" /><figcaption>Figure 2. Real part of the first harmonic contribution to the relative permittivity (a), phase-angle (b), real part (c) and imaginary part (d) of the third harmonic contribution as a function of the AC measuring field. Measurement has been made by decreasing the frequency, i.e. from 100 kHz to 31 Hz, to avoid misinterpretation due to cycling effect discussed afterward. Dashed line corresponds to the measurement at 100 kHz made after all the other measurements, to verify that the evolution of the shape of the curve is not related to cycling but a true frequency effect.</figcaption></figure>

measured device, which implies a decreasing current value when the measurement frequency decreases. When the measuring field increases, the relative permittivity (Fig. 2a) increases, corresponding to irreversible contribution from the motion of polar boundaries (domain walls pinning/unpinning, polar cluster boundaries, or phase boundaries) [20, 24, 33]. When the frequency change, a very similar shape of the curve is observed. The only difference consists in a decrease of the permittivity, when the frequency increases, the slope versus electric field staying constant.
The influence of the frequency on the phase-angle of the third harmonic is visible in Fig. 2b. For frequency between 10 kHz and 100 kHz, $\delta_3$ goes from -180°, divergent like contribution, to -20°, saturation like contribution, and passing by -90°. This kind of response is common for relaxor material [8, 24, 32, 33].

For a frequency of 3 kHz, the limit values at low and high fields are the same, but the transition is very steep. For frequency below 1 kHz, the shape of the curve is different. Instead of increasing after 8 kV cm⁻¹, the phase-angle decreases and does not pass through -90° but to -270°/+90° which is characteristic of a pinching of the hysteresis loop [31].

The pinching of the hysteresis loop can correspond to (i) an electric field-induced transition from a paraelectric to a ferroelectric phase, as in BaTiO₃ [75], (1-x)Bi₁/₂Na₁/₂TiO₃ - x BaTiO₃ [76], (ii) a fresh state in PbZr₀.₅₈Ti₀.₄₂O₃ [31] or (iii) an aged state and oxygen vacancies migration in Cu-doped in BaTiO₃ [77], Ce-doped Ba(Ti₀.₉₉Mn₀.₀₁)O₃ [78]. Here, the AC field is low and thus the sample aging and charge migration are not expected. Another possibility is the presence of PNRs that can be oriented by the application of an AC electric field during the measurement [24]. In the present case, the pinching is well visible even below the transition temperature (T_m = 330 K here) and is more visible for low frequencies. This frequency sensitivity supports the idea that the observed pinching signature ($\delta_3$ ≃ -270° or +90°) may be related to the PNRs since their response is very frequency dependent [16-18].

To verify that the evolution of the $\delta_3$ curve is not related to cycling but a true frequency effect, a measurement at 100 kHz has been performed after all the measurements presented in Fig. 2. The shape of the phase-angle 3rd harmonic is globally the same as the one presented in Fig. 2b, confirming the frequency influence on the pinching effect. The difference is attributed to a slight cycling effect, detailed later on.

It is also possible to analyze the third harmonic response by studying it in-phase and out-of-phase, respectively real and imaginary parts (Figs. 2c and 2d). When the AC field increases, the sign of the real part of the third harmonic contribution changes, when the AC field becomes sufficiently high to overcome the quenched state of the PNR[16, 79]. When the frequency increases, the real part in the range [3 kV cm⁻¹, 20 kV cm⁻¹] significantly decreases in absolute value, diminution of 50%. This stronger frequency decay for the third harmonic than for the first one is common for relaxors in which the nonlinearities have for origin PNRs [16-18].

For relaxors in ceramic form, the imaginary part is very small compared to the real part
[8]. In the present case the imaginary part is in the same order of magnitude as the real part, which suggests that the total third harmonic response is a sum of different contributions i.e. PNRs and domain boundaries contributions. For $f = 31$ kHz, the imaginary part is negative and its absolute value increases linearly between $4 \text{ kV cm}^{-1}$ and $14 \text{ kV cm}^{-1}$, which corresponds to conventional Rayleigh response (i.e irreversible domain wall motion contribution). The almost constant value and close to zero for fields below $4 \text{ kV cm}^{-1}$ corresponds to the low field region in which the polar boundaries contribution is mostly reversible [22, 80]. Above $24 \text{ kV cm}^{-1}$, the imaginary part decreases, which reveals that the domain wall motion contribution no longer follows the Rayleigh law, in agreement with the saturation observed on the first harmonic permittivity (Fig. 2a). When the frequency decreases, the imaginary part at AC fields higher than $25 \text{ kV cm}^{-1}$ progressively goes from negative values to positive values. If only domain wall motion would be present, the imaginary part should stay negative. The progressive change of sign indicates that another contribution to the permittivity (PNRs) progressively takes place and induces the pinching visible for low frequencies. For low AC field (below $10 \text{ kV cm}^{-1}$), the evolution with frequency is almost not visible.

Fitting of the $ε'_r(E_{AC})$ curves (Fig. 2a) allows the extraction of the polar boundaries contributions to the permittivity as function of the frequency, reported in Supplementary Material [60] (see also references [80–89] therein and Fig. S2).

## B. Influence of the number of measurements (cycling)

Fig. 3 shows the first and third harmonics contribution to the permittivity as a function of the measuring electric field, for successive measurements. The relative permittivity exhibits a very small increase after the first measurement and for the following measurements, the increase is almost imperceptible (Fig. 3a). By difference, the third harmonic response shows a clear evolution as a function of cycling. For the first 5 measurements, the phase-angle of the third harmonic exhibits the pinching response which has been seen in the previous study concerning the frequency effect (Fig. 3b). After the 6th measurement, the pinching signature ($\delta_3 = -270°$) is no more visible and the phase-angle response corresponds more to what is obtained for relaxor ferroelectrics ($-180° \rightarrow -90° \rightarrow 0°$). When increasing the number of measurements, the phase response shape slightly evolves: the transition is less steep and the field at which $\delta_3 = -90°$ decreases. The phase-angle of the third-harmonic
<figure><img src="image_4.png" /><figcaption>Figure 3. Real part of the permittivity (a), the phase-angle (b), the real part (c) and the imaginary part (d) of the third harmonics, as a function of the measurement electric field, for different number of measurements, for $f = 1$ kHz.</figcaption></figure>

at low field exhibits a plateau between 4 kV cm⁻¹ and 20 kV cm⁻¹, which value increases from −180° to −145°. This later value is close to what is found for the ferroelectric-relaxor (RE-FE) category given in ref [8].

The existence of a pinched hysteresis loop upon the first field cycles results from the following dipole orientation sequence. At zero field, dipoles are mostly randomly oriented (Fig. 4a), leading to a small but non-zero polarization. As the electric field increases in magnitude, dipoles orient in the same direction (Fig. 4b, like for a conventional ferroelectric) and a large polarization is present. When the field magnitude is decreased, large backswitching occurs (Fig. 4a) and contrary to conventional ferroelectric, a strong decrease of the polarization is visible. This leads to low value of the polarization around zero field, and thus
<figure><img src="image_5.png" /><figcaption>Figure 4. Illustration of the PNR configuration for a fresh sample (a), for a sample under field (b) and for a cycled sample (c).</figcaption></figure>

the pinching/constricted loop [31, 78]. However, when repeating the measurement multiple
times, less and less backswitching occurs. Most PNRs remain oriented in the same direction
and thus applying again the AC field does not lead to a large re-orientation i.e. the pinching
is fading upon cycling (Fig. 4c).

The effect of the number of measurements performed on the real and imaginary parts of
the third harmonic contribution, visible in Figs. 3c and 3d, is very similar to the frequency
effect (Figs. 2c and 2d). When repeating the measurement, the real part at low AC field
decreases, in absolute value, and the field for which the sign changes is lower. This confirms
that after multiple measurements, the applied AC electric field leads to a lower response
since PNRs have been oriented in the same direction. The value at high AC field slightly
increases, corresponding again to a saturation, (δ₃ ≃ 0°) of the permittivity, which is more
visible on the first harmonic after a large number of measurement (Fig. 3a). The evolution
of the pinching when measuring multiple times the material indicates that even at fields
lower than the coercive field the domain/phase structure can be influenced.

Fitting of the $\epsilon'_r(E_{AC})$ curves (Fig. 3a) allows the extraction of the polar boundaries
contributions to the permittivity as function of the cycle number, reported in Supplementary
information [60] (see Fig. S3 and refs. [25, 29, 81, 90–92] therein).

C. First Order Reversal Curves and unipolar Rayleigh measurements

FORC analysis has been performed to estimate local internal field $E_i$ / coercive field $E_c$
distribution and understand the effect of successive measurements on the material response.
$P(E)/I(E)$ loops used for the FORC extraction are visible in Supplementary material [60]
(see Figs. S5 and S6 and references [24, 33, 77, 93] therein.). The distribution given by
first symmetric FORC measurement is presented in Fig. 5a. Close to the axis $E = E_r$, at $(E, E_r) = (-30 \text{ kV cm}^{-1}, -30 \text{ kV cm}^{-1})$ (noted R), a peak is visible corresponding to a large reversible contribution, typical from soft ferroelectric/relaxors [94]. Moreover, it exhibits two peaks (noted P1 and P2), similarly to what is observed for unpoled PZT [73] or pristine state $Hf_{0.5}Zr_{0.5}O_2$ [95–97], and these two peaks are typical of the pinched hysteresis loop [56]. The peak positions are almost symmetric $(E, E_r) = (69 \text{ kV cm}^{-1}, -17 \text{ kV cm}^{-1})$ and $(6.2 \text{ kV cm}^{-1}, -66 \text{ kV cm}^{-1})$. These two peaks in the FORC distribution corresponds to the current peaks visible in Fig. S6(a) at low field. The peak separation is not visible for large field measurements due to the orientation of the PNR, similar to what is found for the cycling study Figs. 3 and 4.

When repeating the symmetric FORC measurement a second time (Fig. 5b), the distribution still presents two peaks at the same positions but their intensity is greatly reduced, especially the peak in the upper right corner (at $6.2 \text{ kV cm}^{-1}, -66 \text{ kV cm}^{-1}$). The hysteron decomposition obtained for the second measurement is thus closer to what is expected for a soft ferroelectric [94], a maximum of hysteron density close to the axis $E_r = E$. Fig. 5e shows the difference between the first and the second symmetric FORC measurements. In addition to the large decrease of the main peaks, the hysteron density close to the origin increases.

To quantify the contribution of each peak and of the reversible part to the total polarization, FORC distribution has been integrated using the following expression: [43]

$$P = \iint p(E, E_r)dEdE_r + \int p_{rev}(E)dE \quad (9)$$

where $p(E, E_r)$ and $p_{rev}(E)$ are the Preisach distribution for irreversible ($E > E_r$) and reversible hysterons ($E = E_r$), respectively. The results of the integration are given in Table 5f. For the first and second measurements, the reversible contribution is 8.4 and 16 times larger than the sum of the 2 peaks related to PNR's, respectively. This explains the reason why the pinching is not visible on the $P(E)$ loops (Fig. S5), enhancing the interest of the third harmonic measurement and the FORC to reveal a pinching contribution to the polarization.

In addition, FORC measurements have been done by starting from a saturated state (Figs. 5c and 5d). In both cases, the distribution presents only one peak, very close to the $E_r = E$ axis, which corresponds to what is expected for a soft ferroelectric [94]. Depending
<figure><img src="image_6.png" /><figcaption>Figure 5. FORC distribution for the first (a) and second (b) symmetric measurements and for positive saturation (c) and negative saturation (d) measurements. Difference between the symmetric measurements (e). Solid white line in (e) indicates the frontier between positive and negative values. Black dashed lines correspond to the area which is probe for symmetric AC measurements. Red lines in (e) indicates the area probe for unipolar impedance measurement. (f) Contribution to the polarization of the reversible contribution (R) and the two peaks (P1) and (P2).</figcaption></figure>

on the initial state (positive or negative saturation), a slight shift of the peak is visible. This peak shift indicates an imprint is present [98], which may arise from a different domain structure for positive and negative bias field. The peak position is very close to the one subsists while doing bipolar FORC measurement. The presence of only one peak shows that the application of a large electric field remove the pinching response of the material and the response correspond to a conventional soft ferroelectric.

A third symmetric FORC measurement has been made after the two starting from saturation. No significant changes compared to the second symmetric FORC is visible, the distribution still exhibits two peaks, but with a magnitude lower than for the first measurement (given in supplementary material [60], Fig. S7).

For an ideal material following the Rayleigh law, a homogeneous distribution of hysteron is
<figure><img src="image_7.png" /></figure>

<figure><img src="image_8.png" /></figure>

Figure 6. Preisach plane areas probed when doing conventional AC measurement (a) and using unipolar AC signal (b). Cyan area corresponds to the portion of the plane which is probed using positive signal and orange to the negative. Real part of the permittivity (c) and phase-angle of the third harmonic (d) measured using the unipolar AC signal given Eq. (12) (colored curves) and using conventional impedance spectroscopy (black lines). Gray dash-dotted line in (d) corresponds to phase angle for a homogeneous hysteron density.

present at low fields [99]. When doing conventional AC measurement, the part of the FORC distribution that is probed is symmetric according to the axis $E_r = -E$, and corresponds to
the triangle $\Delta^{sym}$ defined by the vertices $(-E_{AC}; -E_{AC})$, $(+E_{AC}; -E_{AC})$ and $(+E_{AC}; +E_{AC})$ as shown Fig. 6a. The increase of the permittivity assuming a homogeneous irreversible hysteres distribution may be calculated: [43, 94, 99]

$$ \Delta\epsilon_r = \frac{\rho A^{sym}}{2\epsilon_0 E_{AC}} \qquad (10) $$

$A^{sym}$ corresponds to the area of the triangle $\Delta^{sym}$ and $\rho$ to the hysteron density in the homogeneous part. This gives the linear increases of the permittivity expected by the Rayleigh law:

$$ \Delta\epsilon_r = \frac{\rho E_{AC}}{\epsilon_0} = \alpha^{sym} E_{AC} \qquad (11) $$

$\alpha^{sym}$ refers to the Rayleigh coefficient in the case of a symmetrical AC measurement.

To probe the asymmetry of the FORC distribution, we propose to use instead a unipolar waveform for the applied electric field which consists in:

$$ E(t) = \pm \frac{E_{AC}}{2} (1 + \sin(2\pi ft + \phi)) \qquad (12) $$

The sign $\pm$ denotes if we want to probe the positive hysteres or the negative hysteres (Fig. 6b). This signal has the properties to be always positive (or negative) and has a maximum magnitude of $E_{AC}$. It can be generated directly using the lock-in amplifier by adding to the generated signal the demodulated magnitude. Permittivity is obtained by using the same formula as conventional AC measurement. One can note the $\frac{1}{2}$ factor to obtain the same peak value of the signal.

Probing only a part of the Preisach plane is possible using a recently proposed extension of the FORC, namely the Unipolar Reversal Curve (URC) [100]. Fundamental differences between URC and the proposed unipolar impedance measurement exist since they are complementary. URC is a large signal method, and thus intends to explore the switching dynamic of the material. The proposed unipolar impedance measurement is a small signal method and thus aims to explore non switching dynamic of the material.

Since the proposed method is based on impedance measurement, long integration time can be used to increase the signal-to-noise ratio, allowing high sensitivity of the measurement and higher order harmonics are also extracted (magnitude and phase). This averaging is not straightforward for URC due its large signal nature which can induce fatigue of the sample. Proposed unipolar impedance measurement should not induce large fatigue effects, due to its small signal nature. Nevertheless, for the proposed method, the excitation signal needs
to be sinusoidal, contrary to URC where delay between excitation pulses may be added, allowing to studied time dependant imprint effect.

The portion of hysteron which are switched correspond to $\Delta^+$ and $\Delta^-$ for positive and negative unipolar signal, respectively. Due to symmetry considerations, they are the same and can be calculate using the following formulas:

$$ \Delta\varepsilon_r^- = \frac{\rho A^-}{2\varepsilon_0 \frac{E_{AC}}{2}} \quad \text{and} \quad \Delta\varepsilon_r^+ = \frac{\rho A^+}{2\varepsilon_0 \frac{E_{AC}}{2}} \qquad (13) $$

Which gives:

$$ \Delta\varepsilon_r^+ = \frac{\rho E_{AC}}{2\varepsilon_0} = \alpha^+ E_{AC} \quad \text{and} \quad \Delta\varepsilon_r^- = \frac{\rho E_{AC}}{2\varepsilon_0} = \alpha^- E_{AC} \qquad (14) $$

$\alpha^+$ and $\alpha^-$ refer to the Rayleigh coefficients in the case of a positive and negative unipolar measurement, respectively. With this new methodology, it would possible to apply the Rayleigh law with only a positive (or negative) signal. In the case of a homogeneous distribution of hysteron, the Rayleigh coefficient for the unipolar signal is the half of the one for bipolar signal, i.e. $\alpha^+ = \alpha^- = \frac{\alpha^{sym}}{2}$.

Fig. 6c shows the real part of the permittivity measured using the unipolar AC signal, before and after FORC measurement. For all the measurements, the permittivity increases when the AC field increases, similarly to what is obtained for conventional AC measurement. This increase is due to the irreversible polar boundaries contribution and corresponds to the (almost) homogeneous hysteron distribution close to the origin of the Preisach plane.

Rayleigh coefficients have been extracted for all the measurement conditions and are reported in Table I. The difference between positive and negative unipolar measurements indicates that the hysteron density is not the same for the positive and negative parts, which is consistent with the FORC measurement presented in Fig. 5a. According to Eqs. (11) and (14) the increases should be two times lower for unipolar measurements than for bipolar measurements. In the present case, the increase for negative unipolar and bipolar signals are very similar (values of $\alpha$ close to each other), whereas the increase is lower for unipolar positive than for negative. This reveals that the hysteron density is not homogeneous and the use of both unipolar and bipolar signals permits to measure this. After the FORC measurement, the irreversible contribution, represented by the $\alpha$ parameter, is slightly higher (for symmetric and negative unipolar measurements), similar with what has been observed for AC-poling [90–92]. The Rayleigh coefficient for positive unipolar measurement decreases,
<table><caption>Table I. Rayleigh coefficient extraction using Eq. (2).</caption><thead><tr><th>Condition</th><th>Measurement type</th></tr></thead><tbody><tr><td></td><td>Conventional AC Unipolar negative Unipolar positive</td></tr><tr><td>Before FORC</td><td>7.0 ± 0.4 cm kV<sup>-1</sup> 7.3 ± 0.1 cm kV<sup>-1</sup> 5.0 ± 0.8 cm kV<sup>-1</sup></td></tr><tr><td>After FORC</td><td>8.4 ± 0.1 cm kV<sup>-1</sup> 7.6 ± 0.1 cm kV<sup>-1</sup> 3.4 ± 0.2 cm kV<sup>-1</sup></td></tr></tbody></table>

which corresponds well to what is seen on the FORC density difference (Fig. 5e), a slight decrease of the hysteron density in the triangle $\Delta^+$.

The phase-angle of the third harmonic is also measured and is reported Fig. 6d. Before FORC, the symmetric measurement exhibits the same shape as seen in Fig. 3b (one can note the pinching is not visible due to lower maximum field to avoid cycling influence). For all the others conditions the same global evolution is visible. For low fields (around 5 kV cm$^{-1}$), $\delta_3$ is between -180° and -135°, corresponding to a mainly reversible contribution, and for high fields (above 20 kV cm$^{-1}$), the limit value is between -90° and 0°, corresponding to a mix between irreversible domain wall motion contribution and saturation. The main difference is the limit value which is closer to the ideal value of -90° for unipolar measurement, indicating areas $\Delta^+$ and $\Delta^-$ are more homogeneous than $\Delta^{sym}$.

## IV. CONCLUSIONS

In this article, the dielectric non-linearities in 0.5 (Ba$_{0.7}$Ca$_{0.3}$TiO$_3$) - 0.5 (BaZr$_{0.2}$Ti$_{0.8}$O$_3$) thin film are studied using impedance spectroscopy and harmonic measurements, as a function of the AC measuring field, at different frequencies and upon cycling. The measurements reveal a phase angle of the third harmonic close to -270°/+90°, corresponding to a pinching of the hysteresis loop, especially for low frequency. This strong frequency dependency of the pinching response (i.e. $\delta_3$ close to -270°/+90°) indicate that the pinching is induced by the presence of polar nanoregions (PNRs), whose responses are known to be strongly frequency-dependent.

When repeating the measurement, the PNR contribution changes as indicated by the evolution of the third harmonic from a pinched to a conventional relaxor response. This shows that strong changes in the PNR configuration can be induced even for low AC fields.

First Order Reversal Curve (FORC) measurements confirm the presence of a pinched
hysteresis loop by the presence of two peaks in the Preisach distribution. When repeating
the FORC measurement a second time, the distribution drastically changes towards to a
soft ferroelectric/relaxor. This Preisach distribution presents an asymmetry in the low field
range and we propose a new method to probe this asymmetry, called unipolar impedance
measurement. This method allows to observe only a specific part of the Preisach plane. We
shown the Rayleigh law is still valid even when measuring the dielectric permittivity with
only positive (or only negative) electric field. In addition to the obvious interest of this
new methodology from a characterization point of view, numerous applications use unipolar
signals, in which case the variation of the dielectric permittivity may need to be taken into
account.

SUPPLEMENTARY MATERIAL

Analysis of the different contributions to the permittivity (lattice, reversible and irre-
versible polar boundary motion contributions) as a function of the frequency and the mea-
surement number, *P*(*E*) and *I*(*E*) loops for the FORC extraction and comparison between
the second and third FORC measurements are given in supplementary material.

DATA AVAILABILITY

The data that support the findings of this study are available from the corresponding
author upon reasonable request.

ACKNOWLEDGMENTS

This work has been performed with the means of the CERTeM (microelectronics techno-
logical research and development center) of French region Centre Val de Loire. This work
was funded through the project MAPS in the program ARD+ CERTeM 5.0 by the Région
Centre Val de Loire co-funded by the European Union (ERC, DYNAMHEAT, N°101077402)
and by the French State in the frame of the France 2030 program through the French Re-
search Agency (ANR-22-PEXD-0018). Views and opinions expressed are however those of
the authors only and do not necessarily reflect those of the European Union or the Euro-
pean Research Council. Neither the European Union nor the granting authority can be held
responsible for them.

CONFLICT OF INTEREST

The authors declare no competing financial interest.

[1] A. Jayakrishnan, J. Silva, K. Kamakshi, D. Dastan, V. Annapureddy, M. Pereira, and K. Sekhar, Are lead-free relaxor ferroelectric materials the most promising candidates for energy storage capacitors?, *Progress in Materials Science* **132**, 101046 (2023).

[2] V. Veerapandiyan, F. Benes, T. Gindel, and M. Deluca, Strategies to improve the energy storage properties of perovskite lead-free relaxor ferroelectrics: A review, *Materials* **13**, 5742 (2020).

[3] H. Palneedi, M. Peddigari, G.-T. Hwang, D.-Y. Jeong, and J. Ryu, High-performance dielectric ceramic films for energy storage capacitors: Progress and outlook, *Advanced Functional Materials* **28**, 1803665 (2018).

[4] Prateek, V. K. Thakur, and R. K. Gupta, Recent Progress on Ferroelectric Polymer-Based Nanocomposites for High Energy Density Capacitors: Synthesis, Dielectric Properties, and Future Aspects, *Chemical Reviews* **116**, 4260 (2016).

[5] E. V. Colla, N. K. Yushin, and D. Viehland, Dielectric properties of (PMN)(1-x)(PT)x single crystals for various electrical and thermal histories, *Journal of Applied Physics* **83**, 3298 (1998), publisher: American Institute of Physics.

[6] R. Cowley, S. Gvasaliya, S. Lushnikov, B. Roessli, and G. Rotaru, Relaxing with Relaxors: A Review of Relaxor Ferroelectrics, *Advances in Physics - ADVAN PHYS* **60**, 229 (2011).

[7] A. A. Bokov and Z.-G. Ye, Recent progress in relaxor ferroelectrics with perovskite structure, *Journal of Materials Science* **41**, 31 (2006).

[8] M. Otoničar, A. Bradeško, L. Fulanović, T. Kos, H. Uršič, A. Benčan, M. J. Cabral, A. Henriques, J. L. Jones, L. Riemer, D. Damjanovic, G. Dražić, B. Malič, and T. Rojac, Connecting the Multiscale Structure with Macroscopic Response of Relaxor Ferroelectrics, *Advanced Functional Materials* **30**, 2006823 (2020).
[9] F. Li, S. Zhang, T. Yang, Z. Xu, N. Zhang, G. Liu, J. Wang, J. Wang, Z. Cheng, Z.-G. Ye, J. Luo, T. R. Shrout, and L.-Q. Chen, The origin of ultrahigh piezoelectricity in relaxor-ferroelectric solid solution crystals, Nature Communications 7, 10.1038/ncomms13807 (2016).

[10] F. Li, S. Zhang, D. Damjanovic, L.-Q. Chen, and T. R. Shrout, Local Structural Heterogeneity and Electromechanical Responses of Ferroelectrics: Learning from Relaxor Ferroelectrics, Advanced Functional Materials 28, 1801504 (2018).

[11] F. Li, S. Zhang, Z. Xu, and L.-Q. Chen, The contributions of polar nanoregions to the dielectric and piezoelectric responses in domain-engineered relaxor-PbTiO 3 crystals, Advanced Functional Materials 27, 1700310 (2017).

[12] D. Damjanovic, Contributions to the piezoelectric effect in ferroelectric single crystals and ceramics, Journal of the American Ceramic Society 88, 2663 (2005).

[13] L. M. Garten, M. Burch, A. S. Gupta, R. Haislmaier, V. Gopalan, E. C. Dickey, and S. Trolier-McKinstry, Relaxor ferroelectric behavior in barium strontium titanate, Journal of the American Ceramic Society 99, 1645 (2016).

[14] L. M. Garten, P. Lam, D. Harris, J.-P. Maria, and S. Trolier-McKinstry, Residual ferroelectricity in barium strontium titanate thin film tunable dielectrics, Journal of Applied Physics 116, 044104 (2014).

[15] E. K. H. Salje, M. A. Carpenter, G. F. Nataf, G. Picht, K. Webber, J. Weerasinghe, S. Lisenkov, and L. Bellaiche, Elastic excitations in BaTiO₃ single crystals and ceramics: Mobile domain boundaries and polar nanoregions observed by resonant ultrasonic spectroscopy, Physical Review B 87, 10.1103/physrevb.87.014106 (2013).

[16] J. Dec, S. Miga, W. Kleemann, and B. Dkhil, Nonlinear dielectric properties of PMN relaxor crystals within landau-ginzburg-devonshire approximation, Ferroelectrics 363, 141 (2008).

[17] A. E. Glazounov and A. K. Tagantsev, Phenomenological model of dynamic nonlinear response of relaxor ferroelectrics, Physical Review Letters 85, 2192 (2000), publisher: American Physical Society.

[18] J. Dec, W. Kleemann, S. Miga, C. Filipic, A. Levstik, R. Pirc, T. Granzow, and R. Pankrath, Probing polar nanoregions in Sr0.61Ba0.39Nb2O6 via second-harmonic dielectric response, Physical Review B 68, 092105 (2003), publisher: American Physical Society.

[19] D. V. Taylor and D. Damjanovic, Evidence of domain wall contribution to the dielectric permittivity in PZT thin films at sub-switching fields, Journal of Applied Physics 82, 1973
(1997).

[20] D. Damjanovic, Ferroelectric, dielectric and piezoelectric properties of ferroelectric thin films and ceramics, *Reports on Progress in Physics* **61**, 1267 (1998).

[21] D. A. Hall, Rayleigh behaviour and the threshold field in ferroelectric ceramics, *Ferroelectrics* **223**, 319 (1999).

[22] T. Schenk, M. Hoffmann, M. Pešić, M. H. Park, C. Richter, U. Schroeder, and T. Mikolajick, Physical Approach to Ferroelectric Impedance Spectroscopy: The Rayleigh Element, *Physical Review Applied* **10**, 064004 (2018).

[23] D. V. Taylor and D. Damjanovic, Domain wall pinning contribution to the nonlinear dielectric permittivity in Pb(Zr, Ti)O₃ thin films, *Applied Physics Letters* **73**, 2045 (1998).

[24] K. Nadaud, G. F. Nataf, N. Jaber, M. Bah, B. Negulescu, P. Andreazza, P. Birnal, and J. Wolfman, Subcoercive field dielectric response of 0.5 (Ba₀.₇Ca₀.₃TiO₃)–0.5 (BaZr₀.₂Ti₀.₈O₃) thin film: Peculiar third harmonic signature of phase transitions and residual ferroelectricity, *Applied Physics Letters* **124**, 042901 (2024).

[25] C. Borderon, R. Renoud, M. Ragheb, and H. W. Gundel, Description of the low field nonlinear dielectric properties of ferroelectric and multiferroic materials, *Applied Physics Letters* **98**, 112903 (2011).

[26] J. Bai, J. Yang, Y. Zhang, W. Bai, Z. Lv, K. Tang, J. Sun, X. Meng, X. Tang, and J. Chu, The ac sub-coercive-field dielectric resonances of (Pb, Sr)TiO₃ films at low temperature, *Ceramics International* **43**, S516 (2017).

[27] O. Boser, Statistical theory of hysteresis in ferroelectric materials, *Journal of Applied Physics* **62**, 1344 (1987).

[28] K. Nadaud, M. Sadl, M. Bah, F. Levassort, and H. Ursic, Effect of thermal annealing on dielectric and ferroelectric properties of aerosol-deposited 0.65 Pb(Mg₁/₃Nb₂/₃)O₃–0.35 PbTiO₃ thick films, *Applied Physics Letters* **120**, 112902 (2022).

[29] C. Borderon, A. E. Brunier, K. Nadaud, R. Renoud, M. Alexe, and H. W. Gundel, Domain wall motion in Pb(Zr₀.₂₀Ti₀.₈₀)O₃ epitaxial thin films, *Scientific Reports* **7**, 3444 (2017).

[30] L. M. Riemer, L. Jin, H. Uršič, M. Otonicar, T. Rojac, and D. Damjanovic, Dielectric and electro-mechanic nonlinearities in perovskite oxide ferroelectrics, relaxors, and relaxor ferro-
electrics, Journal of Applied Physics **129**, 054101 (2021).

[31] M. I. Morozov and D. Damjanovic, Hardening-softening transition in Fe-doped Pb(Zr,Ti)O3 ceramics and evolution of the third harmonic of the polarization response, Journal of Applied Physics **104**, 034107 (2008).

[32] M. Otoničar, M. Dragomir, and T. Rojac, Dynamics of domain walls in ferroelectrics and relaxors, Journal of the American Ceramic Society **105**, 10.1111/jace.18623 (2022).

[33] S. Hashemizadeh and D. Damjanovic, Nonlinear dynamics of polar regions in paraelectric phase of (Ba$_{1-x}$Sr$_x$)TiO$_3$ ceramics, Applied Physics Letters **110**, 192905 (2017).

[34] M. Morozov, D. Damjanovic, and N. Setter, The nonlinearity and subswitching hysteresis in hard and soft PZT, Journal of the European Ceramic Society **25**, 2483 (2005).

[35] S. Miga, J. Dec, and W. Cao, Time evolution of polar regions in (Sr$_{0.65}$Ba$_{0.35}$)$_2$NaNb$_5$O$_{15}$ single crystal, Phase Transitions **79**, 513 (2006).

[36] S. Miga, J. Dec, and W. Kleemann, Computer-controlled susceptometer for investigating the linear and nonlinear dielectric response, Review of Scientific Instruments **78**, 033902 (2007).

[37] S. Miga and J. Dec, Non-linear Dielectric Response of Ferroelectric and Relaxor Materials, Ferroelectrics **367**, 223 (2008).

[38] W. Kleemann, S. Miga, J. Dec, and J. Zhai, Crossover from ferroelectric to relaxor and cluster glass in BaTi1-xZrxO3 (x = 0.25-0.35) studied by non-linear permittivity, Applied Physics Letters **102**, 232907 (2013).

[39] S. Shetty, J. Kim, L. W. Martin, and S. Trolier-McKinstry, Non-linearity in engineered lead magnesium niobate (PbMg1/3Nb2/3O3) thin films, Journal of Applied Physics **128**, 194102 (2020), publisher: American Institute of Physics.

[40] R. Pirc, R. Blinc, and Z. Kutnjak, Nonlinear dielectric response of relaxor ferroelectrics, Physical Review B **65**, 214101 (2002), publisher: American Physical Society.

[41] R. Pirc, Z. Kutnjak, and N. Novak, Compressible spherical dipolar glass model of relaxor ferroelectrics, Journal of Applied Physics **112**, 114122 (2012).

[42] R. Pirc and R. Blinc, Spherical random-bond-random-field model of relaxor ferroelectrics, Physical Review B **60**, 13470 (1999).

[43] I. Fujii, E. Hong, and S. Trolier-Mckinstry, Thickness Dependence of Dielectric Nonlinearity of Lead Zirconate Titanate Films, Ultrasonics, Ferroelectrics, and Frequency Control, IEEE Transactions on **57**, 1717 (2010).
[44] L. Mitoseriu, L. Stoleriu, A. Stancu, C. Galassi, and V. Buscaglia, First order reversal curves diagrams for describing ferroelectric switching characteristics, *Processing and Application of Ceramics* **3**, 3 (2009).

[45] D. Ricinschi, A. Stancu, L. Mitoseriu, P. Postolache, and M. Okuyama, First order reversal curves diagrams applied for the ferroelectric systems, *Journal of Optoelectronics and Advanced Materials* **6**, 623 (2004).

[46] R. Pike, First-order reversal-curve diagrams and reversible magnetization, *Phys. Rev. B* **68**, 1 (2003).

[47] R. J. Harrison and J. M. Feinberg, FORCinel: An improved algorithm for calculating first-order reversal curve distributions using locally weighted regression smoothing, *Geochemistry, Geophysics, Geosystems* **9**, 10.1029/2008GC001987 (2008).

[48] L. Mitoseriu, C. E. Ciomaga, V. Buscaglia, L. Stoleriu, and D. Piazza, Hysteresis and tunability characteristics of Ba(Zr,Ti)O₃ ceramics described by first order reversal curves diagrams, *Journal of the European Ceramic Society* **27**, 3723 (2007).

[49] W. Zhu, I. Fujii, W. Ren, and S. Trolier-McKinstry, Influence of Mn doping on domain wall motion in Pb(Zr₀.₅₂Ti₀.₄₈)O₃ films, *Journal of Applied Physics* **109**, 064105 (2011).

[50] L. Stoleriu, A. Stancu, L. Mitoseriu, D. Piazza, and C. Galassi, Analysis of switching properties of porous ferroelectric ceramics by means of first-order reversal curve diagrams, *Phys. Rev. B* **74**, 10.1103/physrevb.74.174107 (2006).

[51] M. Hoffmann, T. Schenk, M. Pešić, U. Schroeder, and T. Mikolajick, Insights into antiferroelectrics from first-order reversal curves, *Applied Physics Letters* **111**, 10.1063/1.5003612 (2017).

[52] A. P. Roberts, Q. Liu, C. J. Rowan, L. Chang, C. Carvallo, J. Torrent, and C.-S. Horng, Characterization of hematite (α-Fe₂O₃), goethite (α-FeOOH), greigite (Fe₃S₄), and pyrrhotite (Fe₇S₈) using first-order reversal curve diagrams, *Journal of Geophysical Research: Solid Earth* **111**, n/a (2006).

[53] A. R. Muxworthy, Assessing the ability of first-order reversal curve (FORC) diagrams to un-ravel complex magnetic signals, *Journal of Geophysical Research* **110**, 10.1029/2004jb003195 (2005).

[54] S. Saremi, R. Xu, F. I. Allen, J. Maher, J. C. Agar, R. Gao, P. Hosemann, and L. W. Martin, Local control of defects and switching properties in ferroelectric thin films, *Physical Review*
Materials **2**, 084414 (2018).

[55] K. Nadaud, C. Borderon, R. Renoud, M. Bah, S. Ginestar, and H. W. Gundel, Evidence of residual ferroelectric contribution in antiferroelectric lead-zirconate thin films by first-order reversal curves, *Applied Physics Letters* **118**, 042902 (2021).

[56] G. Robert, D. Damjanovic, and N. Setter, Preisach modeling of ferroelectric pinched loops, *Applied Physics Letters* **77**, 4413 (2000).

[57] K. Nadaud, G. F. Nataf, N. Jaber, B. Negulescu, F. Giovannelli, P. Andreazza, P. Birnal, and J. Wolfman, Enhancement of piezoelectric properties in a narrow cerium doping range of Ba$_{(1-x)}$Ca$_x$Ti$_{(1-y)}$Zr$_y$O$_3$ evidenced by high throughput experiment, *ACS Applied Electronic Materials* **6**, 7392 (2024).

[58] W. Liu and X. Ren, Large Piezoelectric Effect in Pb-Free Ceramics, *Physical Review Letters* **103**, 257602 (2009).

[59] D. S. Keeble, F. Benabdallah, P. A. Thomas, M. Maglione, and J. Kreisel, Revised structural phase diagram of (Ba$_{0.7}$Ca$_{0.3}$TiO$_{3}$) – (BaZr$_{0.2}$Ti$_{0.8}$O$_{3}$), *Applied Physics Letters* **102**, 092903 (2013).

[60] See supplemental material at [link] analysis of the different contributions to the permittivity (lattice, reversible and irreversible polar boundary motion contributions) as a function of the frequency and the measurement number, *P(E)* and *I(E)* loops for the force extraction and comparison between the second and third FORC measurements are given in supplementary material.

[61] B. Ma, Z. Hu, S. Liu, S. Tong, M. Narayanan, R. E. Koritala, and U. Balachandran, Temperature-dependent dielectric nonlinearity of relaxor ferroelectric Pb$_{0.92}$La$_{0.08}$Zr$_{0.52}$Ti$_{0.48}$O$_{3}$ thin films, *Applied Physics Letters* **102**, 202901 (2013).

[62] D. Damjanovic, A. Biancoli, L. Batooli, A. Vahabzadeh, and J. Trodahl, Elastic, dielectric, and piezoelectric anomalies and Raman spectroscopy of 0.5 (Ba$_{0.7}$Ca$_{0.3}$TiO$_{3}$) – 0.5 (BaZr$_{0.2}$Ti$_{0.8}$O$_{3}$), *Applied Physics Letters* **100**, 192907 (2012).

[63] F. Feng and Y. Yan, Large electrostrictive effect in Mn-doped BCZT ferroelectric ceramics, *Ceramics International* **45**, 21315 (2019).

[64] M. Gao, W. Ge, X. Li, H. Yuan, C. Liu, H. Zhao, Y. Ma, and Y. Chang, Enhanced Dielectric and Energy Storage Properties in Fe-Doped BCZT Ferroelectric Ceramics, *physica status solidi* (a) **217**, 2000253 (2020).
[65] X. Liu, Y. Dai, X. Pei, and W. Chen, Giant room-temperature electrocaloric effect within wide temperature span in Sn-doped (Ba<sub>0.85</sub>Ca<sub>0.15</sub>)(Zr<sub>0.1</sub>Ti<sub>0.9</sub>)O<sub>3</sub> lead-free thin films, *Ceramics International* **49**, 1846 (2023).

[66] X. Xu, G. Yuan, B. Xu, J. Yin, J. Liu, and Z. Liu, Influence of the strain on dielectric and ferroelectric properties of 0.5 (Ba<sub>0.7</sub>Ca<sub>0.3</sub>TiO<sub>3</sub>)–0.5 (BaZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub>), *Journal of the American Ceramic Society* **98**, 2823 (2015).

[67] A. B. Haugen, J. S. Forrester, D. Damjanovic, B. Li, K. J. Bowman, and J. L. Jones, Structure and phase transitions in 0.5 (Ba<sub>0.7</sub>Ca<sub>0.3</sub>TiO<sub>3</sub>)–0.5 (BaZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub>) from -100 °c to 150 °c, *Journal of Applied Physics* **113**, 10.1063/1.4772741 (2013).

[68] F. Benabdallah, A. Simon, H. Khemakhem, C. Elissalde, and M. Maglione, Linking large piezoelectric coefficients to highly flexible polarization of lead free BaTiO<sub>3</sub>–CaTiO<sub>3</sub>–BaZrO<sub>3</sub> ceramics, *Journal of Applied Physics* **109**, 10.1063/1.3599854 (2011).

[69] C. Bhardwaj, B. S. S. Daniel, and D. Kaur, Pulsed laser deposition and characterization of highly tunable (1-x)(Ba<sub>0.7</sub>Ca<sub>0.3</sub>TiO<sub>3</sub>)-x (BaZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub>) thin films grown on LaNiO<sub>3</sub>/Si substrate, *Journal of Physics and Chemistry of Solids* **74**, 94 (2013).

[70] Y. Lin, N. Qin, G. Wu, T. Sa, and D. Bao, Dielectric relaxor behaviors and tunability of (1-x)(Ba<sub>0.7</sub>Ca<sub>0.3</sub>TiO<sub>3</sub>)-x (BaZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub>) thin films fabricated by sol-gel method, *Applied Physics A* **109**, 743 (2012).

[71] K. Nadaud, C. Borderon, R. Renoud, M. Bah, S. Ginestar, and H. W. Gundel, Dielectric, piezoelectric and electrostrictive properties of antiferroelectric lead-zirconate thin films, *Journal of Alloys and Compounds* **914**, 165340 (2022).

[72] K. Nadaud, C. Borderon, R. Renoud, M. Bah, S. Ginestar, and H. W. Gundel, Metastable and field-induced ferroelectric response in antiferroelectric lead zirconate thin film studied by the hyperbolic law and third harmonic response, *Journal of Applied Physics* **133**, 174102 (2023).

[73] L. Cima, E. Laboure, and P. Muralt, Characterization and model of ferroelectrics based on experimental Preisach density, *Review of Scientific Instruments* **73**, 3546 (2002).

[74] K. Nadaud, C. Borderon, R. Renoud, M. Bah, S. Ginestar, and H. W. Gundel, Study of the long time relaxation of the weak ferroelectricity in PbZrO<sub>3</sub> antiferroelectric thin film using Positive Up Negative Down and First Order Reversal Curves measurements, *Thin Solid Films* **773**, 139817 (2023).
[75] W. J. Merz, Double hysteresis loop of BaTiO₃ at the curie point, *Physical Review 91*, 513 (1953).

[76] E. Sapper, N. Novak, W. Jo, T. Granzow, and J. Rödel, Electric-field-temperature phase diagram of the ferroelectric relaxor system (1-x)Bi₁/₂Na₁/₂TiO₃ - xBaTiO₃ doped with man- ganese, *Journal of Applied Physics* **115**, 194104 (2014).

[77] Z. Li, J. Yu, S. Hao, and P.-E. Janolin, Enhancing properties of lead-free ferroelectric BaTiO₃ through doping, *Journal of the European Ceramic Society* **42**, 4693 (2022).

[78] Y. Zhao, J. Wang, L. Zhang, C. Wang, and S. Liu, Aging rate of cerium doped Ba(Ti₀.₉₉Mn₀.₀₁)O₃, *Ceramics International* The 10th Asian Meeting on Electroceramics (AMEC-10), **43**, S70 (2017).

[79] W. Kleemann, S. Miga, and J. Dec, Non-Linear Dielectric Response of Relaxor Ferroelectrics Within Ginzburg-Landau-Devonshire Approximation, *Ferroelectrics* **354**, 13 (2007).

[80] N. Bassiri-Gharb, I. Fujii, E. Hong, S. Trolier-McKinstry, D. V. Taylor, and D. Damjanovic, Domain wall contributions to the properties of piezoelectric thin films, *Journal of Electroceramics* **19**, 49 (2007).

[81] N. B. Gharb and S. Trolier-McKinstry, Dielectric nonlinearity of Pb(Yb₁/₂Nb₁/₂)O₃/PbTiO₃ thin films with {100} and {111} crystallographic orientation, *Journal of Applied Physics* **97**, 064106 (2005).

[82] N. B. Gharb, S. Trolier-McKinstry, and D. Damjanovic, Piezoelectric nonlinearity in ferroelectric thin films, *Journal of Applied Physics* **100**, 044107 (2006).

[83] M. T. Becker, C. J. Burkhardt, R. Kleiner, and D. Koelle, Impedance spectroscopy of ferroelectrics: The domain wall pinning element, *Journal of Applied Physics* **132**, 044104 (2022), publisher: American Institute of Physics.

[84] E. V. Colla, E. L. Furman, S. M. Gupta, N. K. Yushin, and D. Viehland, Dependence of dielectric relaxation on ac drive in [Pb(Mg₁/₃Nb₂/₃)O₃](₁-x) - (PbTiO₃)ₓ single crystals, *Journal of Applied Physics* **85**, 1693 (1999).

[85] D. Damjanovic, S. Bharadwaja, and N. Setter, Toward a unified description of nonlinearity and frequency dispersion of piezoelectric and dielectric responses in pb(zr, ti)O₃, *Materials Science and Engineering: B* **120**, 170 (2005).

[86] A. K. Jonscher, Dielectric relaxation in solids, *Journal of Physics D: Applied Physics* **32**, R57 (1999).
[87] K. Nadaud, C. Borderon, R. Renoud, and H. W. Gundel, Effect of manganese doping of BaSrTiO₃ on diffusion and domain wall pinning, *Journal of Applied Physics* **117**, 084104 (2015).

[88] M. D. Coulibaly, C. Borderon, R. Renoud, and H. W. Gundel, Effect of ferroelectric domain walls on the dielectric properties of PbZrO₃ thin films, *Applied Physics Letters* **117**, 142905 (2020).

[89] A. Jonscher, *Dielectric relaxation in solids* (Chelsea Dielectrics Press Ltd, 1983).

[90] C. Qiu, J. Liu, F. Li, and Z. Xu, Thickness dependence of dielectric and piezoelectric properties for alternating current electric-field-poled relaxor-PbTiO₃ crystals, *Journal of Applied Physics* **125**, 014102 (2019).

[91] C. He, T. Karaki, X. Yang, Y. J. Yamashita, Y. Sun, and X. Long, Dielectric and piezoelectric properties of Pb[(Mg₁/₃Nb₂/₃)₀.₅₂(Yb₁/₂Nb₁/₂)₀.₁₅Ti₀.₃₃]O₃ single-crystal rectangular plate and beam mode transducers poled by alternate current poling, *Japanese Journal of Applied Physics* **58**, SLLD06 (2019), publisher: IOP Publishing.

[92] G. Tian, F. Liu, J. Du, L. Zhao, X. Qi, M. Zhao, and L. Zheng, Simultaneous enhancement of piezoelectricity and mechanical quality factor in relaxor-ferroelectric crystals via defect engineering, *Applied Physics Letters* **122**, 102901 (2023).

[93] X. Cao, H. Tian, C. Hu, F. Huang, Y. Wang, X. Sun, and Z. Zhou, Defect dipole evolution and its impact on the ferroelectric properties of Fe-doped KTN single crystals, *Journal of the American Ceramic Society* **102**, 3117 (2019).

[94] I. Fujii, S. Trolier-Mckinstry, and C. Nies, Effect of grain size on dielectric nonlinearity in model BaTiO₃ -based multilayer ceramic capacitors, *Journal of American Ceramic Society* **199**, 194 (2011).

[95] F. P. G. Fengler, M. Pešić, S. Starschich, T. Schneller, C. Künneth, U. Böttger, H. Mulaosmanovic, T. Schenk, M. H. Park, R. Nigon, P. Muralt, T. Mikolajick, and U. Schroeder, Domain Pinning: Comparison of Hafnia and PZT Based Ferroelectrics, *Advanced Electronic Materials* **3**, 1600505 (2017).

[96] F. P. G. Fengler, M. Hoffmann, S. Slesazeck, T. Mikolajick, and U. Schroeder, On the relationship between field cycling and imprint in ferroelectric Hf₀.₅Zr₀.₅O₂, *Journal of Applied Physics* **123**, 204101 (2018).
[97] T. Schenk, M. Hoffmann, J. Ocker, M. Pešić, T. Mikolajick, and U. Schroeder, Complex internal bias fields in ferroelectric hafnium oxide, *ACS Applied Materials and Interfaces* **7**, 20224 (2015).

[98] P. Chojecki, G. Walters, Z. Forrester, and T. Nishida, Preisach modeling of imprint on hafnium zirconium oxide ferroelectric capacitors, *Journal of Applied Physics* **130**, 094102 (2021).

[99] G. Robert, D. Damjanovic, N. Setter, and A. V. Turik, Preisach modeling of piezoelectric nonlinearity in ferroelectric ceramics, *Journal of Applied Physics* **89**, 5067 (2001).

[100] S. Vecchi, F. M. Puglisi, P. Appelt, R. Guido, X. Wang, S. Slesazeck, T. Mikolajick, and S. Lancaster, Evaluation of Imprint and Multi-Level Dynamics in Ferroelectric Capacitors, *Advanced Electronic Materials* **n/a**, 2400204 (2024).