Periodically driven small polarons

TL;DR

This paper investigates how strong external electric fields influence small lattice polarons, revealing that resonance with ion oscillations can liberate polarons and exponentially increase tunneling frequency.

Contribution

It introduces a direct time-domain integration of the Schrödinger equation to study polaron dynamics under external fields, highlighting resonance effects.

Findings

Resonance with ion oscillations liberates polarons.

Tunneling frequency increases exponentially at resonance.

External fields modulate carrier density and ion motion.

Abstract

Small lattice polarons driven by strong external electric fields are considered. A time-dependent Schroedinger equation is integrated directly in time domain. The field agitates ions both directly and through modulation of carrier density. It is found that when the field is in resonance with local ion oscillations, the polaron is liberated from its self-induced trap and the tunneling frequency increases exponentially.

Equations

$$H$$

$$H_c$$

$$H_i$$

$$H_{ci}$$

$$H_{cf}$$

$$H_{if}$$

$$\tau =\Omega t,\zeta =\frac{y}{y_0},\omega =\frac{{\omega }_0}{\Omega },j=\frac{J}{\hslash \Omega }.$$

$${\alpha }_c=\frac{|e|\mathcal{E}b}{\hslash \Omega },{\alpha }_i=\frac{|e|\mathcal{E}}{\hslash \Omega }\sqrt{\frac{\hslash }{M\Omega }}.$$

$$\lambda =\frac{g^2}{2M{\Omega }^{2}J}.$$

$$z_{ci}=\frac{z_c}{|z_c|}\cdot \frac{z_i}{|z_i|}.$$

$$g=z_{ci}\sqrt{2M{\Omega }^{2}J\lambda }.$$

$$i\frac{\partial}{\partial\tau}\left(\begin{array}[]{c}\psi_{1}\\ \psi_{2}\end{array}\right)=\left(\begin{array}[]{cc}\hat{h}_{1\zeta}&-j\\ -j&\hat{h}_{2\zeta}\end{array}\right)\left(\begin{array}[]{c}\psi_{1}\\ \psi_{2}\end{array}\right),$$

$${\hat{h}}_{1\zeta ,2\zeta }$$

$${\int }_{-\infty }^{\infty }\left\{|{\psi }_1(\tau ,\zeta ){|}^2+|{\psi }_2(\tau ,\zeta ){|}^2\right\}d\zeta =1.$$

$$e^{-i({\hat{h}}_{\zeta }+{\hat{h}}_j)△\tau }=\frac{\left(1-\frac{i}{2}{\hat{h}}_{\zeta }△\tau \right)\left(1-\frac{i}{2}{\hat{h}}_j△\tau \right)}{\left(1+\frac{i}{2}{\hat{h}}_{\zeta }△\tau \right)\left(1+\frac{i}{2}{\hat{h}}_j△\tau \right)}+O(△{\tau }^3),$$

$$\hat{h}_{\zeta}=\left(\begin{array}[]{cc}\hat{h}_{1\zeta}&0\\ 0&\hat{h}_{2\zeta}\end{array}\right),\hskip 14.22636pt\hat{h}_{j}=\left(\begin{array}[]{cc}0&-j\\ -j&0\end{array}\right).$$

$$\displaystyle\left(1+\frac{i}{2}\hat{h}_{\zeta}\triangle\tau\right)\left(1+\frac{i}{2}\hat{h}_{j}\triangle\tau\right)\left\{\begin{array}[]{c}\psi_{1}(\tau_{k+1},\zeta)\\ \psi_{2}(\tau_{k+1},\zeta)\end{array}\right\}=$$

$$\displaystyle\left(1-\frac{i}{2}\hat{h}_{\zeta}\triangle\tau\right)\left(1-\frac{i}{2}\hat{h}_{j}\triangle\tau\right)\left\{\begin{array}[]{c}\psi_{1}(\tau_{k},\zeta)\\ \psi_{2}(\tau_{k},\zeta)\end{array}\right\}.$$

$$\left(1+\frac{i}{2}\hat{h}_{\zeta}\triangle\tau\right)\left\{\begin{array}[]{c}\psi^{\ast}_{1}(\zeta_{l})\\ \psi^{\ast}_{2}(\zeta_{l})\end{array}\right\}=\left\{\begin{array}[]{c}b^{l}_{1}\\ b^{l}_{2}\end{array}\right\}.$$

$$\left(1+\frac{i}{2}\hat{h}_{j}\triangle\tau\right)\left\{\begin{array}[]{c}\psi_{1}(\tau_{k+1},\zeta_{l})\\ \psi_{2}(\tau_{k+1},\zeta_{l})\end{array}\right\}=\left\{\begin{array}[]{c}\psi^{\ast}_{1}(\zeta_{l})\\ \psi^{\ast}_{2}(\zeta_{l})\end{array}\right\}.$$

$$w(\zeta )=\frac{1}{2}{\zeta }^2-\sqrt{j^2+2\lambda j\cdot {\zeta }^2}.$$

$${\zeta }_0=\pm \sqrt{\frac{2j}{\lambda }\left({\lambda }^2-{\lambda }_{\mathrm{cr}}^2\right)}.$$

$$\tilde{\Omega }=\Omega \sqrt{1-\frac{{\lambda }_{\mathrm{cr}}^2}{{\lambda }^2}}\equiv \Omega \cdot \tilde{\omega }.$$

