Activation of nominally silent domain wall-localized phonons from GHz to THz frequencies

TL;DR

This paper reveals that ferroelectric domain walls host phonons from GHz to THz frequencies, explaining their broad GHz conductivity signature and enabling new nanophononic device designs.

Contribution

It uncovers the broad frequency range of DW phonons and explains the activation of silent modes due to wall tilting, advancing understanding of DW nanophononics.

Findings

DW phonons disperse from GHz to THz frequencies

Activation of silent DW modes is due to wall tilting and asymmetry

Simulations predict scanning impedance microscopy as a probe

Abstract

Ferroelectric domain walls (DWs) are nanoscale topological defects that can be easily tailored to create nanoscale devices. Their excitations, recently discovered to be responsible for DW GHz conductivity, hold promise for faster signal transmission and processing speed compared to the existing technology. Here we find that DW phonons disperse from GHz to THz frequencies, thus explaining the origin of the surprisingly broad GHz signature in DW conductivity. Puzzling activation of nominally silent DW sliding modes in BiFeO3 is traced back to DW tilting and resulting asymmetry in wall-localized phonons. The obtained phonon spectra and selection rules are used to simulate scanning impedance microscopy, emerging as a powerful probe in nanophononics. The results will guide experimental discovery of the predicted phonon branches and design of DW-based nanodevices.

Equations

$$m_i{\ddot{x}}_i+\gamma {\dot{x}}_i+m_i{\omega }_i^2{x}_i=\int dr\delta {\vec{P}}_i(\vec{r})\cdot \vec{E}(\vec{r})e^{i\omega t},$$

$$I(\vec{r},\omega )={\omega }^2\gamma \sum _i\frac{{\left[\int d{r}^{\prime }\vec{E}(\overrightarrow{r^{\prime }}-\vec{r})\cdot \delta {\vec{P}}_i(\overrightarrow{r^{\prime }})\right]}^2}{m_i^2({\omega }_i^2-{\omega }^2{)}^2+{\omega }^2{\gamma }^2}.$$

$$f=f_L+f_G+f_c+f_q+f_{fl},$$

$$f_L={\alpha }_i{P}_i^2+\frac{1}{2}{\alpha }_{ij}{P}_i^2{P}_j^2+{\alpha }_{ijk}{P}_i^2{P}_j^2{P}_k^2,$$

$$f_G=\frac{1}{2}{G}_{ijkl}{\partial }_j{P}_i{\partial }_l{P}_k,$$

$$f_c=\frac{1}{2}{ϵ}_{ij}{C}_{ijkl}{ϵ}_{kl},$$

$$f_q=-\frac{1}{2}{ϵ}_{ij}{q}_{ijkl}{P}_k{P}_l,$$

$$f_{fl}=\frac{1}{2}{f}_{ijkl}(ϵij{\partial }_l{P}_k-{\partial }_l{ϵ}_{ij}{P}_k)$$

$$\overrightarrow{d^{\prime }}=\vec{d}+{\vec{\omega }}^{\prime }\times \vec{d}+{\widehat{ϵ}}^{\prime }\vec{d},$$

$${\vec{d}}^{\prime }-{\vec{d}}^{\prime \prime }=({\vec{\omega }}^{\prime }-{\vec{\omega }}^{\prime \prime })\times \vec{d}+({\widehat{ϵ}}^{\prime }-{\widehat{ϵ}}^{\prime \prime })\vec{d}=0.$$

$$d_i({ϵ}_{ij}^{\prime }-{ϵ}_{ij}^{\prime \prime })d_j=0$$

$${ϵ}_{ij}=\frac{1}{2}{Q}_{ijkl}{P}_k{P}_l,$$

Full text

Ferroelectric domain walls (DWs) are nanoscale topological defects that can be easily tailored to create nanoscale devices Scott (2007); Tagantsev et al. (2010). Their excitations, recently discovered to be responsible for DW GHz conductivity, hold promise for faster signal transmission and processing speed compared to the existing technology Lee et al. (2003); Scott et al. (2012); Wu et al. (2017). Here we find that DW phonons disperse from GHz to THz frequencies, thus explaining the origin of the surprisingly broad GHz signature in DW conductivity Wu et al. (2017). Puzzling activation of nominally silent DW sliding modes in BiFeO3 Huang et al. is traced back to DW tilting and resulting asymmetry in wall-localized phonons. The obtained phonon spectra and selection rules are used to simulate scanning impedance microscopy, emerging as a powerful probe in nanophononics. The results will guide experimental discovery of the predicted phonon branches and design of DW-based nanodevices.

