Suppression of vortex – antivortex structures by anti-trimer point defects in hexagonal manganites
Cite as: J. Appl. Phys.127, 194106 (2020);doi: 10.1063/5.0001649
View Online Export Citation CrossMark
Submitted: 20 January 2020 · Accepted: 4 May 2020 · Published Online: 20 May 2020
H. L. Lin,1K. L. Yang,1,a) P. Z. Chen,1 G. Z. Zhou,1C. F. Li,1S. H. Zheng,1 L. Lin,1Z. B. Yan,1X. P. Jiang,2 and J.-M. Liu1,3
AFFILIATIONS
1Laboratory of Solid State Microstructures and Department of Physics, Nanjing University, Nanjing 210093, China
2School of Materials Science, Jingdezhen Ceramic Institute, Jingdezhen 333403, China
3Institute for Advanced Materials, South China Normal University, Guangzhou 510036, China
a)Author to whom correspondence should be addressed:[email protected]
ABSTRACT
The topologically protected vortex–antivortex (V–AV) domain structure in ferroelectric hexagonal manganites has been highly concerned recently, but its stability against intrinsic defects remains to be understood, given the claim that a topological structure would be robust against defects and other perturbations. In fact, it is also known that the V–AV structure is sensitive to the sample quality, and such a well- developed structure is hardly observed in thin films and defective single crystals. In this work, we investigate the influence of anti-trimer point defects on the stability of the V–AV domain structure by employing the phase-field simulation based on the Landau–Devonshire phe- nomenological theory. It is revealed that the characteristic V–AV structure essentially relies on the anti-trimer point defects under consider- ation. These defects lower the trimerization transition temperature on one hand and produce pinning effect on the vortex cores/walls on the other hand. However, the V–AV structure does remain robust if the anti-trimer magnitude of these defects is relatively weak but will be eventually destroyed if the anti-trimer magnitude is strong.
Published under license by AIP Publishing.https://doi.org/10.1063/5.0001649
I. INTRODUCTION
In recent years, topological structures in real space are becoming concerned in condensed matters and materials sci- ences. This is partially due to the breakthroughs in quantum topological phenomena in the momentum space such as topolog- ical insulators and Weyl semimetals, although they may not have much connection except the topology terminology.1–4 For real space topology, a number of examples can be found and here we address those phenomena in ferroic systems.5–7Besides the vortex structure and Skyrmions observed in magnetic systems, emergent topologically protected and cloverleaf (vortex–antivortex, V–AV) ferroelectric domain structures observed in hexagonal manganites (h-RMnO3, R = Sc, Y, In, Dy to Lu) have also been receiving attention.8–13The impact of such V–AV structures on our under- standing of ferroelectric physics and their application potentials are under intense investigations in the last several years.
We may highlight the underlying physics for addressing our main motivation. Structurally, h-RMnO3consists of corner-sharing MnO5trigonal bipyramid layers and R3+layers in theabplane. At the high temperature (T), an h-RMnO3has a centrosymmetric and paraelectric (PE)P63/mmcsymmetry, which is lowered into the fer- roelectric (FE)P63cmsymmetry in the low-Tregion.14,15The struc- tural origin for such a transition is the trimerization due to the mismatch between MnO5layers and R3+ layers and thus coherent tilting of the MnO5trigonal bipyramids.16Therefore, the structural trimerization is the core ingredient for developing such topological domain structures and considered as the primary order parameter.
As the primary order parameter, the trimerization needs two param- etersQandfto characterize.Qis the trimerization amplitude and fis the trimerization angle.17It is also noted that the MnO5trigonal bipyramid tilting induces the coherent displacement of R3+ions with two third displacing upward and one third downward or vice versa.
The unequal population of upward/downward-displacing R3+ions leads to the ferroelectric polarization (P) along thecaxis.8So, the ferroelectric polarization appears as the secondary order parameter in h-RMnO3.
Given the structural trimerization scenario, one immediately sees the three equivalent structural antiphase states (denoted by α,β, andγ) and two equivalent FE states (labeled by + or−for the upward or downwardPalong the caxis), leading to the six energetically equivalent domains (α+,α−,β+,β−,γ+, andγ−). The six interlocked structural antiphase and FE domains emerge from one point and form a vortex or antivortex. Degenerately, there are only two kinds of cycling sequences of the six domains around a central point, which are called a vortex (α+,β−,γ+,α−, β+, andγ−) and an antivortex (α+,γ−,β+,α−,γ+, andβ−), respec- tively, and the neighboring vortex and antivortex always form one V–AV pair.8 The whole space is thus occupied by infinite numbers of such V–AV pairs, constituting a topologically pro- tected V–AV structure in real space.
Nevertheless, one major issue here is to check whether such a V–AV structure is topologically protected or not. This issue has been addressed from various aspects in the literature. On one hand, it was experimentally confirmed that such a V–AV structure is indeed topologically robust against various stimuli such as tempera- ture variation, stress, and electric field, although the domain mor- phology does change seriously.18–20On the other hand, it was also found that it is hard to observe such a V–AV structure in some single crystals or thin films that may not have sufficiently high microstructural quality, hinting that such a structure is sensitive to crystal defects as an intrinsic ingredient.21,22In fact, it was once a mystery that there has been no intensive work reported in subse- quence to the initial studies from Cheong’s group, and one of the reasons was ascribed to the high sensitivity of such a V–AV struc- ture to microstructural quality of crystals under investigation.
Subsequently, quite a few earlier works checked the role of chemical doping as a kind of chemical defects. It was revealed that oxygen interstitials, indium vacancies, and indium-oxygen vacancy pairs in h-InMnO3 favor to stabilize the ferroelectricP63cm phase, while oxygen vacancies favor the non-polarP-3c1 phase.23On the other hand, a substitution of Mn by Ga was found to suppress the trime- rization and thus theP63cmphase.23,24In such a sense, it would be of high interest to test the stability of the ferroelectric V–AV struc- ture against perturbations in local trimerization tendency, to be induced by, e.g., vacancies, interstitial atoms, proper chemical sub- stitution of Mn, and so on.
In this work, we intend to investigate this issue within the Landau–Devonshire phenomenological framework, which was pro- posed for the h-RMnO3systems by Artyukhinet al.25As the first step without adding the complexity, we do not intend to induce direct perturbations to trimerization angle f and ferroelectric polarizationP. Our strategy is to introduce local negative perturba- tion toQwhich is equivalent to a suppression of the local trimeri- zation amplitude. Such a perturbation is viewed as an anti-trimer point defect with a non-zero magnitude. Subsequently, we employ the phase-field simulation to track the domain structure evolution in lattices with spatially randomly distributed defects. The defect density and anti-trimer magnitude will be carefully discussed in terms of the V–AV domain structure evolution. This work thus
represents the first and useful approach to simulate the stability of topological V–AV domain structures against such point defects.
The remaining part of this paper is organized as follows. In Sec.II, the used theoretical model addressing the anti-trimer point defects and the procedure of phase field simulations are described.
The simulated domain structures in lattices with anti-trimer point defects will be presented in Sec.III and the underlying physics of the anti-trimer point defect is discussed. A brief conclusion will be given in Sec.IV.
II. MODEL AND PHASE FIELD SIMULATION A. Landau–Devonshire model
The Landau theory for h-RMnO3has been well-established in terms of the variations of trimerization and polarization. The total free energy was given by Artyukhinet al.,25
Ftot¼ððð
V
(fbulkþfgrad)dV, (1)
where fbulk is the bulk free energy density and expressed in the powers ofQ,f, andP,
fbulk¼a(T) 2 Q2þb
4Q4þc 6Q6þc0
6Q6cos 6f gQ3Pcos 3fþg0
2Q2P2þaP
2 P2, a(T)¼a0(TTC)/TC,
(2)
wherea(T),b,c,c0,g,g0, andaPare the bulk energy coefficients,a0
is a positive constant, andTCis the critical temperature. Termfgrad
is the gradient energy density and can be written as
fgrad¼1 2
X
l¼x,y,z
slQ @Q
@l
2
þQ2 @Q
@l
" 2#
þslP @P
@l
( 2)
þ1 2txP @2P
@x2þ@2P
@y2
2
, (3)
whereslQ,slP, andtxPare the elastic stiffness coefficients which were discussed in the literature.26
The values of order parameters Q, f, and P in equilibrium can be obtained by minimizing termfbulkwith respect toQ,f, and P. Along this line,Qis vanishing atT≥TCand generally increases with decreasing T. The trimerization angle f takes six possible values: 0, ±p/3, ±2π/3,π. The secondary order parameterPcan be expressed as a function ofQandf:P=gQ3cos 3f/(aP+g0Q2). The six possiblefvalues correspond to the six types of domains emerg- ing in h-RMnO3. It is the coexistence of the six types of domains that adds domain wall (gradient) energy to the total free energy.
B. Anti-trimer point defects
Now, we consider how to characterize the point defects to be considered here. Given a defect site in a lattice, we assume that this defect is of anti-trimer functionality, and the trimerization ampli- tude Qwill be modulated by the defect, while the trimerization phase angle is not affected directly for simplicity on one hand and for the fact that this angle prefers to take one of the discrete values nπ/3 (nis an integer∈[0, 6]) and a transition from one angle to another would need to overcome a large energy barrier. Generally, such a modulation can be either suppression or enhancement of the trimerization amplitude. Here, we only consider the case where the defect will suppress the trimerization, following the experimen- tal facts, and we call this defect as anti-trimer point defect.
To characterize such an effect, a straightforward approach is to impose a perturbing term to the bulk free energy. In details, one looks at the bulk Landau free energy density defined in Eq.(2)and assumes that the randomly distributed defect sites have different coefficients from the perfect sites. For an h-RMnO3 system, the high-symmetry paraelectricP63/mmcphase is stabilized asa(T) > 0, while the ferroelectricP63cmphase is stabilized asa(T) < 0. Here, without losing the generality, we only assume coefficienta(T) to be site-dependent asa(T) plays the most important role of all.
Similar to an earlier scheme in the literature,27we re-writea (T) as
a(ri,T)¼a(T)þAtc0(ri),c0¼ 1 for a defect site, 0 for a non-defect site,
(4)
wherec0(ri) defines the defect status for sitei, andAtmay be called the anti-trimer magnitude, which is used to characterize the strength for suppressing trimerization since the trimerization amplitude is associated with the value ofa(ri,T).28A positive At implies a suppression of trimerization by such defects.
Given an isolated defect site, one plots the bulk energy densityfbulkas a function of Qfor different anti-trimer magni- tudesAt, as shown inFig. 1. For each curve, positionQemarking the minimalfbulkdefines the trimerization amplitude in equilib- rium. It is seen thatQe gradually decreases with increasingAt
and vanishes when a≥0, i.e., this site is non-trimerized. For a reference, one may estimate the anti-trimer strength of a point defect in the unit of energy. As shown inFig. 1, the maximal energy variation for a fully trimerized unit cell (Q=Qe∼0.96 Å) and a non-trimerized unit cell (Q= 0) is∼0.6 eV. It implies that an anti-trimer point defect needs to lower the system energy by
∼0.6 eV. This energy scale seems to be accessible by chemical substitution and lattice distortion. For example, for typical ferro- electrics such as BiFeO3and BaTiO3, the ferroelectric double-well potential barrier is ∼0.5 eV in height.29,30 Therefore, the non- polar substitution of a polar lattice unit cell would be sufficient for fully suppressing the local trimerization in h-RMnO3.
C. Phase field simulations
We track the domain structure evolution in h-RMnO3employ- ing the phase-field method. Following the standard procedure in the earlier literature,31 the phase-field simulations are performed by
solving the time-dependent Ginzburg–Landau (TDGL) equations, δQx
δt ¼ LQ δF δQx
, δQy
δt ¼ LQ δF δQy
, δP
δt ¼ LPδF δP, (5)
where LQ andLPare the kinetic coefficients related to the domain wall mobility, and the Cartesian coordinates QxandQyare used to replace the polar coordinates Qandfto measure the trimerization withQx=Q⋅cosfandQy=Q⋅sinf.32
The TDGL equations are solved by the semi-implicit spectral method and the simulations are carried out on a two-dimensional (2D) square lattice with periodic boundary conditions applied along the two directions.33The system size is 512Δx× 512Δy× 1Δz withΔx=Δy=Δz= 3 Å. It should be mentioned that such a simu- lation would be sufficient for the present investigations since the inter-layer interaction is not the dominant factor for illustrating the consequence of anti-trimer point defects.26
We impose a set of point defects into the lattice and the defect density is denoted asC,and they are distributed randomly into the lattice. The simulations start at a temperature Tsta= 1400 K >TC= 1260 K. The initial values of Qx, Qy, and P are chosen randomly within the range of [−1.0, 1.0]. Then, the lattice is cooled down step by step with a step differenceΔT= 1.0 K until the ending temperature Tend= 300 K. At any arbitraryT, the FIG. 1. Trimerization amplitudeQdependence of the bulk free energy density fbulkfor different anti-trimer magnitudesAt.Qemarking the minimalfbulkdefines the trimerization amplitude in equilibrium.
lattice evolves by 10 000 steps and then the data are saved for output. This implies that the outputed domain structure at eachTis just a snapshotted state toward the equilibrium state, and thus the whole simulation sequence is a replica of the continuous cooling process as done experimentally.34 All the parameters used in our simulations are listed inTable I, taking h-YMnO3as a reference.25
III. RESULTS AND DISCUSSION A. Effect on transition temperature
The exact nature of phase transition from the high-symmetry P63/mmcphase to the polarP63cmphase in h-RMnO3 has been a subject of much debate. Now, it is widely accepted that the symmetry changes in a single step by undergoing the so-called trimerization- type structural transition as described above.11Here, it is interesting to know how the anti-trimer point defects affect the transition.
In order to check the effect, one plots the variations of the averaged order parameters, <Q> and <|P|>, over the whole lattice against T, as shown in Fig. 2. Here, the point at which <Q>
changes from zero to nonzero defines the transition temperature Ttr from the high-symmetry P63/mmc phase to the polar P63cm phase. It is seen thatTtrdrops asAtincreases, as shown inFig. 2 (a). As a summary, the dependences ofTtronCfor differentAtare investigated and shown inFig. 2(c). This dependence can be well fitted by the linear relationship. The biggerAtis, the larger the neg- ative slope is. Polarization <|P|> shows the similar behaviors as
<Q> does and the transition temperature forPcoincides with that
for Q, as shown in Fig. 2(b). The reduced <Q> and <|P|> with increasing At is also identified, in accordance with the nature of anti-trimer point defects.
It is worth to point out that earlier experiments on the Al- and Ga-doped h-YMnO3 at the Mn site did show the similar results.24 Likewise, the observed Ttr shows a negative linear dependence on the doping level in both cases. Furthermore, the negative slope in the Al-doping case is larger than that in the Ga-doping case, suggesting that the Al-doping allows stronger anti-trimer effect than the Ga-doping does. More discussion will be given in Sec.III C.
B. Pinning effect
Besides the effect on the transition temperature, the anti-trimer point defects also change the V–AV domain structure remarkably.
The V–AV domain patterns at differentT upon a series of defect densityC have been simulated and one set of the simulated results are presented in Fig. 3, where the snapshotted domain structures mapped by the P-contour at severalTare plotted. Certainly, several features deserve highlighting.
First, as shown in Figs. 3(a)–3(c) for the case of no point defect (C= 0), the well-known V–AV domain structure is devel- oped at T= 1100 K. Upon decreasingT, the domain growth pro- ceeds in accompanying with the annihilation of V–AV pairs, leading to a larger domain size. The local polarization inside the domains increases with decreasing T. It is noted that the domain growth process goes fast at the high-Tregion and slows down when TABLE I.Parameters used in the present calculations.
a0(eV Å−2) b(eV Å−4) c(eV Å−6) c0(eV Å−6) aP(eV Å−2) g(eV Å−4) g0(eV Å−4)
2.626 3.375 0.117 0.108 0.866 1.945 9.931
sxQ(eV) szQ(eV) sxP(eV) szP(eV) txP(eV Å2) At(eV Å−2) TC(K)
5.14 15.4 −8.88 52.7 73.56 1∼30 1260 K
FIG. 2.Temperature dependence of the averages ofQ(a) and |P| (b) for different anti-trimer magnitudesAt. (c) Transition temperatureTtras a function of defect density Cfor different anti-trimer magnitudesAt. Solid lines in (c) represent linear fits.
FIG. 3.Snapshotted domain structures mapped byP-contour during a continuous cooling process. (a)−(c) Cooling process in a defect-free lattice. (d)−(f ) Cooling process in a defective lattice with the anti-trimer magnitudeAt= 3 and defect densityC= 0.05. (g) and (h) TemperatureTdependence of the vortex densitynvfor different anti-trimer magnitudesAtand defect densitiesC. (i) The domain coarsening dynamics. Red and blue circles correspond to the results with and without the anti-trimer defects, respectively. Red line is fitted based on the formulaAlog(A)∼p+qtr.
the lattice goes down to the low-Trange. For a quantitative charac- terization of the domain growth process, the variation of vortex density againstTis drawn and shown using red squares inFig. 3 (g). It seems that the domain growth process is kinetically frozen gradually with decreasing T. However, this is not the case. The domain coarsening dynamics is usually characterized by plotting the average domain areaA∼1/nvas a function of simulation steps, i.e., time t. To study the coarsening dynamics under variational temperature, the variation ofAagainsttis counted in a lattice with larger system size 2048Δx× 2048Δy× 1Δzand drawn in Fig. 3(i).
The temperature linearly falls with increasingtduring the domain coarsening process. As the red line shows inFig. 3(i),Aas a func- tion of t is well fit by the formula AlogA∼p+qtr with r∼1.00, indicating that the power law with a logarithmic correction for the domain coarsening dynamics in 2D vortex domains still works well even against varying temperature.35,36 In other words, domain coarsening dynamics is independent with temperature.
Second, we come to look at the domain structure with anti- trimer point defects. As an example, we discuss the results for At= 3 and C= 0.05, as shown inFigs. 3(d)–3(f ). It is seen that the domain growth process seems to be kinetically terminated since the domain patterns at the three different T’s (1100 K, 700 K, and 300 K) are nearly the same, although the local polari- zation inside the domains still increases with decreasingT. Most of the domain walls and vortex cores seem to be fixed and do not move. Thus, there will be almost no change in the vortex density and domain size, either.
Apparently, it is the pinning effect from the anti-trimer point defects that should be responsible for the terminated domain growth. In order to quantitatively characterize the pinning effect, we also evaluate the variations ofnvagainstTfor differentAtand C, as plotted in Figs. 3(g) and 3(h). It is seen that the pinning effect becomes more serious whenAtis larger andCis higher. The evaluatednv(T) becomes nearly constant asAt= 3.0 andC= 0.05, indicating that the V–AV pair annihilation is seriously suppressed due to the pinning effect from the point defects and the domain coarsening is almost terminated. Let us check the domain coarsen- ing dynamics in the domain structure with the anti-trimer defects.
As the blue circles show inFig. 3(i),Aas a function oftdoes not follow the power law with a logarithmic correction any longer. The reason underlying is yet to study.
C. Instability of V–AV structure
In the above discussion on the pinning effect from the anti- trimer point defects, the chosen anti-trimer magnitudeAt and defect density C are relatively small. Indeed, the topological V–AV domain structure remains robust against such inhomoge- neities. Hereafter, it would be interesting to check the cases of largeAtand highC. Some of the simulated results are presented inFig. 4for defect densityC up toC= 0.2 and anti-trimer mag- nitudeAtup toAt= 12.
Several distinct features deserve for more detailed discussion here. First, the influence of defect densityCon the domain struc- ture is similar for all the cases, and the pinning effect delays seri- ously the domain coarsening. Obviously, the domain structure at C= 0.2 is several times smaller in the characteristic scale than that
at C= 0. Second, the anti-trimer magnitude At has remarkable influence on the domain structure. For a small At (At= 1), the V–AV pattern remains robust against the defects even ifCis up to 0.2, as shown in Fig. 4(d). All the vortex/antivortex cores follow the Z6-symmetry. For a large At (e.g.,At= 3 and 12), the V–AV pattern can still survive as long as Cis low enough, as shown in Figs. 4(e)–4(g) and 4(i). Nevertheless, a slightly high defect density, such as C= 0.1 and 0.2, most likely destroys the V–AV structures in some regions, as shown in Fig. 4(h)forAt= 3.0 and C= 0.2, and Figs. 4(k) and 4(l) for At= 12 and C> 0.05. A zoom-in imaging of the black square area inFig. 4(h), as amplified inFig. 4(m), clearly indicates the broken Z6-symmetry. The small circle area shows a vortex-core, while the small square area indi- cates a vortex core whose Z6-symmetry is broken. When At
becomes even larger, e.g.,At= 12, the V–AV structure is seriously damaged and there is nearly no well-identified Z6-type vortex or antivortex in Fig. 4(l) whereC= 0.2. While the V–AV structure can be well retained for C= 0.05 atAt= 3, as shown inFig. 4(f ), parameterAtas large asAt= 12 is sufficient to damage the V–AV structure forC= 0.05, as shown inFig. 4( j). A zoom-in imaging of the wine square area produces the pattern shown in Fig. 4(n), where a damaged vortex core is marked by a small square. The damaged V–AV structure can be seen more clearly inFigs. 4(k) and 4(l)forC= 0.1 and 0.2, where more vortices/antivortices are damaged. In short, the so-called topological V–AV domain struc- ture would lose its stability against strong perturbations such as strong anti-trimer point defects here.
Now, we present a brief discussion on the underlying physics for the V–AV structure instability against the anti-trimer point defects. The physics can be understood from the free energy sce- nario. In Section II B, the free energy related issue for an isolated point defect is discussed. However, for a practical domain structure, the interaction between a defect site and its neighbors, i.e., the gra- dient energy, is non-negligible and must be considered. As shown in Fig. 1, for example, the free energy density fbulk at At= 3.0 reaches the minimal ifQat the defect site equals toQe= 0. The gra- dient energyfgradbetween the defect site and its nearest-neighbors can be as large as 0.6 eV, implying that trimerization amplitudeQ at the defect site has to increase in order to reducefgrad. As a result, fbulkwill increase. Nevertheless, the gain offbulkcan be completely compensated by the reduction of fgradto ensure that the total free energy is the lowest.
For better illustration of the free energy variations around a defect site, it is useful to look at the order parameter profiles across the defect site, noting that this site is inside a domain and far away from any wall or core. The data are obtained by setting a sufficiently lowCvalue so that the inter-defect interaction can be excluded. For convenience, we set the defect site as the zero point and then extract parameters QandPprofiles from a horizontal line across the defect site. The profiles at At= 1, 3, 12, and 30 (from left to right) for several temperatures (T) are plotted in Fig. 5, where the Q-profiles are plotted in the top row and the P-profiles in the bottom row. It is seen that any defect withAt> 0 always favors the suppression of order parameters Q and P.
WhenAtis small (At= 1), one sees only the weak suppression of Qat the defect site. The suppression ofPas the secondary order parameter seems to be weaker. As At is increased up to 3, the
FIG. 4. The simulated domain structure mapped byf-contour for different anti-trimer magnitudesAtand defect densitiesC. (a)−(d)At= 1; (e)−(h)At= 3; (i)−(l)At= 12;
(m)−(p) zoom-in images from the color-box areas in (h), ( j), (k), and (l), respectively. The intact vortex is marked with a white circle, while the damaged vortex is marked with a white square.
suppression effect becomes stronger but the V–AV domain struc- ture can be retained to some extent. According to Eq. (6), the terma(ri, T) is positive at a defect site as At= 3, implying that the high-symmetryP63/mmcphase (Q= 0) should be stabilized.
However, for the practical domain structure,Qat the defect site may be nonzero due to the effect from the gradient energy term.
In this case, sufficient high defect density can suppress the trime- rization completely and destroy vortex, as shown inFig. 4(m).
ForAt= 12, the case of very strong anti-trimer defects, more and more vortices/antivortices become destabilized and eventually the whole V–AV domain structure will be melt-away if C is suffi- ciently high.
Here, it should be mentioned that the influence of the anti- trimer point defects is of short-range. As illustrated inFig. 5, apart from the defect site, theQandPprofiles at a position roughly 5 Å away from the defect site becomeAt-independent even ifAtis as high as 30. It implies that the spatial radius for the anti-trimer effect from point defect is∼5 Å. Considering that the lattice cons- tant for h-RMnO3is∼6 Å, one understands that a point defect has remarkable impact only on its nearest neighbors. It is necessary to mention that previous works on the interstitial and Ti doping done with DFT show the similar result. But in a different form, they automatically show the influence range by the relation between atomic displacement and distance from a defect.37,38This result is reasonable since the spatial variations ofQand Pare determined
by the gradient energy term (fgrad) which only takes account of the nearest-neighbor sites.
D. Discussion
It is understood that the structural trimerization (i.e., tilting of MnO5 polyhedral and subsequent buckling of R layers) in h-RMnO3 results from the lattice mismatch between the atomic layer consisting of small R ions and layer consisting of large MnO5 polyhedral units. That is to say that the R ionic radius is too small to form a hexagonal structure fitting to the MnO5layers. Therefore, the lattice has to distort to fill up the space around R ions.16This scenario for the trimerization is strongly supported by sufficient experimental evidences.34,39 The most important evidence is that the Curie pointTcin h-RMnO3increases with decreasing R ionic radius. The smaller the R ionic radius is, the stronger the mismatch is. Given the large mismatch, the high-symmetry P63/mmc phase can be stabilized only at highT, while this phase may not survive at the low-Tregion unless the mismatch is relatively weak or vanish- ing. For example, h-InMnO3 adopts the ferroelectric FE P63cm phase at room temperature, while h-InGaO3 adopts the high- symmetryP63/mmcphase at room temperature due to that the Ga3
+ radius (∼0.76 Å) is smaller than that of Mn3+ (∼0.785 Å).
Besides, that the Al and Ga doping on the Mn site in hexagonal YMnO3mentioned above lead to reduced transition temperature is
FIG. 5.The spatial profiles of order parametersQ(top row) andP(bottom row) across a single defect site which is inside a domain far away from any domain wall for different temperaturesTat the anti-trimer magnitudeAt= 1, 3, 12, and 30. The defect site is taken as the zero point.
also consistent with the mismatch scenario. Because Al3+(∼0.675 Å) radius is also smaller than that of Mn3+.
Regarding the effect of the point defects in this work, there is a similar relationship betweenTcfor the trimerization-type struc- tural phase transition and anti-trimer magnitude At. Here, it should be mentioned that Tc is the Curie point for an ideal h-RMnO3lattice. It is known that the condition a(T) = 0 ora(ri, T) = 0 determines the value ofTc. For a defect site, the transition pointTtrbecomes Tc(1−At/a0). Apparently,Ttrat the defect site declines with increasingAtuntilAt=a0, beyond which this defect unit adopts theP63/mmcphase. Regarding the dependence of Ttr on At, it is seen that At actually inversely scales the mismatch between the R layers and MnO5 layers, i.e., the larger At is, the weaker the mismatch is. In fact, Al3+is smaller than Ga3+, corre- sponding to a weaker mismatch. This is the reason why the Al-doping suppresses the transition more seriously.
Experimentally, the Ca substitution of Er and Ti substitution of Mn in h-ErMnO3 were performed to optimize the electronic properties inside domains and at domain walls without destroying the characteristic sixfold vortex domain patterns.38,40 However, these experiments were done at very low substitution levels which never exceed 0.01. Our calculations on the cases of largeAt(∼12) show a similar result, although there is a quantitative difference in substitution level. In fact, the Ca2+radius (∼1.14 Å) is larger than the Er3+ radius (∼1.03 Å), while the Ti4+ radius (∼0.745 Å) is smaller than the Mn3+radius (∼0.785 Å). The both cases fit well the mismatch scenario we have discussed here, i.e., the trimerization originates from the mismatch between small R layers and large MnO5layers.
IV. CONCLUSION
In summary, we used the phase-field method to evolve the domain structure in h-RMnO3in lattices with spatially randomly distributed anti-trimer point defects. Our calculations have revealed that the characteristic sixfold V–AV domain structure in h-RMnO3remains robust, when the anti-trimer magnitude is rel- atively weak. But these anti-trimer point defects can low the trimerization-type transition temperature and produce a pinning effect on the vortex cores/domain walls. On the other hand, some vortices begin to be destroyed at a very low defect density if the anti-trimer magnitude is strong. Also, as defect density increases further, the whole V–AV domain structure will be destroyed. Our calculations are very consistent with the results in a series of doping experiments on h-RMnO3. We argue that the anti-trimer magnitude physically inversely scales the mismatch between the R layers and MnO5 layers in h-RMnO3. The consistency has been well explained based on this mismatch scenario. This present work for the first time provides a comprehensive understanding of the response of the topological V–AV domain structure against defects and may inspire new ideas for tuning functional properties of V–AV domain structures by introducing different kinds of defects as needed.
ACKNOWLEDGMENTS
This work was financially supported by the National Key Research Program of China (Grant No. 2016YFA0300101) and
the National Natural Science Foundation of China (NNSFC) (Grant Nos. 11874031, 11834002, 11774106, 51431006, 51721001, and 11974167).
DATA AVAILABILITY
The data that support the findings of this study are available from the corresponding author upon reasonable request.
REFERENCES
1L. Fu, C. L. Kane, and E. J. Mele,“Topological insulators in three dimensions,”
Phys. Rev. Lett.98, 106803 (2007).
2Y. L. Chen, J. G. Analytis, J. H. Chu, Z. K. Liu, S. K. Mo, X. L. Qi, H. J. Zhang, D. H. Lu, X. Dai, Z. Fang, S. C. Zhang, I. R. Fisher, Z. Hussain, and Z. X. Shen,
“Experimental realization of a three-dimensional topological insulator, Bi2Te3,”
Science325, 178 (2009).
3Y. Du, F. Tang, D. Wang, L. Sheng, E.-J. Kan, C.-G. Duan, S. Y. Savrasov, and X. Wan,“Cate: A new topological node-line and dirac semimetal,”npj Quant.
Mater.2, 3 (2017).
4X. Wan, A. M. Turner, A. Vishwanath, and S. Y. Savrasov,“Topological semi- metal and Fermi-arc surface states in the electronic structure of pyrochlore iri- dates,”Phys. Rev. B83, 205101 (2011).
5X. Z. Yu, Y. Onose, N. Kanazawa, J. H. Park, J. H. Han, Y. Matsui, N. Nagaosa, and Y. Tokura,“Real-space observation of a two-dimensional skyrmion crystal,” Nature465, 901 (2010).
6Y. Li, Y. Jin, X. Lu, J.-C. Yang, Y.-H. Chu, F. Huang, J. Zhu, and S.-W. Cheong,
“Rewritable ferroelectric vortex pairs in BiFeO3,” npj Quant. Mater. 2, 43 (2017).
7Z. Li, Y. Wang, G. Tian, P. Li, L. Zhao, F. Zhang, J. Yao, H. Fan, X. Song, D. Chen, Z. Fan, M. Qin, M. Zeng, Z. Zhang, X. Lu, S. Hu, C. Lei, Q. Zhu, J. Li, X. Gao, and J.-M. Liu,“High-density array of ferroelectric nanodots with robust and reversibly switchable topological domain states,” Sci. Adv. 3, e1700919 (2017).
8T. Choi, Y. Horibe, H. T. Yi, Y. J. Choi, W. Wu, and S. W. Cheong,“Insulating interlocked ferroelectric and structural antiphase domain walls in multiferroic YMnO3,”Nat. Mater.9, 253 (2010).
9S. C. Chae, Y. Horibe, D. Y. Jeong, S. Rodan, N. Lee, and S. W. Cheong,“Self-organization, condensation, and annihilation of topological vortices and antivortices in a multiferroic,” Proc. Natl. Acad. Sci. U.S.A.107, 21366 (2010).
10S.-Z. Lin, X. Wang, Y. Kamiya, G.-W. Chern, F. Fan, D. Fan, B. Casas, Y. Liu, V. Kiryukhin, W. H. Zurek, C. D. Batista, and S.-W. Cheong, “Topological defects as relics of emergent continuous symmetry and Higgs condensation of disorder in ferroelectrics,”Nat. Phys.10, 970 (2014).
11M. Lilienblum, T. Lottermoser, S. Manz, S. M. Selbach, A. Cano, and M. Fiebig,“Ferroelectricity in the multiferroic hexagonal manganites,”Nat. Phys.
11, 1070 (2015).
12Q. N. Meier, M. Lilienblum, S. M. Griffin, K. Conder, E. Pomjakushina, Z. Yan, E. Bourret, D. Meier, F. Lichtenberg, E. K. H. Salje, N. A. Spaldin, M. Fiebig, and A. Cano,“Global formation of topological defects in the multifer- roic hexagonal manganites,”Phys. Rev. X7, 041014 (2017).
13S. M. Griffin, M. Lilienblum, K. T. Delaney, Y. Kumagai, M. Fiebig, and N. A. Spaldin,“Scaling behavior and beyond equilibrium in the hexagonal man- ganites,”Phys. Rev. X 2, 041022 (2012).
14C. J. Fennie and K. M. Rabe,“Ferroelectric transition in YMnO3from first principles,”Phys. Rev. B72, 100103 (2005).
15H. L. Yakel, Jr., W. C. Koehler, E. F. Bertaut, and E. F. Forrat,“On the crystal structure of the manganese(III) trioxides of the heavy lanthanides and yttrium,”
Acta Crystallogr.16, 957 (1963).
16T. Katsufuji, M. Masaki, A. Machida, M. Moritomo, K. Kato, E. Nishibori, M. Takata, M. Sakata, K. Ohoyama, K. Kitazawa, and H. Takagi, “Crystal
structure and magnetic properties of hexagonal RMnO3(R = Y, Lu, and Sc) and the effect of doping,”Phys. Rev. B66, 134434 (2002).
17M. E. Holtz, K. Shapovalov, J. A. Mundy, C. S. Chang, Z. Yan, E. Bourret, D. A. Muller, D. Meier, and A. Cano,“Topological defects in hexago- nal manganites: Inner structure and emergent electrostatics,” Nano Lett. 17, 5883 (2017).
18X. Wang, M. Mostovoy, M. G. Han, Y. Horibe, T. Aoki, Y. Zhu, and S. W. Cheong,“Unfolding of vortices into topological stripes in a multiferroic material,”Phys. Rev. Lett.112, 247601 (2014).
19M.-G. Han, Y. Zhu, L. Wu, T. Aoki, V. Volkov, X. Wang, S. C. Chae, Y. S. Oh, and S.-W. Cheong, “Ferroelectric switching dynamics of topological vortex domains in a hexagonal manganite,”Adv. Mater.25, 2415 (2013).
20K. L. Yang, Y. Zhang, S. H. Zheng, L. Lin, Z. B. Yan, J. M. Liu, and S. W. Cheong,“Electric field driven evolution of topological domain structure in hexagonal manganites,”Phys. Rev. B96, 144103 (2017).
21H. Pang, F. Zhang, M. Zeng, X. Gao, M. Qin, X. Lu, J. Gao, J. Dai, and Q. Li,
“Preparation of epitaxial hexagonal YMnO3thin films and observation of ferro- electric vortex domains,”npj Quant. Mater.1, 16015 (2016).
22J. Nordlander, M. Campanini, M. D. Rossell, R. Erni, Q. N. Meier, A. Cano, N. A. Spaldin, M. Fiebig, and M. Trassin,“The ultrathin limit of improper ferro- electricity,”Nat. Commun.10, 5591 (2019).
23S. M. Griffin, M. Reidulff, S. M. Selbach, and N. A. Spaldin,“Defect chemistry as a crystal structure design parameter: Intrinsic point defects and Ga substitu- tion in InMnO3,”Chem. Mater.29, 2425 (2017).
24H. Sim, H. Kim, K. Park, M. Lilienblum, M. Fiebig, and J.-G. Park,“Doping effects on the ferroelectric transition of multiferroic Y(Mn, Al/Ca)O3,” Phys.
Rev. B98, 085132 (2018).
25S. Artyukhin, K. T. Delaney, N. A. Spaldin, and M. Mostovoy,“Landau theory of topological defects in multiferroic hexagonal manganites,”Nat. Mater.13, 42 (2014).
26K. L. Yang, Y. Zhang, S. H. Zheng, L. Lin, Z. B. Yan, J. M. Liu, and S. W. Cheong,“Spatial anisotropy of topological domain structure in hexagonal manganites,”Phys. Rev. B95, 024114 (2017).
27J. M. Liu, X. Wang, H. L. W. Chan, and C. L. Choy,“Monte Carlo simulation of the dielectric susceptibility of Ginzburg-Landau mode relaxors,”Phys. Rev. B 69, 094114 (2004).
28C. X. Zhang, K. L. Yang, P. Jia, H. L. Lin, C. F. Li, L. Lin, Z. B. Yan, and J. M. Liu,“Effects of temperature and electric field on order parameters in ferro- electric hexagonal manganites,”J. Appl. Phys.123, 094102 (2018).
29Y. Sui, C. Xin, X. Zhang, Y. Wang, Y. Wang, X. Wang, Z. Liu, B. Li, and X. Liu, “Enhancement of multiferroic in BiFeO3 by Co doping,” J. Alloys Compd.645, 78 (2015).
30Y.-W. Fang, H.-C. Ding, W.-Y. Tong, W.-J. Zhu, X. Shen, S.-J. Gong, X.-G. Wan, and C.-G. Duan,“First-principles studies of multiferroic and magne- toelectric materials,”Sci. Bull.60, 156 (2015).
31L.-Q. Chen, “Phase-field models for microstructure evolution,” Annu. Rev.
Mater. Res.32, 113 (2002).
32P. Tolédano and J.-C. Tolédano,“Order-parameter symmetries for improper ferroelectric nonferroelastic transitions,”Phys. Rev. B14, 3097 (1976).
33L. Q. Chen and J. Shen, “Applications of semi-implicit Fourier-spectral method to phase field equations,”Comput. Phys. Commun.108, 147 (1998).
34S. C. Chae, N. Lee, Y. Horibe, M. Tanimura, S. Mori, B. Gao, S. Carr, and S. W. Cheong, “Direct observation of the proliferation of ferroelectric loop domains and vortex-antivortex pairs,”Phys. Rev. Lett.108, 167603 (2012).
35A. J. Bray,“Theory of phase-ordering kinetics,”Adv. Phys.43, 357 (1994).
36F. Xue, N. Wang, X. Wang, Y. Ji, S.-W. Cheong, and L.-Q. Chen,“Topological dynamics of vortex-line networks in hexagonal manganites,” Phys. Rev. B97, 020101 (2018).
37S. H. Skjærvø, E. T. Wefring, S. K. Nesdal, N. H. Gaukås, G. H. Olsen, J. Glaum, T. Tybell, and S. M. Selbach,“Interstitial oxygen as a source of p-type conductivity in hexagonal manganites,”Nat. Commun.7, 13745 (2016).
38T. S. Holstad, D. M. Evans, A. Ruff, D. R. Småbråten, J. Schaab, C. Tzschaschel, Z. Yan, E. Bourret, S. M. Selbach, S. Krohns, and D. Meier,
“Electronic bulk and domain wall properties in B-site doped hexagonal ErMnO3,”Phys. Rev. B97, 085143 (2018).
39H. Sim, J. Jeong, H. Kim, S. W. Cheong, and J.-G. Park,“Studies on the high- temperature ferroelectric transition of multiferroic hexagonal manganite RMnO3,”J. Phys. Condens. Matter30, 105601 (2018).
40E. Hassanpour, V. Wegmayr, J. Schaab, Z. Yan, E. Bourret, T. Lottermoser, M. Fiebig, and D. Meier,“Robustness of magnetic and electric domains against charge carrier doping in multiferroic hexagonal ErMnO3,” New J. Phys. 18, 043015 (2016).