$${\psi }_{10,20}={\kappa }_{1,2}\cdot \frac{{\tilde{\omega }}^{1/4}}{{\pi }^{1/4}}{e}^{-\frac{1}{2}\tilde{\omega }(\zeta -{\zeta }_0{)}^2}.$$

$$\frac{{\kappa }_2}{{\kappa }_1}=2\lambda -2\sqrt{{\lambda }^2-{\lambda }_{\mathrm{cr}}^2},$$

$${\psi }_{0\beta }=e^{-\beta ({\hat{h}}_{\zeta }+{\hat{h}}_j)}{\psi }_0,$$

$$e^{-(△\beta )({\hat{h}}_{\zeta }+{\hat{h}}_j)}\approx$$

$$\left(1-\frac{1}{2}{\hat{h}}_j△\beta \right)\left(1-{\hat{h}}_{\zeta }△\beta \right)\left(1-\frac{1}{2}{\hat{h}}_j△\beta \right).$$

$$P_{1,2}(\tau )={\int }_{-\infty }^{\infty }d\zeta {\left|{\psi }_{1,2}(\tau ,\zeta )\right|}^2,$$

$$P_{(L,R)}(\tau )={\int }_{(-\infty ,0)}^{(0,\infty )}\left\{|{\psi }_1(\tau ,\zeta ){|}^2+|{\psi }_2(\tau ,\zeta ){|}^2\right\}d\zeta .$$

$$E(\tau )$$

$$i\frac{\partial \varphi (\zeta ,\tau )}{\partial \tau }=\left\{-\frac{1}{2}\frac{{\partial }^2}{\partial {\zeta }^2}+\frac{1}{2}{\zeta }^2-z_i{\alpha }_i\zeta \sin (\omega \tau )\right\}\varphi (\zeta ,\tau ),$$

$$i\frac{\partial {\chi }_1}{\partial \tau }$$

$$i\frac{\partial {\chi }_2}{\partial \tau }$$

Full text

I

Introduction

Polaron is a carrier in a condensed-matter system that deforms surrounding atoms or ions from their equilibrium positions. A stable polaron moves through the system together with the deformation. The deformation impedes polaron motion, increases its effective mass and reduces mobility. In polar solids with low carrier density, the carrier-ion interaction derives from an unscreened Coulomb interaction which is strong at short distances. Polaron transport properties sharply, often exponentially, depend on ion positions, and are a sensitive probe of the ion subsystem. Traditionally, ion motion has been controlled through global variables such as temperature or pressure, or through local modifications such as isotope substitutions. Recently, however, direct excitation of ions with lasers became possible. Hu2014 Ions are driven into resonance by strong laser fields and acquire oscillation amplitudes as large as several percent of a lattice constant. Hu2014 Such displacements significantly reduce or even completely eliminate the deformation caused by the interaction with a carrier. As a result, the carrier is “liberated” from the potential trap created by the ions and its transport changes dramatically. External driving of ions becomes therefore a powerful tool to probe intrinsic polaron properties.

Polaron physics has a rich history that dates back to Landau, Landau1933 Pekar, Pekar1946 ; Landau1948 and Fröhlich. Froehlich1954 The development of the field is well documented in several books and reviews. Pekar1951 ; Kuper1963 ; Appel1968 ; Devreese1972 ; Levinson1973 ; Firsov1975 ; Devreese1984 ; Alexandrov1996 ; Devreese1996 ; Alexandrov2007 ; Devreese2009 ; Emin2013 In addition, for decades the polaron has served as a testing ground for novel analytical and numerical techniques including perturbative expansions, Roeseler1968 ; Miyake1975 ; Smondyrev1986 ; Stephan1996 ; Alexandrov2003 path integrals, Feynman1955 ; Feynman1962 ; Verbist1991 Monte Carlo methods, Hirsch1982 ; Fradkin1983 ; DeRaedt1982 ; Alexandrou1990 ; Alexandrou1992 ; McKenzie1996 ; Kornilovitch1997 ; Kornilovitch1998 ; Prokofev1998 ; Mishchenko2000 ; Mishchenko2003 ; Macridin2004 ; Hague2007 cluster diagonalization, Ranninger1992 ; Kabanov1994 ; Marsigio1995 ; Wellein1996 ; Wellein1997 variational, Feynman1955 ; Lee1953 ; Romero1998 ; Bonca1999 ; Bonca2000 ; Chakrabarti2007 ; Chakraborty2011 ; Chandler2016 and other Tjablikov1951 ; Sewell1958 ; Holstein1959a ; Holstein1959b ; Friedman1963 ; Eagles1963 ; Lang1963 ; Emin1969 ; Emin1976 ; Kabanov1993 ; Freericks1993 ; Ciuchi1997 ; Zhang1998 ; Jeckelmann1998 ; Hohenadler2004 ; Berciu2006 ; Goodvin2006 ; Fratini2006 methods. However, there has hardly been any work on externally driven polarons or bipolarons. Most theoretical treatments of polaron-light interactions are based on the linear response theory, Feynman1962 ; Lang1963 ; Eagles1963 ; Friedman1963 ; Devreese1972 ; Mishchenko2003 ; Fratini2006 ; Devreese2009 and as such assume weak coupling between a charge carrier and an external electric field. The field is treated as a perturbation that does not change polaron states but only causes transitions between them. In this work, we are interested in the opposite limit of strong coupling when the field directly alters the local state of the polaron. Arguably, this is a much harder mathematical problem, and this is perhaps one reason why it has not been widely addressed before. If fact, we are not aware of any publication on this subject prior to our own recent paper. Kornilovitch2016

