Flexoelectricity-driven giant polarization in (Bi, Na)TiO_3-based ferroelectric thin films

Domain: physics
Published: 2025-11-26
Editorial outcome: accept

Abstract
This study demonstrates the use of oxygen vacancy-induced planar defects to significantly enhance electrical polarization through a local flexoelectric effect. By introducing an appropriate level of aliovalent dopants, numerous local planar defects are induced in (Bi0.5Na0.5)TiO3-based thin films. These defects, identified as oxygen-deficient structures through direct visualization of oxygen atoms and oxygen vacancies using integrated differential phase-contrast microscopy, result in the formation of head-to-head domain structures. Geometric phase analysis confirms that these structures exhibit a substantial local strain gradient of up to 109 m-1, contributing significantly to the flexoelectric polarization. Consequently, a giant maximum polarization (Pm) of 161 μC cm-2 under 750 kV cm-1 and a remanent polarization Pr = 115 μC cm-2 along with a coercive field of 250 kV cm-1 are achieved, allowing these (Bi0.5Na0.5)TiO3-based thin films to be used in low-power electronic applications. Crucially, the Pm and Pr of the thin films can be sustained at 133 and 98 μC cm-2, respectively, at 230 °C. Additionally, they exhibit exceptional high-temperature fatigue endurance, with Pm and Pr demonstrating a negligible reduction of less than 9% after 107 cycles under 750 kV cm-1 at 230 °C. These values surpass those previously reported for oxide perovskite thin films at elevated temperatures, demonstrating potential applications of our thin films in high-temperature environments. Our findings offer promising avenues for advancing the application fields of ferroelectric thin films.

Introduction
Large polarization is consistently pursued and anticipated in ferroelectric thin films due to its ability to facilitate superior energy storage capacity, improved electro-mechanical coupling efficiency, and elevated dielectric constants1, which significantly impact the functional properties and performance of various devices, including solid-state cooling systems, electrostatic energy storage units, non-volatile random-access memory (FeRAM), and ultrasound medical imaging technologies2. In conventional thin films possessing perovskite structures, enhanced polarization is typically attained through employing highly lattice-mismatched substrates. For example, BiFeO3 (BFO) thin films grown on LaAlO3 substrates experiences significant compressive epitaxial strain (~4.5%), stabilizing a maximum polarization of ~130 μC cm-2 3–5. Another viable route to improve polarization involves utilizing interphase strain in nanocomposite-like films to generate supertetragonal (super-T) structures. A super-tetragonal phase is used to describe the unique lattice structure in perovskite oxide thin films that exhibits a markedly large c/a ratio, e.g., 1.10-1.25 6–8, often resulting in an enhanced polarization. Such super-T structures arise from mismatches in c-lattice parameter between the thin film material and a secondary phase, such as PbO in PbTiO3 thin films6, BaO in BaTiO3 thin films7, and β-Bi2O3 in (Bi0.5Na0.5)TiO3-based thin films8, precipitated due to local chemical inhomogeneities, leading to unusually large polarizations exceeding 100 μC cm-2.

Despite the efficacy of these established methodologies, their limitations are evident. The improvement in polarization achieved through compressive strain from the substrate is finite. As the thickness of thin films increases, the strain imposed by the substrate diminishes, causing the film to gradually revert to their bulk lattice parameters and properties3–5,9. Similarly, when leveraging interphase strain from a secondary phase, a high coercive field (Ec), typically exceeding 2000 kV cm-1, is often necessary for switching domains. The requirement for such substantial power consumption to achieve heightened polarization raises serious concerns regarding Joule heating and catastrophic dielectric breakdown6,7, which are particularly severe at elevated temperatures and also pose challenges for the supporting insulation system, potentially restricting their use in highly integrated miniaturized devices. Balancing the operating field and polarization, i.e., achieving both a low coercive field and large polarization simultaneously, remains a challenge for ferroelectric thin films, especially in the development of beyond-complementary metal oxide semiconductor devices.

In this work, we demonstrate the feasibility of enhancing polarization in relatively thick films under a low driving field by leveraging the flexoelectric effect. The incorporation of Fe cations into 0.94(Bi0.5Na0.5)TiO3-0.06BaTiO3 (BNBT) thin films facilitate the creation of oxygen vacancies, which accumulate in tail-to-tail (T-T) domain walls as charged defects by creating oxygen-deficient superstructures (planar defects). Between the T-T domain walls, head-to-head (H-H) domains walls emerge, generating a giant local strain gradient that can reach a level of up to 109 m-1, contributing significantly to flexoelectric polarization. Consequently, we achieve a giant maximum polarization (Pm) of 161 μC cm-2 under 750 kV cm-1 and a remanent polarization Pr = 115 μC cm-2 along with a coercive field of 250 kV cm-1, demonstrating superior polarization compared with previously reported lead-free ferroelectrics under a comparable electric field.

Results
Room-temperature ferroelectric hysteresis (P–E) loops of Fe-doped BNBT (Fe-BNBT) thin films under selected electric fields are presented in Fig. 1a. The maximum polarization values observed are 161, 118, and 95 μC cm-2 for thin films grown on [001], [110], and [111]-oriented SrTiO3 substrates, respectively, consistent with the empirical relationship of polarization along different crystallographic orientations of a typical tetragonal lattice10, unequivocally demonstrating the tetragonal structure of our Fe-BNBT thin films. To quantify the Pr, the remanent P–E loops were measured under 750 kV cm-1, as illustrated in Fig. 1b. A significant Pr of 115 μC cm−2 was achieved for the [001]-oriented Fe-BNBT thin films. Figure 1c shows the out-of-plane X-ray diffraction (XRD) θ-2θ scan results of Fe-BNBT thin films, indicating that the growth orientation of the films aligns with their respective substrates. Given the tetragonal structure of our Fe-BNBT thin films and their optimal ferroelectric properties when grown on a [001]-oriented substrate, we will concentrate on [001]-oriented films to understand the origin of the large polarization. The XRD reciprocal space map (XRD-RSM) measurements were performed to ascertain the in-plane lattice parameters of Fe-BNBT thin films. As illustrated in Fig. 1(d), diffraction spots referring to LSMO, STO, and Fe-BNBT layers are evident. The identical Qx values of these diffraction spots in the RSM of Fe-BTBT suggest that both Fe-BNBT and LSMO layers are fully strained on the STO substrate. The lattice parameters a and c are calculated as 3.905 and 3.919 Å, respectively. Importantly, neither super-T having elongated c lattice parameter were observed nor was substantial in-plane strain (-0.5%) detected in the XRD and RSM results, suggesting the presence of some other potential contributions to the observed large polarization.

Orbital hybridization plays an essential role for the emergence of ferroelectricity in the thin films of complex oxides11. To understand the electronic structure and orbital hybridization, Fe-BNBT thin films were examined using X-ray absorption spectroscopy (XAS). Figure 2(a) illustrates the Fe K-edge XAS spectra of Fe-BNBT thin films. The absorption edge energy is determined to be 7126.14 eV, consistent with that of Fe2O3, indicating the valent state of Fe is 3+12. Additionally, the pre-edge peak A1 and post-edge A2 and A3 are also observed, attributed to the electric quadrupole-forbidden transition from the O 1s level to Fe 3d, O 2p band to Fe 3d orbit, and O 1s to Fe 4p dipole-allowed transition, respectively. This confirms that the Fe3+ cations enter into the B sites by replacing Ti4+13. Fig. 2(b) illustrates the O K-edge XAS spectra of Fe-BNBT thin films, revealing signals originating from the hybridization of the O 2p and unoccupied Ti/Fe 3d orbitals. The diminished peak intensity ratios of eg and t2g signals (<1) imply a weakening of the Ti-O and Fe-O bonds in the hybridization of O 2p and Ti/Fe 3d orbitals, suggesting the possible formation of oxygen vacancies within the Fe-BNBT thin films14.

The Ti L-edge XAS spectra of Fe-BNBT are illustrated in Fig. 2c to investigate the symmetry of the TiO6 octahedron. L3 and L2 edges are generated by the spin-orbit interaction in the spectra, and each edge contains t2g and eg peaks. The L3 eg state is oriented toward ligand anions, making it more responsive to deviations from the Ti octahedron symmetry15. As shown in the inset of Fig. 2(c), the observed splitting of L3 eg is associated with the non-centrosymmetric TiO6 octahedron, affirming the tetragonal structures of Fe-BNBT films16. Notably, the relatively higher intensity of the peak localized at low photon energy (yellow) compared to the high-photon energy peak (pink) is connected with the localized character of the eg states, demonstrating the existence of defective structures in Fe-doped BNBT thin films16. These distortions are most likely induced by oxygen vacancies and associated planar defects, stemming from the replacement of Ti by Fe cations. Such structural and chemical disorder reduces the local symmetry and modifies the crystal field environment, leading to orbital energy shifts and increased localization of the eg electrons16. Moreover, based on the configuration interaction cluster model calculation, an increased ratio of the peak intensity of the eg and t2g peaks (Ieg/t2g) usually indicates enhanced tetragonal distortion17. The Ieg/t2g value for Fe-BNBT films is calculated as 1.1, suggesting a high level of local strain in Fe-BNBT films16.

X-ray Photoelectron Spectroscopy (XPS) was employed to estimate the level of oxygen vacancy in Fe-BNBT thin films. As shown in Fig. 2(d), the deconvolution of the O 1s peaks suggests the presence of the lattice oxygen OL (blue peak), the lattice oxygen connected with the lower valence cation after vacancy creation OV (green peak), and adsorbed surface water and/or hydroxyl molecules OC (purple peak)18,19. The observation of OV XPS peak echoes a strong correlation between the defective structure observed via XAS and the presence of oxygen vacancies.

To clarify the critical role of oxygen vacancies in producing large polarization in our Fe-BNBT thin films, high-angle annular dark-field (HAADF) scanning transmission electron microscopy (STEM) measurements were conducted. The thickness of the films was confirmed to be ~550 nm through low-magnification HAADF-STEM images (Supplementary Fig. 1). Notably, numerous two-dimensional and oriented segments along [101] and [-101] (referred to as bright streaks) are observed in the Fe-BNBT thin films (Fig. 3(a)). Polarization mapping was performed to examine the impact of bright streaks on the domain structures of Fe-BNBT thin films. The length and direction of the arrows represent the strength and orientation of the polarization, respectively. Significantly, T-T domain walls appear near the bright streak regions in our Fe-BNBT thin films (Fig. 3(b)). Density functional theory (DFT) calculations indicate that the aggregation of oxygen vacancies can give rise to charged T-T domain walls20. The formation of oxygen vacancies in Fe-doped BNBT can be described using the defect chemistry notation: Fe2O3 ! 2TiO2 2Fe'Ti + VO·· + 3OXO. To further investigate the structural features of these bright streaks, integrated differential phase-contrast (iDPC) imaging was employed. This advanced imaging technique can enhance the image contrast of weakly absorbing, low-atomic-number elements, enabling the simultaneous determination of heavy and light element positions21. The iDPC-STEM result in Fig. 3(d) effectively distinguishes A/B-site cations and oxygen atoms, revealing that the bright streaks are planar defects enriched with oxygen vacancies, as delineated by the yellow dotted rectangle.

To provide a more intuitive representation of the oxygen vacancies within the bright streaks, the intensity of oxygen atom signals from two regions—denoted by green and yellow dotted rectangles—was analyzed (Fig. 3(e)). Compared to the surrounding areas (green), the bright streaks (yellow) exhibit a significantly lower contrast, as indicated by blue arrows, underscoring the presence of a nonstoichiometric Fe-BNBT lattice in the vicinity of the streaks. These findings provide compelling evidence that oxygen vacancies are localized within the bright streaks, serving as charged defects, contributing to the formation of T-T domain walls in our Fe-BNBT thin films.

We conducted first-principles calculations based on density functional theory (DFT, Experimental Section) to gain further insight into the onset of the nonstoichiometry (planar defects) in Fe-BNBT. Given that XRD and XAS results confirm the tetragonal structure of our BNBT thin films, we adopted a tetragonal phase with a space group of P4/mbm22 for the DFT calculation, as illustrated in Supplementary Fig. 2a. To investigate the structural stability of BNBT, the formation energy of BNBT unit cells with various dopants, i.e., Mn and Fe cations, were simulated. The formation energies of BNBT unit cells with oxygen vacancies and without oxygen vacancies are given in Supplementary Fig. 2(b). In comparison to Mn-doped BNBT, whose formation energy of the pristine structure is situated between the BNBT unit cells exhibiting oxygen vacancy-related distortions, defect-free Fe-doped BNBT has the highest formation energy. This indicates that the creation of oxygen vacancy favors a more stable Fe-BNBT structure, which is in line with our experimental observations.

Given the high density of planar defects in the films, an intriguing question arises regarding how the domain structures transform between two adjacent T-T domain walls. To address this question, we conducted polarization mapping over a large area (55×20 nm2), as shown in Fig. 4a. Interestingly, H-H domains were found to form between adjacent pairs of bright streaks. The polarization mapping reveals enhanced polarization strength near the H-H domain walls. Figure 4(c) depicts the enlarged polarization maps of Fe-BNBT thin films, with areas 1 and 2 corresponding to the regions outlined by red dashed boxes in Fig. 4(a), which highlights distinct disparities. In area 1, BNBT unit cells adjacent to H-H domain walls exhibit substantial polarization (represented by longer arrows), with polarization diminishing progressively in unit cells further from the domain walls. Conversely, the polarization in area 2, which features H-T domain structures, appears weak and uniformly distributed. To analyze the mechanism behind the enhanced polarization observed near H-H domain walls, as opposed to areas with T-T and H-T domain structures, geometric phase analysis (GPA) is utilized to map the strain distribution in Fe-BNBT thin films, using the strain level of stoichiometric BNBT unit cells (area 2) as a reference. Figure 4(b) displays the strain mapping results of Fe-BNBT thin films generated by GPA. In contrast to the negligible strain levels found in areas with T-T and H-T domain structures, a significant local strain of ~10% is observed near the H-H domain walls, which is believed to be responsible for the increased polarization observed in area 1. Moreover, a strain gradient oriented along the white arrow is also obtained. Figure 4(d) illustrates the level of strain and strain gradient as a function of distance along the direction of the white arrow in Fig. 4(b). Notably, Fe-BNBT thin films exhibit a substantial local strain gradient, reaching up to 109 m-1.

Strain gradients inherently give rise to flexoelectric polarization. A previous study has reported that a large strain gradient on the order of 106 m-1 was found in LaFeO3 thin films grown on (001) LaAlO3 single crystal substrates23, giving rise to a flexoelectric polarization of 9.27 μC cm−2. The strain-gradient-induced electrical polarization (flexoelectric polarization, PF) can be determined by24:

PF = μijkl ∂ϵij / ∂xk (1)

where μijkl is the fourth-rank flexoelectric coefficient and ∂ϵij / ∂xk is the strain gradient. Thus, the polarization for ferroelectrics involving flexoelectric effect can be expressed by:

P = PD + PF = χij Ej + μijkl ∂ϵij / ∂xk (2)

where χij is the dielectric susceptibility of the ferroelectrics, Ej the applied electric field. PD = χij Ej describes the polarization stemming from the dielectric response. To assess the contribution of flexoelectricity to the observed polarization of Fe-BNBT thin films, the determination of the flexoelectric coefficient is necessary. Unfortunately, no commercially available techniques currently exist that enable a standardized process for measuring the flexoelectric coefficient in thin films with compelling reliability and accuracy24. Previous reports revealed that the flexoelectric coefficient of perovskite-structured single crystals was on the order of magnitude of 10-1 to 101 nC m-1, e.g., |μ11|, |μ12|, and |μ44| are 0.15, 5.5, and 1.9 nC m-1 for BaTiO3 single crystal and 0.2, 7 and 5.8 nC m-1 for SrTiO3 single crystal, respectively24. By incorporating this range of flexoelectric coefficient into Equation (1), the flexoelectric polarization of our Fe-BNBT thin films can be up to 100 μC cm-2. We also computed the polarization of typical BNBT with a tetragonal structure (c/a = 1.004). The calculated electrical polarization is 60.88 μC cm-2, which closely aligns with our experimental findings, i.e., 160.88 (60.88 μC cm-2 from the dielectric response + 100 μC cm-2 from the flexoelectric effect) vs 161 μC cm-2. Considering the high density of such local “planar defects” and the associated H-H domain structures in Fe-BNBT thin films, the exceptionally high local strain gradient is believed to play a pivotal role in generating the overall giant polarization.

The phase-field simulation was employed to verify the contribution of the flexoelectric effect to the improved polarization. Details of the phase-field simulations are provided in the Experimental Section. A typical ferroelectric with a tetragonal structure was established as the initial state in a BNBT single crystal. Bound charges were introduced to simulate the oxygen vacancies concentrated at bright streaks (planar defects). Compared to typical BNBT with H-T domain configurations (Fig. 5a, g), T-T domain walls are created upon the introduction of local charges, with H-H domain structures forming between adjacent T-T domain walls (Fig. 5b, h), consistent with observations from STEM. Additionally, strain mapping analysis, presented in Fig. 5c, d, revealed elevated strain levels at the H-H domain walls, while typical BNBT exhibited relatively low and uniformly distributed strain, confirming the presence of strain gradient in films. This finding aligns with the simulated P-E loops, where the bound charge-modified Fe-BNBT demonstrates enhanced Pm and Pr, compared to the typical BNBT, while their Ec remains comparable (Fig. 5e, f). This supports the role of the flexoelectric effect in producing the substantial polarization observed in our BNBT thin films. Furthermore, the T-T and H-H domains exhibit remarkable switchable behavior under electric field stimulation. These domains can be reoriented in response to the applied electric field, and upon the removal of the field, only a small fraction of the domains revert to their initial configuration, while the majority maintain their alignment (Supplementary Fig. 3). This behavior accounts for the exceptional Pr observed in the Fe-BNBT thin film.

To assess the thermal stability of the ferroelectricity in Fe-BNBT thin films, temperature-dependent P–E loops were measured, as illustrated in Fig. 6a–c. Remarkably, the Pm and Pr of Fe-BNBT thin films can be maintained at 133 and 98 μC cm-2 at 230 °C, respectively, indicating that temperature has minimal impact on the strength of the flexoelectric effect. To the best of our knowledge, these values surpass those previously reported for oxide perovskite thin films at elevated temperatures, highlighting their potential for applications in high-temperature environments, particularly in advanced FeRAM systems for aerospace and defense technologies. A charging and discharging frequency of 50 kHz was employed to conduct the fatigue measurement of P–E loops under an AC electric field of 750 kV cm-1, both at room temperature and at 230 °C. Figure 6d, e present the P-E loops of 550 nm-thick Fe-BNBT thin films with respect to the electric field cycles, indicating the variation in Pm and Pr. At room temperature, a reduction in Pm and Pr is found to be 14% and 23%, respectively, compared to the respective values before cycling, indicating a reasonably good endurance to cycling stress. Notably, Fe-BNBT thin films exhibited exceptional fatigue endurance at 230 °C, with Pm and Pr showing only marginal decreases to 122 and 90 μC cm−2, respectively, representing less than 9% variation. This highlights the extraordinary thermal stability of the ferroelectric properties in the flexoelectricity-enhanced Fe-BNBT thin film, emphasizing their suitability for high-temperature applications. The superior fatigue endurance at elevated temperatures, compared to room temperature, maybe attributed to the reduced impact of the pinning effect at higher temperatures25, which further reinforces the feasibility of Fe-BNBT thin films for deployment in extreme thermal environments.

Table 1 compares the Pm, Pr, and Ec of the best-performing ferroelectric thin films to date. BFO-based thin films exhibit relatively low Ec and can simultaneously maintain high values of Pm and Pr, however, polarizations of the biaxial strain-modified BFO thin films are prone to temperature, as even a slight increase in temperature can lead to substantial electric leakage, severely degrading their performance at elevated temperatures. Likewise, although incredibly large remanent polarization of 100 and 236 μC cm−2 were reported in BaTiO3 and PbTiO3 thin films, respectively, the intensely high coercive field, i.e., >2000 kV cm-1, trigger grave concerns over electric breakdown, rapid degradation, and Joule heating, which greatly restricts their practical applications. In contrast, our Fe-BNBT thin films demonstrate excellent thermal stability, and exhibit a significantly lower Ec of 250 kV cm-1 and a giant Pr of 115 μC cm−2. This performance surpasses that of the majority of state-of-the-art ferroelectrics under a comparable strength of driving electric fields (Emax = 750 kV cm-1).

Methods
Experimental

Sample Preparation. Undoped BNBT, along with Mn, La, Sm, Eu, Co, and Fe-doped BNBT ceramic targets for thin film deposition, were prepared by the conventional sintering process. The nominal doping level in our BNBT-based thin films is based on the composition of the ceramic target used for Pulsed Laser Deposition (PLD). 10 mol% excess of Bi2O3 was added to the target for compensating the loss of Bi during the sintering and subsequent laser ablation process. BNBT thin films with various thicknesses and Fe contents were grown on La0.7Sr0.3MnO3 (LSMO)-electroded (001) SrTiO3 (STO) single crystal substrates (a = 0.3905 nm, cubic) by PLD (PASCAL Combi, Japan). A KrF excimer laser (λ = 248 nm) with an energy density of ~1.45 J cm-2 and frequency of 5 Hz was adopted for the deposition. The deposition temperature and oxygen partial pressure are 800 °C and 200 mTorr, respectively. Gold top electrodes with a 200 μm diameter were deposited through a shadow mask by a sputter coater (Leica EM ACE600) prior to the electrical measurements.

Characterizations. Temperature-dependent ferroelectric hysteresis (P-E) loops were measured by the Radiant ferroelectric workstation connected with a program-controlled chamber (Linkam, HFS600E-PB4). The room-temperature X-ray diffraction (XRD) patterns and X-ray diffraction reciprocal space mapping (XRD-RSM) were obtained by Smartlab (Rigaku, 9kW rotating anode thin-film XRD). The XPS measurements were carried out by an X-ray photoelectron spectrometer (Thermo Escalab 250Xi). The orbital hybridization of thin films was evaluated at the Australian Synchrotron and at the Cornell High Energy Synchrotron Source.

Specimens for STEM observations were prepared by focused ion beam (FEI Helios G4 PFIB and Zeiss Auriga dual beam systems). Cross-sectional lamellae were lifted out from the thin films and mounted to copper grids by platinum. Then the lamellae were thinned to the thickness allowing the electron-transparency. HAADF-STEM images were acquired using a HAADF detector in an aberration-corrected FEI-Themis Z operated at 300 kV. The images were obtained with a semi-convergence angle of 17.9 mrad and a collection angle range of 48-200 mrad.

Details of density-functional theory calculations

The formation energy of Mn- and Fe-doped BNBT structures (referred to here as “pristine”) is calculated by substituting a Ti atom in the unit cell with the respective dopant. Starting with this pristine configuration, which has a P4/mbm symmetry, different configurations are generated by placing a single oxygen vacancy (VO) in the unit cell, using the configuration tools available in Pymatgen26. Each resulting structure is fully relaxed, allowing atomic positions and lattice parameters to adjust, using the Perdew-Burke-Ernzerhof (PBE)27 functional within the VASP28 code, along with the standard projector augmented-wave (PAW)29 method. A +U correction is applied to properly model the localized d-electronic states of the transition metals. The optimization process proceeds until the residual forces are below 10-4 eV Å-1. An energy cutoff of 520 eV is set for the PAW pseudopotential, and a 4×4×4 Monkhorst-Pack k-point grid is used consistently in all simulations. Spin polarization is activated in all calculations. The formation energy for each optimized structure is then calculated as the difference between its total energy and the stoichiometric sum of the chemical potentials of its constituent elements.

The energetic, structural, and electric polarization properties of Fe-BNBT thin films were calculated with first-principles methods based on density functional theory (DFT). Due to its relatively small presence and to improve computational affordability, Ba ions were neglected in our system modeling thus we ended up simulating (Bi0.5Na0.5)TiO3 thim films. The PBE variant of the generalized gradient approximation to DFT27 was used as it is implemented in the VASP package28. The “projector augmented wave” method was employed to represent the ionic cores30, and the following electrons were considered as valence: Bi 6s, 6p, and 5d; Na 3s, 2p, and 2s; Ti 3p, 4s, and 3d; and O 2s and 2p. Wavefunctions were represented in a plane-wave basis truncated at 650 eV. We used a 40-atom simulation cell in which the usual ferroelectric and anti-ferrodistortive distortions occurring in oxide compounds could be reproduced31–34 as well as a chemical disorder in the cation positions. For integrations within the first Brillouin zone, we adopted a Gamma-centered k-point grid of 8x8x8. Geometry relaxations were performed by using a conjugate-gradient algorithm that optimized the volume and shape of the simulation cell as well as the atomic positions. The imposed tolerance on the atomic forces was 0.005 eV Å-1. By using these parameters we obtained total energies that converged to within 0.5 meV per formula unit (f.u.). The electric polarization, P, was computed with the Born-effective charges approach32,35,36 using the formula:

