Domain wall mobility and roughening in doped ferroelectric hexagonal manganites

The macroscopic performance of ferroelectric and piezoelectric devices depends strongly on domain wall dynamics. It is clear that structural defects, such as vacancies, interstitials, and dopants codetermine the dynamics, but the microscopic understanding of the wall-defect interactions is still at an early stage. Hexagonal manganites are among of the most intensively studied systems with respect to static domain wall properties and thus are ideal model materials for studying domain wall mobility in the presence of defects. Here we study the mobility of domain walls in the hexagonal manganites and how it is affected by cation dopants using density functional theory calculations. The results are correlated with scanning probe microscopy measurements on single crystals, to conﬁrm an increasing domain wall roughness for the dopants we predict to pin the walls. The pinning originates from elastic strain ﬁelds around the walls interacting with the local crystal perturbations surrounding a dopant. The pinning strength is correlated with the local change in order parameter amplitude caused by the dopant. As a computationally friendly alternative to large supercell calculations, we demonstrate that domain wall pinning can be predicted from the dopants’ effect on the free-energy landscape of polarization switching. This approach allows to directly probe the effect of defects on domain wall mobility using a fraction of the computational cost, opening the door to detailed modeling and understanding of the critical pinning process of domain walls.


I. INTRODUCTION
Since the demonstration of enhanced conductivity at ferroelectric domain walls about a decade ago [1][2][3], much effort has been devoted to understanding the intrinsic domain wall properties [4] as well as developing their functionality for nanoscale device concepts [5,6].A challenging aspect of understanding domain walls is to characterize the local stoichiometry and the influence of point defects on the measured properties in order to establish compositionproperty relationships and disentangle extrinsic from intrinsic effects.Ferroelectric domain walls are structurally subnanometer wide regions that separate two symmetry equivalent domain states [7,8] that are known to accumulate defects, and these small variations in the stoichiometry of transition metal oxides are known to induce large changes in their properties [9].Density functional theory (DFT) calculations have been highly successful in rationalizing the impact of stoichiometric variations on the properties of bulk ferroelectrics.Explicit calculations on domain walls, however, require impedingly large supercells.Increasing computational power and the concerted application of high-resolution transmission electron microscopy investigations has led to rapidly growing understanding of the structure, composition, and properties of domain walls [1,[10][11][12][13][14][15][16].
Point defects, particularly oxygen vacancies, have been known for many decades to impede the macroscopic properties of ferroelectrics and piezoelectrics [17,18].This has been attributed to point defects pinning the domain walls, inhibiting their movement during polarization reversal, hence yielding pinched ferroelectric hysteresis loops and "hard" ferroelectrics [19].Despite the established empirical knowledge that point defects strongly affect the mobility of ferroelectric domain walls, the microscopic understanding of point defect segregation and pinning is very limited.
The multiferroic hexagonal manganites [20,21] are ideal model systems for studying interactions between point defects and domain walls, because the improper ferroelectricity is giving rise to stable charged head-to-head and tail-to-tail walls [4,22] and neutral side-by-side walls, that is, all types of 180 • domain walls.Such walls have been studied intensively, and their basic structural and electronic properties are well understood [10,16].Thus, it is now possible to go beyond as-grown domain walls and include interactions with point defects like oxygen vacancies [11], oxygen interstitials [6], and cation dopant atoms and to study how such point defects affect the domain wall dynamics.
The prototypical hexagonal manganite YMnO 3 is an improper ferroelectric [20] with a T C of 1250 K, where it transitions from paraelectric to ferroelectric, with the space group symmetry changing from P6 3 /mmc to P6 3 cm [23][24][25].The primary order parameter of the improper ferroelectric transition, the K 3 trimerization mode, corresponds to a tilt of the MnO 5 trigonal bipyramids, and an up-down-down corrugation of the Y layers, Fig. 1.Polarization emerges due to inharmonic coupling between the polar − 2 mode and the K 3 mode [26] and the Y corrugation is strongly coupled to the MnO 5 tilting through the strong, ionic Y-O bonds [25].Two polarization directions [27,28] and three possible trimerization centers [29] give rise to six ferroelectric domains labeled α ± , β ± , and γ ± , respectively, each uniquely described by a phase and an amplitude Q.
The improper ferroelectricity results in stable charged domain walls with enhanced or suppressed conductivity [4,5,30], coexisting with the more common neutral domain walls.The ferroelectric domains and domain walls are found to be robust for both donor and acceptor doping of the A site [31][32][33], as well as for donor doping of the B site [34].However, neither the binding energy nor the effect of dopants on the mobility of the domain walls is known.
Here, we use DFT calculations to study the binding energy and pinning of domain walls in YMnO 3 in the presence of different cation dopants (Zr 4+ , Ca 2+ , and Ti 4+ ), complemented by experimental data recorded by piezoelectric force microscopy (PFM) on doped and undoped ErMnO 3 single crystals.A microscopic origin for domain wall pinning and roughening in the hexagonal manganites is suggested by self-consistently solving the free-energy expansion.This approach is proposed as a computationally friendly route towards a fundamental understanding of domain pinning and roughening phenomena in ferroelectrics.

A. Experimental
The experimental PFM data were recorded on single crystals of hexagonal ErMnO 3 [34], Er 1-x Zr x O 3 (x = 0.010) [35], Er 1-x Ca x O 3 (x = 0.010) [32], and ErMn 1-x Ti x O 3 (x = 0.002) [34], grown using the pressurized floating zone method [36].All data were measured on as-grown, cut, oriented, and electromechanically polished, disk shaped samples, with the polarization oriented in-plane (thickness ≈500 μm).The PFM measurements were performed in ambient conditions on an NT-MDT NTEGRA Prima SPM in conjunction with Stanford Research 830R lock-in amplifiers.Commercially obtained MikroMasch NSC35:HQ cantilevers with hard diamondlike carbon coated tips were used, applying a probing signal with an rms amplitude of U rms = 5 V and a frequency of ω = 40 kHz.

A. Domain structures of undoped and doped hexagonal manganites
YMnO 3 is isostructural with ErMnO 3 , which is known to display qualitatively equivalent ferroelectric properties, domain structures, and domain walls, and, importantly, systematic microscopy studies on doped ErMnO 3 are available in literature [4,29,31,34].Representative PFM images of the domain structure in a-axis oriented ErMnO 3 single crystals (undoped and doped) with in-plane polarization are shown in Fig. 2. The domain walls of undoped ErMnO 3 in Fig. 2(a) are smooth and curved, giving rise to both neutral (side-by-side) as well as charged (head-to-head and tail-to-tail) walls, in line with previous reports [4,29].Domain wall roughening is observed with 1% Zr 4+ doping on the Er 3+ site [Fig.2(b)], both for charged and neutral walls, indicating that the walls may be pinned by Zr 4+ dopants.Interestingly, 1% Ca 2+ doping on the Er 3+ site [Fig.2(c)] or 0.2% Ti 4+ doping on the Mn 3+ site [Fig.2(d)] do not give rise to detectable domain wall roughening compared to undoped ErMnO 3 [Fig.2(a)].The differences in the domain sizes are determined by the cooling rate, and not by doping, as previously explained in the literature [46][47][48].The PFM images in Fig. 2 are scaled so that they contain comparable number of domain walls, FIG. 2. Piezoelectric force microscopy (PFM) in-plane contrast on (a) undoped ErMnO 3 [34], (b) Er 1-x Zr x MnO 3 (x = 0.010) [35], (c) Er 1-x Ca x MnO 3 (x = 0.010) [32], and (d) ErMn 1-x Ti x MnO 3 (x = 0.002) [34] (100)-oriented single crystals.The polarization directions are shown by white arrows in (a), where the green domains have a polarization direction to the left in the figure and blue domains have a polarization direction to the right.The PFM images are scaled so that they contain comparable number of domain walls.
since different domain wall densities, i.e., domain sizes, give different length scales for which the domain walls curve from vertex to vertex, and thus different apparent roughness.A comparison of the PFM images with the same absolute scaling is provided in Fig. S1 of the Supplemental Material [49].
From simple electrostatic considerations, the donor dopants Zr 4+ and Ti 4+ are expected to promote pinning of charged tail-to-tail walls, because they can contribute to screen the electrostatic field.In contrast, the head-to-head domain walls are expected to be pinned by the acceptor dopant Ca 2+ .It is then surprising that only Zr 4+ is observed to cause domain wall roughening occurring for all walls, on the nanometer length scale resolved in our PFM data.The unexpected differences in domain wall roughening between differently doped ErMnO 3 single crystals motivate our DFT calculations of binding energies between dopants and domain walls, as well as domain wall migration energy barriers in the absence and presence of dopants.

B. Dopants and local crystal structure distortions
To develop a microscopic understanding of the experimental observations, we start by investigating how the dopants locally perturb the crystal lattice using DFT calculations.The crystal structure of YMnO 3 is shown in Fig. 1.Since the bonds between yttrium (Y) and the planar oxygen are critical to understand the ferroelectricity [25] and the local crystal structure around domain walls and point defects [10,11,50], we emphasize the changes in the YO 7 polyhedra marked in Fig. 1(b).The resulting Y corrugation z Y , A-O4 bond lengths (A = Y 3+ , Zr 4+ , Ca 2+ ), average A-O bond lengths, and AO 7 polyhedral volumes in DFT relaxed 30 atom unit cells of YMnO 3 , Y 1-x Zr x MnO 3 , Y 1-x Ca x MnO 3 , and YMn 1-x Ti x O 3 (x = 1/6) are summarized in Table I.Compared to undoped YMnO 3 , the Y layer corrugation z Y is enhanced by Zr 4+ , subtly reduced by Ti 4+ , and strongly reduced by Ca 2+ .The AO 7 polyhedral volumes respond oppositely to these dopants, with a strong contraction and expansion found for Zr 4+ and Ca 2+ , respectively, following naturally from their Shannon radii [51]: r Y 3+ (VII) = 0.96 Å, r Zr 4+ (VII) = 0.78 Å, and r Ca 2+ (VII) = 1.06 Å. Ti 4+ dopants on the manganese (Mn) sublattice only subtly affect the AO 7 polyhedral volume, and Ti 4+ (0.51 Å) and Mn 3+ (0.58 Å) are more similar in size than Y 3+ and Zr 4+ /Ca 2+ .
The A-O4 bond lengths follow the same trend with dopants as the AO 7 polyhedral volumes.We discuss the importance of these bond lengths to domain wall mobility and roughening further below.A decrease (increase) in the A-O bond lengths correspond to an increase (decrease) in the order parameter Q and the accompanying polarization [26,28,52] (see also Sec.III E).To quantify the effect of doping on the polarization, we calculate the bulk polarization for the DFT relaxed 30 atom unit cells by a simple point charge model using the formal charges.We find bulk polarizations of 7.04, 9.89, 4.38, and 5.88 μC c m -2 for undoped, Zr-doped, Ca-doped, and Ti-doped YMnO 3 , respectively, as summarized in Table I.We note that while the dopant conctration in these calculations are artificially high, they qualitatively show the effect of the different dopants on the local polarization.Electronic densities of states (DOS) for the different dopants are given in the Supplemental Material (Fig. S2 [49]).

C. Domain wall-dopant binding energy
Next, we investigate the crystal structure across two neutral domain walls in the supercell illustrated in Fig. 3(a), starting with the same domain wall structure as reported in Ref. [16].From the Y1-O3 and Y2-O4 bond lengths between Y and planar oxygen atoms mapped across the supercell in Fig. 3(b), two important observations can be made: (i) these bond lengths are significantly different in the vicinity of the domain wall compared to the bulk, and (ii) the Y1 and Y2 walls are not symmetric, consistent with the YMnO 3 crystal structure with two symmetry-inequivalent Y positions in the space group P6 3 cm.
For the segregation of domain walls toward aliovalent point defects, there are two possible driving forces, that is, electrostatic fields between the walls and local strain fields in the vicinity of the walls.For neutral domain walls, the electrostatic potential profile across the 1 × 6 × 1 supercell is flat (Fig. S3(a) [49]) and only the local strain fields are expected to cause pinning of neutral walls by dopants.The variation in Y-O bond lengths across the domain walls, shown in Fig. 3(b), is a visualization of these local strain fields, which arise from discontinuities in the polarization and trigonal bipyramid tilt patterns across the walls.The Y2-O4 bonds are found to be shorter at the Y2-terminated wall (2.35 Å) compared to in bulk (2.41 Å), while oppositely the Y1-O3 bonds are found to be elongated at the Y1-terminated wall (2.32 Å) compared to in bulk (2.29 Å).In contrast, the Mn-O p and Mn-O a bond lengths (Fig. S3(c) [49]) show only small variations of ∼5 mÅ across the supercell.Thus, the local strain fields in the vicinity of neutral domain walls mainly affect the local structure of the YO 7 polyhedra.
The tendency for pinning of domain walls by defects, such as substitutional dopants, can be quantified by the binding energy or segregation enthalpy profiles for the different dopants.This segregation enthalpy profile of a defect can be expressed as ( Here, E (x) is the calculated total energy for placing a defect at a distance x from the domain wall, and E bulk is the total energy for placing the defect within the center of a domain of the supercell, which corresponds to the bulk material.A negative E seg (x) at a domain wall means that a defect costs less energy here, indicating possible pinning of a moving domain wall.The calculated segregation enthalpies for Zr 4+ , Ca 2+ , and Ti 4+ across 1 × 6 × 1 supercells with neutral domain wall pairs are given in Figs.3(c)-3(e).
Tendencies for segregation of the Y2-terminated wall towards Zr 4+ and segregation of the Y1-terminated wall towards Ca 2+ are observed, with calculated segregation enthalpies of −0.03 eV for both dopants.These values are comparable to those previously calculated for oxygen interstitials [6], but significantly smaller in magnitude than the positive segregation enthalpy previously found for oxygen vacancies [11].As noted above, the Zr-O bonds in Zr-doped YMnO 3 are found to be shorter (average 2.17 Å) than the Y-O bonds (average 2.30 Å) in undoped YMnO 3 , while the Ca-O bonds are found to be longer (average 2.37 Å).Hence, Y2-terminated walls with intrinsically shorter Y-O bonds are expected to show affinity for Zr 4+ dopants, while Y1-terminated walls with longer Y-O bonds are expected to bind to Ca 2+ dopants, consistent with our calculated segregation enthalpy profiles in Figs.3(c)-3(d).The rationalization is that it costs the system less energy to have a point defect and a planar defect like a domain wall together if these defects induce similar local structural distortions; the volume of perturbed bulk material is minimized by such segregation.No strong trends for domain wall segregation with Ti 4+ doping are observed [Fig.3(e)], as expected from the unperturbed Mn-O bond lengths across the walls (Fig. S3(c) [49]).Note that the total energy for Zr 4+ at the Y2-site in bulk is ∼0.1 eV higher than at the Y1-site, and opposite for Ca 2+ .This stems from shorter Y1-O4 bonds compared to longer Y2-O4 bonds, which implies that Zr 4+ prefers shorter bonds to O at the Y1 site while Ca 2+ prefers longer bonds at the Y2 site.
For the charged domain walls, simple electrostatic considerations imply that donor and acceptor dopants should pin tail-to-tail and head-to-head walls, respectively.This is also what we find from static mapping of the domain wall segregation enthalpy, as shown in Fig. S4 [49].Attempts to relax supercells with charged domain walls and aliovalent dopants were, however, unsuccessful as the walls migrate to the pinning dopant during relaxation (Supplement note 1 [49]).Hence, the charged domain walls will not be addressed further.
To summarize, we can rationalize the segregation energy, or binding energy, between dopants and domain walls from similarities between bond lengths at the walls and surrounding a dopant in bulk.When a domain wall is pinned by a dopant causing similar lattice distortions, the volume of unperturbed material increases and the total energy of the system decreases.

D. Domain wall mobility and dopants
Having established the static binding energies between neutral domain walls and aliovalent dopants, we next investigate the effect of dopants on the domain wall mobility by calculating the migration energy barriers.The calculated migration energy barrier and corresponding crystal structure of undoped YMnO 3 during migration of equidistant parallel domain walls are shown in Figs.4(a)-4(d).A migration energy barrier of 0.24 eV (27.61 mJ/m 2 ) is found, which is in agreement with Ref. [16].The domain walls broaden during migration towards the saddle point, where the domain wall adopts the P 3c1 space group symmetry [labeled 2/c model in Fig. 4(a)] as further elaborated in Ref. [16].The most pronounced change in the local crystal structure during the migration occurs in the Y layer, where Y2 [blue in Figs.4(a)-4(c)] displaces along the polar c axis, going through the highsymmetry position and regaining their initial Y2 symmetry after the migration.A subtle rotation of the MnO 5 trigonal bipyramids in the ab plane also occurs during migration, corresponding to a change in the order parameter phase of 30 • at the saddle point and a total change in of 60 • at the final point with negligible displacements of Mn.
Before we investigate the coupling between domain wall mobility and doping, we first need to define the two possible domain wall migration paths across the dopants, illustrated in Figs.4(e) and 4(f).In path (1), Fig. 4(e), the dopants are substituted in the Y/Mn columns marked with blue spheres.The substituted Y2 site, or the corresponding Y2 site above the Ti 4+ dopant, displaces along c during the migration, passing through the centrosymmetric position and regains the Y2 symmetry in the final structure.The sign of the relative displacement compared to the high-symmetry position, c, changes through this migration path, while the magnitude is preserved.In migration path (2) in Fig. 4(f), the substituted Y2 site, or the corresponding Y2 site above Ti 4+ , switches symmetry from Y2 to Y1 (or vice versa).The sign of c is, however, retained and only a subtle change in the magnitude of c is observed, corresponding to the inherent difference in the relative Y1 and Y2 displacements.For both migration paths ( 1) and ( 2) the Mn site, with and without Ti 4+ , shows no discernable displacements.The calculated domain wall migration energy barriers for the for the two distinct paths (1) and ( 2) when migrating in bulk (I → II) and across a dopant (II → III) are shown in Figs.4(g) and 4(h), respectively, discussed below.The calculated domain wall migration energy barriers are summarized in Table II.For migration path (1), presented in Fig. 4(g), all three doped systems show a bulk domain wall migration energy barrier of 0.24 eV, (I → II), equal to that of undoped YMnO 3 described above.The three doped systems show, however, significantly different behavior when the walls are migrating across the dopants, (II → III).The domain wall migration energy barrier is increased to 0.47 eV when the wall is migrating across Zr 4+ , which indicates significant domain wall pinning with Zr 4+ doping.Oppositely, the migration energy barrier is reduced to 0.14 eV for domain wall migration across Ca 2+ , suggesting enhanced domain wall mobility with Ca 2+ doping.Finally, the domain wall migration barrier shows only a subtle decrease to 0.21 eV when migrating across Ti 4+ and, hence, no changes in the domain wall mobility are expected for Ti 4+ doping.This means that for migration path (1) we expect Zr 4+ dopants to pin the domain walls, enhancing the domain wall roughness; Ca 2+ to promote domain wall movement, smoothening the domain walls; and no significant change for T 4+ domain wall roughness compared to undoped system.Next, we consider migration path (2), Fig. 4(g).In comparison to the results for path (1), described above, we see no apparent changes in the calculated migration energy barriers when migrating across the dopants, as compared to migrating in defect-free bulk.The relative energy differences for configurations II and III with Zr 4+ or Ca 2+ in Fig. 4(g) stems from the aforementioned site preferences for the two dopants as the corresponding lattice site changes symmetry from Y2 to Y1 when going from II to III.Note that we here only show the results for migrating the Y2-terminated wall across the dopants.Similar results are obtained for the Y1-terminated wall (Fig. S5 [49]).As such, when the domain walls follow this migration path, we do not expect to see any effect on roughness from Zr 4+ , Ca 2+ , or Ti 4+ dopants.
To summarize this part, the domain wall migration energy barrier only changes when a wall migrates across a dopant which induces local changes in Y-O bond lengths, meaning that the domain walls are expected to roughen for Zr 4+ doping if the domain walls move via path (1), consistent with the domain roughness seen in Fig. 2.

E. Point defects and the free-energy landscape of polarization switching
To provide a computationally friendly alternative to explicit domain wall migration calculations using large supercells, we now investigate the coupling between point defects and the free-energy landscape of the ferroelectric domain state, i.e., polarization switching.The free-energy expansion with respect to the phase , amplitude Q, and polarization P for the hexagonal manganites was first described in Ref. [28], where they arrived at the characteristic sombrero potential with six minima along the brim of the hat, corresponding to the six ferroelectric domain states.As a simple approximation, the sombrero potential can be extrapolated from the static energy landscape of the K 3 + 2 mode amplitude from paraelectric (space group P6 3 /mmc) to ferroelectric (P6 3 cm), and the domain switching energy barrier between each domain state with = 60 • .This approximated sombrero potential should then capture the microscopic nature of both the energy lowering of the domain state and the domain wall migration energy barrier.
The static energy landscape of the K 3 + 2 mode as a function of the relative mode amplitude in a stoichiometric 30 atom unit cell of YMnO 3 is shown in Fig. 5(a), where we assume a relative mode amplitude of 1.0 corresponding to a DFT relaxed ferroelectric unit cell.The ferroelectric domain state is 0.82 eV lower in energy than the paraelectric structure.The domain switching energy barrier between two unit cells with a change in of = 60 • is calculated to be 0.15 eV, Fig. 5(b).From the results in Figs.5(a) and 5(b), we extrapolate the sombrero potential surface in Fig. 5(c).
Here we assume that each domain is located at its corresponding (Q cos , Q sin ) value, with linear trajectories between each domain and to the top of the hat.This resulting DFT approximated sombrero potential is in good agreement with the reported free-energy Landscape from Landau modeling [28] and experimental observations [53].
Now that we have a self-consistent approach to solve the free-energy expansion, we expand the perturbations of the sombrero potential to include doping.A comparison of the K 3 + 2 mode amplitude energy landscape, the domain switching energetics, and the resulting approximated sombrero potential surfaces for the three doping schemes is given in Figs.5(d)-5(f).The Zr-doped unit cell shows a K 3 + 2 mode amplitude energy lowering equal to 2.14 eV, as compared to 0.82 eV in stoichiometric YMnO 3 , Fig. 5(a).5(a)] and a domain switching energy of 0.18 eV [Fig.5(b)], the latter is comparable to stoichiometric YMnO 3 .For Ti 4+ doping, the sombrero potential in Fig. 5(f) is also similar to that of the stoichiometric system, indicating negligible perturbations of the free-energy landscape and, hence, a similar robustness of the ferroelectric domain state.
Note that for the doped systems in Fig. 5(a), the dopants are added to a preconverged high-symmetric stoichiometric unit cell.The doped high-symmetric cells are not relaxed in order to avoid symmetry breaking; any energy contributions from changes in bond lengths in the high-symmetry phase upon doping are not included.Thus, for the doped systems, we expect the high-symmetric energies and the energy lowering of the ferroelectric states to be lower than what is reported in Fig. 5(a).
These results indicate that the ferroelectric domain state is sensitive to doping of the Y sublattice, and close to insensitive to doping of the Mn sublattice.This is due to the fact that the ferroelectricity in the hexagonal manganites is mainly governed by the rigidity of the YO 7 polyhedra, which is again governed by the strength of the Y-O bonds [25].The resulting Sombrero potentials in Figs.5(c)-5(f) can hence be reasoned from changes in the local chemical environments of the AO 7 polyhedra (A = Y , Zr, Ca) with doping illustrated in Fig. 2.During the domain switching, the A cation will displace along the polar c axis, where the ferroelectric displacement is accompanied by breaking and forming of A-O4 bonds.The robustness of the ferroelectric domain state will thus be governed by the A-O4 bond strength, as well as the rigidity of the AO 7 polyhedra.As Zr 4+ doping gives the shortest A-O4 bond, during domain switching the breaking of this strong bond requires more energy than for the other systems.In addition, with the shortest average Zr-O bonds and corresponding lowest polyhedral volume, the ZrO 7 polyhedra should in general be the most rigid, indicating a robust ferroelectric domain state.In contrast, Ca 2+ doping gives the longest Ca-O4 bond, which requires the least energy to break during domain switching.Ca 2+ doping also leads to the longest average Ca-O bonds and corresponding lowest polyhedral volume, which gives the least rigid AO 7 polyhedra.The latter indicates a weak ferroelectric domain state.Ti 4+ doping, on the other hand, does not give any significant changes in the Y-O bonding in the YO 7 polyhedra.The microscopic origin for any perturbation in the free-energy landscape with doping in the hexagonal manganites can thus be explained by changes in chemical bonding in the YO 7 polyhedra, which governs the robustness of the ferroelectric domain state.
To summarize this section, dopants which increase the amplitude of the Landau free-energy landscape for ferroelectric distortions are also predicted to pin domain walls, while dopants which reduce this amplitude are not.A dopant which locally increases (reduces) the ferroelectric distortions will (not) pin domain walls.Calculating the free-energy landscape is thus a computationally efficient alternative to explicit supercell calculations with domain walls.
Note that in our explicit domain wall migration calculations, the domain walls migrate as equidistant walls, while a nucleation and growth mechanism would be more realistic [54][55][56][57][58].This is not easily investigated explicitly by DFT calculations, but may become feasible with large-scale DFT approaches [59,60].However, we believe that even though nucleation and growth mechanisms are not captured in these calculations, the energetics of the domain wall propagation front should still be described well.Hence, any predictions from our simplified DFT approach should provide a further understanding of domain wall pinning phenomena observed experimentally.We have demonstrated the experimental relevance of our approach in a recent study on h-LuFeO 3 [61].
Here our DFT calculations predicted that oxygen vacancies will enhance the domain wall mobility, similar to the case with Ca 2+ doping described above, which was confirmed experimentally on h-LuFeO 3 thin films with different thermoatmospheric histories.
It should also be noted that while the negative domain wall segregation enthalpies for Zr 4+ and Ca 2+ (Fig. 3) may suggest that both these dopants should pin the domain walls, we observe no clear correlation between segregation enthalpies and migration energy barriers.In contrast, we observe a strong correlation between the migration energy barrier and the local order parameter amplitude Q.Hence, domain wall pinning and domain roughening phenomena in hexagonal manganites can be understood from domain wall migration energy barriers rather than segregation enthalpies.

IV. CONCLUSION
In conclusion, the domain wall mobility, inferred from DFT calculated migration energy barriers, is affected by the presence of cation dopants and can explain differences in domain wall roughening for different dopants.The microscopic origin for domain wall pinning in hexagonal manganites can be explained by dopant-induced changes in the Y-O bonding, which is crucial for stabilizing ferroelectricity in YMnO 3 .Zr 4+ doping gives shorter and thus stronger Zr-O bonds compared to the Y-O bonds in the stoichiometric material, acting as pinning centers for domain walls and promoting domain wall roughening.In contrast, Ca 2+ and Ti 4+ doping gives weaker or similar Ca-O or Y-O bonds as in pure YMnO 3 and will thus not act as pinning centers, and thus not induce domain wall roughening.
Understanding how dopants and other defects cause domain wall pinning and domain roughening is imperative to control the performance of ferroelectric and piezoelectric devices.Our results imply that the domain wall pinning of stationary cation defects can, in principle, be predicted from any structural changes that results in deeper energy wells in the free-energy landscape.We suggest that our first-principles approach is universally applicable to all ferroelectrics where the microscopic origins of energy lowering are known.

FIG. 1 .
FIG. 1. Crystal structure of YMnO 3 represented by (a) MnO 5 trigonal bipyramids, and (b) YO 7 polyhedra.The different Wyckoff positions of Y, Mn, and O in the hexagonal P6 3 cm space group are labeled accordingly.The Y corrugation z Y is defined as the difference in the z component of the two Y positions Y1 and Y2, as illustrated by the horizontal black lines in (b).The corresponding Y 2 O 7 polyhedra for which the bond lengths and polyhedral volumes reported in Table I are calculated is marked in (b).

FIG. 3 .
FIG. 3. (a) Structurally relaxed 1 × 6 × 1 stoichiometric YMnO 3 supercell with a neutral domain wall pair.(b) Calculated Y-O p bond lengths across the supercell, illustrating the strain fields in the vicinity of the walls.Mapped domain wall segregation enthalpy profiles for substitution of the Y1 site and Y2 site with Zr 4+ or Ca 2+ are shown in (c) and (d), respectively.The segregation enthalpy profile for substitution of Mn with Ti 4+ is shown in (d).

FIG. 4 .
FIG. 4. (a)-(c) The local crystal structure evolution and (d) the calculated migration energy barrier for migrating a neutral domain wall pair in a 1 × 6 × 1 YMnO 3 supercell.Here, for simplicity, only the local crystal structure evolution of the Y2-terminated wall is shown.The initial (a) and final (c) structures correspond to the sharp model, and the saddle point (b) corresponds to the 2/c model [16].The Y that goes through the high-symmetry position during the migration is marked in blue.Vertical lines in (a)-(c) correspond to the position of the walls during the migration.Illustration of domain wall migration paths (1) (e) and path (2) (f) and (g), (h) corresponding calculated migration energy barriers from bulk to the vicinity of a dopant (I → II) and across a dopant (II → III) in 1 × 6 × 1 supercells of YMnO 3 , Y 1-x Zr x MnO 3 , Y 1-x Ca x MnO 3 , and YMn 1-x Ti x O 3 (x = 1/36).

FIG. 5 .
FIG. 5. Calculated (a) static energy landscape of the K 3 + 2 mode amplitude and (b) energy barrier for domain switching with = 60 • , in 30 atom unit cells of YMnO 3 , Y 1-x Zr x MnO 3 , Y 1-x Ca x MnO 3 , and YMn 1-x Ti x O 3 with x = 1/6.(c)-(f) The corresponding extrapolated sombrero potentials from the static K 3 + 2 energy profiles in (a) and the domain switching energetics in (b), illustrating differences in the energy landscape with respect to doping.The arrows in the inset Sombrero potentials in (a) and (b) illustrate the corresponding energy paths for the K 3 + − 2 mode energy lowering and the domain switching, respectively.

TABLE II .
(2)culated domain wall migration energy barriers in bulk, and across the dopants for the two migration paths (1) [DW(1)] and path(2)[DW(2)].The average values in both direction of the asymmetric minimum energy paths for the domain wall migration across the dopants in path (2) in Fig.4(h) are reported.