The basic physics of driven polarons can be understood on a two-site model system. In the absence of both ion displacement and electric field, Fig. 1(a), two atomic levels are in resonance, leading to carrier transfer with a bare amplitude $J$. An external electric field, Fig. 1(b), periodically drives the levels on and off resonance, which on average reduces the tunneling amplitude to $J_{1}<J$. When an ion is included in the system, Fig. 1(c), its displacement drives the levels off resonance on a more sustained basis. Tunneling is possible only with simultaneous reversal of the deformation, and as such is exponentially suppressed. This constitutes polaron formation. However, a sufficiently strong external field can drive the ion in resonance, see Fig. 1(d). That symmetrizes the potential so that the two levels come into resonance periodically. As a result, Kornilovitch2016 the tunneling amplitude increases exponentially compared to the undriven case, $J_{3}\gg J_{2}$, but of course remains smaller than bare amplitude $J$.

The goal of the present work is to extend the treatment of Ref. [Kornilovitch2016, ] beyond the adiabatic approximation. Additionally, indirect ion-field interaction through modulation of carrier density is also included in the picture. In general, the direct and indirect interactions compete with each other and their balance should be carefully analyzed. Only single polaron case is studied while the bipolaron case is left for future work.

II

Two-site driven polaron model

We consider one carrier of charge $z_{c}$ moving between two tight-binding sites $|1\rangle$ and $|2\rangle$ with a hopping integral $J$. The distance between the sites is $2b$. One ion of mass $M$, frequency $\Omega$, charge $z_{i}$, and dynamic coordinate $y$ is positioned symmetrically between the sites. The symmetric position and polarization of ion displacement (parallel to the line between the sites) is inspired by the layered structure of high-temperature superconductors. Note that this arrangement is different from the Holstein model of molecular crystals, Holstein1959a where each ion interacts with one lattice site only. The full Hamiltonian consists of five terms: free carrier $H_{c}$, free ion $H_{i}$, carrier-ion $H_{ci}$, carrier-field $H_{cf}$, and ion-field $H_{if}$:

[TABLE]

Here ${\cal E}$ and $\omega_{0}$ are the field’s amplitude and frequency, and $g$ is the force between the ion and the carrier on either site. The force is taken to be independent of $y$, i.e., carrier-ion interaction is treated in the linear approximation. Note that Eqs. (1)-(6) are invariant under simultaneous site inversion $|1\rangle\leftrightarrow|2\rangle$, coordinate inversion $y\rightarrow-y$ and time inversion $t\rightarrow-t$.

Because of a large number of physical parameters, it is useful to transform to dimensionless variables. It is convenient to measure everything in oscillator units. We choose $\Omega^{-1}$, $\hbar\Omega$, and $y_{0}=\sqrt{\hbar/(M\Omega)}$ as the units of time, energy, and length, respectively. Accordingly, dimensionless time $\tau$, coordinate $\zeta$, frequency $\omega$, and hopping integral $j$ are introduced as follows:

[TABLE]

Dimensionless coupling constants between the field and the carrier and between the field and the ion follow from the transformation:

[TABLE]

To estimate $\alpha_{c}$ and $\alpha_{i}$, one can use values reported by Hu et al Hu2014 in their studies of dynamically stabilized superconductivity in YBa2Cu3O6.5: ${\cal E}=3.0$ MV cm*-1* and $\hbar\Omega=83$ meV. Using the distance between copper-oxygen bilayers (8.2 Å) for $2b$ and the mass and charge of the oxygen ion, one obtains $\alpha_{c}\approx 1.5$ and $\alpha_{i}\approx 0.02$. Despite the fact that $\alpha_{i}\ll\alpha_{c}$, effects of the two coupling types are comparable, as will be detailed below.

The carrier-ion coupling constant $\lambda$ is now discussed. Commonly, Alexandrov1996 $\lambda$ is defined as the ratio of polaron shift $E_{p}=g^{2}/(2M\Omega^{2})$ (polaron energy in the atomic limit $J=0$) to half of bare bandwidth $D$ (i.e. the lowest energy of a free carrier). The half bandwidth in a two-site system is $D_{2}=J$. Accordingly, $\lambda$ is defined here as

[TABLE]

Expression (9) is quadratic in $g$, so the same $\lambda$ describes both attraction and repulsion. To distinguish between the two possibilities, a sign variable $z_{ci}$ is introduced

[TABLE]

Accordingly, $g$ in Eq. (4) can be written as:

[TABLE]

For the sake of visual clarity and keeping in mind the physics of YB2C3O6.5, hereafter we will focus on the case of attraction, i.e. positive $z_{c}=+1$ (a hole inside a copper-oxygen plane) and negative $z_{i}=-2$ (an apical oxygen ion), which implies $g<0$. With such a choice, the on-site energy of $|1\rangle$ decreases for negative $y$, see Eq. (4).

The full wave function is a two-element array $\{\psi_{1};\psi_{2}\}$, both elements being functions of time $\tau$ and ion displacement $\zeta$. Collecting all the definitions, the Schrödinger equation reads

[TABLE]

where

[TABLE]

Normalization is chosen to be

[TABLE]

III

Numerical method