The seminal work of Seidel et al. Seidel et al. (2009) demonstrated that DC conductivity is higher at DWs in BiFeO3 (BFO) than in the bulk, enabling signal transmission along the walls. This inspired a new paradigm for the design of DW-based nanoelectronic devices Catalan et al. (2012); Matsubara et al. (2015); McQuaid et al. (2017); Sharma et al. (2017); Huang and Cheong (2017); Mundy et al. (2017); Turner et al. (2018); Akamatsu et al. (2018); Schaab et al. (2018), leading to demonstration of DW-based diode and transistor operating in kHz range Schaab et al. (2018); McQuaid et al. (2017). The push for higher frequencies, used in modern computers, has led to a discovery of giant increase of the effective DW conductivity at GHz frequencies in recent microwave microscopy experiments Tselev et al. (2016); Wu et al. (2017); Prosandeev et al. (2018). The conductivity results from excitation of soft DW-localized phonon modes, that are at the heart of switching, microwave dielectric loss and dielectric constant enhancement in ferroelectrics Wu et al. (2017); Hlinka et al. (2017). They correspond to oscillations of the DW plane and can be excited by an AC electric field that favors one of the domains during its half period (Fig. 1(a)). Here we show that the frequency of these phonons can go from GHz all the way up to THz, going beyond the frequency range of modern surface acoustic wave-based cell phone transducers.

The results explain why the frequency range of the SMIM anomaly is so wide and reveal DW-localized bands in the phonon spectra of ferroelectrics. We also address the peculiar selection rules for DW excitation. When the ferroelectric polarization component along the driving field is the same in the two domains across the wall, the DW vibration should not be excited (Fig. 1(b)), as seen at 180∘ walls on the $ac$ face of YMnO3 Wu et al. (2017). However, this rule is violated in recent experiments on rhombohedral BFO, where polarization components along the field across the wall do not differ (Fig. 1(c)), while the DW mode is still excited Huang et al. .

The results could inspire the exploration of DW-based nanosystems in THz frequency range, and lead to novel DW-based phononic nanodevices.

BFO is a rare room-temperature multiferroic with the largest spontaneous polarization among single-phase compounds known to date, and may be soon used in nano-devices Manipatruni et al. (2018). Three types of ferroelectric DWs in the rhomohedral phase are distinguished by the number of polarization components being reversed Marton et al. (2010). The most puzzling data is on R71∘ DWs, at which the $P_{2}$ component is reversed, and the polarization rotates by [math] across the wall, as seen in Fig. 1(c). R71 DWs are often seen in BFO films grown on a $[001]$ oriented substrate Ziegler et al. (2013); Domingo et al. (2017), such as SrTiO3 and DyScO3. The sliding mode of this wall should not be excited by the field normal to the (100) surface, since neither domain is favored by the field, as seen in Fig. 1(c). However, the impedance signal is still observed experimentally Huang et al. . The important distinction between [math] DWs in YMnO3 and R71 DWs in BFO is their orientation with respect to the surface. YMnO3 DWs are mostly normal to the surface, while BFO R71 DWs are tilted away from the surface normal. In that case the combination of the surface and DW orientations breaks the symmetry, as seen in Fig. 1(c). This symmetry breaking plays a key role in their response to electric fields, even when the simple selection rule, summarised in Fig. 1(a,b), suggests no response.

Results

We use the Ginzburg-Landau-Devonshire (GLD) model for DWs in BFO to capture the energetics of interacting ferroelectric polarization and strain (see Methods for details). The DW is tilted 45∘ away from the surface normal and bends near the surface to reduce the DW area, as shown in Fig. 2. Long-range strain textures emanate from the bent segments. $P_{3}$ decrease in the blue area near the surface is driven by the compression $\epsilon_{33}$ and the tensile strain $\epsilon_{11}$ via electrostriction.

