A Predictive Theory for Domain Walls in Oxide Ferroelectrics Based on Interatomic Interactions and its Implications for Collective Material Properties

Domain walls separating regions of ferroelectric material with polarization oriented in different directions are crucial for applications of ferroelectrics. Rational design of ferroelectric materials requires the development of a theory describing how compositional and environmental changes affect domain walls. To model domain wall systems, a discrete microscopic Landau–Ginzburg–Devonshire (dmLGD) approach with A‐ and B‐site cation displacements serving as order parameters is developed. Application of dmLGD to the classic BaTiO3, KNbO3, and PbTiO3 ferroelectrics shows that A–B cation repulsion is the key interaction that couples the polarization in neighboring unit cells of the material. dmLGD decomposition of the total energy of the system into the contributions of the individual cations and their interactions enables the prediction of different properties for a wide range of ferroelectric perovskites based on the results obtained for BaTiO3, KNbO3, and PbTiO3 only. It is found that the information necessary to estimate the structure and energy of domain‐wall “defects” can be extracted from single‐domain 5‐atom first‐principles calculations, and that “defect‐like” domain walls offer a simple model system that sheds light on the relative stabilities of the ferroelectric, antiferroelectric, and paraelectric bulk phases. The dmLGD approach provides a general theoretical framework for understanding and designing ferroelectric perovskite oxides.


While
LGD is a powerful tool for understanding material behavior by compactly describing the system in terms of the order parameter (polarization P for ferroelectric materials), the compactness comes at the price of providing only a phenomenological description of the material that makes LGD fundamentally incapable of predicting changes in the properties due to the changes in composition. Furthermore, to consider the effect of spatial variation in polarization (i.e., the formation of domain structure and domain walls), gradient parameters g ij are used to characterize the energy cost of P variation. These parameters must be derived from scarce experimental data or from first-principles calculations on large cells and are usually assumed to weakly vary with changes in the system such as strain, temperature, and composition. To the best of our knowledge, there is no model of interatomic interactions available for estimating g ij and for evaluating whether or not they can be expected to be constant, for example, with temperature or strain variation. Considering the wide use of LGD theory for modeling of ferroelectrics, this lack of physical understanding of the origin of g ij is an important fundamental problem. In the study of domain-wall kinetics, LGD theory cannot account for the difference in the energies between the A-and B-site centered 180° tetragonal-phase domain walls that give rise to the Peierls potential [16] (see Supporting Information), since such a difference would require two different values of the g 44 parameter for the same material which is clearly unphysical. Additionally, LGD calculations cannot be used to estimate barriers for domain nucleation because an arbitrarily small step can be chosen to obtain arbitrarily small nucleation energy.
At the other extreme, DFT calculations can obtain the electronic wave function of the complex system and provide a precise description of the desired properties for essentially any composition and environmental conditions. The complexity of the wave function(s), however, leads to the difficulty in interpreting the DFT results and makes DFT calculations the in silico equivalent of experimental measurements that can accurately evaluate the necessary properties but do not, in and of, themselves provide insight or guidance for materials design. Furthermore, for complex multicomponent solid solutions, the computational cost of DFT evaluation of the many possible domain-wall configurations rapidly becomes prohibitive so that to the best of our knowledge, DFT studies of ferroelectric oxide domain walls have been limited to classic simple materials (e.g., BaTiO 3 , PbTiO 3 , LiNbO 3 , and BiFeO 3 [16][17][18][19][20][21][22][23][24][25][26][27] ) with the exception of one study of two ordered PbZr 0.5 Ti 0.5 O 3 phases. [28] The effective Hamiltonian method is a DFT-based extension of the LGD approach that can accurately treat the effects of temperature and other perturbations such as strain and electrical boundary conditions. [29][30][31] However, a separate first-principles-based derivation for the effective Hamiltonian parameters must be carried out for each composition, making it unsuitable for rapid exploration of compositionalphase space. Therefore, for a fundamental understanding of ferroelectric domain walls and to guide the exploration of the vast possible design space of ferroelectric materials, an analytical predictive theory that describes domain-wall structure and energies in terms of interatomic interactions and relates the properties of the bulk material to the properties of the DW-defect state is necessary.
Here, we use first-principles calculations for ferroelectric bulk unit cells and domain walls to reveal the limitations of the standard LGD approach. We then develop an extension of the LGD theory based on the use of atomic displacements as discrete microscopic order parameters and the interatomic interactions as the terms coupling the order parameters; we call this approach discrete microscopic LGD (dmLGD) theory. Based on the dmLGD approach, we identify A-B repulsion as the interaction that predominantly mediates the coupling of polarization in neighboring unit cells and predicts domain wall energies based on 5-atom, single-domain calculations. Using i) oxygen vacancy pinning of domain walls, ii) thermal fluctuations of the domain walls in the Peierls potential, iii) coercive fields in single-crystal ferroelectric solid solutions, iv) relative stabilities of ferroelectric and antiferroelectric phases and v) relative stabilities of ferroelectric and paraelectric phases as examples, we then demonstrate that the accurate decomposition of the energy of the system into contributions of individual interactions and lattice sites enables predictions of a wide range of properties in different ferroelectrics and the understanding of the effect of compositional variations on material properties based on the results obtained from the simple BaTiO 3 , KNbO 3, and PbTiO 3 systems only. Our results reveal that there exists a bidirectional relationship between the single-domain material (with uniform polarization) and its DW "defects" (Figure 1). On the one hand, while even the most simple models (i.e., 5-atom unit cell, single domain) contain the necessary information about the interactions to estimate the structure and energy for the DW "defects," the DW "defects" themselves contain information for predicting the relative stabilities of the ferroelectric, antiferroelectric, and paraelectric single-domain bulk phases. Thus, this work provides an analytical theoretical framework for understanding and designing ferroelectric perovskite oxides.