Equation (12) does not admit analytical solutions and has to be integrated numerically. One approach is based on Floquet theory and spectral expansion of $\psi(\tau)$ in a Floquet basis. Grossmann1991 ; Grifoni1998 Here we adopt an alternative approach of direct integration in time domain. The particular differencing scheme is based on the Cayley form of the finite-time evolution operator Goldberg1967 ; Press1992

[TABLE]

where

[TABLE]

This representation is unitary and accurate to order $(\triangle\tau)^{2}$ even for non-commuting operators $\hat{h}_{\zeta}$ and $\hat{h}_{j}$, which is the case here. Equation (15) leads to an unconditionally stable implicit time-stepping rule Crank1947

[TABLE]

where $\tau_{k}$ is the $k$-th time step. Factorization of the evolution operator into $\zeta$ and $j$ parts enables sequential inversion: first along $\zeta$ axis and then along the site index. The time-stepping procedure is as follows. (i) Both wave function components are discretized along $\zeta$ dimension: $\psi_{1,2}(\tau_{k},\zeta)\rightarrow\psi_{1,2}(\tau_{k},\zeta_{l})$. (ii) Starting with $\psi_{1,2}(\tau_{k},\zeta_{l})$, the right hand side of Eq. (22) is computed. First, $j$ part is applied which amounts to a $(2\times 2)$ matrix multiplication for each $l$. Second, $\zeta$ part is applied, which is performed independently for $\psi_{1}$ and $\psi_{2}$ [$\hat{h}_{\zeta}$ in Eq. (16) is diagonal]. The resulting right hand side is a two-row matrix $\{b^{l}_{1};b^{l}_{2}\}$. (iii) The left hand side of Eq. (22) is inverted, again in two consecutive steps. First, an intermediate function $\psi^{\ast}$ is found from the following relation

[TABLE]

Since $\hat{h}_{\zeta}$ is diagonal, $\psi^{\ast}_{1}$ and $\psi^{\ast}_{2}$ are computed independently using the method of Ref. [Goldberg1967, ]. (iv) Finally, $\psi(\tau_{k+1},\zeta_{l})$ is found by solving

[TABLE]

Since operator $\hat{h}_{j}$ does not mix spatial coordinates, solving the last equation amounts to inverting a $(2\times 2)$ matrix for each $l$. (v) The above sequence is repeated for the entire time interval of interest. Most results presented in this paper were obtained with $\triangle\tau=0.001$ and $\triangle\zeta=0.1$. Thus, time evolution over a total time of $\tau_{\rm max}=10^{4}$ requires $10^{7}$ time steps.

Choice of initial conditions $\psi_{10,20}$ is now discussed. The primary quantity of interest in this work is polaron tunneling frequency between the sites. Therefore, we seek an initial state that is localized at one of the two sites. Adiabatic approximation provides a good starting point. Remove the second derivative and assume $\alpha_{c}=\alpha_{e}=0$ in Eq. (13). Then fix $\zeta$ and solve the two-site one-carrier problem to obtain the polaron adiabatic potential

[TABLE]

At $\lambda>\lambda_{\rm cr}=\frac{1}{2}$, $w(\zeta)$ develops two symmetric minima at

[TABLE]

Near $\zeta_{0}$, the potential is quadratic with renormalized frequency (in physical units)

[TABLE]

Accordingly, a good starting wave function describing the polaron in one of its minima is the ground state of an oscillator with frequency $\tilde{\omega}$ shifted by $\zeta_{0}$

[TABLE]

Weights $\kappa_{1,2}$ follow from solving the two-site problem

[TABLE]

and from the normalization condition $|\kappa_{1}|^{2}+|\kappa_{2}|^{2}=1$.

The initial wave function can be further improved by applying the projection operator

[TABLE]

where $\beta$ is a positive dimensionless number. Calculation of Eq. (30) is also performed in a step-like fashion, by splitting total imaginary time $\beta$ into small steps $\triangle\beta\ll 1$ and replacing the short-time evolution operator with a $(\triangle\beta)^{2}$-accurate representation DeRaedt1987

[TABLE]

Application of all the factors is done in the same way as in the case of real-time evolution, as described above. Since the projection operator does not conserve normalization, $\psi_{0}$ is normalized after every time step $\triangle\beta$. Most results presented below were obtained with $\beta=2.0$ and $\triangle\beta=0.001$.

For the purposes of this paper, a full wave function contains too much information. To visualize tunneling rate, it is more convenient to consider an integral probability for the carrier to occupy either site 1 or site 2

[TABLE]

or probability for the ion to be on the left (L) or on the right (R) from the symmetry point $\zeta=0$

[TABLE]

Another quantity of interest is instantaneous polaron energy without ion-field and carrier-field contributions

[TABLE]

It can be shown by standard means that $E(\tau)$ is real for all times $\tau$.

IV

Limit cases

In this section, some limit cases of general model (12) are reviewed. Consider first zero carrier-ion interaction, $\lambda=0$. That splits the system in two independent subsystems: (i) a free ion interacting with an external field, and (ii) a free carrier hopping between two sites and interacting with an external field. Solution to Eq. (12) is sought in a factorized form: $\psi_{1,2}(\zeta,\tau)=\phi(\zeta,\tau)\cdot\chi_{1,2}(\tau)$. Equation (12) decouples into one equation for $\phi$:

[TABLE]

and two coupled equations for $\chi_{1,2}$

[TABLE]