The vibrational modes of the system with DWs are computed from the linearized equations of motion within a discretized Ginzburg-Landau model, and shown in Fig. 3(a,b) (see the Methods section for details). The DW breaks translational invariance and therefore the modes are not characterized by a well-defined quasimomentum. The traditional phonon band structure is then replaced by a spectral function, as seen in Fig. 3(a), where the intensity at frequency $\omega$ and wavevector $k$ indicates the content of plane waves with that wavevector in the eigenmodes at that frequency. For a translationally invariant system, sharp peaks would then appear in the spectral function ${A(\omega,\vec{k})=\sum_{\nu}\delta(\omega\pm\omega_{\nu})|\langle\delta P_{2,\nu}|\exp(i\vec{k}\cdot\vec{r})\rangle|^{2}}$. The long-ranged strain texture of the wall, shown in Fig. 2, is responsible for strong mixing between the DW sliding modes and acoustic phonons.

Fig. 3(a,b) shows the content of plane waves with $P_{2}\sim e^{i\vec{k}\vec{r}}$, while Fig. 3(d) shows the phonon dispersion along the wall. At low frequencies acoustic phonons are observed (their intensity is divided by 10 in Fig. 3(d) to make the DW branch visible). They correspond to strain modulations and mix with $P_{2}$ modes due to electrostriction, ${f_{q}=-\frac{1}{2}\epsilon_{ij}q_{ijkl}P_{k}P_{l}}$. The V-shaped low frequency branch, marked with (1-3) in Fig. 3(a), extends from around 10 GHz at $\Gamma$-point all the way to the bulk polar phonons, and corresponds to the DW sliding and wobbling modes, illustrated in Fig. 3(e-g). The phonons in this branch disperse along the wall (along $k_{\parallel}$), but are localized in the perpendicular direction, and therefore their Fourier components extend to $k_{\perp}\sim\pm\pi/\lambda$, where $\lambda$ stands for the DW width. The intensity of the DW-localized branch is lower than that of the bulk phonons due to the low volume fraction occupied by DWs. The higher frequency modes are the bulk polar phonons.

Fig. 3(e-i) show $\delta P_{2}$ profiles corresponding to some of the lowest frequency DW sliding modes(e,f,g) and breathing modes (h,i). The nodeless mode shown in the upper panel of Fig. 3(e) corresponds to the $P_{2}$ increase at the wall during half-period of the oscillation, therefore adding the DW area to the positive domain. During another half-period the negative domain grows. Therefore this mode can be thought of as DW sliding, depicted schematically in the lower panels of Fig. 3(e). A higher energy mode with a node in that branch, shown in the upper panel of Fig. 3(f), corresponds to DW wobbling. The DW breathing modes, schematically shown in Fig. 3(h,i), are found below the band of bulk polar modes, and are indicated with markers (4,5) in Fig. 3(a).

Now we move to the origin of the SMIM signal. In the experiments an AC electric field is applied between the tip and the back electrode, as shown in Fig. 1, and the current is measured. The phonons that have a non-zero energy in the oscillating field of the tip, ${-\int dr\vec{\delta P}(\vec{r})\cdot\vec{E}(\vec{r})e^{i\omega t}}$, are excited and give rise to a displacement current component in phase with the field. The corresponding loss, ${\mathrm{Re}\int dr\vec{j}(r)\cdot\vec{E}(r)}$ is measured in addition to the Ohmic losses due to itinerant electrons. The amplitude $x_{i}$ of a phonon mode $i$ with an eigenfrequency $\omega_{i}$ is governed by the equation of motion

[TABLE]

where the mode is characterized by its effective mass $m_{i}$, damping $\gamma_{i}$ and the spatial polarization profile $\delta P_{i}(r)$. The oscillating driving force on the right hand side is due to the electric field $\vec{E}(\vec{r})e^{i\omega t}$ of the SMIM tip. Looking for the solution in the form $x_{i}=x_{i0}e^{i\omega t}$, we obtain the loss power at the tip position $\vec{r}$,

[TABLE]

The frequency dependence is characterized by a Lorentzian. The integral in the numerator is negligible for phonons whose spatial oscillation period is much smaller that the length scale of the electric field inhomogeneity, roughly determined by the tip radius $R$, therefore the phonons with wavevectors $k\lessapprox R^{-1}$ are excited. Fig. 3(c) shows the simulated signal for a BFO sample containing a R71 DW. An asymmetric peak at the DW and a weak signature due to the long-range strain features are visible.

Discussion