P = 1/V Σj Z*ijk × uk

where V is the volume of the simulation cell, index j runs over all the atoms in the simulation cell, indexes i and k represent Cartesian directions, uk is the displacement vector of the j-th atom as referred to as a non-polar reference phase, and Z*ijk the Born effective charge tensor calculated for a non-polar reference state. The Born effective charge tensor was estimated with first-principles density functional perturbation theory28.

Details of phase-field simulation

The phase-field simulation of BNBT is based on the Ginzburg-Landau-Devonshire theory. By solving the phase-field equation, we can obtain the temporal and spatial evolution of the system’s microscopic structure. In this study, the domain structure uses spontaneous polarization (P) as the order parameter. The energy minimization principle is applied to solve both the time-dependent Ginzburg-Landau equation and the mechanical equilibrium equation37–39:

∂Pi / ∂t = -L δF / δPi

∂(σij(r,t)) / ∂xj = 0

here, L is a kinetic coefficient related to domain wall mobility, F is the total free energy of the system, δF / δPj is the thermodynamic driving force for polarization evolution, σij is the stress tensor, and r and t represent the spatial coordinate and time, respectively. The total free energy of a bulk system can be defined as follows:

Ftotal = Fbulk(P) + Fgrad(P) + Felas(P) + Felec(P,E)

where, the total free energy Ftotal is composed of the bulk energy Fbulk(P), gradient energy Fgrad(P), elastic energy Felas(P), and electrical Felec(P,E), where E is the external electric field. The bulk free energy is expressed as a sixth-order polynomial expansion, as follows40:

Fbulk = α1(P1^2+P2^2+P3^2) + α11(P1^4+P2^4+P3^4) + α12(P1^2P2^2+P2^2P3^2+P1^2P3^2) + α111(P1^6+P2^6+P3^6) + α112(P1^4(P2^2+P3^2) + P2^4(P1^2+P3^2) + P3^4(P1^2+P2^2)) + α123(P1^2P2^2P3^2)

where, α, αij, αijk represent the Landau parameters, which can be obtained experimentally. The elastic energy density is given by:

Felas = 1/2 Cijkl(εij - ε0ij)(εkl - ε0kl)

here, Cijkl is the elastic stiffness tensor, εij is the total strain, and ε0ij is the spontaneous strain during the phase transformation. The spontaneous strain is related to the polarization through the electrostrictive coefficients, ε0ij = Qijkl Pk Pl, where Qijkl is the electrostrictive coefficient. The anisotropy factor is defined as A = C44 / C11. The gradient energy density is expressed as:

Fgrad = Gijkl Pi,j Pk,l

here, Gijkl is the gradient coefficient, and Pi,j denotes the derivative of the ith component of the polarization vector Pi with respect to the jth coordinate. In three-dimensional space, i,j=1,2,3. The electrostatic field energy felec in this system can be expressed as:

felec = -EiPi - 1/2 ε0 κij Ei Ej

where, Ei represent the electric field component, ε0 is the vacuum permittivity, and κij is the dielectric constant.