Equation (35) admits an explicit solution by means of Husimi substitution Husimi1953 ; Kerner1958 that maps a forced harmonic oscillator onto a free oscillator. For example, the lowest ($m=0$) Floquet state has the form

[TABLE]

with

[TABLE]

[TABLE]

[TABLE]

This explicit formula can be used to validate the numerical method described in the preceding section. In case of resonant excitation, $\omega=1$, no stable periodic solution is possible. Instead, the amplitude of wave function oscillations $\eta(\tau)$ grows linearly with time. If the oscillator is in its ground state when the field is turned on, the time-dependent solution still has the form of Eq. (38) but with $\varepsilon_{0}=\frac{1}{2}$ and new functions $\eta(\tau)$ and $S(\tau)$:

[TABLE]

[TABLE]

Equations (36) and (37) and their analogues in quantum optics have been studied in several papers Bloch1940 ; Shirley1965 ; Aravind1984 ; Agarwal1994 ; Raghavan1996 but no analytical solution has been reported. (Somewhat surprisingly, the more general problem of an infinite chain in a periodically varying electric field is exactly solvable, providing analytical description of dynamic localization. Dunlap1986 ) To gain some insight, we show several numerical solutions in Fig. 2. At low frequencies, $\omega\ll j$, the on-site energies change slowly compared to inter-site hopping. When the levels are out of resonance, the carrier gets confined to one of the sites and the frequency of oscillations increases, as can be seen in the left panels of Fig. 2. At high frequencies, $\omega\gg j$, fast oscillations of on-site energies average to zero and the dynamics approaches that of a free carrier with an intersite tunneling amplitude $j$, see Fig. 2 on the right. At intermediate frequencies, carrier dynamics is complex, as illustrated in Fig. 2, center.

In the absence of an external field, $\alpha_{i}=\alpha_{c}=0$, the general problem, Eq. (12), reduces to a free two-site polaron. Note that although the present model involves an intersite ion that couples to both sites simultaneously, the model can be mapped to the two-site Holstein polaron that has been studied in detail. Ranninger1992 ; Kabanov1994 ; Holstein1959b A relevant quantity is the energy split of the lowest level pair, as it defines polaron “mass” and tunneling frequency in the undriven case [$f_{\rm tun}=\triangle\epsilon_{12}/(2\pi)$]. In this paper, the ground state energy and the first excited state energy are computed by discretizing polaron Hamiltonian of Eq. (12) in $\zeta$ space and directly diagonalizing a resulting matrix, see Fig. 3. Tunneling frequencies thus obtained are later used as a reference for the driven case.

There are two well understood limits of the (undriven) two-site polaron. In the adiabatic limit, Holstein1959b ; Emin1976 ; Kabanov1993 ; Kabanov1994 $j\gg 1$, the two-function Schrödinger equation, Eq. (12), is reduced to a one-function equation with Hamiltonian

[TABLE]

where $w(\zeta)$ is the adiabatic potential, Eq. (25). Eigenvalues of $\hat{h}_{\rm ad}$ can be found numerically by discretizing Eq. (44) using standard rules. Level splitting can also be estimated with the instanton technique, Coleman1977 ; Rajaraman1982 ; Kleinert2004 leading to the following result: Kornilovitch2016 ; Kabanov1994

[TABLE]

[TABLE]

[TABLE]

where $\lambda_{\rm cr}=\frac{1}{2}$. Figure 3 compares Eq. (45) with numerical diagonalization of adiabatic Hamiltonian $\hat{h}_{\rm ad}$ and of full Hamiltonian, Eq. (12).

In the opposite, anti-adiabatic limit, $j\ll 1$, a fast ion instantaneously follows carrier transitions. Level splitting is derived from the sudden approximation Tjablikov1951 ; Sewell1958 ; Holstein1959a

[TABLE]

We now turn to the full time-dependent problem, Eq. (12), with nonzero $\alpha_{i}$ and $\alpha_{c}$. We begin with the adiabatic limit $j\gg 1$.

V

Two-site driven polaron: Adiabatic case

Effects of direct ion-field interaction, $\alpha_{i}\neq 0$, $\alpha_{c}=0$, are considered first. Figure 4 compares time evolution of the driven and undriven polaron for $j=10$ and $\lambda=0.9$. Period of undriven oscillations is $T_{0}=7080$ which is consistent with the lowest doublet split of $\triangle\epsilon_{12}=8.81\cdot 10^{-4}$ [formula (45) is accurate to 5.7%]. In contrast, the driven case shows much faster oscillations with a period $T\approx 410$, i.e. a decrease of 17.3 times. Notice that the driven $P_{L}(\tau)$ is not a clean sinusoid but rather a complex function with multiple overtones.

A possible mechanism behind faster oscillations is now discussed. The argument is based on the property of quantum oscillators that the centroid of a driven wave function follows a classical equation of motion. Husimi1953 Consider an ion near the bottom of the adiabatic potential, as shown in Fig. 5. Frequency of small-amplitude oscillations is $\tilde{\omega}=\sqrt{1-(\lambda_{\rm cr}/\lambda)^{2}}$, see Eq. (27). The classical equation of motion at resonance reads

[TABLE]

An explicit solution with arbitrary starting point $x_{0}$ and zero initial velocity is

[TABLE]

The oscillation-averaged total energy is

[TABLE]