Figure 3(c) shows that SMIM conductivity is much higher at DW than in the bulk due to low-frequency DW sliding and wobbling modes. Their dispersion spans the whole GHz range and extends to the bulk polar phonon band, usually positioned at THz frequency. This explains why the GHz microwave conductivity at DWs in hexagonal manganites does not show a narrow peak in frequency domain, but rather rises monotonically towards higher frequencies. The proximity of DW phonons to acoustic branches also leads to phonon scattering and affects thermal conductivity Royo et al. (2017). The low frequency of sliding modes is also behind the giant enhancement of the static dielectric permittivity Hlinka et al. (2017). In addition, soft DW phonons are responsible for mechanical softening phonomena, recently observed in several ferroelectrics Stefani et al. .

In materials with strong strain-polarization interactions, DWs are aligned into regular lattices by strain fields Artyukhin et al. (2013); Wang et al. (2014), and spectral features similar to the ones in Fig. 3(a,b) are expected. Alternatively, when DWs are randomly oriented, the spectral weigth of DW wobbling and breathing modes spreads along different directions, but the phonon density of states is expected to retain its shape, with the peak at GHz frequencies. The directionality of DWs and their populations in the sample can then be inferred from the orientation and intensity of low-energy branches, as seen in Fig. S2.

The recent observation of nominally silent DW modes in BFO Huang et al. being activated is a surprising evidence that DW type and orientation control local AC conductivity. Possible DW orientations have been classified using symmetry and compatibility analyses Janovec (1976); Fousek and Janovec (1969). In essence, the distances between ionic planes in the two domains must match along the wall. This is only possible for particular DW orientations, for which the strain components in the DW plane match at the wall (see Eq. S3). The energy penalty for violating this condition scales linearly with the wall area, therefore fixing the DW orientation in the bulk. However, in thin films surface strain relaxation allows the DW to deviate from its optimal bulk orientation.

In an infinite bulk sample with a DW, a diagonal mirror plane $M_{dw}(-1,0,1)$, shown by the dashed line in Fig. 1(c), is a symmetry operation. This symmetry imposes the requirements: $\epsilon_{11}=\epsilon_{33},\epsilon_{12}=\epsilon_{23}$. The surface of the thin film and the domain wall orientation together break this mirror symmetry $m_{dw}$ if they are not orthogonal. The stress at the tilted domain wall $\sigma_{ij}=\partial f/\partial\epsilon_{ij}$ then acts on geometrically different regions from the left and right sides of the wall. The resulting asymmetric strain induces asymmetric polarization across the wall via electrostriction Eq. 7. This phenomenon is observed in Fig. 4, where we simulated thin films with different surface orientations with respect to the crystallographic directions (different ways of cutting a film out of a slab). It is seen that the DW orients perpendicular to the $[P_{1},0,P_{3}]$ direction, and a small bending near the surface is observed in Fig. 4(a) Ishibashi and Salje (2002); Ishibashi et al. (2005); Conti and Weikard (2004); Li et al. (2018). In very thin films, the wall deviates from the $[101]$ orientation towards the surface normal. Fig. 4(b) shows the polarization component normal to the surface. When the wall is tilted at 45∘, the asymmetry of the polarization $P_{3}$ is evident in Fig. 2(a), especially near the surface where the polarization deviates the most from the bulk value. This asymmetry comes from the asymmetric strain, plotted in Fig. 2(b,c). Even though the DW polarization profile is very narrow, the accompanying strain texture, emanating from the areas of DW bending near the surface Ishibashi and Salje (2002); Ishibashi et al. (2005); Conti and Weikard (2004); Li et al. (2018), extends surprisingly far Salje (1993); Salje et al. (2016), as seen in Fig. 2(b,c).

The polarization asymmetry eventually results in the asymmetry in the phonon polarization profile of the DW sliding mode, as seen in the lower panel of Fig. 3(c). When the polarization $P_{3}$ interacts with the external electric field, normal to the surface, the electrostatic energy $-\int dr\vec{E(r)}\cdot\vec{P}(r-r_{DW})$ results in a force on the wall due to $\vec{E}\cdot\frac{\partial\vec{P}(r-r_{DW})}{\partial x_{1,DW}}\neq 0$ and excitation of the DW vibration by the SMIM field.