The parameters are provided as follows: α1 = 4.124×105(T-T0), α11 = -5×108(c-c0), α12 = 8.0×108(c-c0) + 2α11, α111 = 1.294×109(c-c0), α112 = 1.5×109(c-c0) + 3α111, α113 = 7.5×108(c-c0) + 6α111, Q11 = 0.1, Q12 = -0.034, Q44 = 0.029, s11 = 9.1×10-12, s12 = -3.2×10-10, s44 = 8.2×10-10, c = 0.06. These parameters are all in SI units. The compliance coefficients are represented by ‘s’, the Curie temperature is denoted by ‘T0’, and the doping level and critical doping level are represented by ‘c’ and ‘c0’, respectively. The simulation scale is 128 dx × 128 dy × 12 dz, with grid scales dx, dy, and dz set to 1 nm. The thickness of the substrate and air layers is set to 2 nm. The semi-implicit Fourier-spectral method is used to solve the equations.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "The manuscript \"Flexoelectricity-driven giant polarization in (Bi, Na)TiO3-based ferroelectric thin films\" by Sun and coworkers investgates the flexoelectric polarization of Fe-doped (Bi0.5, Na0.5)TiO3-based thin films. The authors identify oxygen-deficient structures as the local source of large polarization values and remanent polarization with a reasonable low coercive field. The authors further report a reduction of less than 9% after 10^7 cycles under 750 kV cm^-1 at 230 °C for the remanent and maximum polarization.\n\nThe origin of the large polarization values seems to be oxygen-vacancy-induced planar defects contributing probably to a large extend, to the stabilization of H-H domain structures necessary to reach and maintain these large values.\nOverall, an interesting manuscript but I have a couple of questions and comments.\n\nConcerning the film growth:\nWhat is the Fe-doping level and how did you verify it? XPS?\nHow did you arrive at this doping level and judged it to be optimal?\nYou are using quite a large fluence, how large is the area of your laser beam on the target?\nAlso, one should not use the term laser MBE because it is misleading even if the company selling your PLD system calls it this way. You still do PLD with the control capabilies originally available in MBE systems. In addition, the kinetic energy of your arriving species is and will be very different to thermally evaporated species, this is what MBE is about. The only \"correctly\" labelled laserbased MBE system the referee is aware of is the one developed at the Max-Planck Institut in Stuttgart and in the mean time commercially available. There CO2 lasers are used, not excimers operating in the UV.\n\nXRD: what do the authors mean with elongated lattice parameters? I am also a little puzzled about the lattice values you mention. LSMO growth with a tensile strain on STO, whereas BNBT should grow more or less strain free and the lattice constant is clearly larger than STO? Any explanation why? Looking at your RSMs, your film(s) must contain a lot of defects (your rocking curve must be very broad), but you have the same for the STO position. Is this a film contribution to the substrate or a substrate property because you induce defects in STO?\n\nXAS: the referee assumes, since it is not written, that the authors have recorded (and report) TEY data. The depth of information is only, at best, a few nm otherwise you probe surface/surface-near electronic properties. Identifying absorped water and hydroxyle molecules with TEY, fair enough, since I assume, a non-vacuum transfere to the XAS beamline has been done. How come that you generalize a very surface-near feature to the entire film structure? This cannot be correct. If you would like to do this, than show the corresponding TFY data and explain where the water in your films (d>500nm) is coming from. TFY data are more representative to bulk properties, you comment only on surface properties and deduce bulk properties?? What is the justification for it?\n\nTEM: As far as I can tell, and maybe I overlooked it, you forgot to state the obvious. The oxygen deficient structures you report are oriented along [110] and [1-10]. Is this correct? From a growth point of view, one would naivly expect a random distribution of defects. Is thermodynamics really the only explanation for it? I am not so certain. It has the look of twinning on a nano-scale.\n\nFor your calculations, you do not discuss the influence of growth properties like diffusion (length) of species. The saturation of species on a surface will for sure influence the species diffusion and you use a very large fluence meaning potentially a large number of species arriving at the substrate with a very wide kinetic energy distribution.\n\nYou use Fe-doped BNBT and observed these oxygen deficient structures supporting the H-H and T-T domains. What makes the Fe-doping so special? Is the density of these oxygen deficient structures doping level dependent? Would it also work without Fe-doping? From this point of view, the manuscript is very incomplete.\n\nIf an in-build electric field makes those structures \"visible\" or better recognizible as shown in Fig. 4 SI, where does this field comes from while imaging?\n\nAlso the word \"twilt\", I am not sure what it means in this context. I can see what you are using it for, but this is not the same. I know it in the context of (fish)bait or folds in a fabric. Looking at the structures you show, is there a more appropriate word for it? Those are 1D, oriented segments as far as I can see. Figure 3, caption: contract difference should probably read contrast difference\n\nSince you have a lot of these 1D defects, what is the consequence for leakage properties of these structures? Some of the PE-loops you show, albeit at 230°C, you could label as charge affected.\n\nWhat you show is potentially technically of interest, how reproducible is it? You do not discuss this point.\n\nI am also wondering, is the large polarization substrate dependent? STO is know to induce FE into other materials due to the very large dielectric constant, and if strained. How would the FE properties of Fe-doped BNBT look like if you grow it on LaAlO3 and ignore that this subtrate is twinned but nicely lattice matched. Or LSAT as an alterantive. Would you observe the oxygen deficient structures if you grow BNBT directly on STO?", "reviewer_2": "The study reported by Y. Sun et al. is focused on how planar defects coupled with oxygen vacancies can enhance electrical polarization in (Bi,Na)TiO3-based (BNT) thin films through the flexoelectric effect. By introducing iron dopants in 0.94(Bi0.5,Na0.5)TiO3-0.06BaTiO3 (BNBT) thin films, the researchers achieved a remarkable maximum polarization of 161 μC/cm² and a remnant polarization of 115 μC/cm² under a relatively low coercive field. With the help of HRTEM, they attributed this to the formation of charged domain walls and demonstrated that head-to-head (H-H) domain walls exhibit strain and strain gradients. The authors propose that local strain gradients enhance the thin-film ferroelectric properties through the flexoelectric effect, as rationalized by DFT and phase-field modeling.\n\nWhile this manuscript presents a quite comprehensive work with interesting findings, particularly the identification of strained regions in the films, the arguments provided by the authors are insufficient to explain the observed enhanced polarization. Firstly, the authors’ interpretation is directly related only to the vicinity of the defects associated with the H-H domain walls, and hence, to some, relatively small, fraction of the material volume. Therefore, if the defects do contribute to the increased polarization, these contributions should be weighted by their relative volumes. Secondly, after polarization switching, H-H domain walls should become T-T domain walls and acquire different properties because of the change of their charge sign. In turn, the oxygen vacancies should follow the change of the wall type quickly enough for the proposed mechanism to work. These aspects are not mentioned in the manuscript. Therefore, the authors’ picture that the observed effects are entirely a result of the flexoelectric effect is not substantiated. In fact, BNBT exhibits relaxor properties, and the used composition is close to a morphotropic phase boundary (MPB), where small variations in doping can significantly impact polarization (see, e.g., https://doi.org/10.1039/C7CP03530F). While the flexoelectric effect may contribute to the observed properties, it is likely only one of several contributing factors, and its contribution may be negligible.\n\nAdditionally, HRTEM measurements were conducted under “zero-field” conditions, while an applied electric field can induce a rhombohedral-to-tetragonal phase transition at the MPB. This raises questions about whether the observed effects would remain under applied fields. Furthermore, flexoelectric coupling at defects such as grain boundaries and domain walls is currently well known. However, it remains generally unclear whether enhanced polarization leads to larger strain or vice versa.\n\nSpecific Comments:\n• \"Super-tetragonal\" (Line 73) – This term should be introduced (defined and explained) in Introduction.\n• \"Facilitate the creation of oxygen vacancies, which accumulate in tail-to-tail (T-T) domain configurations\" (Line 95) -> Should it be “T-T domain walls”?\n• “Between the adjacent T-T domain walls, H-H domains emerge\"(Line 97). Please reformulate more carefully. Should it be “H-H domain walls”? The word “adjacent” is not correctly used in this sentence. It may be omitted.\n• The charge compensation of the T-T walls is due to accumulation of oxygen vacancies at the walls. How the charge compensation of the H-H walls is achieved? It is worth illustrating the oxygen vacancy formation with equations of defect chemistry.\n• The value of the ratio Aeg/t2g of 3.51 looks very large compared to the typical values of 1–1.5. This may suggest that the background has not been accurately subtracted. Please comment.\n• The statement “connected with the localized character of the eg states, demonstrating the existence of defective structures in Fe-doped BNBT thin films” (Line 159) – needs more explanation.\n• What is meant by \"bright twill\"? “Twill” is not a technical term but rather describes perception of the pattern in the image by the observer using some analogy. Therefore, it needs to be descriptively introduced before the use in the image analysis. The text is confusing as now and does not assist readers in following the authors’ logic. For example, are the twills 1D or 2D objects? To facilitate the readers’ work, it looks beneficial to move Figure S1 into the main text and use it to explain the meaning of the word “twills”.\n• The “twill” areas apparently possess polarization. Would it be more appropriate to understand the bright stripes as narrow, platelet-like, ferroelastic domains oriented along the normal to the image plane? The statement \"bright twills function as T-T domain walls in Fe-BNBT thin films\" (Line 192) is unclear. There appears to be a non-trivial polarization distribution within and outside the “bright twill” region. Are the boundaries between this region and neighboring domains uncharged?\n• The claim that \"T-T configurations are more energetically favorable for generating charged domain walls\" (Line 195) is puzzling, and generally, incorrect. T-T domain walls are charged, and charged domain walls are energetically unfavorable, requiring stabilization, which is kinetic in nature. For example, T-T domain walls carry negative polarization charge, which can be compensated by positively charged oxygen vacancies.\n• In Figure 3b, where is the yellow region located? Is it at the boundary or inside the “twill”? Could you provide a clearer depiction in the figure? Additionally, could you compare multiple line profiles from different locations (middle, boundary) inside the bright twill and the rest of the film, similar to panels (c) and (d)? The phrase \"intensity of oxygen atom signals\" is unclear. Please explain how you distinguish oxygen from other elements.\n• The statement \"H-H domain walls form between adjacent bright twills, acting as transition zones for T-T domain structures\" (Line 140) is unclear. A wall is a narrow transition zone between domains, but what is meant here is a relatively extended zone. What are “T-T domain structures”? In Figure 4, there is a region with head-to-head polarization. Inside this region, there is also in-plane polarization, which does not appear to be an artifact. Can this region be a domain?\n• The difference between regions like 1 and 2 in Figure 4 is due to different polarization orientations. The formation of polarization clusters with random orientations is common in relaxors.\n• The phase-field modeling suggests domains extending through the entire film thickness, but this is not observed experimentally. Would this discrepancy significantly affect the estimated polarization values?\n\nGrammar:\n--\"indicates enhanced tetragonally distortion\"\n--\"the area intensity of the eg and t2g peaks\"\n\nIn summary, this manuscript, while reporting interesting and potentially important experimental outcomes, provides no evidence that the achieved high values of polarization can be attributed solely to the flexoelectric effect. Given that this assumption is central to the paper and that the whole concept may need a revision, the manuscript cannot be recommended for publication in Nature Communications. The authors may consider publication in a more specialized journal.", "reviewer_3": "This study significantly enhanced the polarization of (Bi,Na)TiO₃-based ferroelectric thin films by introducing a defect structure with oxygen vacancies through Fe doping and leveraging the flexoelectric effect. The formation of H-H domain structures and a high strain gradient was confirmed via iDPC-STEM analysis, resulting in a maximum polarization of 161 μC cm⁻² and a remanent polarization of 115 μC cm⁻². Furthermore, the films exhibited excellent thermal stability, maintaining their performance even at 230°C, highlighting their potential for high-temperature electronic device applications.\n\nHowever, several revisions are required for publication.\n1. The core finding of this study is that the polarization of (Bi,Na)TiO₃-based ferroelectric thin films was maximized by leveraging the flexoelectric effect induced by oxygen vacancies. While the flexoelectric effect is identified as the primary cause of the polarization enhancement, the flexoelectric coefficients were estimated using values from existing literature rather than direct experimental measurement. To improve accuracy, the flexoelectric coefficient should be measured using AFM techniques, or the induced polarization should be further analyzed through XRD and Raman spectroscopy.\n2. This study demonstrates that the flexoelectric effect induces a significant increase in polarization; however, a quantitative separation of ferroelectric and flexoelectric polarization was not performed. To clarify the contribution of the flexoelectric effect, a more detailed P-E loop analysis should be conducted, including a comparison with conditions where the flexoelectric effect is not considered.\n3. This study mentions that oxygen vacancies play an important role in inducing the flexoelectric effect, but it lacks an explanation of how the polarization characteristics change with varying oxygen vacancy densities. It is necessary to analyze whether the flexoelectric effect increases as the number of oxygen vacancies increases or whether it saturates beyond a certain threshold.\n4. Since the flexoelectric effect is determined by the strain gradient, it is likely that the effect varies with film thickness. Although the experiment was primarily conducted on a 550 nm thick film, further research is needed to understand how the flexoelectric effect changes in thinner or thicker films.\n5. Although this study claims that Pm = 133 μC cm⁻² and Pr = 98 μC cm⁻² are maintained at 230°C, the specific reason for the high-temperature stability is not explained. Generally, problems such as oxygen vacancy movement, increased leakage current, and domain wall movement occur at high temperatures, and this study should explain how these factors were suppressed.\n6. The current study focused on Fe doping, but it is necessary to compare and analyze with previous studies how the flexoelectric effect changes when doping with other elements (Mn, Co, V, etc.).\n7. In general, typos (emperirical, etc.) and awkward expressions should be corrected.", "reviewer_4": "I co-reviewed this manuscript with one of the reviewers who provided the listed reports. This is part of the Nature Communications initiative to facilitate training in peer review and to provide appropriate recognition for Early Career Researchers who co-review manuscripts."}, "author_response": "Re: NCOMMS-25-12919, Title: Flexoelectricity-driven giant polarization in (Bi,\nNa)TiO -based ferroelectric thin films\nXun Geng, Ji Zhang, Jack Yang, Claudio Cazorla, Changqing Guo, Shery L. Y. Chang,\nXiaojie Lou, Houbing Huang, Zibin Chen, Shujun Zhang and Danyang Wang\nWe sincerely appreciate the referees’ insightful comments and valuable\nsuggestions on our manuscript entitled “Flexoelectricity-driven giant polarization in\n(Bi, Na)TiO -based ferroelectric thin films” (NCOMMS-25-12919A). We have\n3\nmeticulously revised the manuscript, with all modifications clearly highlighted. The\nclarifications and revisions made in accordance with the referees’ comments are\nprovided below.\nResponse to the comments of Referee #1\nThe manuscript \"Flexoelectricity-driven giant polarization in (Bi, Na)TiO3-\nbased ferroelectric thin films\" by Sun and coworkers investgates the flexoelectric\npolarization of Fe-doped (Bi0.5, Na0.5)TiO3-based thin films. The authors identify\noxygen-deficient structures as the local source of large polarization values and\nremanent polarization with a reasonable low coercive field. The authors further report\na reduction of less than 9% after 10^7 cycles under 750 kV cm^-1 at 230 °C for the\nremanent and maximum polarization.\nThe origin of the large polarization values seems to be oxygen-vacancy-induced\nplanar defects contributing probably to a large extend, to the stabilization of H-H\n1\ndomain structures necessary to reach and maintain these large values.\nOverall, an interesting manuscript but I have a couple of questions and comments.\nResponse: We sincerely thank the reviewer for the positive comments on our work.\nComments (1): Concerning the film growth: What is the Fe-doping level and how did\nyou verify it? XPS? How did you arrive at this doping level and judged it to be optimal?\nYou are using quite a large fluence, how large is the area of your laser beam on the\ntarget? Also, one should not use the term laser MBE because it is misleading even if\nthe company selling your PLD system calls it this way. You still do PLD with the control\ncapabilies originally available in MBE systems. In addition, the kinetic energy of your\narriving species is and will be very different to thermally evaporated species, this is\nwhat MBE is about. The only \"correctly\" labelled laser based MBE system the referee\nis aware of is the one developed at the Max-Planck Institut in Stuttgart and in the mean\ntime commercially available. There CO2 lasers are used, not excimers operating in the\nUV.\nResponse: The Fe-doping level in our BNBT thin films is 2 mol%, which is determined\nby the composition of the ceramic targets for PLD deposition. One of the hallmark\nadvantages of PLD lies in its ability to reproduce the stoichiometry of the multi-\ncomponent target material in the resulting thin films.1-4 This attribute is particularly\nbeneficial in applications requiring precise control of chemical compositions, making\nPLD a preferred technique for the growth of complex oxide thin films. Therefore, it is\nreasonable to assume that the deposited BNBT thin films exhibit a comparable Fe\nconcentration to that in the PLD target, e.g., 2 mol%. To address this query from the\nreviewer, we explicitly stated the doping level in the revised manuscript by adding the\nfollowing sentence in the Methods section:\n2\n“The Fe-doping level in our BNBT thin films is 2 mol%, which is determined by the\ncomposition of the ceramic targets for Pulsed Laser Deposition (PLD)”\nTo understand the effect of Fe doping level on the ferroelectric properties of\nBNBT-based thin films, we fabricated a series of thin films (350-nm-thick for a quicker\nscreening of doping level) with various Fe contents, i.e., 0.5, 1.0, 1.5, 2.0, 2.5, and 3.0\nmol%, by PLD. As shown in Figure R1, both the maximum and remanent polarization\nare optimal when the Fe doping level is 2 mol%. Notably, while the best ferroelectric\nperformance is observed at 2 mol% doping, a significant increase in polarization only\noccurs when the Fe content increases from 0.5 to 1 mol%, beyond which the growth\nrate of polarization greatly slows down till a Fe content of 2 mol% (Figure R1(h)). A\nhigher Fe doping level (2.5 and 3.0 mol%) incurs a drop in polarizations, as shown in\nFigure R1. This observation suggests that the Fe content already reaches an optimal\nlevel at 2.0 mol%. We have added these experimental results with detailed discussions\nto the Supplementary Information in the revised manuscript:\n“To understand the effect of Fe doping level on the ferroelectric properties of BNBT-\nbased thin films, we fabricated a series of thin films (350-nm-thick for a quicker\nscreening of doping level) with various Fe contents, i.e., 0.5, 1.0, 1.5, 2.0, 2.5, and 3.0\nmol% by PLD. As shown in Supplementary Figure 12, both the maximum and remanent\npolarization are optimal when the Fe doping level is 2 mol%.”\nloops of 350 nm-thick BNBT thin films with (a) 0, (b) 0.5, (c) 1.0, (d) 1.5, (e) 2.0, (f) 2.5, and\n(g) 3.0 mol% Fe doping levels measured under selected electric fields. (h) Evolution of\n4\nmaximum and remanent polarizations of the 350 nm-thick BNBT thin films as a function of Fe\ndoping level.\nBecause of a recent upgrade to the optical components in our PLD system, we\nhave re-measured the laser spot size, which was determined to be ~5.5 mm2.\nAccordingly, the laser fluence has been recalculated as ~1.45 J cm-2, which is\ncomparable to the values reported in previous studies.5, 6 In the original manuscript, the\nlaser fluence was estimated using an outdated spot size value provided by the lab\nmanager. This laser fluence has now been corrected in the revised manuscript.\nWe also greatly appreciate the reviewer’s clarification on the difference between\nPLD and MBE. We replaced the term “Laser MBE” with “PLD” in the revised\nmanuscript.\nComments (2): XRD: what do the authors mean with elongated lattice parameters? I\nam also a little puzzled about the lattice values you mention. LSMO growth with a\ntensile strain on STO, whereas BNBT should grow more or less strain free and the\nlattice constant is clearly larger than STO? Any explanation why? Looking at your\nRSMs, your film(s) must contain a lot of defects (your rocking curve must be very broad),\nbut you have the same for the STO position. Is this a film contribution to the substrate\nor a substrate property because you induce defects in STO?\nResponse: The term “elongated lattice parameter” in our manuscript particularly refers\nto the pronounced increase in c-axis lattice parameter, which often results in a \"super\ntetragonal\" structure in perovskite oxide thin films.7-9 The super tetragonal structure\nexhibits a very large c/a ratio, e.g., 1.238 in PbTiO thin films7, 1.12 in BaTiO thin\n3 3\nfilms8, 1.23 in BiFeO thin films10, and 1.25 in BNBT-based thin films9. Nevertheless,\n3\nthe measured c/a ratio of our Fe-BNBT thin films is ~1.004 (with a = 0.3905 nm and c\n5\n= 0.3919 nm, as stated in the manuscript), suggesting that a super tetragonal structure\nor a substantially elongated c-axis does not exist.\nTo avoid any confusion, we changed the term \"elongated lattice parameter\" to \"super-\nT having elongated c lattice parameter\" in the revised manuscript.\nRSM results clearly indicate that both LSMO and BNBT layers remain fully\nstrained on the STO substrates, as all the diffraction spots align to a specific Q position.\nx\nUnder appropriate growth conditions, perovskite oxide films (such as LSMO and\nBNBT) tend to adopt the in-plane lattice parameters of the substrate below a critical\nthickness without relaxation towards their respective bulk state. The bulk crystal\nstructures and lattice parameters of LSMO, BNBT, and the STO substrate, along with\nthe lattice parameters of the LSMO and BNBT layers epitaxially grown on the STO\nsubstrate, are summarized in Table R1. It is evident that both LSMO and BNBT layers\nexperience tensile in-plane strain when grown on STO, which is expected to result in a\nshortened out-of-plane lattice parameter relative to their respective bulk values.\nHowever, in addition to the epitaxial strain, oxygen vacancies can also lead to lattice\nexpansion in perovskite oxides.11 It is known that the PLD-produced oxide thin films\nare always oxygen deficient,3 as noted by the reviewer regarding the presence of defects\nin our films. Since the in-plane lattice of our thin film is constrained by the STO\nsubstrate, the lattice expansion primarily occurs along the c-axis, resulting in an\nelongated c lattice parameter (0.3919 nm in our thin films vs. 0.3903 nm in bulk BNBT).\nTable R1. Summary of lattice parameters for bulk LSMO, BNBT, and the STO substrate used\nin this study, as well as those of the LSMO and BNBT layers epitaxially grown on the STO\n(001) substrate.\nStructure\nLSMO (bulk)12 Pseudocubic 0.3876 0.3876 0.3876\nBNBT (bulk)13 Tetragonal 0.3885 0.3885 0.3903\nSTO substrate (single crystal)14 Cubic 0.3905 0.3905 0.3905\nLSMO layer (~15 nm thick) Tetragonal 0.3905 0.3905 0.3861\nBNBT layer (550 nm thick) Tetragonal 0.3905 0.3905 0.3919\nWe measured the rocking curves of our BNBT-based thin films, both undoped\nand doped with 2% mol of Mn and Fe, as shown in Figure R2. Although the Fe-BNBT\nthin films exhibit a larger full width at half maximum (FWHM) of 0.029° in contrast to\nthat of the Mn-doped counterpart, such FWHM remains still quite low compared to\nthose reported for typical perovskite oxide thin films,15, 16 suggesting the minimal\nimpact of the Fe-doping induced defects on the epitaxial quality of the films.\nThe identical in-plane lattice parameters observed among STO, Fe-BNBT, and\nLSMO do not result from the “strain” in the substrate imposed by the films or defects\ngenerated within the STO single crystal. Instead, this observation reflects the dominant\nrole of the single-crystal substrates in dictating the in-plane lattice parameters of the\nepitaxially grown overlayers. Because both bulk LSMO and BNBT have lattice\nparameters quite similar to those of STO, they adopt the in-plane lattice of the substrate\nduring epitaxial growth, as confirmed by the XRD-RSM patterns.\nundoped and doped with 2% mol of Mn and Fe. The rocking curve of the STO single crystal\nsubstrate is also shown for comparison purposes.\nWe hope the above clarifications have adequately addressed the reviewer’s\nconcerns.\nComments (3): XAS: the referee assumes, since it is not written, that the authors have\nrecorded (and report) TEY data. The depth of information is only, at best, a few nm\notherwise you probe surface/surface-near electronic properties. Identifying absorped\nwater and hydroxyle molecules with TEY, fair enough, since I assume, a non-vacuum\ntransfere to the XAS beamline has been done. How come that you generalize a very\nsurface-near feature to the entire film structure? This cannot be correct. If you would\n8\nlike to do this, than show the corresponding TFY data and explain where the water in\nyour films (d>500nm) is coming from. TFY data are more representative to bulk\nproperties, you comment only on surface properties and deduce bulk properties?? What\nis the justification for it?\nResponse: We sincerely appreciate the reviewer’s insightful suggestion. The Fe K-edge\nXAS data presented in our original study were acquired using the TFY mode. To further\naddress the reviewer’s concern, we conducted the XAS measurements for the O K-edge\nand Ti L-edge using the TFY mode. As illustrated in Figure R3, the TFY spectra are\nconsistent with the TEY spectra, confirming that the conclusions based on TEY data in\nthe original manuscript remain valid.\nAdditionally, no spectral features indicative of the presence of adsorbed water or\nhydroxyl species were detected in our XAS spectra.\nWe have clarified in the revised manuscript that the XAS measurements were\nconducted using the TFY mode. The revised manuscript has also been updated to\nincorporate the TFY-based XAS results, as shown in Figure R4.\nin total electron yield (TEY) mode. The (c) O K-edge and (d) Ti L-edge XAS spectra of Fe-\nBNBT thin films acquired in total fluorescence yield (TFY) mode. The insets in (b) and (d) are\nthe results of peak fitting of the L e peak of Fe BNBT thin films.\n3 g\nL-edge XAS spectra of Fe-BNBT thin films acquired in total fluorescence yield (TFY) mode.\nThe inset in (a) refers to the full spectra of the Fe K-edge. The inset in (c) is the result of peak\nfitting of the L e peak of Fe-BNBT thin films. The configuration interaction cluster model\n3 g\ncalculation is a computational technique used in quantum chemistry and solid-state physics to\nstudy the electronic structure of atomic systems, which is particularly valuable in systems with\nstrong electron correlation. (d) XPS spectra of the O 1s electronic level for Fe-BNBT thin films.\nBlue, green, and purple peak is ascribed to the lattice oxygen, lattice oxygen connected with\nthe lower valence cation after vacancy creation, and the adsorbed surface water and/or hydroxyl\nmolecules, respectively.\nComments (4): TEM: As far as I can tell, and maybe I overlooked it, you forgot to state\nthe obvious. The oxygen deficient structures you report are oriented along [110] and\n11\n[1-10]. Is this correct? From a growth point of view, one would naivly expect a random\ndistribution of defects. Is thermodynamics really the only explanation for it? I am not\nso certain. It has the look of twinning on a nano-scale.\nResponse: We thank the reviewer for flagging up this critical observation. It is true that\nthe oxygen-deficient structures line up in [101] and [-101] directions. We added the\nfollowing statements to the revised manuscript:\n\"Notably, numerous two-dimensional and oriented segments along [101] and [-101]\n(referred to as “bright streaks”) are observed in the Fe-BNBT thin films.\"\nAlthough one might expect a random distribution of oxygen vacancies\nthroughout the film due to thermodynamic effects, the actual scenario is more complex.\nAt low concentrations, oxygen vacancies initially appear as discrete point defects that\nare randomly distributed within the lattice.17 As the concentration of Vo increases, these\nvacancies begin to aggregate into finite-length chains. With further increase in the\nconcentration of Vo, the chains extend further and align along the [101]\ncrystallographic direction, forming an ordered structure in which layers of tetrahedrally\ncoordinated cations alternate with those of octahedrally coordinated cations, as\ndemonstrated in Figure R5.17, 18\nFigure Redacted\nFigure R5. Schematic illustration depicting the evolution of an oxygen vacancy chain along\nthe [101] direction, originating from an isolated vacancy and progressing with increasing\noxygen vacancy concentration. 17, 18\nTwinning in perovskites typically refers to the occurrence of two or more regions\n(called twin domains) within a crystal that are related by a specific crystallographic\nsymmetry operation—such as a mirror plane, rotation, or inversion—but are not related\n12\nby simple translation. These domains share a coherent interface known as the twin\nboundary.19-22 This structural characteristic is schematically illustrated in Figure R6.\nBased on our experimental observations and the nature of the identified oxygen-\ndeficient regions, it is more appropriate for these features to be described as planar\ndefects rather than conventional twinning structures.\nFigure redacted\nFigure R6. (a) Schematic drawing of a coherent (111) twin in BaTiO in the [110] projection.\n3\n(b) Experimental lattice images of the (111) twin boundary in BaTiO as seen in the [110]\n3\nprojection.22\nComments (5): For your calculations, you do not discuss the influence of growth\nproperties like diffusion (length) of species. The saturation of species on a surface will\nfor sure influence the species diffusion and you use a very large fluence meaning\npotentially a large number of species arriving at the substrate with a very wide kinetic\nenergy distribution.\nResponse: As discussed in our response to Comment (1), the laser fluence is\nrecalculated to be ~1.45 mm2, which is consistent with those reported in previous\nstudies.5, 6\nWe agree with the reviewer that the laser fluence undoubtedly affects diffusion\nbehaviours of the species on the substrate, thereby influencing the structural and\nfunctional properties of the resulting films. As a critical processing parameter, laser\nfluence can be strategically adjusted to manipulate the microstructural and functional\nproperties of thin films. For instance, increasing the laser fluence from ~1.3 to 7.6 J cm-\n2 has been reported to enhance the c/a ratio from 1.06 to 1.10 and increase the remanent\npolarization from 70 to 100 µC cm-2 in PbZr Ti O films, attributed to the generation\n0.2 0.8 3\n13\nof charged defects induced by the elevated laser fluence during pulsed laser\ndeposition.23\nComments (6): You use Fe-doped BNBT and observed these oxygen deficient structures\nsupporting the H-H and T-T domains. What makes the Fe-doping so special? Is the\ndensity of these oxygen deficient structures doping level dependent? Would it also work\nwithout Fe-doping? From this point of view, the manuscript is very incomplete.\nResponse: BNBT films with various Fe contents, i.e., 0.5, 1.0, 1.5, 2.0, 2.5, and 3.0\nmol%, were grown on LSMO-electroded STO (001) single crystal substrates using PLD\nunder the same deposition conditions. Both maximum and remanent polarization are\npeaked at 2.0 mol% Fe doping, as detailed in the response to Comment (1). X-ray\nphotoelectron spectroscopy (XPS) was employed to estimate the level of oxygen\nvacancy in Fe-BNBT thin films, as shown in Figure R7. The relative level of oxygen\ndefects ( ) in the sample could be estimated by the:\n𝑂𝑉%\n𝐴𝑂𝑉\n𝑂𝑉% =\n𝐴𝑂𝐿 +𝐴𝑂𝑉 +𝐴𝑂𝐶\nwhere , and are the areas enclosed under O , O and O peak,\nL V, C\nrespecti𝐴v 𝑂 el 𝐿 y. 𝐴T 𝑂 h 𝑉 e XPS𝐴 𝑂 re 𝐶 sults reveal a progressive increase in oxygen vacancy\nconcentration with the increase in Fe doping level. Nevertheless, an excessive amount\nof oxygen vacancy (in 2.5 and 3.0 mol% Fe doped samples) will deteriorate the\nferroelectric properties, as shown in Figure R1, owing to the more pronounced domain-\nwall pinning effect24 and increased leakage current25.\n14\nelectronic level for Fe-BNBT thin films with (a) undoped and (b) 0.5, (c) 1.0, (d) 1.5, (e) 2.0,\n(f) 2.5, and (g) 3.0 mol% Fe doping levels. Blue (O ), green (O ), and purple (O ) peak is\nL V C\nascribed to the lattice oxygen, lattice oxygen connected with the lower valence cation after\nvacancy creation, and the adsorbed surface water and/or hydroxyl molecules, respectively. (h)\nEvolution of the O level in BNBT-based thin films as a function of Fe doping level.\nV\nThe undoped BNBT thin films exhibited a maximum polarization of 62 cm-2,\nwhereas 2.0 mol% Fe-doped films showed a significantly enhanced P of 161 µC cm-\nm\n2. HAADF STEM images confirm the absence of planar defects in the undoped film, as\nshown in Figure R8(d), suggesting that the oxygen-deficient structures originate from\nthe introduction of Fe dopants and the critical role of the oxygen-deficient structures in\nboosting the polarization.\nloops of (a) 2.0 mol% Fe-doped BNBT thin films and (b) undoped BNBT thin films. HAADF\nSTEM images of (c) 2.0 mol% Fe-doped BNBT thin films and (d) undoped BNBT thin films.\nThese results unequivocally suggest that Fe doping plays a decisive role in\nenhancing the polarization of BNBT-based thin films. We attribute the giant\npolarization of Fe-doped BNBT thin films to several fundamental mechanisms as\nfollows:\n1) Fe3+ (r = 58.5 pm) preferentially replaces Ti4+ (r = 60.5 pm) at the B-site of\nFe3+ Ti4+\nthe perovskite lattice, as an acceptor, leading to the formation of oxygen vacancies to\ncompensate for the charge imbalance. These vacancies tend to organize along specific\n17\ncrystallographic planes,17, 18 as detailed in the response to Comment (4). In contrast,\nother commonly employed dopants in BNBT lattice such as La3+ (r = 136 pm), Sm3+\nLa3+\n(r = 124 pm), and Eu3+ (r = 129 pm) substitute for A-site cations (r = 139 pm,\nSm3+ Eu3+ Na+\nr = 120 pm, r = 135 pm) as donors, due to their similar ionic radii to that of the\nBi3+ Ba2+\nA-site cations, which do not promote the formation of oxygen vacancies.26\n2) The multivalent nature of Mn (Mn2+, Mn3+, Mn4+) in principle allows them to enter\ndifferent lattice sites of BNBT. Considering the similarity of the ionic radii between\nMn cations (r ~ 0.83 Å, r ~ 0.65 Å, r ~ 0.53 Å) and the respective cations\nMn2+ Mn3+ Mn4+\nin the host lattice, as well as the global instable index indices calculated in our previous\nwork27, 28, Mn2+ prefers to reside in the A-site, and Mn3+ and Mn4+ are inclined to\nsubstitute Ti4+ at the B-site. It should be noted that the presence of Mn2+ often impedes\nthe formation of oxygen vacancies.29 More importantly, Mn dopants tend to introduce\nlocalized point defects rather than extended planar faults.30\n3) Co3+ (r = 55 pm) can also substitute for Ti4+ at the B-site and is reported to be\nCo3+\nmore effective than Fe in promoting the creation of oxygen vacancies.31 However, this\nincreased oxygen vacancy concentration is often accompanied by an increased\nelectrical conductivity, which may incur leakage current pathways and deteriorate\nferroelectric performance, as shown in Figure R19.\nWe added these experimental results and the corresponding discussions to the\nSupplementary Information of the revised manuscript:\n“X-ray photoelectron spectroscopy (XPS) was employed to estimate the level of oxygen\nvacancy in Fe-BNBT thin films, as shown in Supplementary Figure 13. The relative\nlevel of oxygen defects ( ) in the sample could be estimated by the:\n𝑂𝑉%\n18\n𝐴𝑂𝑉\n𝑂𝑉% =\n𝐴𝑂𝐿 +𝐴𝑂𝑉 +𝐴𝑂𝐶\nwhere , and are the areas enclosed under O , O and O peak,\nL V, C\nrespectively. 𝐴T 𝑂 h 𝐿 e 𝐴X 𝑂 P 𝑉 S res𝐴u 𝑂 lt 𝐶 s reveal a progressive increase in oxygen vacancy\nconcentration with the increase in Fe doping level. Nevertheless, an excessive amount\nof oxygen vacancy (in 2.5 and 3.0 mol% Fe doped samples) will deteriorate the\nferroelectric properties, as shown in Supplementary Figure 12, owing to the more\npronounced domain-wall pinning effect and increased leakage current.\nThe undoped BNBT thin films exhibited a maximum polarization of 62.1 cm-2,\nwhereas 2.0 mol% Fe-doped films showed a significantly enhanced P of 161 µC cm-\nm\n2. HAADF STEM images confirm the absence of planar defects in the undoped film, as\nshown in Supplementary Figure 14, suggesting that the oxygen-deficient structures\noriginate from the introduction of Fe dopants and the critical role of the oxygen-\ndeficient structures in boosting the polarization.\nThese results unequivocally suggest that Fe doping plays a decisive role in\nenhancing the polarization of BNBT-based thin films. We attribute the giant\npolarization of Fe-doped BNBT thin films to several fundamental mechanisms as\nfollows:\n1) Fe3+ (r = 64.5 pm) preferentially replaces Ti4+ (r = 60.5 pm) at the B-site of\nFe3+ Ti4+\nthe perovskite lattice, as an acceptor, leading to the formation of oxygen vacancies to\ncompensate for the charge imbalance. These vacancies tend to organize along specific\ncrystallographic planes, as detailed in the response to Comment (4). In contrast, other\ncommonly employed dopants in BNBT lattice such as La3+ (r = 136 pm), Sm3+ (r\nLa3+ Sm3+\n= 124 pm), and Eu3+ (r = 129 pm) substitute for A-site cations (r = 139 pm, r\nEu3+ Na+ Bi3+\n19\n= 120 pm, r = 135 pm) as donors, due to their similar ionic radii to that of the A-site\nBa2+\ncations, which do not promote the formation of oxygen vacancies.\n2) The multivalent nature of Mn (Mn2+, Mn3+, Mn4+) in principle allows them to enter\ndifferent lattice sites of BNBT. Considering the similarity of the ionic radii between\nMn cations (r ~ 0.83 Å, r ~ 0.65 Å, r ~ 0.53 Å) and the respective cations\nMn2+ Mn3+ Mn4+\nin the host lattice, as well as the global instable index indices calculated in our previous\nwork, Mn2+ prefers to reside in the A-site, and Mn3+ and Mn4+ are inclined to substitute\nTi4+ at the B-site. It should be noted that the presence of Mn2+ often impedes the\nformation of oxygen vacancies. More importantly, Mn dopants tend to introduce\nlocalized point defects rather than extended planar faults.\n3) Co3+ (r = 55 pm) can also substitute for Ti4+ at the B-site and is reported to be\nCo3+\nmore effective than Fe in promoting the creation of oxygen vacancies. However, this\nincreased oxygen vacancy concentration is often accompanied by an increased\nelectrical conductivity, which may incur leakage current pathways and deteriorate\nferroelectric performance, as shown in Supplementary Figure 15.”\nWe hope this response sufficiently clarifies the unique role of Fe doping in\nstabilizing oxygen-deficient structures and enhancing ferroelectric properties of BNBT-\nbased thin films.\nComments (7): If an in-build electric field makes those structures \"visible\" or better\nrecognizible as shown in Fig. 4 SI, where does this field comes from while imaging?\nResponse: Supplementary Figure 4 in the original manuscript is phase-field simulation\nresults, in which the in-built electric field is not specified. In our phase-field simulations,\nthe bound charges are intentionally introduced at domain walls to mimic the\nexperimentally observed domain structures. The in-built electric field in the films may\n20\noriginate from localized charged defects accumulation. Specifically, induced oxygen\nvacancies are preferentially aligned along the [101] crystallographic direction,17, 18,\nwhich facilitates the formation of T-T domain walls.32 The oxygen vacancies are\nelectropositive in nature, thereby contributing to the development of the built-in electric\nfield.\nComments (8): Also the word \"twilt\", I am not sure what it means in this context. I can\nsee what you are using it for, but this is not the same. I know it in the context of (fish)bait\nor folds in a fabric. Looking at the structures you show, is there a more appropriate\nword for it? Those are 1D, oriented segments as far as I can see. Figure 3, caption:\ncontract difference should probably read contrast difference.\nResponse: We appreciate the reviewer's comment on the critical terms we used. Based\non the morphology and orientation of these structures, we have replaced \"bright twills\"\nwith \"bright streaks\" throughout the manuscript to better reflect their two-dimensional\nand oriented nature.\nAdditionally, we changed \"contract difference\" to \"contrast difference\" in the\ncaption of Figure 3.\nComments (9): Since you have a lot of these 1D defects, what is the consequence for\nleakage properties of these structures? Some of the PE-loops you show, albeit at 230°C,\nyou could label as charge affected.\nResponse: While oxygen vacancies and planar defects may potentially incur high\nleakage currents, our Fe-BNBT thin films render a relatively low leakage density over\nthe entire measuring temperature range, as shown in Figure R9. Even at 250°C, the\nleakage current density remains as low as ~10-4 A cm-2 under an applied field of 750\n21\nkV cm-1, which is considered minimal for typical ferroelectric thin films possessing a\nperovskite structure.33-35 This suggests that the induced defects do not significantly\ncompromise the insulating properties of our thin films.\nFigure R9 (Supplementary Figure 7 in the revised manuscript). Leakage current density as\na function of the applied electric fields of Fe-BNBT thin films at various temperatures.\nIn dielectric thin films, the charge transport is often dominated by bulk-limited\nconduction mechanisms, including ohmic conduction, space-charge-limited conduction\n(SCLC), and Poole-Frenkel emission (PF), or interface-limited Schottky emission\n(SE).36 The space-charge-limited current obeys quadratic electric field dependence in\nthe case of discrete traps,37\n2\n𝐽𝑆𝐶𝐿𝐶 = 9𝜇ɛ0ɛ𝑟𝐸 /8𝑑\nwhere μ is the carrier mobility and d the film thickness. The current density across a\nSchottky barrier is37\n1\n∗ 2 22 2\n𝐽𝑆 = 𝐴 𝑇 exp− [(𝜙𝑏 − 𝛽𝑠𝐸 )/𝑘𝑇]\nwhere , is the Richardson constant, T the temperature, q the\n3 1/2 ∗\nelectro𝛽n 𝑠 c=ha(r𝑞ge,/ 4𝜋ɛ 0 thɛe 𝑟 )Scho𝐴ttky barrier height, k the Boltzmann constant, the\npermittivity of fre𝜙e 𝑏 space, and the dielectric constant of the thin film. The ɛP 0 ool-\nFrenkel conductivity is given\nbyɛ37𝑟\n1\n2\nwhere 𝐽𝑃𝐹 =, 𝜎0 i𝐸s ea xspam−p [l(e𝐸-d 1 e−pe n𝛽d 𝑃 e 𝐹 n𝐸t z)e/r𝑘o𝑇-f]ield conductivity, and\n3 1/2\nis the tr𝛽a 𝑃 p 𝐹 io=ni(z𝑞ati/o𝜋nɛ e 0 nɛe 𝑟 r)gy. T𝜎h 0 e functional relationship of SE and PF is similar so th𝐸a 1 t\nthe plots of ln(J T-2) versus E0.5 and ln(J V-1) versus E0.5 for SE and PF, respectively,\nwill yield linear relationships.\nWe plotted and fitted the leakage current data based on the above-mentioned\ntransport mechanisms, as shown in Figure R10. At room temperature, the leakage\ncurrent is governed by PF under the electric field < 150 kV cm-1, given the relative\ndielectric constant extracted from the slope of the linear fits to PF is 118.8, which lies\nwithin the physically reasonable range between the optical dielectric constant (6.71)38\nand quasi-static dielectric constant (500- 800)39 typical for BNT-based thin films.\nApart from the low field region (< 150 kV cm-1), the Schottky emission\nmechanism becomes dominant, supported by the extracted dielectric constant of 52.8,\nwhich is also physically plausible. SCLC can be ruled out in the entire range of\nmeasuring temperature and electric field, given slope of the LnJ-LnE curves evidently\ndeviates from 2, as shown in Figure R10(a). At 250 °C, the leakage current follows the\nPF mechanism across the entire range of the measuring electric field.\ncurrent mechanisms in Fe-BNBT thin films based on (a) SCLC, (b) SE, and (c) PF emissions.\nRegarding the P–E loops measured at elevated temperatures, e.g., 230 °C, we\nacknowledge that thermal activation of charges will lead to a rise in leakage current\ndensity, potentially influencing the measured polarizations. However, the well-defined\nP-E loops with decent fatigue resistance at elevated temperatures (Figure 6e in the\noriginal manuscript), along with the clearly resolved temperature-dependent remanent\nP-E loops (Figure 6b in the original manuscript), suggest that the effect of thermally\nactivated charges on the polarization response remains modest. The relatively low\nleakage current density can be attributed to the trapping of oxygen vacancies and free\ncharges at T-T and H-H domain walls, respectively. These domain walls act as potential\nbarriers that impede charge transport, thereby suppressing the formation and growth of\nelectrical trees and effectively reducing leakage currents, particularly at high\ntemperatures.40\nWe added these experimental results and the corresponding discussion to the\nSupplementary Information of the revised manuscript.\n“The leakage current density vs electric field curves at various temperatures of our Fe-\nBNBT thin films were shown in Supplementary Figure 7. The leakage current density\nremains as low as ~10-4 A cm-2 even at elevated temperatures, which is considered\n24\nminimal for typical ferroelectric thin films possessing a perovskite structure.7-9 The\nrelatively low leakage current density can be attributed to the trapping of oxygen\nvacancies and free charges at T-T and H-H domain walls, respectively. These domain\nwalls act as potential barriers that impede charge transport, thereby suppressing the\nformation and growth of electrical trees and effectively reducing leakage currents,\nparticularly at high temperatures.10\nIn dielectric thin films, the charge transport is often dominated by bulk-limited\nconduction mechanisms, including ohmic conduction, space-charge-limited conduction\n(SCLC), and Poole-Frenkel emission (PF), or interface-limited Schottky emission\n(SE).11 The space-charge-limited current obeys quadratic electric field dependence in\nthe case of discrete traps,12\n2\n𝐽𝑆𝐶𝐿𝐶 = 9𝜇ɛ0ɛ𝑟𝐸 /8𝑑\nwhere μ is the carrier mobility and d the film thickness. The current density across a\nSchottky barrier is12\n1\n∗ 2 2\n𝐽𝑆 = 𝐴 𝑇 exp− [(𝜙𝑏 − 𝛽𝑠𝐸 )/𝑘𝑇]\nwhere , is the Richardson constant, T the temperature, q\n3 1/2 ∗\nthe electron c𝛽h 𝑠 a=rge(,𝑞 /4 t𝜋heɛ 0 Sɛc 𝑟 h)ottky𝐴 barrier height, k the Boltzmann constant, the\npermittivity of free s𝜙p 𝑏 ace, and the dielectric constant of the thin film. The ɛP 0 ool-\nFrenkel conductivity is given byɛ 𝑟\n1\n2\n𝐽𝑃𝐹 = 𝜎0𝐸exp− [(𝐸1− 𝛽𝑃𝐹𝐸 )/𝑘𝑇]\nwhere , is a sample-dependent zero-field conductivity, and\n3 1/2\nis the tr𝛽a 𝑃 p 𝐹 io=ni(z𝑞ati/o𝜋nɛ e 0 nɛe 𝑟 r)gy. T𝜎h 0 e functional relationship of SE and PF is similar, so th𝐸a 1 t\n25\nthe plots of ln(J T-2) versus E0.5 and ln(J V-1) versus E0.5 for SE and PF, respectively,\nwill yield linear relationships.\nWe plotted and fitted the leakage current data based on the above-mentioned\ntransport mechanisms, as shown in Supplementary Figure 8. At room temperature, the\nleakage current is governed by PF under the electric field < 150 kV cm-1, given the\nrelative dielectric constant extracted from the slope of the linear fits to PF is 118.8,\nwhich lies within the physically reasonable range between the optical dielectric constant\n(6.71) and quasi-static dielectric constant (500- 800) typical for BNT-based thin films.\nApart from the low field region (< 150 kV cm-1), the Schottky emission\nmechanism becomes dominant, supported by the extracted dielectric constant of 52.8,\nwhich is also physically plausible. SCLC can be ruled out in the entire range of\nmeasuring temperature and electric field, given slope of the LnJ-LnE curves evidently\ndeviates from 2, as shown in Supplementary Figure 8(a). At 250 °C, the leakage current\nfollows the PF mechanism across the entire range of the measuring electric field.”\nComments (10): What you show is potentially technically of interest, how reproducible\nis it? You do not discuss this point.\nResponse: We confirm the giant polarization in our Fe-doped BNBT thin films is\nhighly reproducible, as demonstrated by the consistent P-E loops (Figure R11) obtained\nfrom five samples prepared under optimal processing conditions.\nTo further strengthen the manuscript, we have added a brief discussion on the\nreproducibility of our results in the revised manuscript:\n“To confirm the reproducibility of the giant polarization, we fabricated multiple Fe-\ndoped BNBT thin films under optimized processing conditions. The P–E loops of five\nconsistent ferroelectric properties observed across all samples clearly demonstrate the\nexcellent reproducibility of the giant polarization in our Fe-doped BNBT thin films.”\nFigure R11 (Supplementary Figure 4 in the revised manuscript). (a)-(e)Room-temperature\nP–E loops obtained from five Fe-BNBT thin films prepared under optimal processing\nconditions. (f) Comparison of the P–E loops of five Fe-BNBT thin films under 750 kV cm-1.\nComments (11): I am also wondering, is the large polarization substrate dependent?\nSTO is know to induce FE into other materials due to the very large dielectric constant,\nand if strained. How would the FE properties of Fe-doped BNBT look like if you grow\nit on LaAlO3 and ignore that this subtrate is twinned but nicely lattice matched. Or\nLSAT as an alterantive. Would you observe the oxygen deficient structures if you grow\nBNBT directly on STO?\nResponse: To further investigate the effect of substrate on the polarization of the\nresulting films, we fabricated 550-nm-thick 2 mol% Fe-doped and Mn-doped BNBT\nthin films on LSMO-electroded (LaAlO ) (Sr AlTaO ) (LSAT) (001) single crystal\n3 0.3 2 6 0.7\nsubstrates (cubic, a = 3.868 Å) using PLD under the identical deposition conditions to\n27\nthose grown on STO substrates. As shown in Figure R12(a) and (b), although the\nferroelectric polarizations of both Fe-doped and Mn-doped BNBT films grown on\nLSAT are slightly lower compared to their counterparts grown on STO, the Fe-BNBT\nfilm still exhibits a notably high P of 145 µC cm-2 and a P of 96.7 µC cm-2. These\nm r\nvalues remain substantially higher than those of the Mn-BNBT films on LSAT (P\nm\n=49.1 µC cm-2, P =19.4 µC cm-2). These results suggest that the large polarization of\nr\nour Fe-BNBT thin films is not strongly dependent on the substrates.\n2 mol% Fe-doped BNBT thin films were also directly deposited on Nb-doped\nSTO (Nb:STO) (001) single crystal substrates. As shown in Figure R12(c), the\nconsiderably large polarizations (P = 110 µC cm-2 and P =157 µC cm-2) are highly\nm r\ncomparable to those of samples grown on LSMO-electroded STO (001) substrates,\nstrongly underpinning the presence of the oxygen-deficient structures.\nloops of (a) Fe-BNBT thin films and (b) Mn-BNBT thin films deposited on LSAT (001) single\ncrystal substrates, and (c) Fe-BNBT thin films deposited on Nb:STO (001) substrates. The inset\nin (b) is the P-E loops of Mn-BNBT deposited on STO (001) substrate. (d) Comparison of the\nP–E loops of Fe-BNBT thin films deposited on different substrates.\nWe added these experimental results and the corresponding discussion to the\nSupplementary Information of the revised manuscript.\n“To further investigate the effect of substrate on the polarization of the resulting films,\nwe fabricated 550-nm-thick 2 mol% Fe-doped and Mn-doped BNBT thin films on\nLSMO-electroded (LaAlO ) (Sr AlTaO ) (LSAT) (001) single crystal substrates\n3 0.3 2 6 0.7\n(cubic, a = 3.868 Å) using PLD under the identical deposition conditions to those grown\non STO substrates. As shown in Supplementary Figure 11(a) and (b), although the\nferroelectric polarizations of both Fe-doped and Mn-doped BNBT films grown on\nLSAT are slightly lower compared to their counterparts grown on STO, the Fe-BNBT\nfilm still exhibits a notably high P of 145 µC cm-2 and a P of 96.7 µC cm-2. These\nm r\nvalues remain substantially higher than those of the Mn-BNBT films on LSAT (P\nm\n=49.1 µC cm-2, P =19.4 µC cm-2). These results suggest that the large polarization of\nr\nour Fe-BNBT thin films is not strongly dependent on the substrates.\n2 mol% Fe-doped BNBT thin films were also directly deposited on Nb-doped\nSTO (Nb:STO) (001) single crystal substrates. As shown in Supplementary Figure\n11(c), the considerably large polarizations (P = 110 µC cm-2 and P =157 µC cm-2) are\nm r\nhighly comparable to those of samples grown on LSMO-electroded STO (001)\nsubstrates, strongly underpinning the presence of the oxygen deficient structures.”\nThe study reported by Y. Sun et al. is focused on how planar defects coupled with\noxygen vacancies can enhance electrical polarization in (Bi,Na)TiO3-based (BNT) thin\nfilms through the flexoelectric effect. By introducing iron dopants in\n0.94(Bi0.5,Na0.5)TiO3-0.06BaTiO3 (BNBT) thin films, the researchers achieved a\nremarkable maximum polarization of 161 μC/cm² and a remnant polarization of 115\nμC/cm² under a relatively low coercive field. With the help of HRTEM, they attributed\nthis to the formation of charged domain walls and demonstrated that head-to-head (H-\nH) domain walls exhibit strain and strain gradients. The authors propose that local\nstrain gradients enhance the thin-film ferroelectric properties through the flexoelectric\neffect, as rationalized by DFT and phase-field modeling.\nWhile this manuscript presents a quite comprehensive work with interesting findings,\nparticularly the identification of strained regions in the films, the arguments provided\nby the authors are insufficient to explain the observed enhanced polarization.\nResponse: We sincerely appreciate the reviewer’s insightful feedback.\nIn addition to addressing the specific concerns the reviewer raised, we have also\nsupplemented the experimental results pertaining to the effects of dopant species,\ndoping level, substrates, and film thickness on the flexoelectric-enhanced polarization\nof BNBT-based thin films in the revised manuscript. We hope our additional results\nand revisions not only substantially strengthen our arguments in the manuscript but also\nimprove the clarity, depth, and scientific rigor of our conclusions.\nComments (1): Firstly, the authors’ interpretation is directly related only to the vicinity\nof the defects associated with the H-H domain walls, and hence, to some, relatively\nsmall, fraction of the material volume. Therefore, if the defects do contribute to the\n30\nincreased polarization, these contributions should be weighted by their relative\nvolumes.\nResponse: It is important to note that although the bright streaks observed in the STEM\nimages exhibit a two-dimensional contrast, i.e., linear features in the projection, the\nassociated planar defects extend beyond these streaks into the adjacent lattice, spanning\nover approximately 5 unit cells. This broader spatial distribution is confirmed by the\nline profiles of oxygen atom signal intensities taken from the eight sequentially adjacent\nregions, delineated by blue to green dashed rectangles, as shown in Figure R13.\nFigure R13 (Supplementary Figure 9 in the revised manuscript). (a) HAADF-STEM image\nof Fe-BNBT thin films. (b) STEM-iDPC image of the region outlined by the red dashed square\nin (a). The yellow dashed line indicates the position of the bright streak, while the series of\ncolored dashed rectangles represents eight sequentially adjacent regions extending outward\nfrom the bright streak. (c) Line profiles of oxygen atom signal intensities extracted from the\nregions marked by the blue to green dashed rectangles in (b), illustrating a progressive decrease\n31\nin oxygen vacancy concentration from the bright streak toward the surrounding lattice. The blue\narrow mark the positions of oxygen vacancies.\nBased on this analysis, the volume fraction of oxygen vacancy-associated planar\ndefects can be estimated from the STEM images (Figure R14) using the line intercept\nmethod, which is a common approach to evaluate the dislocation density in alloys.41, 42\nThe volume fraction is calculated using the following equation:\n𝑉𝑑𝑒𝑓𝑒𝑐𝑡\n𝑜\n𝐴𝑑𝑒𝑓𝑒𝑐𝑡 𝐿𝑠𝑡𝑟𝑒𝑎𝑘 ×𝐴𝑢𝑛𝑖𝑡 𝑐𝑒𝑙𝑙 ×8×sin45\n𝑉𝑑𝑒𝑓𝑒𝑐𝑡 = =\n𝐴𝑓𝑖𝑙𝑚 𝐴𝑓𝑖𝑙𝑚\nwhere and represent the total projected area of the planar defects and the\ncross-s𝐴ec 𝑑 t 𝑒 i 𝑓 o 𝑒 n 𝑐 a 𝑡 l area𝐴 𝑓 o 𝑖𝑙 f 𝑚 the thin film, respectively. denotes the cumulative\nlength of the bright streaks, and corresponds 𝐿to 𝑠𝑡 𝑟 th 𝑒𝑎 e 𝑘 projected area of a single\nFe-BNBT unit cell. The volume f𝐴ra 𝑢 c 𝑛 t 𝑖 i 𝑡 o 𝑐 n 𝑒 𝑙 o 𝑙 f the planar defects in the Fe-BNBT thin films\nis estimated to be ~16%.\nof Fe-BNBT thin films, with orange-shaded regions indicating the presence of planar defects.\nAs noted in the manuscript, the flexoelectric coefficient for perovskite-structured\nsingle crystals typically ranges from 10-1 to 101 nC m-1. Within this range, the estimated\nflexoelectric polarization is up to 1000 μC cm-2. Considering the calculated density of\nplanar defects above, the expected contribution of the flexoelectric effect to the\npolarization in our Fe-BNBT films amounts to ~160 μC cm-2. The experimentally\nmeasured flexoelectric polarization in Fe-BNBT thin films is determined to be ~100.2\nμC cm-2, obtained by subtracting the polarization of typical BNBT (60.9 μC cm-2) thin\nfilms from that of Fe-BNBT thin films (161.1 μC cm-2). This result strongly underpins\nour conclusion that the flexoelectric effect, driven by oxygen vacancy-induced planar\ndefects, is a dominant contributor to the observed giant polarization in our Fe-BNBT\nthin films.\nWe added these experimental results and the corresponding discussion to the\nSupplementary Information in the revised manuscript. We hope the above clarifications\nhave effectively addressed the reviewer’s concern.\n“Although the bright streaks observed in the STEM images exhibit a one-dimensional\ncontrast, i.e., linear features in the projection, the associated planar defects extend\nbeyond these streaks into the adjacent lattice, spanning over approximately 5 unit cells.\nThis broader spatial distribution is confirmed by the line profiles of oxygen atom signal\nintensities taken from the eight contiguous regions, delineated by blue to green dashed\nrectangles, as shown in Supplementary Figure 9. Based on this analysis, the volume\nfraction of oxygen vacancy-associated planar defects can be estimated from the STEM\nimages (Supplementary Figure 10) using the line intercept method, which is a common\n33\napproach to evaluate the dislocation density in alloys.13,14 The volume fraction\nIt is calculated using the following equation:\n𝑉𝑑𝑒𝑓𝑒𝑐𝑡\n𝐴𝑑𝑒𝑓𝑒𝑐𝑡 𝐿𝑠𝑡𝑟𝑒𝑎𝑘 ×𝐴𝑢𝑛𝑖𝑡 𝑐𝑒𝑙𝑙 ×8×sin45°\n𝑉𝑑𝑒𝑓𝑒𝑐𝑡 = =\n𝐴𝑓𝑖𝑙𝑚 𝐴𝑓𝑖𝑙𝑚\nwhere and represent the total projected area of the planar defects and the\ncross-s𝐴ec 𝑑 t 𝑒 i 𝑓 o 𝑒 n 𝑐 a 𝑡 l area𝐴 𝑓 o 𝑖𝑙 f 𝑚 the thin film, respectively. denotes the cumulative\nlength of the bright streaks, and corresponds 𝐿to 𝑠𝑡 𝑟 th 𝑒𝑎 e 𝑘 projected area of a single\nFe-BNBT unit cell. The volume f𝐴ra 𝑢 c 𝑛 t 𝑖 i 𝑡 o 𝑐 n 𝑒 𝑙 o 𝑙 f the planar defects in the Fe-BNBT thin films\nis estimated to be ~16%.”\nComments (2): Secondly, after polarization switching, H-H domain walls should\nbecome T-T domain walls and acquire different properties because of the change of\ntheir charge sign. In turn, the oxygen vacancies should follow the change of the wall\ntype quickly enough for the proposed mechanism to work. These aspects are not\nmentioned in the manuscript. Therefore, the authors’ picture that the observed effects\nare entirely a result of the flexoelectric effect is not substantiated.\nResponse: While H-H and T-T domain walls are conceptually linked, their behavior\nunder an applied electric field is more complex than a simple one-to-one transformation.\nRather than directly converting from H-H to T-T upon polarization switching, domain\nwalls predominantly migrate in the direction of the applied electric field,43, 44 as\ndemonstrated by the results of phase-field simulations in the manuscript (Figure R15;\nor Supplementary Figure 4 in the original manuscript).\nresults depicting the polarization distribution of bound charge-modified Fe-BNBT in various\nstates: (a) initial configuration, (b) under the application of an upward-directed electric field\n(E ), (c) following the removal of the upward electric field, (d) under the application of a\nup\ndownward-directed electric field (E ), and (e) following the removal of the downward\ndown\nelectric field.\nNotably, the emergence of H-H and T-T domain walls is intrinsically linked to\nthe alignment of oxygen vacancies along specific crystallographic directions. 17, 18, 32\nWhile the domain walls can reorient in response to an external electric field, the ordered\ndistribution of oxygen vacancies, which manifest as planar defects, remains largely\nunaltered. This stability arises from the elevated activation energy barriers for the\nmigration of oxygen vacancies, particularly in the immediate proximity of dislocations,\ncompared to those in the bulk lattice.45 As a result, the local strain gradients and the\nassociated flexoelectric effects are preserved under the applied electric fields. This\nscenario is also evidenced by our fatigue measurement presented in the original\nmanuscript (Figure 6), which reveal a sustained high P of 139 μC cm-2 after 107 cycles,\nm\n35\nstrongly indicating a persistent flexoelectric effect in the Fe-BNBT thin films under the\napplied electric field.\nWe added the corresponding discussion to the Supplementary Information in the\nrevised manuscript. We hope the above clarifications have effectively addressed the\nreviewer’s concern.\n“The emergence of H-H and T-T domain walls is intrinsically linked to the\nalignment of oxygen vacancies along specific crystallographic directions. While the\ndomain walls can reorient in response to an external electric field, the ordered\ndistribution of oxygen vacancies, which manifest as planar defects, remains largely\nunaltered. This stability arises from the elevated activation energy barriers for the\nmigration of oxygen vacancy, particularly in the immediate proximity of dislocations,\ncompared to those in the bulk lattice. As a result, the local strain gradients and the\nassociated flexoelectric effects are preserved under the applied electric fields.”\nComments (3): In fact, BNBT exhibits relaxor properties, and the used composition is\nclose to a morphotropic phase boundary (MPB), where small variations in doping can\nsignificantly impact polarization (see, e.g., https://doi.org/10.1039/C7CP03530F).\nWhile the flexoelectric effect may contribute to the observed properties, it is likely\nonly one of several contributing factors, and its contribution may be negligible.\nResponse: We acknowledge the reviewer’s comment regarding the composition of our\nFe-BNBT being near the morphotropic phase boundary (MPB), which typically leads\nto the coexistence of rhombohedral (R) and tetragonal (T) phases in ceramics. However,\nthe lattice structures of epitaxial thin films may differ from their bulk counterparts. In\nepitaxial thin films, the lattice structure can be strongly influenced by the biaxial strain\nimposed by the single-crystal substrate46 or secondary phase7, 8.\na predominantly tetragonal structure through RSM (Figure 1(d) in the original\nmanuscript), XAS (Figure 2 in the original manuscript) analyses, and most notably, the\npolarization measurements of differently oriented Fe-BNBT thin films (Figure 1(a) in\nthe original manuscript).\nIt is known that the ferroelectric properties of both tetragonal and rhombohedral\nstructures are strongly dependent on the crystallographic orientation.47-49 In the\nrhombohedral phase, the polarization along different crystallographic directions\ntypically follows the ratio P :P :P = 1: : .49 Whilst in a tetragonal\n(001) (110) (111)\nstructure, the polarization magnitudes obey a relation√ 2of P√3 :P :P = : :1, as\n(001) (110) (111)\nillustrated in Figure R16(b).48 For our Fe-BNBT thin films, the measured p√o3la√ri2zation\nvalues along the (001), (110), and (111) orientations are 161, 118, and 95 μC cm-2,\nrespectively, as shown in Figure R16(a) (Figure 1(a) of the original manuscript). These\nvalues closely match the relation of : :1, indicating a tetragonal structure.\n√3 √2\nFigure Redacted\nFigure R16. (a) The room-temperature P–E loops of [001], [110], and [111]-oriented Fe-BNBT\nthin films under selected electric fields. (b) Schematic illustration of domain configurations in\ntetragonal Fe-BNBT thin films subjected to electric fields applied along the (001), (110), and\n(111) crystallographic directions, respectively.48\nAdditionally, we assessed the relaxor properties of our thin films by calculating\nthe relaxor diffuseness factor, γ, using the modified Curie-Weiss law50:\n𝛾\n1 1 (𝑇−𝑇𝑚)\n− =\n𝜀𝑟 𝜀𝑟,𝑚 𝐶\n37\nwhere is dielectric permittivity at the temperature T, the maximum dielectric\npermitt𝜀i 𝑟 vity at the temperature , and C the Curie const𝜀a 𝑟 n , t 𝑚 . In typical ferroelectrics,\nγ is expected to be 1, while γ = 2𝑇 𝑚 suggests an ideal relaxor ferroelectric. Our calculation\nyielded a γ value of 1.013, as shown in Figure R17, which unambiguously suggests a\ntypical ferroelectric behavior of our Fe-BNBT thin films, rather than a relaxor at a\ncomposition close to MPB as suggested by the reviewer.\nFigure R17 (Supplementary Figure 6 in the revised manuscript). Fitting of the relaxor\ndiffuseness factor of Fe-BNBT films based on the temperature-dependent dielectric data.\nIn comparison, the work by Zhou et al. (https://doi.org/10.1039/C7CP03530F)\ndiscusses polycrystalline BNBT thin films (BT content ranging from 0–16%) prepared\nby a sol-gel method. It may not be appropriate to compare the polarization behavior of\na PLD-produced epitaxial thin film to that of polycrystalline relaxor thin films (near\nMPB compositions as suggested in the work by Zhou et al) anyhow, even the nominal\ncompositions are almost identical. Furthermore, Zhou et al.\n38\n(https://doi.org/10.1039/C7CP03530F) did not examine the role of dopants in\nmodifying the ferroelectric behavior of BNBT thin films.\nWe also studied the effect of dopant species on the polarization of resultant\nepitaxial BNBT-based thin films. To that end, 550 nm-thick BNBT-based thin films\ndoped with 2.0 mol% of Mn, La, Sm, Eu, and Co, along with undoped BNBT thin films\nwere deposited on STO (001) single crystal substrate by PLD under optimal processing\nconditions. As shown in Figure R18, the maximum polarization of the undoped, Mn-,\nLa-, Sm-, Eu, and Co-doped BNBT films are 62, 73, 68, 73, 74, and 97 µC cm-2,\nrespectively. In stark contrast, Fe-doped BNBT films exhibit a considerably higher\nmaximum polarization of ~161 µC cm-2, which cannot be solely ascribed to the MPB\ncomposition or doping effect. In the revised manuscript, we’ve provided both\nexperimental and simulation evidence to support our findings that oxygen defects\ninduced flexoelectric effect play a critical role in boosting the polarization of Fe-BNBT\nthin films.\nTable R2 summarizes the reported ferroelectric properties of typical doped\nBNBT-based thin films with various dopants produced using different methods. To date,\nwe’re not aware of any dopants other than Fe that can enhance the polarization of BNBT\nthin films having an MPB composition, regardless of polycrystalline or epitaxial, to\nabove 100 µC cm-2. This fact suggests that another contributing factor, namely the\nflexoelectric effect as explained in detail in this work, is NOT negligible. It is evident\nthat the contribution of the flexoelectric effect to the giant polarization in our Fe-BNBT\nthin films is nontrivial.\n39\nloops of (a) undoped BNBT thin film and BNBT-based thin films doped with 2.0 mol% of (b)\nMn, (c) La, (d) Sm, (e) Eu, (f) Co, and (g) Fe, respectively. (h) Comparison of P–E loops of\nBNBT thin films with different dopants.\nTable R2. Comparative analysis of ferroelectric properties in BNBT thin films doped with\nvarious elements and synthesized via different fabrication techniques.\nSamples Pr Pm Ec Em Substrate Method Thickness Ref.\nµC cm-2 µC cm-2 kV cm-1 kV cm-1 nm\nUndoped 19 62.1 135 750 (001) SrTiO PLD 550 This\nUndoped 18 73 100 750 Pt/TiO/SiO/Si Chemical 200 51\n2 2\nBNBT (MPB) solution\ndeposition\nUndoped 16 48 115 1780 (111) Chemical 680 52\nBNBT (MPB) Pt/Ti/SiO/Si solution\n2\ndeposition\nUndoped 12 23 180 800 Pt/TiO/SiO/(1 Spin-coating 125 53\n2 2\nBNBT (MPB) 00)Si\n0.96(Bi ,Na ) 12 38 100 650 Pt(111)/Ti/SiO Sol-gel 200 Noted by\n0.5 0.5 2\nTiO- /Si(100) reviewer\n3\n0.04BaTiO\nLa- and Ce- 29.5 70 74 325 (001) PLD 400 48\ncodoped BNBT (LaAlO)TiO(\n3 3\nSrAlTaO)\n2 60.35\nRb-doped 28.9 50 200 600 Pt(111)/Ti/SiO Metal-organic - 54\ndeposition\nMn-doped 25 48 150 650 (001) LaAlO PLD 162 55\nLi-doped 23.9 60 124 400 Pt(111)/Ti/SiO Metal-organic - 56\ndeposition\nFe-doped 115 161 250 750 (001) SrTiO PLD 550 This\nCo-doped 51 97 284 750 (001) SrTiO PLD 550 This\nEu-doped 41 74 165 750 (001) SrTiO PLD 550 This\nSm-doped 26 73 100 750 (001) SrTiO PLD 550 This\nMn-doped 19 73 120 750 (001) SrTiO PLD 550 This\nLa-doped 17 68 85 750 (001) SrTiO PLD 550 This\nWe have integrated these experimental results and the corresponding discussion\ninto the Supporting Information. We hope this clarification effectively addresses the\nconcern.\n“We assessed the relaxor properties of our thin films by calculating the relaxor\ndiffuseness factor, γ, using the modified Curie-Weiss law6:\n𝛾\n1 1 (𝑇−𝑇𝑚)\n− =\n𝜀𝑟 𝜀𝑟,𝑚 𝐶\nwhere is dielectric permittivity at the temperature T, the maximum dielectric\npermitt𝜀i 𝑟 vity at the temperature , and C the Curie const𝜀a 𝑟 n , t 𝑚 . In typical ferroelectrics,\nγ is expected to be 1, while γ = 2𝑇 𝑚 suggests an ideal relaxor ferroelectric. Our calculation\nyielded a γ value of 1.013, as shown in Supplementary Figure 6, which unambiguously\nsuggests a typical ferroelectric behaviour of our Fe-BNBT thin films, rather than a\nrelaxor at a composition close to MPB.”\n“To confirm the unique role of Fe doping in enhancing the polarization, 550 nm-thick\nBNBT thin films doped with 2.0 mol% of a few other commonly adopted dopants in\nperovskite ferroelectrics, including Mn, La, Sm, Eu, and Co, as well as undoped BNBT\nthin films were grown on LSMO-electroded STO (001) single crystal substrates using\nPLD under the same deposition conditions as that for Fe-doped counterparts. Room-\ntemperature P-E loops of these BNBT-based thin films are presented in Supplementary\nFigure 15.\nThe maximum polarization P of undoped, Mn-, La-, Sm-, Eu-, and Co-doped\nm\nBNBT thin films are 62, 73, 68, 73, 74, and 97 µC cm-2, respectively. In contrast, Fe-\ndoped BNBT thin films exhibit a drastically higher P of 161 µC cm-2.”\nm\nComments (4): Additionally, HRTEM measurements were conducted under “zero-field”\nconditions, while an applied electric field can induce a rhombohedral-to-tetragonal\nphase transition at the MPB. This raises questions about whether the observed effects\nwould remain under applied fields. Furthermore, flexoelectric coupling at defects such\n42\nas grain boundaries and domain walls is currently well known. However, it remains\ngenerally unclear whether enhanced polarization leads to larger strain or vice versa.\nResponse: As discussed in our manuscript (and in the response to Comment (3)), we\nhave confirmed that our Fe-BNBT thin films predominantly exhibit a tetragonal\nstructure under zero-field conditions. Therefore, a field-induced rhombohedral-to-\ntetragonal phase transition in the MPB compositions of bulks57, 58 is not anticipated in\nour Fe-BNBT thin films. In addition, T-T and H-H domain walls will reorient to the\ndirection of the applied electric field, rather than undergoing a simple one-to-one\ndomain transformation, as detailed in our original manuscript (and in the response to\nComment (2)). The emergence of H-H and T-T domain walls is intrinsically linked to\nthe alignment of oxygen vacancies along specific crystallographic directions.17, 18, 32\nWhile the domain walls can reorient in response to an external electric field, the ordered\ndistribution of oxygen vacancies, which manifest as planar defects, remains largely\nunaltered. This stability arises from the elevated activation energy barriers for the\nmigration of oxygen vacancies, particularly in the immediate proximity of dislocations,\ncompared to those in the bulk lattice.45 As a result, the local strain gradients and the\nassociated flexoelectric effects are preserved under the applied electric fields. The P-E\nand remanent P-E loops results (Figure 1(a), 1(b), 6(a), and 6(b) in the original\nmanuscript clearly demonstrate that the local strain gradients and the associated\nflexoelectric effects remain under applied electric fields.\nRegarding the correlation between polarization and strain, we infer that the\nreviewer refers to the electrostrain instead of the epitaxial strain. In general, field-induced\nphase transitions in ferroelectrics with MPB composition or relaxor ferroelectrics may\nlead to both enhanced polarization and electrostrain.59-61 Nonetheless, it is clear that\nan enhanced polarization does not necessarily result in a large electrostrain, e.g., in a\nferroelectric material where polarization switching is dominated by 180o domains which\ndo not contribute to electrostrain,62 nor does a large electrostrain always originate from\na high polarization, as detailed in Table R3.\nWhile polarization and electrostrain are often correlated as domain switching is\ninvolved in both processes, distinct physical mechanisms can contribute to the\nmeasured polarization and electrostrain responses. Polarization typically arises from\ndomain reorientation and lattice distortion7. In contrast, the mechanisms underlying\nelectrostrain can be more complex. For example, when mechanisms such as\nelectrostriction63, hopping of oxygen vacancies near ceramic surfaces64, or enhanced\nlocal structural heterogeneity65 dominate, a pronounced electrostrain may arise even\nwhen the maximum polarization remains moderate.\nWe thank the reviewer for raising this point. We will investigate the\npiezoresponse of our Fe-BNBT thin films using a similar technique reported by Liu et\nal (Nat. Commun., 15, 10798 (2024)) in the near future.\nTable R3. Summary of the electrostrain and polarization values of some reported perovskite\noxides.\nMaterials Electrostrain P Mechanism behind the Ref.\nm\n(%) (µC cm-2) enhanced electrostrain\n(K0.5,Na0.5)NbO3 ceramics 1.9 10 Hopping of oxygen vacancies 64\nnear ceramic surfaces\n0.92(Bi0.5,Na0.5)TiO3-0.07BaTiO3- 0.36 45 Transition between 66\n0.01Sr(Sn0.5Ta0.4)O3 ceramics ferroelectric and ergodic\nrelaxor states\nSm-doped Pb(Mg1/3Nb2/3)O3- 0.3 30 Enhanced local structural 65\nPbTiO3 single crystals heterogeneity\n0.9(Bi0.5,Na0.5)TiO3-0.1SrTiO3 0.1 61.3 Transition from the nonergodic 59\nsingle crystal to the ergodic relaxor state\nComments (5): \"Super-tetragonal\" (Line 73) – This term should be introduced (defined\nand explained) in Introduction.\nincluded the definition of the term “super-tetragonal” in the revised manuscript for\nbetter clarity as follows:\n“A super-tetragonal phase is used to describe the unique lattice structure in perovskite\noxide thin films that exhibits an extremely large c/a ratio, e.g., 1.10-1.25,6-8 often\nresulting in an enhanced polarization.”\nComments (6): \"Facilitate the creation of oxygen vacancies, which accumulate in tail-\nto-tail (T-T) domain configurations\" (Line 95) -> Should it be “T-T domain walls”?\nResponse: The corrected sentence now reads in the revised manuscript:\n“Facilitate the creation of oxygen vacancies, which accumulate in tail-to-tail (T-T)\ndomain walls.”\nComments (7): “Between the adjacent T-T domain walls, H-H domains emerge\"(Line\n97). Please reformulate more carefully. Should it be “H-H domain walls”? The word\n“adjacent” is not correctly used in this sentence. It may be omitted.\nResponse: We have rephrased the sentence for better clarity. The revised text now\nreads:\n“Between the T-T domain walls, head-to-head (H-H) domain walls emerge.”\nComments (8): The charge compensation of the T-T walls is due to accumulation of\noxygen vacancies at the walls. How the charge compensation of the H-H walls is\nachieved? It is worth illustrating the oxygen vacancy formation with equations of defect\nchemistry.\nResponse: The charge compensation of H-H domain walls is primarily achieved by the\ntrapping of free electrons, which is derived from the excess carriers generated during\nthe formation of oxygen vacancies.67 Specifically, electropositive oxygen vacancies\npreferentially accumulate at T-T domain walls, while free electrons, being negatively\ncharged, become trapped and localized at H-H domain walls.\nThe formation of oxygen vacancies in Fe-doped BNBT can be described using\nthe defect chemistry notation:\n2𝑇𝑖𝑂2\n′ ·· ×\n𝐹𝑒2𝑂3 ⎯⎯ 2𝐹𝑒𝑇𝑖 +𝑉𝑂 +3𝑂𝑂\nWe have added the following explanations to the revised manuscript:\n“The formation of oxygen vacancies in Fe-doped BNBT can be described using the\ndefect chemistry notation: .”\n2𝑇𝑖𝑂2\n′ ·· ×\n𝐹𝑒2𝑂3 ⎯⎯ 2𝐹𝑒𝑇𝑖 +𝑉𝑂 +3𝑂𝑂\nComments (9): The value of the ratio Aeg/t2g of 3.51 looks very large compared to the\ntypical values of 1–1.5. This may suggest that the background has not been accurately\nsubtracted. Please comment.\nResponse: We have carefully re-examined our data and confirmed that the background\nsubtraction was properly carried out during the XAS data processing. The reported\nA value of 3.51 corresponds to the ratio of integrated areas under the e and t\neg/t2g g 2g\npeaks. We also note that significantly lower values of e /t ratios, such as 1.04, 1.55,\ng 2g\nand 1.03, have been reported in other studies;68 however, these typically refer to the\nratio of peak intensity rather than the ratio of integrated area. To avoid further confusion,\nwe also adopt the ratio of peak intensity, which is ~1.1, in the revised manuscript. We\nhope this clarification satisfactorily addresses the reviewer’s concern.\nComments (10): The statement “connected with the localized character of the eg states,\ndemonstrating the existence of defective structures in Fe-doped BNBT thin films” (Line\n159) – needs more explanation.\nspectra, which provide insights into the electronic structure and coordination\nenvironment of Ti ions in our Fe-doped BNBT thin films. The L e region is\n3 g\ndeconvoluted into the low-photon-energy peak (yellow peak) and high-photon-energy\npeak (pink peak), as shown in Figure 2 (c) in the revised manuscript. The deconvolution\nwas performed due to the sensitive asymmetry of the e peak, which is not expected\ng\nwithin an ideal octahedral symmetry.\nFigure R19 (Figure 2(c) in the revised manuscript). The (c) Ti L-edge XAS spectra of Fe-\nBNBT thin films acquired in total fluorescence yield (TFY) mode. The inset in (c) is the result\nof peak fitting of the L e peak of Fe-BNBT thin films.\n3 g\nThe enhanced intensity of the low-energy peak suggests a local splitting of the e\ng\nstates, which arises from the distortions of the TiO octahedra.69 These distortions are\n6\nmost likely induced by oxygen vacancies and associated planar defects. Such structural\nand chemical disorder reduces the local symmetry and modifies the crystal field\n47\nenvironment, leading to orbital energy shifts and increased localization of the e\ng\nelectrons.69 This interpretation is consistent with previously reported XAS studies on\ndefective or distorted perovskite systems, where asymmetry in the e region originates\ng\nfrom the defect-induced orbital localization.69\nWe have added the following explanations to the revised manuscript:\n“These distortions are most likely induced by oxygen vacancies and associated planar\ndefects, stemming from the replacement of Ti by Fe cations. Such structural and\nchemical disorder reduces the local symmetry and modifies the crystal field\nenvironment, leading to orbital energy shifts and increased localization of the e\ng\nelectrons.16”\nComments (11): What is meant by \"bright twill\"? “Twill” is not a technical term but\nrather describes perception of the pattern in the image by the observer using some\nanalogy. Therefore, it needs to be descriptively introduced before the use in the image\nanalysis. The text is confusing as now and does not assist readers in following the\nauthors’ logic. For example, are the twills 1D or 2D objects? To facilitate the readers’\nwork, it looks beneficial to move Figure S1 into the main text and use it to explain the\nmeaning of the word “twills”.\nResponse: We thank the reviewer for raising the issue regarding the use of the term\n\"bright twill.\" To address this, we describe these features as \"two-dimensional and\noriented segments along [101] and [-101] (referred to as “bright streaks”)\" in the revised\nmanuscript. We hope this revision will help better communicate the nature of the\nobserved structures.\nInformation has been integrated into Figure 3 of the revised manuscript to visually\nsupport the description of these features.\nFigure R20 (Figure 3 in the revised Manuscript). (a) HAADF STEM images of Fe-BNBT\nthin films, indicating the high density of bright streaks. (b) The polarization map of Fe-BNBT\n49\nthin films, where the length and direction of the arrows represent the strength and orientation\nof the polarization, respectively. (c) A magnified view illustrating that the bright streaks act as\nT-T domain walls in the Fe-BNBT thin films, extracted from the red-dashed square region in\n(b). The yellow and green dashed rectangles highlight the bright streaks region and its nearby\nareas, respectively. (d) STEM-iDPC image of Fe-BNBT thin films, corresponding to the red-\ndashed square region in (c). The yellow and green dashed rectangles highlight the bright streaks\nregion and its nearby areas, respectively. (e) Oxygen atom signal intensities from the two\nregions marked by the yellow and green dashed rectangles in (c) and (d), demonstrating the\ncontrast difference between the bright streaks and its surrounding areas.\nComments (12): The “twill” areas apparently possess polarization. Would it be more\nappropriate to understand the bright stripes as narrow, platelet-like, ferroelastic\ndomains oriented along the normal to the image plane? The statement \"bright twills\nfunction as T-T domain walls in Fe-BNBT thin films\" (Line 192) is unclear. There\nappears to be a non-trivial polarization distribution within and outside the “bright twill”\nregion. Are the boundaries between this region and neighboring domains uncharged?\nResponse: Ferroelastic domains refer to non-180° domains.70 As shown in Figure R21,\nthe majority of the domains along the bright streaks in our Fe-BNBT thin films are 180°\ndomains. Consequently, we believe that the bright streaks in our Fe-BNBT thin films\nare best described as ferroelectric domain walls rather than ferroelastic domains.\nTo avoid further confusion, the statement “bright twills function as T-T domain\nwalls in Fe-BNBT thin films” has been revised to “T-T domain walls appear near the\nbright streaks in our Fe-BNBT thin films”.\nWe also confirm that the boundary between the bright region and neighboring\ndomains is charged. This is evidenced by the oxygen atom signal intensities at the\nboundary, where the presence of oxygen vacancies corroborates that the boundary is\n50\ncharged (Figure R13 in the response to Comment (1) and Figure R22 in the response to\nComment (14)).\nFigure R21. The polarization map of Fe-BNBT thin films, where the length and direction of\nthe arrows represent the strength and orientation of the polarization, respectively.\nComments (13): The claim that \"T-T configurations are more energetically favorable\nfor generating charged domain walls\" (Line 195) is puzzling, and generally, incorrect.\nT-T domain walls are charged, and charged domain walls are energetically\nunfavorable, requiring stabilization, which is kinetic in nature. For example, T-T\ndomain walls carry negative polarization charge, which can be compensated by\npositively charged oxygen vacancies.\nResponse: We sincerely thank the reviewer for this insightful observation. It is\nimportant to emphasize that T-T domain walls can form through the line-up of oxygen\n51\nvacancies along specific crystallographic direction.32 This vacancy ordering leads to the\nemergence of charged domain walls, whose formation becomes energetically\nfavourable in regions with elevated oxygen vacancy concentrations.\nWe have revised the manuscript to better reflect this understanding. The revised\nstatement now reads:\n“Density functional theory (DFT) calculations indicate that the aggregation of oxygen\nvacancies can give rise to charged T-T domain walls.”\nComments (14): In Figure 3b, where is the yellow region located? Is it at the boundary\nor inside the “twill”? Could you provide a clearer depiction in the figure? Additionally,\ncould you compare multiple line profiles from different locations (middle, boundary)\ninside the bright twill and the rest of the film, similar to panels (c) and (d)? The phrase\n\"intensity of oxygen atom signals\" is unclear. Please explain how you distinguish\noxygen from other elements.\nResponse: In Figure 3(b) of the original manuscript, the region highlighted by the\nyellow dashed box corresponds to the boundary between the bright streak and the\nadjacent unit cell. To improve the clarity, we have updated the figure caption as follows:\n“The yellow dashed box indicates the interface between the bright streak and its\nneighbouring unit cell, whereas the green dashed box denotes a region approximately\nseven unit cells away from the defect zone.”\nOwn to the two-dimensional nature of the bright streaks, manifesting as linear\nfeatures in the STEM projection, oxygen atomic columns are absent from the middle\nportion of these streaks in the imaging plane. To further elucidate the spatial variation\nand persistence of oxygen vacancies, the line profiles of oxygen atom signal intensities\n52\nwere extracted from eight sequentially adjacent regions near the streaks, delineated by\nblue to green dashed rectangles, as shown in Figure R13. These profiles reveal a clear\ngradient in oxygen vacancy concentration, indicating that the associated planar defects\nextend beyond these streaks into the adjacent lattice, spanning over approximately 5\nunit cells. In addition, we analyzed line profiles from two other bright streak regions,\nas presented in Figure R22. These comparative results are in strong agreement with\nthose shown in the manuscript and further substantiate our interpretation regarding the\npresence and distribution of oxygen deficiencies near the planar defects.\nFigure R13 (Supplementary Figure 9 in the revised manuscript). (a) HAADF-STEM image\nof Fe-BNBT thin films. (b) STEM-iDPC image of the region outlined by the red dashed square\nin (a). The yellow dashed line indicates the position of the bright streak, while the series of\ncolored dashed rectangles represent eight sequentially adjacent regions extending outward from\nthe bright streak. (c) Line profiles of oxygen atom signal intensities extracted from the regions\nmarked by the blue to green dashed rectangles in (b), illustrating a progressive decrease in\n53\noxygen vacancy concentration from the core of the bright streak toward the surrounding lattice.\nThe blue arrows mark the positions of oxygen vacancies.\nFigure R22. (a) and (b) HAADF STEM images of Fe-BNBT thin films. (c) and (d) STEM-\niDPC images aligned with the HAADF images in (a) and (b), respectively. (e) and (f) Line\nprofiles of oxygen atom signal intensities extracted from the three regions indicated by the\nblack, light green, and green dashed rectangles in (c) and (d). The blue arrows mark the\npositions associated with oxygen vacancies.\nOxygen atoms were identified by comparing HAADF images with iDPC STEM\nimages. In HAADF mode, atomic contrast is primarily governed by atomic number (Z-\ncontrast), allowing clear visualization of A-site and B-site cations but rendering oxygen\natoms almost invisible due to their low atomic number. iDPC imaging enables the\nsimultaneous detection of both heavy and light elements. In this mode, the oxygen\natomic columns appear as additional low-contrast features along the [001] projection\n54\naxis that are absent in the HAADF images. These comparative analyses allow for the\nreliable identification of oxygen positions and enable a reasonable assessment of\noxygen vacancy distribution.\nComments (15): The statement \"H-H domain walls form between adjacent bright twills,\nacting as transition zones for T-T domain structures\" (Line 140) is unclear. A wall is a\nnarrow transition zone between domains, but what is meant here is a relatively\nextended zone. What are “T-T domain structures”? In Figure 4, there is a region with\nhead-to-head polarization. Inside this region, there is also in-plane polarization, which\ndoes not appear to be an artifact. Can this region be a domain?\nResponse: We sincerely appreciate the reviewer’s insightful comments.\nTo avoid further confusion, we have revised the statement: “H-H domain walls\nwere found to form between adjacent pairs of bright twills, acting as transition zones\nfor T-T domain structures.” to “H-H domains were found to form between adjacent\npairs of bright streaks.”\nRegarding the “in-plane polarization” observed in the region shown in Figure 4\nof the original manuscript, it may be resulted from the presence of a strongly charged\ndomain wall, which induces local polarization rotation.71 It is more appropriate to\ndescribe this as a domain wall rather than a discrete ferroelectric domain, as it spans\nover one unit cell only.\nComments (16): The difference between regions like 1 and 2 in Figure 4 is due to\ndifferent polarization orientations. The formation of polarization clusters with random\norientations is common in relaxors.\nour Fe-BNBT thin film is a typical ferroelectric, instead of a relaxor, and predominantly\nadopts a tetragonal structure.\nThe distinction between regions 1 and 2 in Figure 4 primarily arises from\nsignificant differences in the magnitude of B-site cation displacement, rather than from\nvariations in polarization orientation. In STEM-based polarization mapping, the length\nand direction of the arrows represent the strength and orientation of the polarization,\nrespectively. In region 1, the markedly larger displacements produce well-defined\nvectors, whereas the smaller displacements in region 2 yield vectors that appear only as\n“arrowheads”. This visual discrepancy may cause an impression of relaxor-like\nbehavior.\nTo better reflect the displacement behavior of B-site cations in the areas distant\nfrom planar defects, we performed polarization mapping at multiple planar-defect-free\npositions across the Fe-BNBT thin films, as presented in Figure R23. It is evident that\nthe Fe-BNBT thin films exhibit typical tetragonal domains instead of polarization\nclusters with random orientations.\nComments (17): The phase-field modeling suggests domains extending through the\nentire film thickness, but this is not observed experimentally. Would this discrepancy\nsignificantly affect the estimated polarization values?\nResponse: We would like to clarify that the polarization values obtained from the\nphase-field simulations are employed to qualitatively illustrate the mechanism of\npolarization enhancement induced by local strain gradients, rather than quantifying the\nexact polarization that mirrors experimental results. The primary purpose of these\nsimulations is to conceptually establish the relationship between oxygen vacancy\n57\naccumulation and the formation of T-T domain walls, as well as the link between T-T\nand H-H domain walls and polarization enhancement through the flexoelectric effect.\nTo further address the reviewer’s concern regarding the impact of domain density\nand scale on the resultant polarization, additional phase-field simulations were\nconducted with both reduced density and reduced spatial scale of T-T and H-H domain\nwalls. As shown in Figure R24, even with this lower concentration and smaller length\nof H-H and T-T domain walls, where only part of the regions contain T-T and H-H\ndomain walls while the remainder exhibits typical H-T domain walls, an enhancement\nin polarization is still observed. This result strengthens our conclusion that the\npolarization improvement arises primarily from the localized strain gradients associated\nwith the planar defects without the need for them to be extended throughout the entire\nthickness of the thin films.\nand for BNBT modified by the flexoelectric effect with (b) a high and (c) a low concentration of T-\nT and H-H domain walls. The simulated P-E loops for (d) the pristine BNBT and (e)\nflexoelectric effect-modified BNBT.\nComments (18): Grammar:\n--\"indicates enhanced tetragonally distortion\"\n--\"the area intensity of the eg and t2g peaks\"\nResponse: We thank the reviewer for flagging up the grammatical errors. We have\ncorrected these issues in the revised manuscript:\nThe phrase “indicates enhanced tetragonally distortion” has been revised to\n“indicates enhanced tetragonal distortion”.\nThe phrase “the area intensity of the e and t peaks” has been revised to “the\ng 2g\npeak intensity of the e and t peaks” for better clarity.\ng 2g\nComments (19): In summary, this manuscript, while reporting interesting and\npotentially important experimental outcomes, provides no evidence that the achieved\nhigh values of polarization can be attributed solely to the flexoelectric effect. Given that\nthis assumption is central to the paper and that the whole concept may need a revision,\nthe manuscript cannot be recommended for publication in Nature Communications.\nThe authors may consider publication in a more specialized journal.\nResponse: We sincerely thank the reviewer for acknowledging the significance of our\nexperimental findings and their thoughtful and constructive comments.\nIn addition to addressing the specific comments made by the reviewer above, we\nhave supplemented additional experimental results to clarify the impacts of dopant\n59\nspecies, doping level, substrates, and film thickness on polarization of the resultant\nBNBT-based thin films. Moreover, our modelling based on the phase field simulations\nconfirms qualitatively that the unique domain configurations, i.e., T-T domain walls,\nresulting from the oxygen-deficient planar defects, play a critical role in improving the\npolarization of our Fe-BNBT thin films. In short, we’ve refined our analysis and\nstrengthened our claims through providing compelling evidence, both experimentally\nand theoretically, that the flexoelectric effect induced by the local defects in the Fe-\ndoped BNBT thin films is the key origin of the largely enhanced polarization.\nResponse to the comments of Referee #3\nThis study significantly enhanced the polarization of (Bi,Na)TiO₃-based\nferroelectric thin films by introducing a defect structure with oxygen vacancies through\nFe doping and leveraging the flexoelectric effect. The formation of H-H domain\nstructures and a high strain gradient was confirmed via iDPC-STEM analysis, resulting\nin a maximum polarization of 161 μC cm⁻² and a remanent polarization of 115 μC cm⁻².\nFurthermore, the films exhibited excellent thermal stability, maintaining their\nperformance even at 230°C, highlighting their potential for high-temperature\nelectronic device applications. However, several revisions are required for publication.\nResponse: We greatly appreciate the reviewer’s positive feedback on our work.\nComments (1): The core finding of this study is that the polarization of (Bi,Na)TiO₃-\nbased ferroelectric thin films was maximized by leveraging the flexoelectric effect\ninduced by oxygen vacancies. While the flexoelectric effect is identified as the primary\ncause of the polarization enhancement, the flexoelectric coefficients were estimated\nusing values from existing literature rather than direct experimental measurement. To\nimprove accuracy, the flexoelectric coefficient should be measured using AFM\n60\ntechniques, or the induced polarization should be further analyzed through XRD and\nRaman spectroscopy.\nResponse: We appreciate the significance of quantifying the flexoelectric effect\nthrough direct measurements. Atomic force microscopy (AFM) has been widely\nemployed to evaluate the strain gradient.72 However, we are not aware of any previous\nstudies that determined the flexoelectric coefficient of thin films grown on rigid single-\ncrystal substrates using AFM-based techniques.\nTo the best of our knowledge, high-resolution transmission electron microscopy\n(HRTEM) is the only practical technique that is appropriate for probing the local strain\ngradient at atomic scale in our thin films. This is because the spatial resolution of other\ncommonly employed probes, including laboratory-based XRD (>10 µm)73, synchrotron\nXRD (>0.1 µm)74, and neutron diffraction (>250 µm)75, is orders of magnitude larger\nthan the “feature size” of the localized strain gradient in our Fe-BNBT thin films. XRD\nis a reliable technique for determining the overall lattice symmetries, such as identifying\nthe super-tetragonal phase in ferroelectric films.9 However, in the case of our Fe-BNBT\nthin films, the observed enhancement in polarization primarily arises from the\nflexoelectric effect, which is associated with a pronounced local strain gradient rather\nthan significant changes in the average lattice parameters. Based on the XRD results,\nwe can only conclude that the films exhibit a slight tetragonal distortion (c/a = 1.004),\nas discussed in Figure 1 in the original manuscript. This trivial structural distortion is\ninsufficient to account for the substantially enhanced polarization. In short, XRD\nmeasurements alone cannot uncover the origin of the giant polarization.\nAlthough Raman spectroscopy can, in principle, provide indirect evidence of\nstrain gradients, since strain alters phonon frequencies and can induce local symmetry\n61\nbreaking through mechanical bending or lattice mismatch between thin films and rigid\nsingle-crystal substrates, we are not aware of any previous studies in which Raman has\nbeen employed to quantify ferroelectric or flexoelectric polarization of such films. Akin\nto XRD, Raman spectroscopy primarily captures structural information related to lattice\ndistortions and phase composition. While it is effective for chemical identification and\nphase analysis, it does not directly probe flexoelectric polarization induced by local\nstrain gradients.\nAs reported in prior studies,76, 77 investigation of the flexoelectric effect at the\nnanoscale has predominantly relied on the combination of HRTEM and theoretical\nmodelling, both of which are also employed in our study.\nWe hope these clarifications adequately address the reviewer’s concern.\nComments (2): This study demonstrates that the flexoelectric effect induces a\nsignificant increase in polarization; however, a quantitative separation of ferroelectric\nand flexoelectric polarization was not performed. To clarify the contribution of the\nflexoelectric effect, a more detailed P-E loop analysis should be conducted, including\na comparison with conditions where the flexoelectric effect is not considered.\nResponse: In our original manuscript, we conducted a quantitative analysis to\ndistinguish between ferroelectric and flexoelectric contributions to the overall\npolarization.\nThe P-E loops for undoped BNBT thin films are presented in Figure R8(b)\n(Figure 6 in the Supplementary Information of the original manuscript), revealing a\nmaximum polarization P of 62 μC cm-2. In addition, the STEM observations confirm\nm\nno sensible bright streaks in the undoped BNBT thin films, as shown in Figure R8(d).\nIn contrast, the Fe-doped BNBT thin films exhibit a significantly enhanced P of 161\nm\n62\nμC cm-2, which largely arises from the flexoelectric effect. The contribution of\nflexoelectricity can be roughly separated based on the relation ,78 where\nP represents the total polarization, P refers to the polarizat𝑃io 𝑡𝑜 n 𝑡 o=rig𝑃i 𝐷 na+tin𝑃g 𝐹 from the\ntot D\ndielectric response, and P corresponds to the flexoelectric polarization. Consequently,\nF\nthe ferroelectric polarization of the Fe-BNBT thin films is ~62 μC cm-2, while the\nflexoelectric contribution is ~ 99 μC cm-2 (161 – 62 = 99).\nFigure R8 (Supplementary Figure 14 in the revised manuscript). Room-temperature P–E\nloops of (a) 2.0 mol% Fe-doped BNBT thin films and (b) undoped BNBT thin films. HAADF\nSTEM images of (c) 2.0 mol% Fe-doped BNBT thin films and (d) undoped BNBT thin films.\ninducing the flexoelectric effect, but it lacks an explanation of how the polarization\ncharacteristics change with varying oxygen vacancy densities. It is necessary to analyze\nwhether the flexoelectric effect increases as the number of oxygen vacancies increases\nor whether it saturates beyond a certain threshold.\nResponse: To reveal the dependence of polarization on the oxygen vacancy density,\nwe fabricated a series of thin films (350-nm-thick for a quicker screening) with various\nFe contents, i.e., 0.5, 1.0, 1.5, 2.0, 2.5, and 3.0 mol%, by PLD. As shown in Figure R1,\nboth the maximum and remanent polarization are optimal when the Fe doping level is\n2 mol%. Notably, while the best ferroelectric performance is observed at 2 mol%\ndoping, a significant increase in polarization only occurs when the Fe content increases\nfrom 0.5 to 1 mol%, beyond which the growth rate of polarization greatly slows down\ntill a Fe content of 2 mol% (Figure R1(h)). A higher Fe doping level (2.5 and 3.0 mol%)\nincurs a drop in polarizations, as shown in Figure R1. This observation suggests that\nthe Fe content already reaches an optimal level at 2.0 mol%.\nloops of 350 nm-thick BNBT thin films with (a) 0, (b) 0.5, (c) 1.0, (d) 1.5, (e) 2.0, (f) 2.5, and\n(g) 3.0 mol% Fe doping levels measured under selected electric fields. (h) Evolution of\n65\nmaximum and remanent polarizations of the 350 nm-thick BNBT thin films as a function of Fe\ndoping level.\nX-ray photoelectron spectroscopy (XPS) was employed to estimate the level of\noxygen vacancy in Fe-BNBT thin films, as shown in Figure R7. The relative level of\noxygen defects ( ) in the sample could be estimated by the:\n𝑂𝑉%\n𝐴𝑂𝑉\n𝑂𝑉% =\n𝐴𝑂𝐿 +𝐴𝑂𝑉 +𝐴𝑂𝐶\nwhere , and are the areas enclosed under O , O and O peak,\nL V, C\nrespecti𝐴v 𝑂 el 𝐿 y. 𝐴T 𝑂 h 𝑉 e XPS𝐴 𝑂 re 𝐶 sults reveal a progressive increase in oxygen vacancy\nconcentration with the increase in Fe doping level. Nevertheless, an excessive amount\nof oxygen vacancy (in 2.5 and 3.0 mol% Fe doped samples) will deteriorate the\nferroelectric properties, as shown in Figure R1, owing to the more pronounced domain-\nwall pinning effect24 and increased leakage current25. In summary, these findings\nindicate that the strength of flexoelectric effect in our Fe-BNBT thin is strongly\ndependent on the level of oxygen vacancies. It is evident that the flexoelectric response\nreaches a maximum at an Fe doping concentration of 2.0 mol %, i.e. a moderate level\nof oxygen vacancies as suggested by the XPS results.\n66\nelectronic level for Fe-BNBT thin films with (a) undoped and (b) 0.5, (c) 1.0, (d) 1.5, (e) 2.0,\n(f) 2.5, and (g) 3.0 mol% Fe doping levels. Blue (O ), green (O ), and purple (O ) peak is\nL V C\nascribed to the lattice oxygen, lattice oxygen connected with the lower valence cation after\nvacancy creation, and the adsorbed surface water and/or hydroxyl molecules, respectively. (h)\nEvolution of the Ov level in BNBT-based thin films as a function of Fe doping level.\nComments (4): Since the flexoelectric effect is determined by the strain gradient, it is\nlikely that the effect varies with film thickness. Although the experiment was primarily\nconducted on a 550 nm thick film, further research is needed to understand how the\nflexoelectric effect changes in thinner or thicker films.\nResponse: To understand the influence of film thickness on the flexoelectric effect and\npolarization behavior, 2.0 mol% Fe-doped BNBT thin films with thicknesses of 150,\n350, 750, and 950 nm were also deposited on STO (001) single crystal substrates by\nPLD under the same conditions as those for the 550-nm-thick sample. As shown in\nFigure R26, all the Fe-doped BNBT thin films exhibit well-defined P-E loops at room\ntemperature with very large polarization (both P and P), implying that the\nm r\nflexoelectric effect acts on all the thicknesses.\nP–E loops of Fe-BNBT thin films with a thickness of (a) 150 nm, (b) 350 nm, (c) 550 nm, (d)\n750 nm, and (e) 950 nm. (f) Comparison of P–E loops of Fe-BNBT thin films having different\nthicknesses.\nMulti-dimensional characterizations, including XRD, RSM, and STEM were\nconducted to perceive the structural origin of the thickness dependence of the observed\npolarization. The XRD patterns confirm the preferred (00l) orientation of all the thin\nfilms. Notably, the observation of a broad diffraction peak near 41o on the XRD pattern\nof the 950 nm-thick Fe-BNBT films suggests the presence of super-tetragonal phase,\ncharacterized by a significant c-axis expansion, as confirmed in our previous study9.\nRSM scans around the (103) plane of thinner films (150, 350, and 550 nm) indicate that\nthe Q values of STO, LSMO, and Fe-BNBT diffraction spots are almost identical,\nx\nsignifying that LSMO and Fe-BNBT layers remain fully strained on the STO substrate\n(in-plane direction). However, in the thicker films, i.e., 750 and 950 nm-thick, the Fe-\nBNBT diffraction spots broaden, suggesting a partial strain relaxation as film thickness\nincreases.\nθ-2θ scan results of Fe-BNBT thin films having different thicknesses. (b) XRD-RSM around\n(103) reflections of Fe-BNBT thin films having different thicknesses. Notably, a small bump\naround 41° is perceived in 950 nm-thick FeBNBT thin films, suggesting the creation of a super-\nT structure. The appearance of super-T structures may be attributed to the film being too thick,\npreventing the evaporation of Bismuth from the films during the film growth process.\nLow-magnification STEM images of Fe-BNBT thin films with thicknesses of\n150, 550, and 950 nm, as shown in Figure R28 (a), (b), and (c), respectively, show the\ndifference in the spatial distribution of the defects. We give a closer examination of\nthree representative regions along the thickness of the thin films, i.e. Area 1, near the\nsubstrate; Area 2, in the middle of the film; and Area 3, close to the film surface, as\nmarked by the red boxes in Figure R28 (a)-(c).\nIt is evident that the spatial distribution of planar defects is considerably different\nin the films of different thicknesses:\n1) In the 150 nm-thick film, planar defects are barely seen in Areas 1 and 2, but appear\nsparsely in Area 3.\n2) In the 550 nm-thick film, planar defects remain absent in Area 1 but are loosely\npresent in Area 2, with a noticeable increase in density in Area 3.\n3) In the 950 nm-thick film, planar defects are still absent in Area 1, while they render\na high density in Area 2. However, these defects are barely seen in Area 3, where β-\nBi₂O₃ layers begin to develop9.\nmagnification HAADF-STEM images of Fe-BNBT thin films with thicknesses of (a) 150 nm,\n(b) 550 nm, and (c) 950 nm. (d–f) HAADF-STEM images of the 150 nm-thick Fe-BNBT thin\nfilm, acquired from Areas 3, 2, and 1, respectively. (g–i) HAADF-STEM images of the 550\nnm-thick Fe-BNBT thin film, captured from the corresponding Areas 3, 2, and 1. (j–l) HAADF-\nrespectively.\nIt is known that epitaxial thin films with lattice mismatch relative to their\nsubstrates tend to relieve the epitaxial strains through the creation of\ndislocations/defects when the film thickness exceeds a critical value. The strain\nrelaxation typically proceeds gradually over an extended range of thickness,79 instead\nof an abrupt process at the critical thickness. According to the Matthews-Blakeslee\nmodel,80 which is widely used to estimate the critical thickness at which a coherent\nepitaxial film begins to relax, a smaller lattice mismatch often results in a larger critical\nthickness. For instance, in the case of (La Sr )MnO (LSMO) films with a lattice\n0.7 0.3 3\nmismatch of 0.83% grown on (001) STO substrates via PLD, it has been reported that\nthe film remains fully strained up to 100 nm and that relaxation remains incomplete\neven at a thickness of 432 nm.81 In addition to lattice mismatch, growth conditions,\nsuch as deposition rate, can also significantly influence the relaxation behavior, with\nhigh deposition rate generally delaying the onset of relaxation.82\nGiven the minimal lattice mismatch of 0.51% between bulk BNBT and the STO\nsubstrate, combined with the relatively high deposition frequency of 5 Hz employed for\nour Fe-doped BNBT thin films, it is reasonable to infer that the critical thickness for\nthe onset of strain relaxation substantially exceeds 100 nm. Consequently, the Fe-\nBNBT layers are expected to remain fully strained well beyond this threshold, with\nstrain relaxation remaining incomplete even at thicknesses greater than 432 nm, which\nis consistent with our STEM observations. In regions that are close to the substrate (<\n100 nm to the substrate), the films remain fully strained with a minimal level of planar\ndefects (see Figures R29 (e), (f), (i), (l)). As the film thickness increases (well beyond\n100 nm), partial strain relaxation occurs, accompanied by the emergence of planar\ndefects (Figures R29 (d), (g), (h), (k)). Upon further increase in thickness (i.e.,> ~500\nnm), the BNBT films transition to a fully relaxed state. In addition to the planar defects,\nβ-Bi₂O₃ layers are seen in Figure R26 (j), attributed to chemical inhomogeneity. This\nstructure evolution correlates strongly with the increasing roughness of the growing Fe-\nBNBT surface, compared to the atomically flat STO substrate. On a rougher surface,\nthe kinetic energy of incoming heavy adatoms, e.g., Bi atoms, is insufficient to support\nlong-range surface diffusion. As a result, Bi-rich species tend to accumulate locally,\nleading to the segregation and formation of β-Bi O as a secondary phase.9\n2 3\nThese observations are in good agreement with the strength of the flexoelectric\neffect and resultant polarization of respective Fe-BNBT thin films of different thickness,\ni.e., the higher the density of planar density, the stronger the flexoelectric effect and\noverall polarization, and justify the optimal polarization in the 550 nm-thick sample.\nDespite the reduced density of planar density, the slight increase in polarization for the\n950 nm-thick films, relative to the 750 nm-thick films, is likely attributed to the\nformation of super-tetragonal lattices due to the onset of the secondary phase, namely\nβ-Bi O , which can also enhance polarization as reported before.9\n2 3\nWe added the above experimental results and discussions to the Supporting\nInformation of the revised manuscript. We sincerely hope that we have clarified the\nthickness dependence of the flexoelectric effect and polarization in our Fe-BNBT thin\nfilms.\n“To understand the influence of film thickness on the flexoelectric effect and\npolarization behavior, 2.0 mol% Fe-doped BNBT thin films with thicknesses of 150,\n350, 750, and 950 nm were also deposited on STO (001) single crystal substrates by\nPLD under the same conditions as that for the 550-nm-thick sample. As shown in\nloops at room temperature with very large polarization (both P and P), implying that\nm r\nthe flexoelectric effect acts on all the thicknesses.\nMulti-dimensional characterizations, including XRD, RSM, and STEM were\nconducted to perceive the structural origin of the thickness dependence of the observed\npolarization. The XRD patterns confirm the preferred (00l) orientation of all the thin\nfilms. Notably, the observation of a broad diffraction peak near 41o on the XRD pattern\nof the 950 nm-thick Fe-BNBT films suggests the presence of super-tetragonal phase,\ncharacterized by a significant c-axis expansion, as confirmed in our previous study.\nRSM scans around the (103) plane of thinner films (150, 350, and 550 nm) indicate that\nthe Q values of STO, LSMO, and Fe-BNBT diffraction spots are almost identical,\nx\nsignifying that LSMO and Fe-BNBT layers remain fully strained on the STO substrate\n(in-plane direction). However, in the thicker films, i.e., 750 and 950 nm-thick, the Fe-\nBNBT diffraction spots broaden, suggesting a partial strain relaxation as film thickness\nincreases.\nLow-magnification STEM images of Fe-BNBT thin films with thicknesses of\n150, 550, and 950 nm, as shown in Supplementary Figure 18 (a), (b), and (c),\nrespectively, show the difference in the spatial distribution of the defects. We give a\ncloser examination of three representative regions along the thickness of the thin films,\ni.e. Area 1, near the substrate; Area 2, in the middle of the film; and Area 3, close to the\nfilm surface, as marked by the red boxes in Supplementary Figure 18 (a)-(c).\nIt is evident that the spatial distribution of planar defects is considerably different\nin the films of different thicknesses:\n75\n1) In the 150 nm-thick film, planar defects are barely seen in Areas 1 and 2, but appear\nsparsely in Area 3.\n2) In the 550 nm-thick film, planar defects remain absent in Area 1 but are loosely\npresent in Area 2, with a noticeable increase in density in Area 3.\n3) In the 950 nm-thick film, planar defects are still absent in Area 1, while they render\na high density in Area 2. However, these defects are barely seen in Area 3, where β-\nBi₂O₃ layers begin to develop.\nIt is known that epitaxial thin films with lattice mismatch relative to their\nsubstrates tend to relieve the epitaxial strains through the creation of\ndislocations/defects when the film thickness exceeds a critical value. The strain\nrelaxation typically proceeds gradually over an extended range of thickness, instead of\nan abrupt process at the critical thickness. According to the Matthews-Blakeslee model,\nwhich is widely used to estimate the critical thickness at which a coherent epitaxial film\nbegins to relax, a smaller lattice mismatch often results in a larger critical thickness.\nFor instance, in the case of (La Sr )MnO (LSMO) films with a lattice mismatch of\n0.7 0.3 3\n0.83% grown on (001) STO substrates via PLD, it has been reported that the film\nremains fully strained up to 100 nm and that relaxation remains incomplete even at a\nthickness of 432 nm. In addition to lattice mismatch, growth conditions, such as\ndeposition rate, can also significantly influence the relaxation behavior, with high\ndeposition rates generally delaying the onset of relaxation.\nGiven the minimal lattice mismatch of 0.51% between bulk BNBT and the STO\nsubstrate, combined with the relatively high deposition frequency of 5 Hz employed for\nour Fe-doped BNBT thin films, it is reasonable to infer that the critical thickness for\nthe onset of strain relaxation substantially exceeds 100 nm. Consequently, the Fe-\nBNBT layers are expected to remain fully strained well beyond this threshold, with\nstrain relaxation remaining incomplete even at thicknesses greater than 432 nm, which\nis consistent with our STEM observations. In regions that are close to the substrate (<\n100 nm to the substrate), the films remain fully strained with a minimal level of planar\ndefects (see Supplementary Figure 18 (e), (f), (i), (l)). As the film thickness increases\n(well beyond 100 nm), partial strain relaxation occurs, accompanied by the emergence\nof planar defects (Supplementary Figure 18 (d), (g), (h), (k)). Upon further increase in\nthickness (i.e.,> ~500 nm), the BNBT films transition to a fully relaxed state. In addition\nto the planar defects, β-Bi₂O₃ layers are seen in Supplementary Figure 18 (j), attributed\nto chemical inhomogeneity. This structure evolution correlates strongly with the\nincreasing roughness of the growing Fe-BNBT surface, compared to the atomically flat\nSTO substrate. On a rougher surface, the kinetic energy of incoming heavy adatoms,\ne.g., Bi atoms, is insufficient to support long-range surface diffusion. As a result, Bi-\nrich species tend to accumulate locally, leading to the segregation and formation of β-\nBi O as a secondary phase.\n2 3\nThese observations are in good agreement with the strength of the flexoelectric\neffect and resultant polarization of respective Fe-BNBT thin films of different thickness,\ni.e., the higher the density of planar density, the stronger the flexoelectric effect and\noverall polarization, and justify the optimal polarization in the 550 nm-thick sample.\nDespite the reduced density of planar density, the slight increase in polarization for the\n950 nm-thick films, relative to the 750 nm-thick films, is likely attributed to the\nformation of super-tetragonal lattices due to the onset of the secondary phase, namely\nβ-Bi O , which can also enhance polarization as reported before.”\n2 3\nComments (5): Although this study claims that Pm = 133 μC cm⁻² and Pr = 98 μC cm⁻²\nare maintained at 230°C, the specific reason for the high-temperature stability is not\n77\nexplained. Generally, problems such as oxygen vacancy movement, increased leakage\ncurrent, and domain wall movement occur at high temperatures, and this study should\nexplain how these factors were suppressed.\nResponse: The temperature-dependent leakage current density of our Fe-BNBT thin\nfilms (550-nm thick) is presented in Figure R9. Even at 250 °C, the leakage current\ndensity remains as low as ~10-4 A cm-2 under an applied field of 750 kV cm-1, which is\nconsidered minimal for typical ferroelectric thin films possessing a perovskite\nstructure.33-35\nFigure R9 (Supplementary Figure 7 in the revised manuscript). Leakage current density as\na function of the applied electric fields of Fe-BNBT thin films at various temperatures.\nIn dielectric thin films, the charge transport is often dominated by bulk-limited\nconduction mechanisms, including ohmic conduction, space-charge-limited conduction\n(SCLC), and Poole-Frenkel emission (PF), or interface-limited Schottky emission\n78\n(SE).36 The space-charge-limited current obeys quadratic electric field dependence in\nthe case of discrete traps,37\n2\n𝐽𝑆𝐶𝐿𝐶 = 9𝜇ɛ0ɛ𝑟𝐸 /8𝑑\nwhere μ is the carrier mobility and d the film thickness. The current density across a\nSchottky barrier is37\n1\n∗ 2 2\nwhere 𝐽𝑆 =, 𝐴 𝑇is ethxep R−i c[h(𝜙ar 𝑏 ds−on 𝛽 c 𝑠 o𝐸ns)t/a𝑘n𝑇t, ]T the temperature, q the\n3 1/2 ∗\nelectro𝛽n 𝑠 c=ha(r𝑞ge,/ 4𝜋ɛ 0 thɛe 𝑟 )Scho𝐴ttky barrier height, k the Boltzmann constant, the\npermittivity of fre𝜙e 𝑏 space, and the dielectric constant of the thin film. The ɛP 0 ool-\nFrenkel conductivity is given byɛ 𝑟\n1\n2\n𝐽𝑃𝐹 = 𝜎0𝐸exp− [(𝐸1− 𝛽𝑃𝐹𝐸 )/𝑘𝑇]\nwhere , is a sample-dependent zero-field conductivity, and\n3 1/2\nis the tr𝛽a 𝑃 p 𝐹 io=ni(z𝑞ati/o𝜋nɛ e 0 nɛe 𝑟 r)gy. T𝜎h 0 e functional relationship of SE and PF is similar, so th𝐸a 1 t\nthe plots of ln(J T-2) versus E0.5 and ln(J V-1) versus E0.5 for SE and PF, respectively,\nwill yield linear relationships.\nWe plotted and fitted the leakage current data based on the above-mentioned\ntransport mechanisms, as shown in Figure R10. At room temperature, the leakage\ncurrent is governed by PF under the electric field < 150 kV cm-1, given the relative\ndielectric constant extracted from the slope of the linear fits to PF is 118.8, which lies\nwithin the physically reasonable range between the optical dielectric constant (6.71)38\nand quasi-static dielectric constant (500- 800)39 typical for BNT-based thin films.\nApart from the low field region (< 150 kV cm-1), the Schottky emission\nmechanism becomes dominant, supported by the extracted dielectric constant of 52.8,\n79\nwhich is also physically plausible. SCLC can be ruled out in the entire range of\nmeasuring temperature and electric field, given slope of the LnJ-LnE curves evidently\ndeviates from 2, as shown in Figure R10(a). At 250 °C, the leakage current\nfollows the PF mechanism across the entire range of the measuring electric field.\nFigure R10 (Supplementary Figure 8 in the revised manuscript). Analysis of leakage\ncurrent mechanisms in Fe-BNBT thin films based on (a) SCLC, (b) SE, and (c) PF emissions.\nThe suppression of leakage currents at high temperature in this work can be\nprimarily attributed to the trapping of oxygen vacancies and free charges at T-T and H-\nH domain walls, respectively. These domain walls serve as effective charge sinks,\nimpeding carrier transport and thereby inhibiting the formation and propagation of\nelectrical trees.40 Consequently, this mechanism significantly reduces leakage currents,\nensuring stable ferroelectric performance even at elevated temperatures.\nWe added these experimental results and discussions to the Supporting\nInformation of the revised manuscript.\n“The leakage current density vs electric field curves at various temperatures of our Fe-\nBNBT thin films were shown in Supplementary Figure 7. The leakage current density\nremains as low as ~10-4 A cm-2 even at elevated temperatures, which is considered\n80\nminimal for typical ferroelectric thin films possessing a perovskite structure. The\nrelatively low leakage current density can be attributed to the trapping of oxygen\nvacancies and free charges at T-T and H-H domain walls, respectively. These domain\nwalls act as potential barriers that impede charge transport, thereby suppressing the\nformation and growth of electrical trees and effectively reducing leakage currents,\nparticularly at high temperatures.10\nIn dielectric thin films, the charge transport is often dominated by bulk-limited\nconduction mechanisms, including ohmic conduction, space-charge-limited conduction\n(SCLC), and Poole-Frenkel emission (PF), or interface-limited Schottky emission (SE).\nThe space-charge-limited current obeys quadratic electric field dependence in the case\nof discrete traps,12\n2\n𝐽𝑆𝐶𝐿𝐶 = 9𝜇ɛ0ɛ𝑟𝐸 /8𝑑\nwhere μ is the carrier mobility and d the film thickness. The current density across a\nSchottky barrier is12\n1\n∗ 2 2\n𝐽𝑆 = 𝐴 𝑇 exp− [(𝜙𝑏 − 𝛽𝑠𝐸 )/𝑘𝑇]\nwhere , is the Richardson constant, T the temperature, q the\n3 1/2 ∗\nelectro𝛽n 𝑠 c=ha(r𝑞ge,/ 4𝜋ɛ 0 thɛe 𝑟 )Scho𝐴ttky barrier height, k the Boltzmann constant, the\npermittivity of free s𝜙p 𝑏 ace, and the dielectric constant of the thin film. The ɛP 0 ool-\nFrenkel conductivity is given byɛ 𝑟\n1\n2\n𝐽𝑃𝐹 = 𝜎0𝐸exp− [(𝐸1− 𝛽𝑃𝐹𝐸 )/𝑘𝑇]\nwhere , is a sample-dependent zero-field conductivity, and\n3 1/2\nis the tr𝛽a 𝑃 p 𝐹 io=ni(z𝑞ati/o𝜋nɛ e 0 nɛe 𝑟 r)gy. T𝜎h 0 e functional relationship of SE and PF is similar, so th𝐸a 1 t\n81\nthe plots of ln(J T-2) versus E0.5 and ln(J V-1) versus E0.5 for SE and PF, respectively,\nwill yield linear relationships.\nWe plotted and fitted the leakage current data based on the above-mentioned\ntransport mechanisms, as shown in Supplementary Figure 8. At room temperature, the\nleakage current is governed by PF under the electric field < 150 kV cm-1, given the\nrelative dielectric constant extracted from the slope of the linear fits to PF is 118.8,\nwhich lies within the physically reasonable range between the optical dielectric constant\n(6.71) and quasi-static dielectric constant (500- 800) typical for BNT-based thin films.\nApart from the low field region (< 150 kV cm-1), the Schottky emission\nmechanism becomes dominant, supported by the extracted dielectric constant of 52.8,\nwhich is also physically plausible. SCLC can be ruled out in the entire range of\nmeasuring temperature and electric field, given slope of the LnJ-LnE curves evidently\ndeviates from 2, as shown in Supplementary Figure 8 (a). At 250 °C, the leakage current\nfollows the PF mechanism across the entire range of the measuring electric field.”\nComments (6): The current study focused on Fe doping, but it is necessary to compare\nand analyze with previous studies how the flexoelectric effect changes when doping\nwith other elements (Mn, Co, V, etc.).\nResponse: The reported ferroelectric properties of typical doped BNBT thin films\nfabricated by PLD are summarized in Table 4. Notably, none of the previous studies\nhas reported a polarization of beyond 100 µC cm-2 in BNBT-based thin films.\nIn addition, we also studied the effect of dopant species on the polarization of\nresultant epitaxial BNBT thin films. To that end, 550 nm-thick BNBT-based thin films\ndoped with 2.0 mol% of Mn, La, Sm, Eu, and Co, along with undoped BNBT thin films,\nwere deposited on STO (001) single crystal substrates by PLD under optimal processing\nconditions. As shown in Figure R18, the maximum polarization of the undoped, Mn-,\nLa-, Sm-, Eu, and Co-doped BNBT films is 62, 73, 68, 73, 74, and 97 µC cm-2,\n82\nrespectively. In stark contrast, Fe-doped BNBT films exhibit a considerably higher\nmaximum polarization of ~161 µC cm-2.\nThese results unequivocally suggest that Fe doping plays a decisive role in\nenhancing the polarization of BNBT thin films. We attribute the giant polarization of\nFe-doped BNBT thin films to several fundamental mechanisms as follows:\n1) Fe3+ (r = 58.5 pm) preferentially replaces Ti4+ (r = 60.5 pm) at the B-site of\nFe3+ Ti4+\nthe perovskite lattice, as an acceptor, leading to the formation of oxygen vacancies to\ncompensate for the charge imbalance. These vacancies tend to organize along specific\ncrystallographic planes, such as [101]orientation.17, 18 In contrast, other commonly\nemployed dopants in BNBT lattice such as La3+ (r = 136 pm), Sm3+ (r = 124 pm),\nLa3+ Sm3+\nand Eu3+ (r = 129 pm) substitute for A-site cations (r = 139 pm, r = 120 pm,\nEu3+ Na+ Bi3+\nr = 135 pm) as donors, due to their similar ionic radii to that of the A-site cations,\nBa2+\nwhich do not promote the formation of oxygen vacancies.26\n2) The multivalent nature of Mn (Mn2+, Mn3+, Mn4+) in principle allows them to enter\ndifferent lattice sites of BNBT. Considering the similarity of the ionic radii between\nMn cations (r ~ 0.83 Å, r ~ 0.65 Å, r ~ 0.53 Å) and the respective cations\nMn2+ Mn3+ Mn4+\nin the host lattice, as well as the global instable index indices calculated in our previous\nwork27, 28, Mn2+ prefers to reside in the A-site, and Mn3+ and Mn4+ are inclined to\nsubstitute Ti4+ at the B-site. It should be noted that the presence of Mn2+ often impedes\nthe formation of oxygen vacancies.29 More importantly, Mn dopants tend to introduce\nlocalized point defects rather than extended planar faults.30\n3) Co3+ (r = 55 pm) can also substitute for Ti4+ at the B-site and is reported to be\nCo3+\nmore effective than Fe in promoting the creation of oxygen vacancies.31 However, this\nincreased oxygen vacancy concentration is often accompanied by an increased\n83\nelectrical conductivity, which may incur leakage current pathways and deteriorate\nferroelectric performance, as shown in Figure R18.\nTable R2. Comparative analysis of ferroelectric properties in BNBT thin films doped with\nvarious elements.\nSamples Pr Pm Ec Em Substrate Thickness Ref.\nµC cm-2 µC cm-2 kV cm-1 kV cm-1 nm\nMn-doped BNBT 33 68 200 800 (001) SrTiO 250 27\nBNBT (LaAlO)TiO(SrAlTa\n3 3 2\nO)\n60.35\nMn-doped BNBT 25 48 150 650 (001) LaAlO 162 55\n3\n84\nloops of (a) undoped BNBT thin film and BNBT-based thin films doped with 2.0 mol% of (b)\nMn, (c) La, (d) Sm, (e) Eu, (f) Co, and (g) Fe, respectively. (h) Comparison of P–E loops of\nBNBT thin films with different dopants.\nComments (7): In general, typos (emperirical, etc.) and awkward expressions should\nbe corrected.\nResponse: We have thoroughly checked the manuscript and corrected the\ntypographical errors, such as \"emperirical\", which has been amended to “empirical”,\nas well as improved awkward expressions, including replacing “twills” with the more\nappropriate term of “streaks”.\nResponse to the comments of Referee #4\nI co-reviewed this manuscript with one of the reviewers who provided the listed\nreports. This is part of the Nature Communications initiative to facilitate training in\npeer review and to provide appropriate recognition for Early Career Researchers who\nco-review manuscripts.\nResponse: We truly appreciate the time and effort you’ve invested in reviewing our\nmanuscript and providing constructive feedback. Your insights have been incredibly\nhelpful in improving the quality of our work.\n(1) Shepelin, N. A.; Tehrani, Z. P.; Ohannessian, N.; Schneider, C. W.; Pergolesi, D.; Lippert,\nT.; A Practical Guide to Pulsed Laser Deposition. Chemical Society Reviews 2023, 52 (7), 2294-\n2321.\n(2) Aziz, M. J.; Film Growth Mechanisms in Pulsed Laser Deposition. Applied Physics A 2008,\n93, 579-587.\n(3) Chrisey, D. B.; Hubler, G., Pulsed Laser Deposition. 1994.\n(4) Christen, H. M.; Eres, G.; Recent Advances in Pulsed-Laser Deposition of Complex Oxides.\nJournal of Physics: Condensed Matter 2008, 20 (26), 264005.\n(5) Jiang, Y.; Parsonnet, E.; Qualls, A.; Zhao, W.; Susarla, S.; Pesquera, D.; Dasgupta, A.;\nAcharya, M.; Zhang, H.; Gosavi, T.; Enabling Ultra-Low-Voltage Switching in Batio3. Nat.\nMater. 2022, 21 (7), 779-785.\n(6) Pan, H.; Zhu, M.; Banyas, E.; Alaerts, L.; Acharya, M.; Zhang, H.; Kim, J.; Chen, X.; Huang,\nXuang,\nX\n(7) Zhang, L.; Chen, J.; Fan, L.; Diéguez, O.; Cao, J.; Pan, Z.; Wang, Y.; Wang, J.; Kim, M.;\nDeng, S.; Giant Polarization in Super-Tetragonal Thin Films through Interphase Strain. Science\n2018, 361 (6401), 494-497.\n(8) Wang, Y.; Zhang, L.; Wang, J.; Li, Q.; Wang, H.; Gu, L.; Chen, J.; Deng, J.; Lin, K.; Huang,\nL.; Chemical-Pressure-Modulated Batio3 Thin Films with Large Spontaneous Polarization and\nHigh Curie Temperature. J. Am. Chem. Soc. 2021, 143 (17), 6491-6497.\n(9) Sun, Y.; Zhang, L.; Huang, Q.; Chen, Z.; Wang, D.; Seyfouri, M. M.; Chang, S. L.; Wang,\nY.; Zhang, Q.; Liao, X.; Ultrahigh Energy Storage Density in Glassy Ferroelectric Thin Films\nunder Low Electric Field. Adv. Sci. 2022, 9 (31), 2203926.\n(10) Béa, H.; Dupé, B.; Fusil, S.; Mattana, R.; Jacquet, E.; Warot-Fonrose, B.; Wilhelm, F.;\nRogalev, A.; Petit, S.; Cros, V.; Evidence for Room-Temperature Multiferroicity in a\nCompound with a Giant Axial Ratio. Phys. Rev. Lett. 2009, 102 (21), 217603.\n(11) Marrocchelli, D.; Perry, N. H.; Bishop, S. R.; Understanding Chemical Expansion in\nPerovskite-Structured Oxides. Physical Chemistry Chemical Physics 2015, 17 (15), 10028-\n10039.\n(12) Španková, M.; Štrbík, V.; Chromik, Š.; Zheng, D.; Li, J.; Machajdík, D.; Kobzev, A.;\nPlecenik, T.; Sojková, M.; Characterization of Epitaxial Lsmo Films Grown on Sto Substrates.\nActa Phys. Pol., A 2017, 131 (4), 848-850.\n(13) Zhang, L.; Wang, H.; Wang, D.; Guo, M.; Lou, X.; Wang, D.; A New Strategy for Large\nDynamic Piezoelectric Responses in Lead‐Free Ferroelectrics: The Relaxor/Morphotropic\nPhase Boundary Crossover. Adv. Funct. Mater. 2020, 30 (45), 2004641.\n(14) Arya, M.; Kumar, S.; Hasina, D.; Sen, R.; Ojha, S.; Kumar, V.; Som, T.; Dhar, S.;\nCombining Experimental and Modelling Approaches to Understand the Expansion of Lattice\nParameter of Epitaxial Srti1-Xtaxo3 (X= 0–0.1) Films. Computational Materials Science 2023,\n217, 111917.\n(15) Son, J.; Kim, B. G.; Cho, J.; Thin Film Growth of Epitaxial, Polycrystalline and\nAmorphous Srruo3. Thin Solid Films 2007, 515 (18), 7086-7090.\n(16) Haque, A.; D'Souza, H. J.; Parate, S. K.; Sandilya, R. S.; Raghavan, S.; Nukala, P.;\nHeterogeneous Integration of High Endurance Ferroelectric and Piezoelectric Epitaxial Batio3\nDevices on Si. Adv. Funct. Mater. 2025, 35 (1), 2413515.\n(17) Becerro, A.; McCammon, C.; Langenhorst, F.; Seifert, F.; Angel, R.; Oxygen Vacancy\nOrdering in Catio3–Cafeo2. 5 Perovskites: From Isolated Defects to Infinite Sheets. Phase\nTransitions 1999, 69 (1), 133-146.\n(18) Scott, J.; Dawber, M.; Oxygen-Vacancy Ordering as a Fatigue Mechanism in Perovskite\nFerroelectrics. Appl. Phys. Lett. 2000, 76 (25), 3801.\n(19) Jia, C.; Urban, K.; Mertin, M.; Hoffmann, S.; Waser, R.; The Structure and Formation of\nNanotwins in Batio3 Thin Films. Philos. Mag. A 1998, 77 (4), 923-939.\n(20) Eibl, O.; Pongratz, P.; Skalicky, P.; Schmelz, H.; Formation of (111) Twins in Batio3\nCeramics. J. Am. Ceram. Soc. 1987, 70 (8), C‐195-C‐197.\n87\n(21) Eibl, O.; Pongratz, P.; Skalicky, P.; Crystallography of (111) Twins in Batio3. Philos. Mag.\nB 1988, 57 (4), 521-534.\n(22) Reþnik, A.; Twins in Barium Titanate. Acta Chim. Slov 2001, 48, 1-50.\n(23) Huang, C.; Liao, Z.; Li, M.; Guan, C.; Jin, F.; Ye, M.; Zeng, X.; Zhang, T.; Chen, Z.; Qi,\nY.; A Highly Strained Phase in Pbzr0. 2ti0. 8o3 Films with Enhanced Ferroelectric Properties.\nAdv. Sci. 2021, 8 (8), 2003582.\n(24) Tan, Q.; Li, J.; Viehland, D.; Role of Lower Valent Substituent-Oxygen Vacancy\nComplexes in Polarization Pinning in Potassium-Modified Lead Zirconate Titanate. Appl. Phys.\nLett. 1999, 75 (3), 418-420.\n(25) Singh, S.; Ishiwara, H.; Maruyama, K.; Enhanced Polarization and Reduced Leakage\nCurrent in Bifeo3 Thin Films Fabricated by Chemical Solution Deposition. J. Appl. Phys. 2006,\n100 (6).\n(26) Bendahhou, A.; Chourti, K.; Loutou, M.; El Barkany, S.; Abou-Salama, M.; Impact of\nRare Earth (Re 3+= La 3+, Sm 3+) Substitution in the a Site Perovskite on the Structural, and\nElectrical Properties of Ba (Zr 0.9 Ti 0.1) O 3 Ceramics. RSC Adv. 2022, 12 (18), 10895-10910.\n(27) Lin, Q.; Ding, R.; Li, Q.; Tay, Y. Y.; Wang, D.; Liu, Y.; Huang, Y.; Li, S.; Large\nPiezoelectricity and Ferroelectricity in Mn‐Doped (Bi0. 5na0. 5) Tio3‐Batio3 Thin Film\nPrepared by Pulsed Laser Deposition. J. Am. Ceram. Soc. 2016, 99 (7), 2347-2353.\n(28) Sun, Y.; Yang, J.; Li, S.; Wang, D.; Defect Engineering in Perovskite Oxide Thin Films.\nChem. Commun. 2021, 57 (68), 8402-8420.\n(29) Shuai, Y.; Zhou, S.; Bürger, D.; Reuther, H.; Skorupa, I.; John, V.; Helm, M.; Schmidt,\nH.; Decisive Role of Oxygen Vacancy in Ferroelectric Versus Ferromagnetic Mn-Doped Batio3\nThin Films. J. Appl. Phys. 2011, 109 (8).\n(30) Maier, R.; Pomorski, T.; Lenahan, P.; Randall, C.; Acceptor-Oxygen Vacancy Defect\nDipoles and Fully Coordinated Defect Centers in a Ferroelectric Perovskite Lattice: Electron\nParamagnetic Resonance Analysis of Mn2+ in Single Crystal Batio3. J. Appl. Phys. 2015, 118\n(16).\n(31) Chang, H.; Bjørgum, E.; Mihai, O.; Yang, J.; Lein, H. L.; Grande, T.; Raaen, S.; Zhu, Y.-\nA.; Holmen, A.; Chen, D.; Effects of Oxygen Mobility in La–Fe-Based Perovskites on the\nCatalytic Activity and Selectivity of Methane Oxidation. ACS Catalysis 2020, 10 (6), 3707-\n3719.\n(32) Petralanda, U.; Kruse, M.; Simons, H.; Olsen, T.; Oxygen Vacancies Nucleate Charged\nDomain Walls in Ferroelectrics. Phys. Rev. Lett. 2021, 127 (11), 117601.\n(33) Yan, H.; Song, B.; Zhu, K.; Xu, L.; Shen, B.; Zhai, J.; Achieved High Energy Density and\nExcellent Thermal Stability in (1− X)(Bi0. 5na0. 5) 0.94 Ba0. 06tio3− X Bi (Mg0. 5ti0. 5) O3\nRelaxor Ferroelectric Thin Films. Journal of Materials Science: Materials in Electronics 2021,\n32 (12), 16269-16278.\n(34) Hejazi, M.; Safari, A.; Temperature-Dependent Leakage Current Behavior of Epitaxial\nBi0. 5na0. 5tio3-Based Thin Films Made by Pulsed Laser Deposition. J. Appl. Phys. 2011, 110\n(10).\n(35) Sun, Y.; Du, J.; Jiang, C.; Liang, J.; Geng, X.; Wang, Y.; Wang, D.; Large Negative\nElectrocaloric Response Induced by Nanoscale Phase Transition in (Bi, Na) Tio3-Based Thin\nFilms. Appl. Phys. Lett. 2023, 122 (9).\n(36) Wang, D.; Lin, D.; Kwok, K. W.; Chan, N.; Dai, J.; Li, S.; Chan, H.; Ferroelectric,\nPiezoelectric, and Leakage Current Properties of (K 0.48 Na 0.48 Li 0.04)(Nb 0.775 Ta 0.225)\nO 3 Thin Films Grown by Pulsed Laser Deposition. Appl. Phys. Lett. 2011, 98 (2), 022902.\n(37) Zubko, P.; Jung, D.; Scott, J.; Electrical Characterization of Pbzr 0.4 Ti 0.6 O 3 Capacitors.\nJ. Appl. Phys. 2006, 100 (11), 114113.\n(38) Yang, H.; Suvorova, N.; Jain, M.; Kang, B.; Li, Y.; Hawley, M.; Dowden, P.; DePaula, R.;\nJia, Q.; Lu, C.; Effective Thickness and Dielectric Constant of Interfacial Layers of Pt∕ Bi3.\n15nd0. 85ti3o12∕ Srruo3 Capacitors. Appl. Phys. Lett. 2007, 90 (23).\n(39) Grove, K. M.; Fox, A.; Cann, D. P.; Ko, S. W.; Mardilovich, P.; Gibbons, B. J.; Impact of\nProcess Conditions on Chemical Solution Deposited Bnkt Thin Film Electromechanical\nProperties. SN Applied Sciences 2021, 3, 1-10.\n(40) Ahmad, M.; Bashir, N.; Ahmad, H.; Abd Jamil, A.; Suleiman, A.; An Overview of\nElectrical Tree Growth in Solid Insulating Material with Emphasis of Influencing Factors,\nMathematical Models and Tree Suppression. TELKOMNIKA Indonesian Journal of Electrical\nEngineering 2014, 12 (8), 5827-5846.\n(41) Meng, Y.; Ju, X.; Yang, X.; The Measurement of the Dislocation Density Using Tem.\nMater. Charact. 2021, 175, 111065.\n(42) Ham, R.; The Determination of Dislocation Densities in Thin Films. Philos. Mag. 1961, 6\n(69), 1183-1184.\n(43) Ma, J.; Ma, J.; Zhang, Q.; Peng, R.; Wang, J.; Liu, C.; Wang, M.; Li, N.; Chen, M.; Cheng,\nX.; Controllable Conductive Readout in Self-Assembled, Topologically Confined Ferroelectric\nDomain Walls. Nat. Nanotechnol. 2018, 13 (10), 947-952.\n(44) Lu, H.; Tan, Y.; McConville, J. P.; Ahmadi, Z.; Wang, B.; Conroy, M.; Moore, K.; Bangert,\nU.; Shield, J. E.; Chen, L. Q.; Electrical Tunability of Domain Wall Conductivity in Linbo3\nThin Films. Adv. Mater. 2019, 31 (48), 1902890.\n(45) De Souza, R.; Oxygen Diffusion in Srtio3 and Related Perovskite Oxides. Adv. Funct.\nMater. 2015, 25 (40), 6326-6342.\n(46) Yamada, H.; Garcia, V.; Fusil, S.; Boyn, S.; Marinova, M.; Gloter, A.; Xavier, S.; Grollier,\nJ.; Jacquet, E.; Carrétéro, C.; Giant Electroresistance of Super-Tetragonal Bifeo3-Based\nFerroelectric Tunnel Junctions. ACS nano 2013, 7 (6), 5385-5390.\n(47) Li, X.; Li, J.; Li, Y.; Liu, X.; Yang, S.; Wu, J.; Hou, D.; Zhang, J.; Wu, H.; Zhang, Y.;\nHigh Performance Relaxor Ferroelectric Textured Ceramics for Electrocaloric Refrigeration.\nNat. Commun. 2025, 16 (1), 1-10.\n(48) Ding, R.; Wang, D.; Chu, D.; Li, S.; Crystallographic Orientation Dependence on\nElectrical Properties of (Bi, Na) Tio 3‐Based Thin Films. J. Am. Ceram. Soc. 2013, 96 (11),\n3530-3535.\n(49) Cao, H.; Fang, B.; Xu, H.; Luo, H.; Crystal Orientation Dependence of Dielectric and\nPiezoelectric Properties of Tetragonal Pb (Mg1/3nb2/3) O3–38% Pbtio3 Single Crystal. Mater.\nRes. Bull. 2002, 37 (13), 2135-2143.\n(50) Bokov, A. A.; Ye, Z.-G.; Dielectric Relaxation in Relaxor Ferroelectrics. Journal of\nAdvanced dielectrics 2012, 2 (02), 1241010.\n(51) Kobald, H.; Kobald, A. M.; Panzic, I.; Deluca, M.; Enhanced Energy Storage in Relaxor\n(1-X) Bi0. 5na0. 5tio3-Xbazryti1-Yo3 Thin Films by Morphotropic Phase Boundary\nEngineering. Communications materials 2025, 6 (1), 12.\n(52) Zhao, J.; Niu, G.; Ren, W.; Wang, L.; Zhang, N.; Sun, Y.; Wang, Q.; Shi, P.; Liu, M.;\nZhao, Y.; Polarization Behavior Oflead-Free 0.94 (Bi0. 5na0. 5) Tio3-0.06 Batio3 Thin Films\nwith Enhanced Ferroelectric Properties. J. Eur. Ceram. Soc. 2020, 40 (12), 3928-3935.\n(53) Alonso‐Sanjosé, D.; Jiménez, R.; Bretos, I.; Calzada, M. L.; Lead‐Free Ferroelectric\n(Na1/2bi1/2) Tio3–Batio3 Thin Films in the Morphotropic Phase Boundary Composition:\nSolution Processing and Properties. J. Am. Ceram. Soc. 2009, 92 (10), 2218-2225.\n(54) Acharya, S. K.; Ahn, B.-G.; Jung, C. U.; Koh, J.-H.; Choi, I.-H.; Lee, S.-K.; Effect of Rb\nDoping on Ferroelectric and Piezoelectric Properties of Bi0. 5na0. 5tio3–Batio3 Thin Films. J.\nAlloys Compd. 2014, 603, 248-254.\n(55) Gallegos-Melgar, A.; Espinosa-Arbelaez, D.; Flores-Ruiz, F.; Lahmar, A.; Dellis, J.-L.;\nLemée, N.; Espinoza-Beltran, F.; Muñoz-Saldaña, J.; Ferroelectric Properties of Manganese\nDoped (Bi1/2na1/2) Tio3 and (Bi1/2na1/2) Tio3–Batio3 Epitaxial Thin Films. Appl. Surf. Sci.\n2015, 359, 923-930.\n(56) Acharya, S. K.; Ahn, B. G.; Hyung, J.-H.; Lee, S.-K.; Effect of Li Doping on Ferroelectric\nand Piezoelectric Properties of Ba 0.5 Na 0.5 Tio 3-Batio 3 (Bnt-Bt) Thin Films. J. Korean\nPhys. Soc. 2013, 62, 794-799.\n(57) Daniels, J. E.; Jo, W.; Rödel, J.; Jones, J. L.; Electric-Field-Induced Phase Transformation\nat a Lead-Free Morphotropic Phase Boundary: Case Study in a 93%(Bi0. 5na0. 5) Tio3–7%\nBatio3 Piezoelectric Ceramic. Appl. Phys. Lett. 2009, 95 (3).\n(58) Franzbach, D. J.; Seo, Y.-H.; Studer, A. J.; Zhang, Y.; Glaum, J.; Daniels, J. E.; Koruza,\nJ.; Benčan, A.; Malič, B.; Webber, K. G.; Electric-Field-Induced Phase Transitions in Co-\nDoped Pb (Zr1− Xtix) O3 at the Morphotropic Phase Boundary. Science and Technology of\nAdvanced Materials 2014, 15 (1), 015010.\n(59) Guan, C.; Wei, Y.; Dong, S.; Xia, L.; Jin, C.; Nan, R.; Hu, L.; Gao, L.; Enhanced\nPiezoelectric Performance and Improved Polarization in 0.9 (Na0. 5bi0. 5) Tio3-0.1 Srtio3\nSingle Crystals Via Flux Growth Method. J. Alloys Compd. 2024, 1007, 176365.\n(60) Zhang, L.; He, L.; Zhang, C.; Xu, S.; Li, S.; Yang, Z.; Zhang, Y.; Wu, H.; Wang, D.; Yang,\nS.; Superior Electrostrain with Excellent Thermal Stability by Constructing Glassy\nFerroelectric Crossover. Acta Mater. 2025, 121153.\n(61) Yang, Y.; Liu, C.; Ji, Y.; He, L.; Ren, X.; Designed Morphotropic Relaxor Boundary\nCeramic Exhibiting Large Electrostrain and Negligible Hysteresis. Acta Mater. 2021, 208,\n116720.\n(62) Chaplya, P. M.; Carman, G. P.; Dielectric and Piezoelectric Response of Lead Zirconate–\nLead Titanate at High Electric and Mechanical Loads in Terms of Non-180 Domain Wall\nMotion. J. Appl. Phys. 2001, 90 (10), 5278-5286.\n(63) Chen, K.; Ma, J.; Wu, J.; Wang, X.; Miao, F.; Huang, Y.; Shi, C.; Wu, W.; Wu, B.; Improve\nPiezoelectricity in Batio 3-Based Ceramics with Large Electrostriction Coefficient. Journal of\nMaterials Science: Materials in Electronics 2020, 31, 12292-12300.\n(64) Xu, Z.; Shi, X.; Liu, Y.-X.; Wang, D.; Thong, H.-C.; Jiang, Y.; Sha, Z.; Li, Z.; Yao, F.-Z.;\nCai, X.-X.; High Electrostrain in a Lead-Free Piezoceramic from a Chemopiezoelectric Effect.\nNat. Mater. 2025, 1-9.\n(65) Li, F.; Cabral, M. J.; Xu, B.; Cheng, Z.; Dickey, E. C.; LeBeau, J. M.; Wang, J.; Luo, J.;\nTaylor, S.; Hackenberger, W.; Giant Piezoelectricity of Sm-Doped Pb (Mg1/3nb2/3) O3-Pbtio3\nSingle Crystals. Science 2019, 364 (6437), 264-268.\n(66) Qiu, Q.; Wang, H.; Xu, Y.; Chen, X.; Liu, X.; Zhou, Y.; Lv, Z.; Sun, J.; Zhang, S.; Zhang,\nH.; Achieving High Electrostrain Performance in Bnt‐Based Lead‐Free Piezoelectric Ceramics\nModified by Sr (Sn0. 5ta0. 4) O3. J. Am. Ceram. Soc. 2025, 108 (8), e20595.\n(67) Gunkel, F.; Christensen, D. V.; Chen, Y.; Pryds, N.; Oxygen Vacancies: The (in) Visible\nFriend of Oxide Electronics. Appl. Phys. Lett. 2020, 116 (12).\n(68) Luo, H.; Liu, H.; Huang, H.; Song, Y.; Tucker, M. G.; Sun, Z.; Yao, Y.; Gao, B.; Ren, Y.;\nTang, M.; Achieving Giant Electrostrain of above 1% in (Bi, Na) Tio3-Based Lead-Free\nPiezoelectrics Via Introducing Oxygen-Defect Composition. Sci. Adv. 2023, 9 (5), eade7078.\n(69) Zhang, W.; Feng, Q.; Hosono, E.; Asakura, D.; Miyawaki, J.; Harada, Y.; Tetragonal\nDistortion of a Batio3/Bi0. 5na0. 5tio3 Nanocomposite Responsible for Anomalous\nPiezoelectric and Ferroelectric Behaviors. ACS omega 2020, 5 (36), 22800-22807.\n(70) Nagarajan, V.; Roytburd, A.; Stanishevsky, A.; Prasertchoung, S.; Zhao, T.; Chen, L.-Q.;\nMelngailis, J.; Auciello, O.; Ramesh, R.; Dynamics of Ferroelastic Domains in Ferroelectric\nThin Films. Nat. Mater. 2003, 2 (1), 43-47.\n(71) Moore, K.; O’Connell, E. N.; Griffin, S. M.; Downing, C.; Colfer, L.; Schmidt, M.;\nNicolosi, V.; Bangert, U.; Keeney, L.; Conroy, M.; Charged Domain Wall and Polar Vortex\nTopologies in a Room-Temperature Magnetoelectric Multiferroic Thin Film. ACS Appl. Mater.\nInterfaces 2022, 14 (4), 5525-5536.\n(72) Park, S. M.; Wang, B.; Chen, L.-Q.; Noh, T. W.; Yang, S. M.; Lee, D.; Flexoelectric\nControl of Physical Properties by Atomic Force Microscopy. Appl. Phys. Rev. 2021, 8 (4).\n(73) Epp, J., X-Ray Diffraction (Xrd) Techniques for Materials Characterization. In Materials\nCharacterization Using Nondestructive Evaluation (Nde) Methods, Elsevier: 2016; pp 81-124.\n(74) Harrington, G. F.; Santiso, J.; Back-to-Basics Tutorial: X-Ray Diffraction of Thin Films.\nJ. Electroceram. 2021, 47 (4), 141-163.\n(75) Wu, W.; Stoica, A. D.; Berry, K. D.; Frost, M. J.; Skorpenske, H. D.; An, K.; Pind: High\nSpatial Resolution by Pinhole Neutron Diffraction. Appl. Phys. Lett. 2018, 112 (25).\n(76) Lee, D.; Yoon, A.; Jang, S.; Yoon, J.-G.; Chung, J.-S.; Kim, M.; Scott, J.; Noh, T.; Giant\nFlexoelectric Effect in Ferroelectric Epitaxial Thin Films. Phys. Rev. Lett. 2011, 107 (5),\n057602.\n(77) Chen, H.; Soh, A. K.; Ni, Y.; Phase Field Modeling of Flexoelectric Effects in Ferroelectric\nEpitaxial Thin Films. Acta Mechanica 2014, 225 (4), 1323-1333.\n(78) Zubko, P.; Catalan, G.; Tagantsev, A. K.; Flexoelectric Effect in Solids. Annu. Rev. Mater.\nRes 2013, 43 (1), 387-421.\n(79) Meyer, T. L.; Jiang, L.; Park, S.; Egami, T.; Lee, H. N.; Strain-Relaxation and Critical\nThickness of Epitaxial La1. 85sr0. 15cuo4 Films. APL Mater. 2015, 3 (12).\n(80) Matthews, J. W.; Blakeslee, A. E.; Defects in Epitaxial Multilayers: I. Misfit Dislocations.\nJ. Cryst. Growth 1974, 27, 118-125.\n(81) Maurice, J.-L.; Pailloux, F.; Barthélémy, A.; Durand, O.; Imhoff, D.; Lyonnet, R.; Rocher,\nA.; Contour, J.-P.; Strain Relaxation in the Epitaxy of La2/3sr1/3mno3 Grown by Pulsed-Laser\nDeposition on Srtio3 (001). Philos. Mag. 2003, 83 (28), 3201-3224.\n(82) Kawai, M.; Kan, D.; Isojima, S.; Kurata, H.; Isoda, S.; Kimura, S.; Sakata, O.; Shimakawa,\nY.; Deposition Rate Effect on Critical Thickness of Batio3 Epitaxial Thin Film Grown on Srtio3\n(001). MRS Online Proceedings Library (OPL) 2007, 1034, 1034-K10-04.\nNa)TiO -based ferroelectric thin films\nXun Geng, Ji Zhang, Jack Yang, Claudio Cazorla, Changqing Guo, Shery L. Y. Chang,\nHuizhong Wang, Xiaojie Lou, Houbing Huang, Zibin Chen, Shujun Zhang and Danyang Wang\nWe greatly appreciate the referees’ constructive comments and suggestions on\nour manuscript entitled “Flexoelectricity-driven giant polarization in (Bi, Na)TiO -\n3\nbased ferroelectric thin films” (NCOMMS-25-12919A). The clarifications and changes\nmade in accordance with the referees’ comments are provided below.\nResponse to the comments of Referee #1\nDear Authors, thanks for your very detailed response, and most of what you write is\nfine and acceptable.\nComments (1): On two comments however, I disagree. Your notion and understanding\nof PLD with respect to composition transfere is not correct. One can get a congruent\ntransfere, but as far as I see this issue, this is more the exception rather than the rule.\nYou have to optimize the deposition conditions in order to reach the envisaged\ncomposition like for any other deposition technique. There are a couple of simple\nguidelines which can get you to an almost congruent transfere, but this is not a given.\nThe physics and scattering properties inside these laser-induced plasma plumes are in\nfact rather complex. As for the references you quite, I know the authors. Non of them\nwould sign what you state. The correct statement in your manuscript would be that you\nassume an Fe-concentration in your films as set by the target. Though, it is a nominal\ndoping level, and this is as much as you can go. Still, it does not take anything away\nfrom your data.\nResponse: We sincerely appreciate the reviewer’s insightful comments. We agree that\nachieving congruent transfer in PLD generally requires careful optimization of\ndeposition conditions, and it is challenging to accurately quantify the Fe concentration\nin the films, due to the trace levels of dopants involved. Therefore, the Fe concentration\nin our thin films is more appropriately described as a nominal doping level defined by\nthe target composition. In line with this clarification, we have revised the Methods\nsection of the manuscript as follows:\n“Undoped BNBT, along with Mn, La, Sm, Eu, Co, and Fe-doped BNBT ceramic targets\nfor thin film deposition, were prepared by the conventional sintering process. The\nnominal doping level in our BNBT-based thin films is determined by the composition\nof the ceramic target used for Pulsed Laser Deposition (PLD).”\nComments (2): Also the statement on p6: It is known that PLD-produced oxide thin\nfilms are always oxygen deficient. The word \"always\" is a bold statement, and not really\ncorrect. You can have oxygen deficient films, yes. But this will depend very much on\nyour deposition conditions. Meaning, you can have anything from highly deficient to\nfully stoichiometric. If what you say is correct, how could you grow insulating oxide\nfilms? Say STO since you use STO as a substrate? Or Multiferroics, or other FE\nmaterials? You would have a very difficult time.\nResponse: We sincerely appreciate the reviewer’s insightful comments. We agree that\nthe concentration of oxygen vacancies in PLD-grown thin films can be effectively\ncontrolled by tuning the deposition parameters. Nevertheless, to achieve epitaxial\ngrowth and favourable plume dynamics, relatively low oxygen partial pressures (e.g.,\n10-3 to 10-2 Torr) are often used, resulting in the formation of oxygen vacancies in the\nas-grown thin films.\nWhile near-stoichiometric films can be obtained through rationally controlling the\ndeposition conditions, we would like to clarify that the direct quantification of oxygen\nvacancy concentrations remains extremely challenging, particularly when the vacancies\nare present in trace amounts. We also understand that a few published works1, 2\ndemonstrate a “stoichiometric growth window” for ferroelectric or complex oxide thin\nfilms deposited by PLD. However, this term is more appropriately interpreted as\nsuggesting a regime where oxygen vacancy concentrations are sufficiently low so as\nnot to significantly degrade functional properties, rather than implying the complete\nabsence of oxygen vacancies. For instance, the presence of well-defined ferroelectric\nbehavior does not necessarily confirm perfect stoichiometry or the absence of defects,\nas ferroelectricity is known to tolerate a certain level of oxygen vacancies.3\nIt is also not uncommon that oxygen-deficient perovskite ferroelectric thin films can\nremain insulating due to several factors, including the low concentration of oxygen\nvacancies4, localization of defects5 (defect trapping as in our case), and the formation\nof compensating defects such as defect dipoles or complex6-8, etc.\nAs mentioned in our initial response, our Fe-BNBT thin films exhibit oxygen deficiency,\nas evidenced by the XPS results (Supplementary Figure 13 in the revised manuscript)\nand STEM-iDPC images (Supplementary Figure 9). These localised oxygen vacancies\nfacilitate the formation of T-T domain walls and contribute to the development of the\nbuilt-in electric field. Nevertheless, our thin films are still considered insulating, given\nthe well-defined P-E loops as shown in Supplementary Figure 12.\nComments (3): I am also puzzled about the numbers you quote for FWHM in\nFig. R2. Are you certain that the number for undoped BNBT is correct? It looks wider\nto me. And a FWHM of 0.02° for STO is not good for a single crystal. This would\nindicate either a not so good crystal, or your alignment on the STO Bragg peak could\nhave been better.\nResponse: We sincerely appreciate the reviewer’s insightful comments. We\nrecalculated the FWHM of the rocking curve for the (002) peak of the undoped BNBT\nthin film, as shown in our initial response letter, and determined the correct value to be\n0.058°, rather than the previously stated 0.031°. This correction was made in Figure R1\n(a). We are grateful to the reviewer for bringing this to our attention.\nIn general, a FWHM of ~0.02° for the (002) rocking curve of a SrTiO (001) single\n3\ncrystal is considered indicative of good crystal quality, which is suitable for epitaxial\ngrowth.3 The STO substrates used in our study are commercially sourced from\nShinkosha Co., Ltd. (Japan). The relatively large FWHM value reported in our initial\nresponse letter arises primarily from the XRD instrumental settings, rather than from\nany deficiency in substrate quality. Specifically, the earlier measurement was\nperformed using a divergence slit of 0.6 mm. It is known that a larger slit often results\nin a broader peak due to increased beam divergence. To address this issue, we repeated\nthe rocking curve measurements for both STO substrates and our thin films using a\nnarrower divergence slit of 0.1 mm, while keeping the scan settings (step size, counting\ntime, etc) unchanged. As a result, the FWHM of the STO (002) rocking curve was\ndetermined to be 0.0078°, as shown in Figure R1(b), which closely matches the\nsupplier-provided data (FWHM = 0.0062o for the (200) rocking curve; see\nhttps://www.shinkosha.com/english/). Consistently, the FWHM values for the rocking\ncurves of (002) peaks of our thin films also decrease when measured with the narrower\ndivergence slit, as presented in Figure R1(b). Importantly, these updated FWHM values\ndo not compromise our previous conclusion that Fe-doping-induced defects have only\na minimal impact on the epitaxial quality of the thin films, as stated in our earlier\nresponse letter.\nFigure R1. XRD rocking curves of (002) reflections for BNBT-based thin films,\nincluding undoped and doped with 2% mol of Mn and Fe. The rocking curve of the\nSTO single crystal substrate is also shown for comparison purposes. (a) results provided\nin the previous response letter and (b) newly acquired measurements with a narrower\ndivergence slit.\nResponse to the comments of Referee #2\nThe overall analysis and interpretation presented in the revised versions of the main\ntext of the paper and its supporting information made the paper substantially more\ncomprehensive and clearer than the original manuscript. The newly provided data rule\nout the alternative interpretations suggested by the reviewers. The experimental work\nis of very high quality. Although the claim made in the manuscript title has not been\ndirectly demonstrated and is supported primarily by indirect evidence—through\nexclusion of alternatives and quantitative modeling—the paper may be accepted for\npublication in Nature Communications in its present form.\nResponse: We sincerely appreciate your recommendation to publish our manuscript in\nNature Communications.\nResponse to the comments of Referee #3\nI would like to recommend the publication of this paper in this journal.\nResponse: We sincerely appreciate your recommendation to publish our manuscript in\nNature Communications.\nResponse to the comments of Referee #4\nI co-reviewed this manuscript with one of the reviewers who provided the listed reports.\nThis is part of the Nature Communications initiative to facilitate training in peer review\nand to provide appropriate recognition for Early Career Researchers who co-review\nmanuscripts.\nResponse: We sincerely appreciate your recommendation to publish our manuscript in\nNature Communications.\nReference:\n(1) Liu, G.; Lei, Q.; Xi, X.; Stoichiometry of Srtio3 Films Grown by Pulsed Laser\nDeposition. Appl. Phys. Lett. 2012, 100 (20).\n(2) Golalikhani, M.; Lei, Q.; Wolak, M.; Davidson, B.; Xi, X.; Narrow Growth Window\nfor Stoichiometric, Layer-by-Layer Growth of Laalo3 Thin Films Using Pulsed Laser\nDeposition. J. Cryst. Growth 2016, 443, 50-53.\n(3) Wang, T.; Gong, F.; Ma, X.; Pan, S.; Wei, X.-K.; Kuo, C.; Yoshida, S.; Ku, Y.-C.;\nWang, S.; Yang, Z.; Large Enhancement of Ferroelectric Properties of Perovskite\nOxides Via Nitrogen Incorporation. Sci. Adv. 2025, 11 (2), eads8830.\n(4) You, L.; Chua, N. T.; Yao, K.; Chen, L.; Wang, J.; Influence of Oxygen Pressure\non the Ferroelectric Properties of Epitaxial Bifeo 3 Thin Films by Pulsed Laser\nDeposition. Physical Review B—Condensed Matter and Materials Physics 2009, 80 (2),\n024105.\n(5) Noguchi, Y.; Matsuo, H.; Kitanaka, Y.; Miyayama, M.; Ferroelectrics with a\nControlled Oxygen-Vacancy Distribution by Design. Sci. Rep. 2019, 9 (1), 4225.\n(6) Kim, J.; Saremi, S.; Acharya, M.; Velarde, G.; Parsonnet, E.; Donahue, P.; Qualls,\nA.; Garcia, D.; Martin, L. W.; Ultrahigh Capacitive Energy Density in Ion-Bombarded\nRelaxor Ferroelectric Films. Science 2020, 369 (6499), 81-84.\n(7) Saremi, S.; Kim, J.; Ghosh, A.; Meyers, D.; Martin, L. W.; Defect-Induced (Dis)\nOrder in Relaxor Ferroelectric Thin Films. Phys. Rev. Lett. 2019, 123 (20), 207602.\n(8) Damodaran, A. R.; Breckenfeld, E.; Chen, Z.; Lee, S.; Martin, L. W.; Enhancement\nof Ferroelectric Curie Temperature in Batio3 Films Via Strain‐Induced Defect Dipole\nAlignment. Adv. Mater. 2014, 26 (36), 6341-6347.\nRe: NCOMMS-25-12919B, Title: Flexoelectricity-driven giant polarization in (Bi, Na)TiO -\n3\nbased ferroelectric thin films\nAuthors: Yunlong Sun, Ranming Niu, Zizheng Song, Shiyu Tang, Huizhong Wang, Xun Geng,\nJi Zhang, Jack Yang, Claudio Cazorla, Changqing Guo, Shery L. Y. Chang, Xiaojie Lou,\nHoubing Huang, Zibin Chen, Shujun Zhang and Danyang Wang\nWe greatly appreciate the referee #1’s comments and suggestions on our manuscript\nentitled “Flexoelectricity-driven giant polarization in (Bi, Na)TiO -based ferroelectric thin\n3\nfilms” (NCOMMS-25-12919B). The clarifications and changes made in accordance with the\nreferee#1’s comments are provided below.\nResponse to the comments of Referee #1\nDear Authors, thanks for your revision, and I still have a couple of minor points based on the\nreply of yours to my comments:\nComments (1): You write now that \"The nominal doping level in our BNBT-based thin films is\ndetermined by the composition of the ceramic target used for Pulsed Laser Deposition (PLD)\"\nIf you write \"determined\" you imply that you measure the doping level. Fact remains you didn't.\nThan state something along the line that the stated doping levels correspond to the assumed\nvalues based on the sintering process. You need to change this.\nResponse: We sincerely appreciate the reviewer’s insightful comments. We have replaced the\nterm “determined” with “based on” to avoid implying a direct compositional measurement of\ndopants. Accordingly, we have revised the Methods section of the manuscript as follows:\n“The nominal doping level in our BNBT-based thin films is based on the composition of the\nceramic targets used for Pulsed Laser Deposition (PLD).”\nComments (2): In your second reply, I notice you seemed to have a somewhat strange notion\nabout epitaxial growth in general and PLD in particular (Nevertheless, to achieve epitaxial\ngrowth and favourable plume dynamics, relatively low oxygen partial pressures (e.g., 10-3 to\n10-2 Torr) are often used, resulting in the formation of oxygen vacancies in the as-grown thin\nfilms.) You imply for one epitaxial growth under vacuum conditions cannot be done if I\noveremphasize your wording? This is plain wrong. Or adjust the oxygen content. An oxygen\nbackground can help to adjust the correct oxygen content, but it is far from what you try to tell\nme. Most of the oxygen in an oxide film will originate from the target, and not from the\nbackground.\nResponse: We recognize that epitaxial growth of thin films can occur under both vacuum and\nvarious oxygen partial pressures. Our use of the term \"often\" was intended to acknowledge this\nand was not meant to exclude the vacuum conditions.\nComments (3): What do you mean with a few published work? There are plenty of publication\nbesides the two you quote who measure their composition, some more appropriate than others.\nYou don't even attempt to measure your composition proper and the extend you would need it\nand you critizice others that they do? For your information, normally the oxygen content is\nmeasured, not the complement and the best one can seriously do is with nuclear techniques go\nto 1-2%, this is already excellent. Measuring with XPS is in your case insufficient since you\nrefere from a surface related property to bulk properties. As you may know, surfaces don't\nnecesarely reflect bulk properties and there are plenty of complex oxides where the top couple\nnm are different in oxygen composition than the bulk of the film. This degradation already can\ntake place while the film is still in vacuum. If you want to measure vacancies, than use Raman\nand search if you have around 620 cm^-1 a signal. You better do it with a green laser otherwise\nyou fool yourself since other wavelengt are less sensitive to oxygen vacancies.\nResponse: We did not intend to criticize any prior work. Our intention was merely to note that\noxygen vacancies are commonly present in oxide thin films and their direct quantification is\nchallenging, particularly at low concentrations in films of several hundred nanometres in\nthickness.\nOur XPS results were used in a semi-quantitative manner to illustrate relative near-surface\ntrends in oxygen vacancy levels induced by Fe (acceptor) doping. We did not claim to provide\nbulk quantification.\nFollowing the reviewer’s suggestion, Raman spectra of BNBT-based thin films doped with Fe,\nMn, La, Sm, Eu, and Co, as well as undoped BNBT thin films, were measured with a 532-nm\n(green) laser on an inVia Qontor Raman microscope (Renishaw). However, due to the Raman\nactivity of the SrTiO substrate, its strong signal dominates the spectra, making it extremely\n3\ndifficult to resolve the vibrational features specific to the films. As shown in Figure R1,\nthe pronounced peaks near 620 cm⁻¹ clearly reflect the strong substrate influence, resulting in\nRaman spectra of thin films that closely resemble those of the SrTiO substrate. Consequently,\n3\nthe Raman spectra alone cannot serve as a definitive indicator of oxygen vacancies in our\nBNBT-based thin films.\nIn contrast, iDPC-STEM provides direct, bulk-sensitive imaging of oxygen columns and is a\nmore suitable and widely accepted technique for detecting oxygen vacancies in oxide thin\nfilms.1 Our use of iDPC-STEM imaging was therefore appropriate for identifying localized\noxygen vacancies within the bulk of the thin films.\nFigure R1. Room-temperature Raman spectra of STO substrate, BNBT-based thin films and\n2.0Fe-BNBT ceramic target.\nComments (4): I am glad that you remeasured your STO Bragg peaks and figure out where\nyour obviously not so great values come from. STO from Shinkosha is good and presently the\nbest on the marked in terms of properties. Still, I have seen and measured plenty of Shinkosha-\nSTO which showed values in the range originally given by you. Also STO made by Sinkosha in\nthe long distant past when they have been than the best source for STO. Still, STO is in general\nnot a good single cystal as compared to Si or Scandates. It does not mean, one cannot grow\ngood films on it.\nResponse: We thank the reviewer for sharing their experience with the Shinkosha STO\nsubstrates and for highlighting the challenges of STO substrates may present compared to other\nsingle crystals.\nReferences:\n(1) Zhang, W. D.; Song, Z. Z.; Tang, S. Q.; Wei, J. C.; Cheng, Y.; Li, B.; Chen, S. Y.; Chen, Z.\nB.; Jiang, A. Q. Ultrahigh dielectric permittivity in Hf0. 5Zr0. 5O2 thin-film capacitors.\nNat. Commun. 2025, 16 (1), 2679."}]

Decision letter