Long-term behavior of $\langle\varepsilon(\tau)\rangle$ is quadratic growth. However, short-term details depend on starting position $x_{0}$. For positive $x_{0}$, energy initially decreases and oscillations stop before picking up again. On the basis of this observation, the following mechanism can be proposed, see Fig. 5. Initially, the ion is near the bottom of the left well of the adiabatic potential. The external field drives the ion in resonance and its energy rises. As energy approaches the top of the potential barrier, quantum-mechanical tunneling probability increases exponentially and eventually the ion tunnels under the barrier into the right well. There, the ion finds itself out-of-phase with the field (the in-phase position would be on the right side of the right well) and starts losing energy. Once the ion descends to the bottom of the right well, the process repeats but in the opposite direction.

Reality is more complicated. One complication comes from non-harmonicity of the adiabatic potential. As can be inferred from Eq. (25), the instantaneous frequency decreases with the amplitude of oscillations. As a result, an external field cannot always be in resonance. One can only speak of an “effective” or “average” resonance for the duration of the process. It is clear that the rate of energy absorption or loss is less than for ideal resonance. Another complication is the gradual nature of quantum-mechanical tunneling. Numerical analysis shows that under-barrier tunneling is not a sharp process. Rather, the full wave function is nonzero in both wells at all times, so the “start” or “end” of tunneling is not easy to identify.

The qualitative argument presented above leads to an important conclusion. According to Eq. (50), the rate of amplitude increase is proportional to coupling constant $\alpha_{i}$. If one assumes that tunneling takes place when energy reaches a certain value, and by association when the amplitude reaches a certain value, then the tunneling condition will be reached faster at larger $\alpha_{i}$. This implies that the frequency of driven tunneling is directly proportional to the field strength and to the square root of laser intensity

[TABLE]

To verify this prediction, we show $f_{\rm tun}(\alpha_{i})$ computed for $j=10$ and $\lambda=0.9$ in Fig. 6. $f_{\rm tun}$ has been determined from the location of the largest peak in a Fourier power spectrum of $P_{L}(\tau)$. Relation (52) approximately holds at intermediate $0.03<\alpha_{i}<0.10$. There is a local minimum at weak couplings, $\alpha_{i}\approx 0.01$, which implies that a very weak coupling initially disrupts quantum-mechanical tunneling and decreases the tunneling rate. In the opposite limit of strong coupling, $\alpha_{i}>0.1$, oscillations become highly irregular without a dominant harmonic in the Fourier transform.

Notice that experimentally achievable field strengths, Hu2014 $\alpha_{i}=0.02$, lie near the lower end of the linearity region, according to the present calculation. This suggests that a further increase in the field strength should systematically raise $f_{\rm tun}$. This effect may have implications for superconductivity in the cuprates, as was suggested in Ref. [Kornilovitch2016, ].

In real crystals, local ion frequencies may be shifted by dispersion, lattice imperfections, temperature, and other factors. Therefore it is important to understand the sensitivity of $f_{\rm tun}$ increase to variations of ion frequency $\Omega$. In the present model, $\Omega^{-1}$ is taken to be a unit of time, so one should consider sensitivity to external driving frequency $\omega$ instead. Figure 7 presents $P_{L}(\tau)$ for several values of $\omega$. One can observe that enhanced polaron tunneling persists for detunings of up to several percent of $\Omega$ on both sides of the “optimal” frequency $\omega=0.71$. (Which of course renders the choice of optimal frequency itself somewhat arbitrary.) For larger detunings, $P_{L}(\tau)$ assumes more complex beat-like shapes which makes it difficult to interpret the motion as periodic tunneling between two potential wells. An important conclusion, however, is that polaron delocalization is a robust effect that exists in a finite interval of ion and driving frequencies, and as such should be observable in real materials.

So far, only interaction between the external field and ion has been taken into account. Now we consider the opposite case of carrier-field interaction in the absence of ion-field interaction: $\alpha_{e}\neq 0$ and $\alpha_{i}=0$. Figure 8 compares $P_{L}(\tau)$ for several $\alpha_{e}$. The overall $f_{\rm tun}(\alpha_{e})$ dependence is similar to the $f_{\rm tun}(\alpha_{i})$ behavior: initially, at small $\alpha_{e}$, the frequency decreases relative to the undriven case before increasing at larger $\alpha_{e}$.

A possible underlying mechanism is now discussed, see Fig. 9. In the adiabatic limit, $j\gg 1$, $j\gg\omega$, the carrier wave function equilibrates between the two sites for any instantaneous position of the ion and instantaneous value of the field. The time-dependent field “rocks” the carrier density between the two sites. In turn, the oscillating carrier density pulls the ion with variable strength and excites the latter in its own well if the field frequency is close to the ion frequency. Thus instead of direct excitation by the field, the ion is driven into resonance by oscillating carrier density and nonzero carrier-ion interaction.

Fluctuations of carrier density occur on a background set by the much larger bare kinetic energy $J$. As a result, influence of the external field is effectively reduced. In order to estimate this effect, solve the two-site problem for small $\alpha_{e}\ll j$. In the first order, one obtains for the wave functions

[TABLE]

and for the densities

[TABLE]

where

[TABLE]

(Note that $|\psi_{1}(\tau)|^{2}+|\psi_{2}(\tau)|^{2}=1$.) Thus, a carrier-field interaction plus a carrier-ion interaction is equivalent to an ion-field interaction with some effective coupling constant $\tilde{\alpha}_{i}$. Recasting carrier-ion interaction, Eq. (4), with $\triangle n$ from Eq. (56) into ion-field interaction form, Eq. (6), one obtains