In addition to SMIM, it may be possible to observe the DW-localized modes with inelastic neutron or X-ray scattering, although their intensity may be low due to low volume fraction occupied by DWs. Rapid quenching of a sample across the ferroelectric transition, leading to high DW densities Griffin et al. (2012), or strain-induced generation and alignment of ferroelastic DWs Chae et al. (2012) should increase the intensity of these branches and enable their observation in inelastic neutron scattering experiments and other bulk spectroscopies.

In summary, the simulated phonon spectra reveal the DW-localized phonon branches starting from GHz and extending all the way to THz frequencies. This explains wide frequency range of the conductivity anomaly observed in hexagonal manganites Wu et al. (2017) and other materials. The surprising activation of the nominally silent DW mode at R71∘ DWs in BFO is interpreted in terms of phonon polarization asymmetry due to the interplay of electrostriction and elastic compatibility at a tilted DW. The proposed way to simulate SMIM experiments may be used in emerging second-principles methodologies. Since the low-energy theory used here is rather general, similar phonon spectra must be expected for all ferroic materials, where the order parameter couples strongly to the lattice. We hope this work will motivate the experimental search for DW-localized phonon branches, guide the interpretation of SMIM studies and eventually inspire the development of high-speed DW-based phononic devices, extending surface acoustic wave-based technology to THz DW-based circuits.

References

I Supporting Information

The two polar phases which define neighboring domains have to be mechanically connected by a DW plane described by the normal vector $\vec{n}$. The spontaneous polarizations and deformations in the two domains are denoted as $\{P_{i}^{\prime},\epsilon_{ij}^{\prime}\}$ and $\{P_{i}^{\prime\prime},\epsilon_{ij}^{\prime\prime}\}$ respectively, where $|P_{i}^{\prime}|=|P_{i}^{\prime\prime}|$. Let $\vec{d}$ be a vector in the DW plane within the reference structure. In a polar domain the electrostriction leads to strain, therefore $\vec{d}$ at the boundary of the first domain goes into Landau (2004)

[TABLE]

where the second term on the right-hand side represents a rigid rotation by angle $|\vec{\omega}^{\prime}|$ around vector $\vec{\omega}^{\prime}$. The third term corresponds to a deformation vector due to the strain tensor $\hat{\epsilon}^{\prime}$. An analogous formula holds for $\vec{d}^{\prime\prime}$ in the second domain. To avoid dislocations and cracks between the two domains, these deformed vectors $\vec{d}^{\prime},\vec{d}^{\prime\prime}$ must match at the DW, $\vec{d}^{\prime}=\vec{d}^{\prime\prime}$. Therefore

[TABLE]

Taking the scalar product with $\vec{d}$ eliminates the vector product term leading to:

[TABLE]

Therefore the strain compatibility requires the strain components in the DW plane to be equal in the two domains across the wall.

In the bulk of each domain, the strain tensor $\epsilon_{ij}$ can be obtained from the stress free condition $\sigma_{ij}=\partial f/\partial\epsilon_{ij}$,

[TABLE]

where $Q_{ijkl}=\frac{1}{2}(C^{-1})_{ijmn}\,q_{mnkl}$. We substitute Eq. S4 into Eq. S3, and take into account that only $P_{2}$ changes sign across the R71 DW, to obtain $P_{2}d_{2}(P_{1}d_{1}+P_{3}d_{3})=0$. The two solutions are $d_{2}=0$ and $[P_{1},P_{3}]\perp[d_{1},d_{3}]$, the latter meaning that all the vectors in the wall plane are perpendicular to $[P_{1},0,P_{3}]$. However, in thin films the violation of the constraint does not lead to infinite elastic energy, and allows the strain to relax at the surface. Therefore, the DW may deviate from the lowest energy orientation in the bulk. To confirm this, we fixed $P_{1},P_{3}$ and imposed domain walls in $P_{2}$ with different orientations $P_{2}\sim\tanh\frac{\vec{n}\cdot\vec{r}}{\lambda}$, as shown by black lines in the colored slabs in Fig. S1(a). As seen from Fig. S1(b), the DW energy is minimized when the DW normal $\vec{n}$ is parallel to $[P_{1},0,P_{3}]$. Note that in the bulk the deviations from the optimal DW orientation are penalized with elastic energy scaling linearly with the DW area. Fig. S1(b) shows that the energy penalty indeed comes from the elastic contribution due to strain mismatch between the domains.