Results
To understand the atomistic origins of the ferroelectric domainwall energy σ θ (where θ is the angle between the polarization vectors on the two sides of the domain wall), we study the energies of tetragonal 180° domain walls (σ 180 ) for the classic perovskite ferroelectrics BaTiO 3 (BTO), KNbO 3 (KNO), and PbTiO 3 (PTO) in a variety of strain states (Figure 2). As observed from Figure 2a, for the gradient energy parameter g 44 that specifies the energy cost of P variation across tetragonal 180° domain wall, standard LGD theory based on P as the order parameter (see S1, Supporting Information) obtains variations of approximately 800%, 600%, and 40% for BTO, KNO, and PTO, respectively, for the strain range of 1.2%. Such a large variation in the gradient energy parameter for a small strain is implausible because a small strain leads to small structural changes and should only weakly affect the interatomic interactions in the material. Thus, the standard LGD approach does not provide a physically meaningful picture for the strain variation effect on domain walls even for these simple systems.
As mentioned above, LGD theory also obtains two different values of the g 44 parameter for the same material depending on whether the energy of the A-or B-centered domain wall is used for the calculation of g 44 , which is clearly unphysical www.advmat.de www.advancedsciencenews.com (see S1, Supporting Information). This indicates that it is necessary to go beyond the standard unit-cell-based LGD approach in order to understand the interactions in the material that control σ θ .
In the standard LGD approach that describes the bulk system, a macroscopic observable namely polarization is used as an order parameter to accurately describe the free energy surface and the macroscale phenomena and properties of ferroelectrics. Intuitively, to accurately describe microscopic, angstrom-and nanometer-scale potential surface that governs phenomena such as ferroelectric domain walls, the order parameter should be a microscopic variable. While polarization is a measurable bulk property, local polarization is not a rigorously defined observable. By contrast, atomic displacements are welldefined and measurable microscopic variables that are zero for the high-symmetry cubic structure and are non-zero in the ferroelectric state. Therefore, we choose cation displacements at the A-site (d A ) and at the B-site (d B ) as the order parameters for describing the microscopic potential energy surface of ferroelectric perovskites.
We follow the LGD approach of assuming that the potential energy surface of the ferroelectric perovskite can be expressed in a Taylor-series expansion in the powers of d A and d B around the high-symmetry structure. Here, a four th-order polynomial is used to describe the local energy contribution at the A-site and B-site to the total energy.
Based on previous work, [32] the coupling between the A-and B-site cation displacements is dominated by short-range cationcation repulsion between the A-and B-cations. This interaction decays exponentially with the increasing A-B distance, R A-B . Thus, the coupling energy term U g is given by where ν 0 and k are constants, the index i runs over the A-site nearest neighbor atoms of site B, and R i A-B are the distances between the atom on the B-site and its nearest A-site neighbors. In turn, R i A-B are functions of the off-center displacements of A-site and B-site atoms.
For accurate dmLGD formulation, it is important to distinguish between the different possible definitions of cation displacements. The displacements of the A-and B-site cations can be measured with reference to the center of mass of their O 12 and O 6 oxygen cages; these displacements change the AO and BO bond lengths and bond energies and generate polarization.  [40,43]. Left-hand image: Reproduced with permission. [40] Copyright 2016, Springer Nature. Right-hand image: Reproduced under the terms of the CC-BY Creative Commons Attribution 4.0 International license (https://creativecommons.org/licenses/ by/4.0). [43] Copyright 2021, The Authors, published by Springer Nature.

www.advmat.de www.advancedsciencenews.com
We denote these displacements by d A and d B , respectively. Since d B is defined with respect to the center of the O 6 cages, a shear of the O 6 octahedra even without any shift of the A-and B-cations from the high-symmetry positions will generate non-zero d B and local dipole moments while preserving the symmetric distances between the A-and B-cations ( Figure 2b). Alternatively, it is possible to define the A-and B-cation displacements with reference to the A-and B-site positions in the high-symmetry structure. We denote these displacements as D A and D B . These displacements may not change the AO and BO distances in the case where the oxygen cages displace together with the cations, but will change the A-B distances and the A-B repulsive interactions if D A and D B are different ( Figure 2c). Thus, d A and d B should be used to describe the energetics of the AO and BO interactions, and D A and D B should be used to describe the energetics of the A-B interactions. While the question of which definition of ionic displacement should be used may appear to be highly technical, as we will see it plays a key role in modeling domain energies because it is intimately related to the chemical bonding and interatomic interactions in the material.
For a perovskite with A-and B-cation displacements D A and D B from the high-symmetry A-and B-site positions, the coupling energy term due to A-B cation repulsion can be rewritten as by using the Taylor expansion of the exponential repulsive interaction between the A-and B-site cations (see S2-S4, Supporting Information). The g constant here will not have the same values as the g ij in the LGD theory. However, it can be shown that dmLGD becomes equivalent to LGD theory in the limit of slow spatial variation of polarization and ionic displacements (see Supporting Information). The total energy of the bulk ferroelectric system with respect to the high-symmetry structure is described by and for a system with spatially varying ionic displacements, the total energy is given by where the index i runs over the unit cells of the perovskite in the material and U tot is the total energy difference between the ferroelectric and paraelectric states. Since the driving force for ferroelectric distortion is controlled by changes in the bond lengths in the local A-O 12 and B-O 6 environments, [32] the local  www.advmat.de www.advancedsciencenews.com energy due to the distortion that is represented by the first two terms is controlled by d A and d B , and is parameterized by the material-and strain-state-specific constants a 1 , a 11 , and b 1 . For a bulk system, the parameters a 1 , a 11 , and b 1 can be obtained from fitting the energies of a series of DFT calculations for bulk 5-atom structures with different values of A-and B-cation displacements to Equation 4 as described in Tables S1-S6 (Supporting Information). The calculated a 1 , a 11 , b 1 , and g values of BTO, KNO, and PTO at selected lattice parameters are reported in Table 1 (see Supporting Information for a 1 , a 11 , b 1 , and g values for additional lattice parameter values, Tables S1-S6, Supporting Information).
For bulk systems, dmLGD is equivalent to the standard LGD theory and can be used to make the same predictions regarding the properties of the ferroelectric system at different temperatures, pressure, and applied electric field conditions. Going beyond the capabilities of the standard LGD approach, Equation (5) can be used to accurately the energies of systems with spatially varying polarization such as the domain wall energies σ θ for different θ. This is because the first two terms of Equation (5) give the local energy U loc while the last term represents the gradient-energy term U g because it disfavors variations in the ionic displacements in different unit cells of the material (see Supporting Information). In contrast to the standard LGD approach where a single-domain configuration does not contain any information regarding U g , the more detailed atomistic approach contains information about the energy cost of the polarization gradient U g even at the level of a 5-atom unit cell, single-domain state and therefore can be used to predict the energies of different domain walls. This is a crucial advantage of the discrete microscopic approach.
Furthermore, the discrete nature of the order parameter more accurately describes the potential energy surface that is due to the interactions between discrete atoms. This is illustrated in Figure S6 (Supporting Information) showing the domain wall profiles in the continuous P and discrete d B and d A representations for a system with the preference for off-center displacement at the B-site and the preference to remain at the high-symmetry position for the A-site (e.g., BaTiO 3 ). While the local energy cost of the continuous P representation is large due to the zeroing out of P at the A-site centered domain wall, the local energy cost in dmLGD is low due to the preference of the A-site to remain at the high-symmetry position. Since the true order parameters in ferroelectrics, namely ionic displacements, are in fact discrete and are only approximated in a coarse-grained model by the continuous P order parameter, the minimum grid size for LGD simulations of ferroelectric perovskites must either one half or one lattice constant, leading to a finite energy of nucleation processes.
Formally, the dmLGD theory expressed by Equation (4) describes the potential energy surface of a 1D ferroelectric (e.g., tetragonal BTO or PTO). In a manner similar to that for the standard LGD theory, the dmLGD expression given in Equation (4) can be extended to three dimensions and to include strain effects by using the x-, y-, and z-components of ionic displacements and the strain as the order parameters. Similarly, coupling to additional order parameters such as octahedral tilting and magnetization can also be included in a manner similar to that of the standard LGD theory. [29][30][31] The dmLGD approach can also be extended to describe the effects of compositional changes (A-, and B-site substitution and introduction of O vacancies) through the appropriate changes in the a 1 , a 11 , and b 1 parameters as revealed by DFT calculations. The effects of temperature and pressure can also be included by assuming linear dependence of the a 1 and b 1 parameters on temperature and pressure.
More importantly, as will be illustrated below, even without a formal extension of dmLGD to the description of compositional and environmental variations, the discrete and physically meaningful description of the potential energy surface of the ferroelectric by the dmLGD approach enables predictions of the effects of the changes in composition (e.g., from BaTiO 3 to PbTiO 3 or due to the introduction of O vacancies) on different ferroelectric properties and the relative stabilities of the different phases of the material (e.g., ferroelectric, antiferroelectric, and paraelectric) based on simple structural and energetic criteria.
We now apply the dmLGD approach to the modeling of ferroelectric domain walls. We first calculated 180° domain wall energy (σ 180 ) only allowing the cations to displace while freezing the oxygen ions at their high-symmetry positions (O-fixed) and analyzed the energetics of the system in terms of A-and B-site displacements using Equation (5). We fit a 1 , a 11 , b 1 , and g to the results of DFT calculations on 5-atom unit cells with different A-and B-site cations (see S3, Supporting Information) and using the fit values to predict the energies for the oxygen-fixed 180° DWs. Even though for all of our 5-atom calculations, the system has uniform polarization, we obtained good agreement with the DFT oxygen-fixed σ 180 values for the domain-wall configuration with non-uniform polarization (Figure 2d,e).
Consideration of the dipole-dipole interactions at the 180° domain wall shows that an attractive interaction will be present between the dipoles formed by the titanium (100) displacements on the two opposite sides of the domain wall. This interaction decreases σ 180 and likely accounts for the overestimation of σ 180 by our model. The presence of this dipole-dipole interaction also means that in the absence of A-site (B-site) displacement that couples to the B-site (A-site) displacements in the neighboring cells, the anti-parallel displacement pattern will be preferred. This is demonstrated by the plot of the DFT-calculated energy difference between the antiferroelectric and ferroelectric states with the A-site atoms fixed at highsymmetry positions that shows that the antiferroelectric state is preferred, with the preference increasing with greater c lattice parameter and B-cation displacement (Figure 2f and see S6, Supporting Information). Thus, our results show that it is the A-B repulsion that is responsible for the gradient-energy cost of polarization variation and the appearance of the ferroelectric phase while in the absence of the ability to minimize this repulsion by coupled A-and B-displacement, the usual preference of dipolar systems for antiparallel (antiferroelectric-like) arrangement will be restored.
Comparison of the oxygen-fixed σ 180 and fully relaxed σ 180 values show that oxygen relaxation is quite important for domain-wall energies, with much higher domain-wall energies (by a factor of 1.5-4) obtained for the oxygen-fixed case. Examination of the relaxed DFT domain-wall structure (Figure 2b) shows that the O 6 octahedra at the domain wall shear, substantially decreasing the displacement of the titanium cations at the domain wall relative to the A-site cation. We note that the local energy of the B-site displacement is controlled by d Ti that can be created by the O 6 -shear mode. By contrast, the gradient-energy term in Equation 5 is controlled by D Ti and is not affected by the O 6 -shear mode. Thus, the activation of the O 6 shear simultaneously keeps the local energy U loc low because d Ti is approximately equal to the single-domain d Ti and keeps the gradient U g low due to the small D Ti and small change in the A-B cation distance. The energy cost of O 6 shear is low and thus, avoidance of the trade-off between U g and U loc dramatically lowers the domain-wall energy relative to oxygen-fixed domain wall energies.
Assuming that the energy cost of the shear mode is zero and that D Ti is equal to one half of the single-domain value D A value (see S5, Supporting Information), we evaluated σ based on the a 1 , a 11 , b 1 , and g constants obtained from 5-atom calculations to obtain the lower-bound estimate of σ 180 . Good agreement (Figure 2d,e) is obtained between the model and DFT values with the underestimation accounted for by our neglect of the O 6 -shear-mode cost (which increases with higher strain).
Having elucidated the energetics of σ 180 in terms of A-and B-site displacements, we consider the tetragonal 90° and rhombohedral 71° domain walls (Figure 3 and Figure S7, Supporting Information) that control the switching and E c in bulk ferroelectrics. [33] Figure 3a,e shows the domain wall structure and Figure 3b,d shows the polarization profiles of these domain walls. In contrast to the 180° domain walls, these domain walls exhibit dipole-dipole interactions that are less favorable for the domain wall than for the single domain, giving rise to a contribution to σ 90 and σ 71 that is proportional to P 2 . Interestingly, we find that the A-B cation repulsion energy for the stable oxygencentered 90° domain wall is slightly lower than that in the single-domain structure (Figure 3c, and S7, Supporting Information). By contrast, the A-B centered 90° domain wall exhibits a very high A-B cation repulsion energy and is unstable. Thus, we conclude that σ 90 is essentially due to dipole-dipole interactions and scales approximately as P 2 . For the 71° rhombohedral-phase domain wall, the domain-wall energy σ 71 is also dominated by the dipole-dipole interactions proportional to P 2 and the A-B cation repulsion contribution is slightly positive (Figure 3f, and S7, Supporting Information).
We also consider the application of the dmLGD approach to 71° and 180° conducting domain walls (71°c and 180°c) formed by tail-to-tail and head-to-head arrangements of local dipoles. The polarization profiles for these domain walls are generated following the approach of Petralanda et al. [34] and are shown in   [24] The black, cyan, red spheres represent Pb, Ti, and O atoms. Local polarization direction is shown by the red arrows. www.advmat.de www.advancedsciencenews.com the domain wall energy σ will have local energy and repulsive energy contributions as well as a contribution due to the creation of the bound charge at the domain wall as expressed by tot l oc g c harge where ΔP is the change in the polarization across the domain wall and C 0 and A are constant and area perpendicular to the polarization variation. For different values of the polarization inside the domain, we evaluate the DFT domain wall energies for these domain walls as well as our model local energy and repulsive energy contributions using the a 1 , a 11 , b 1 , and g parameters obtained from 5-atom bulk DFT calculations, with the results shown in Figure 4b,c. We then take the difference between the DFT-calculated σ values and the U loc + U g contributions to the domain wall energy to obtain the U charge contribution to the domain wall energy (σ charge ) for different values of P inside the domain (Figure 4b,c). Plotting the charge σ ) versus P inside the domain (Figure 4d,e), we find that for both 71°c and 180°c walls and for the BTO, KNO and PTO materials, vs P. Thus, our dmLGD approach is also suitable for the modeling of conducting domain walls. [35] The physical understanding of the domain wall energetics in terms of discrete ionic displacements and the interatomic interactions enabled by dmLGD theory makes it a general and broadly useful tool for modeling of ferroelectrics domain walls that can be used to predict various domain-wall related properties of a wide range of ferroelectric perovskite materials based on the findings for the simple BTO, KNO, and PTO systems only. We illustrate this by using the physical understanding obtained by our dmLGD theory for BTO, KNO, and PTO to predict the effects of compositional variations on: i) oxygen vacancy migration energies, ii) Peierls potential and thermal fluctuations of domain walls, iii) coercive fields in ferroelectric solid-solution crystals, iv) relative stabilities of ferroelectric and antiferroelectric phases and v) ferroelectric-to-paraelectric transition temperatures.
Oxygen vacancies have been considered to be a primary cause of domain wall pinning that hinders switching in ferroelectrics and the preference for the vacancy to be located at the 180° domain wall rather than inside the domain was confirmed for PbTiO 3 by DFT calculations [36] but was not studied for other materials. Our dmLGD analysis of the potential energy surface for BTO, KNO, and PTO shows that the b 1 parameter is negative for PTO while it is positive for KNO and BTO, indicating that Ba and K prefer to remain in the high-symmetry position while Pb prefers to displace from the high-symmetry position. This means that the PbO bonds for the Pb at the 180° Pbcentered domain walls are weaker than the PbO bond of Pb atoms in the bulk of the domain, whereas the BaO and KO bonds at the domain wall are stronger than the BaO and KO bonds of the Ba and K atoms inside the domain. Since the introduction of the O vacancy into the material incurs the cost of AO bond breaking, the introduction of the vacancy at the 180° domain wall is favorable in PTO (because weaker bonds are being broken), while for BTO and KNO the vacancy formation energy should be higher at the domain wall. To confirm this prediction, we carried out DFT calculations for BaTiO 3 and PbTiO 3 supercells containing a 180° domain wall and an oxygen vacancy at different positions. As shown in Figure 5, in agreement with our expectations, a strong preference of the vacancy to be located at the 180° domain wall is observed for PbTiO 3 as indicated by the lower vacancy formation energy, while for BaTiO 3 the vacancy formation energy is lowest inside the domain. This suggests that domain wall pinning by vacancies will be either weak or entirely absent in BaTiO 3 and KNbO 3 .
Next, we consider the effect of compositional change in the Ba 1−x Sr x TiO 3 (BST) solid solution on the Peierls potential and thermal oscillations of 180° domain walls. Domain wall motion can be considered to proceed on a Peierls potential formed by the difference between the energies of the A-site and B-site centered domain walls. [37] In particular, domain wall motion can take the form of the random oscillations of the domain wall due to thermal fluctuations. The thermally driven hopping of the 90° domain walls in orthorhombic-phase BST thin films (aa 1 /aa 2 domain variant) with T c close to room temperature has been shown to give rise to resonance effects and spikes in the quality factor (Q) spectrum under an application of an in-plane external bias field. [8] However, such spikes are absent in the absence of the external bias field due to a too-high barrier for thermal domain wall oscillations. By contrast, a system with a low barrier for thermal domain wall oscillations, i.e., a flat Peierls potential can be expected to display thermal domain wall oscillations even in the absence of an external bias field. Here, we use the dmLGD approach to predict that increasing BST content will lead to a rapid flattening of the Peierls potential in tetragonal BST films.
For thin films of tetragonal ferroelectrics, the switching of the out-of-plane polarization is controlled by the movement of 180° domain walls. [33] In BTO, since the energy cost of the Ba-centered 180° domain wall (σ Ba ) is due to the increased Ba-Ti repulsion at the DW, it depends on the difference between the barium and titanium displacements D Ba -D Ti . On www.advmat.de www.advancedsciencenews.com the other hand, the energy of the Ti-centered 180° domain wall (σ Ti ) is due to the energy cost of zeroing out the titanium displacement and depends on d Ti . As seen from Figure 6a, D Ba -D Ti changes more weakly with the changes in the c parameter (and P) than d Ti , suggesting weaker dependence of σ Ba on P compared to σ Ti . In fact, for BTO, both our model and DFT calculations predict that σ Ba scales as P 1.9 while σ Ti scales as P 3 (Figure 6b and S11, Supporting Information).
This difference in the P scaling between σ Ba and σ Ti implies that the application of tensile strain decreases the c/a lattice parameter ratio and P will strongly decrease σ Ti but will more weakly decrease σ Ba . Thus, for a sufficiently small P (due to a small c lattice parameter, see S11, Supporting Information), the Peierls potential given by σ Ti -σ Ba will become flat or close to flat (Figure 6), enabling rapid thermal oscillations of domain walls. Since our analytical model and DFT results show that D A and D B are controlled by the c lattice parameter (see S11, Supporting Information), the flattening of the Peierls poten-tial and thermal domain-wall fluctuations should also occur when the c, P, and U loc of BTO are decreased by mixing with SrTiO 3 to obtain Ba 1−x Sr x TiO 3 (BST) solid solutions. Based on our previous results, [37] the presence of thermal hopping of 180° domain-walls should lead to the appearance of resonance peaks in the Q spectrum of tetragonal c + /c − BST films with higher SrTiO 3 content and a high density of 180° domain walls even in the absence of the external electric field.
To confirm our predictions, we grew x = 0.3 BST films on GdScO 3 (111) substrate because according to the BST phase diagram for this strain and composition the BST film will be found in the tetragonal phase with Tc close to room temperature (see S15-S18, Supporting Information) and examined their Q spectrum. Examination of reciprocal space maps showed that the films are tetragonal and piezoresponse force microscopy (PFM) experiments revealed polarization profiles consistent with the c + /c − domain configuration with a high density of 180° domain walls (Figure 6c). As shown in Figure 6d Figure S18, Supporting Information). The repulsive energy contribution to the σ 180 Ba is given by (for details see Equation S19, Supporting Therefore, slow variation of Δ implies that σ 180 Ba will go down slowly with the c parameter. b) σ 180 Ti and σ 180 Ba as function of c lattice parameter of BTO (DFT and model calculations are considered using LDA functional and Equation (1) ). Δσ = σ Ti-centered − σ Ba-centered is the difference between the Ti-center and Ba-center DW energies of BTO 180° DW that give the height of the Peierls potential. c) Piezoresponse force microscopy map of a written box-in-box pattern. d) Room-temperature microwave frequency-dependent quality (Q) factor as a function of frequency for c + /c − Ba 0.7 Sr 0.3 TiO 3 film (experimental). Resonances are explained by thermal fluctuations of 180° DWs. Q factor of a 100 nm-thick epitaxial Ba 0.7 Sr 0.3 TiO 3 /GdScO 3 (110) interdigitated electrode capacitor, exhibiting zero bias external electric field resonance spikes between 2 and 5 GHz well above the background. The insets show representative out-of-plane (OP) and in-plane (IP) piezoresponse force microscopy (PFM) images for the BST film confirming the presence of the c + /c − DW variant with directions denoted in blue, and the piezoresponse voltage-phase hysteresis further confirming local ferroelectricity in the film. The area of the PFM image is 500 × 500 nm 2 and the color-map corresponds to negative and positive system phase. Q spectrum measurements performed for these BST films without an application of a bias showed high resonance peaks between 2 and 5 GHz. The presence of these resonance peaks indicates the existence of strong thermally driven hopping of 180° domain walls, consistent with the flattening of the Peierls potential predicted by our dmLGD calculations.
For our third example of the application of physical understanding extracted from the dmLGD model, we consider the prediction of the compositional variation of the coercive field E c for bulk rhombohedral ferroelectric perovskite oxides. In this case, the switching of polarization proceeds through the motion of 71° domain walls and the energy required for ferroelectric switching is characterized by the coercive field E c that is a crucial material parameter for many applications. For these materials, the coercive field E c at the temperature T [E c (T)] scales with σ 71 (T) 2 /P(T) [σ 90 (T) 2 /P(T)] [33] and since σ 71 (σ 90 ) is largely due to dipole-dipole interactions and scales with P 2 , we obtain E c (T) ∼ P 3 (T). According to simple Landau theory P(T) = (P 0 2 /T c )(T -T c ) 1/2 , where P 0 is the 0 K polarization. Since our analysis of domain-wall energies showed that σ 71 is dominated by the P 2 contribution and T c is empirically known to scale as P 0 2 , [38][39][40][41][42] we obtain E c (300 K) ∼ (T c -300) 1.5 . To test this prediction, we collected the experimental results for room temperature E c and T c from the literature (References for experimental T c and E c data are given in Section S14, Supporting Information) for 14 different single-crystal ferroelectric perovskite oxide solid solution compositions (see Section S13, Supporting Information) and then plot the room temperature E c values versus the corresponding T c values (Figure 7a). It is found that the experimental E c and T c data show good agreement with the predicted E c (300 K) ∼ (T c -300) 1.5 relationship. We note that this relationship allows an easy evaluation of the intrinsic E c of a material that is difficult to obtain either by direct measurements (due to the presence of imperfections such as grain boundaries, vacancies, and impurities) or by P measurements due to the effects of the electrode-oxide interface on the measured P values. By www.advmat.de www.advancedsciencenews.com Adv. Mater. 2022, 34, 2106021 the vacancy-formation energy and the melting temperature in metals. In both cases, the microscopic features of the material at 0 K (σ 71 and the vacancy-formation energy) measure the energy cost of the configuration that is dominant in the highenergy disordered phase and are therefore correlated with the temperature of the order-disorder phase transition. The relationship between T c and 0 K σ 71 is further discussed in a different paper where we show that it can be used to generalize the previously obtained empirical T c -P 0 2 relationships. [44] The above insights into the relative stabilities of antiferroelectric, ferroelectric and paraelectric phases could not have been obtained by considering the ground state 5-atom structure with fixed polarization but are easily revealed through the consideration of the DW defect states with varying polarization. Thus, our dmLGD model shows that the 180°, 90°, and 71° domain walls provide simple model systems for the elucidation of the physics controlling the properties of bulk single-domain ferroelectric perovskite oxide materials.

Conclusion
We have studied ferroelectric domain walls using a combination of DFT calculations and analytical modeling using the discrete microscopic LGD (dmLGD) theory. We find that A-B cation repulsion is the interaction that predominantly mediates the coupling of polarization in neighboring unit cells and governs the energy of 180° domain walls. Thus, the energy contribution governed by the gradient term is in fact present even in a single-domain state, and single-domain, 5-atom unit cell calculations can be used to analytically predict domain-wall energies that require much larger (at least 40-atom) DFT calculations. We then demonstrate that the understanding of the different contributions to domain-wall energetics obtained from the simple BaTiO 3 , KNbO 3 , and PbTiO 3 systems allows the prediction of the potential energy surface for oxygen vacancy migration, and of the polarization switching and dielectric properties controlled by domain-wall movement for ferroelectric materials with a wide range of compositions. Furthermore, the successful use of the domain wall defect to characterize the interactions in the system that are important for the collective behavior of bulk, single-domain materials, along with the previous similar use of the vacancy energy to predict melting points in metals, suggests that the approach of studying the structure and energy of defect states at 0 K for understanding of collective properties at high temperature is general and applicable to wide variety of physical systems.

Supporting Information
Supporting Information is available from the Wiley Online Library or from the author.