[TABLE]

For $j=10$ and $z_{i}=-2$, the last formula changes from $\tilde{\alpha}_{i}\approx-\frac{\alpha_{e}}{6}$ at $\lambda=\frac{1}{2}$ to $\tilde{\alpha}_{i}\approx-\frac{\alpha_{e}}{36}$ at $\lambda=1$. Using the estimates given in Sec. II, $\alpha_{e}=1.5$ and $\alpha_{i}=0.02$, one concludes that $\tilde{\alpha}_{i}>\alpha_{i}$, that is the indirect light-ion interaction is at least comparable to the direct interaction and may even be dominant. One should note that precise values of $\alpha_{e}$ and $\alpha_{i}$ depend on the local crystal structure and fine details of the dielectric response of a particular solid. These topics are beyond the scope of the present work and for this reason $\alpha_{e}$ and $\alpha_{i}$ will continue to be treated here as phenomenological parameters.

Nonetheless, one general conclusion can be reached by examining Fig. 9: the direct and indirect ion-field interactions are always of opposite signs. For example, in case of attraction (depicted in the figure), the same field pulls the ion in one direction while modifying carrier density to pull the ion in the opposite direction. This can also be seen in Eq. (57) that becomes negative whenever $z_{e}$ and $z_{i}$ are of opposite signs. The main conclusion remains true in case of carrier-ion repulsion.

Thus the two interactions always compete and the combined effect is always smaller than either of the two. This competing nature is confirmed by numerical calculations summarized in Fig. 10. Shown is the polaron tunneling frequency, as determined from the location of the largest peak in a power spectrum of $P_{L}(\tau)$, as a function of $\alpha_{i}$ and $\alpha_{e}$. It is evident from the figure that $f_{\rm tun}$ is in fact a function of $\alpha^{\prime}_{i}=\alpha_{i}-c^{-1}\alpha_{e}$, where $c\approx 22-25$ is a numerical coefficient that is in agreement with Eq. (57). Given that neither $\alpha_{i}$ nor $\alpha_{e}$ is known (unless detailed first-principle calculations of the dielectric response are performed) it is convenient to lump both interactions in one effective ion-field interaction with a new phenomenological coupling constant $\alpha^{\prime}_{i}$. Therefore, $\alpha_{i}=\alpha^{\prime}_{i}$ and $\alpha_{e}=0$ will be assumed for the rest of the paper. Again, we choose to deal with an effective ion-field interaction rather than carrier-field interaction because of easier visualization.

We now return to the main subject of this paper: exponential increase of the tunneling frequency under a driving field. Figure 11 compares $f_{\rm tun}$ of the driven and undriven polarons as a function of $\lambda$. The undriven line is essentially the same that was shown earlier in Fig. 3 but now computed not from the spectrum of a stationary Schrödinger equation but from Fourier analysis of the time-dependent Schrödinger equation, and divided by $2\pi$. The two methods produce numbers that match within less than 1%. The undriven case shows the familiar exponential decrease with $\lambda$ reflecting a deepening potential barrier between two minima. In contrast, $f_{\rm tun}$ of the driven case remains approximately constant. This is understandable, since a small deepening of the potential results in a proportionally small increase of the time the ion needs to gain enough energy to get over the barrier. In other words, in the driven case, $f_{\rm tun}$ is more a function of the field strength rather than details of the potential.

An important consequence is that the ratio between driven and undriven $f_{\rm tun}$ grows exponentially in the deep polaron regime. In the example of Fig. 11, the driven frequency exceeds the undriven one by two orders of magnitude at $\lambda=1.0$ and by three orders at $\lambda=1.1$. Thus field assistance promotes polaron tunneling and increases polaron mobility. This effect has important consequences for transport properties of polaron systems.

VI

Non-adiabatic cases

VI.1

Slow carriers, fast ions: $j\ll 1$

We now turn to the opposite anti-adiabatic case $j\ll 1$. In this regime, ion oscillations are fast and carrier transitions between the sites are rare. When the carrier is confined to one of the two sites, the ion shifts to a new equilibrium position toward the carrier by $(2j\lambda)^{1/2}$. When the carrier tunnels to the second site, the displacement must reverse. Overlap of ion wave functions separated by $2(2j\lambda)^{1/2}$ leads to an exponentially small level splitting given by Eq. (48). Numerical solution confirms this physical picture. Figure 12(a) shows temporal evolution of an undriven polaron for $j=0.1$ and $\lambda=35$. The left-space probability has an ideal sine-wave shape as expected for a closely split level dublet. The numerical time period is in perfect agreement with the analytical formula, as shown in Fig. 12(b).

An external electric field drives the ion in resonance. As the oscillation amplitude rises in accordance with Eq. (50), the overlap integral grows allowing for more frequent transitions. Increased tunneling rate between the two states is evident in Fig. 12(c). Figure 12(d) of the same figure shows the corresponding Fourier spectrum with a prominent peak near 0.00165. Different from the adiabatic case, in the $j\ll 1$ limit the ionic potential remains undistorted and harmonic. As a result, the most effective driving frequency is simply $\omega=1.0$. This is confirmed by numerical results.

Similarly to the adiabatic case, stronger external fields cause more frequent tunneling. This can be seen by comparing the time series in Figs. 12(c) and 12(e). In the latter case, the field is almost twice as strong (0.014 vs. 0.008), which leads to faster oscillations. The overall $f_{\rm tun}(\alpha^{\prime}_{i})$ dependence is roughly linear, as shown in Fig. 12(f). The physical reason is the same as before: the amount of time needed to reach an amplitude at which tunneling takes place “comfortably” is inversely proportional to $\alpha^{\prime}_{i}$ and ${\cal E}$.

VI.2

Intermediate frequencies: $j\sim 1$

As a typical example of intermediate phonon frequencies we consider $j=1.0$. The polaron is expected to behave between the $j\gg 1$ and $j\ll 1$ limits. In particular, an optimal driving frequency should be less than 1.0 but not as much as in the adiabatic regime. Figure 13 shows a sample of numerical results obtained for $\lambda=5.0$. A scan over driving frequencies found that the most prominent tunneling occurs at $\omega=0.99$, i.e., this case is closer to the anti-adiabatic limit. The figure reveals the now familiar physical picture. The undriven case (not shown) corresponds to a deep polaron regime with a large tunneling period of $T\approx 85\,880$. An external force promotes tunneling and increases $f_{\rm tun}$ by 2-to-3 orders of magnitude. $f_{\rm tun}$ increases linearly with the field’s strength in accordance with the physical mechanism described in previous sections.

VII

Summary and discussion

The ability to drive ions in resonance by external laser fields is a powerful technique of experimental solid state physics. In systems where ions are strongly coupled to charge carriers or other degrees of freedom, direct excitation of ions must lead to large and measurable changes in other subsystems. In this paper, we investigated one model system where such effects are particularly strong: the small lattice polaron. By definition, a lattice polaron is a charge carrier that interacts with surrounding ions so strongly that it deforms the lattice and displaces the ions from their equilibrium positions by finite amounts. If those ions are agitated by a different source, in our case by an external laser field, that should change the balance of forces between the ions and the carrier. Specifically, the energy barrier that prevents the carrier from tunneling between lattice sites will be affected. That should lead to exponentially large changes in polaron tunneling rates, which may be detectable through transport or optical response.

In the absence of analytical solutions even for the simplest cases, we employed direct step-by-step integration of the time-dependent Schrödinger equation as the main analysis tool. A two-site polaron model has been chosen as being sufficient to reveal the essential physics. Our findings can be summarized as follows.

(i) Increase of polaron tunneling rate has been found at all ion frequencies $\Omega$ (as compared to the bare transfer integral $J$), as long as the field is appropriately tuned to the ion resonance. The effect exists in a finite range of driving frequencies. In accordance with general polaron theory, an optimal driving frequency has been found to be $\omega=\Omega$ in the anti-adiabatic limit $J\ll\Omega$, and gradually decreasing to $\omega<\Omega$ in the adiabatic limit $J\gg\Omega$.

(ii) Resonant energy build up has been proposed as a common physical mechanism behind enhancement of tunneling. When an oscillator, quantum or classical, is driven in resonance, its amplitude grows linearly with time. As a result, one of the ion’s turning points gets progressively closer to a symmetrical position between two stable states, symmetrizing the polaron potential and reducing its height. After a finite number of oscillations the barrier is reduced so much that the ion tunnels “easily.” (Of course, tunneling is not instantaneous but because of exponential dependence, at a very high level one can assume that tunneling happens as soon as the barrier is reduced to a certain height.) After tunneling, the ion initially loses energy dropping to the potential bottom. After that, the process begins in the opposite direction.

(iii) It follows from the above mechanism that $f_{\rm tun}$ scales linearly with the number of resonant oscillations, i.e., linearly with the potential’s depth. At the same time, frequency of undriven tunneling scales exponentially with the potential’s depth. Thus the ratio of driven to undriven $f_{\rm tun}$ scales exponentially with the coupling constant and as such can reach several orders of magnitude in the deep polaron regime.

(iv) With similar reasoning, one can predict a linear dependence of $f_{\rm tun}$ on the field amplitude ${\cal E}$. (And subsequently, a square root dependence on laser intensity, $f_{\rm tun}\propto\sqrt{I}$.) The rate of amplitude increase is proportional to ${\cal E}$. Therefore, the number of oscillations needed to reach the “easy tunneling” condition is $\propto{\cal E}^{-1}$, from which the stated dependence follows.

(v) In addition to affecting ions, the external field couples directly to carriers. The field modulates carrier density on both sites, which in turn rocks the ion and drives it in resonance even in the absence of direct laser-ion interaction. The two interactions are generally of the same order but they pull the ion in opposite directions, so that the two forces partially compensate each other.

Let us conclude by briefly discussing possible consequences of polaron “undressing”. The present work was motivated by recent experiments on dynamically stabilized superconductivity. Hu2014 We previously proposed Kornilovitch2016 that exponential enhancement of bipolaron tunneling rates between copper-oxygen bilayers could explain the observed increase of the apparent critical temperature in YBCO. This proposal requires additional proof. Given that the observable quantity in Ref. [Hu2014, ] was an optical response, it would make sense to compute the dynamical optical response of a driven (bi)polaron. This calculation is left for future work. More generally, any polaron property that depends on the tunneling probability, for example dc conductivity or Hall coefficient, will be affected by a resonant excitation of ions. A general theory of transport properties of driven polarons seems to be worth developing.

Acknowledgements.

The Author wishes to thank Alexei Voronin for help with computational resources, as well as David Dunlap, David Eagles and Vassilios Kovanis for useful discussions on the physics of driven systems.