Table Of ContentSemiclassical Time Evolution of the Density Matrix
and Tunneling
Joachim Ankerhold and Hermann Grabert
Fakult¨at fu¨r Physik, Albert-Ludwigs-Universit¨at Freiburg, Hermann-Herder-Straße 3, D-79104 Freiburg, Germany
0 (February 1, 2008)
0
0 The time dependent density matrix of a system with potential barrier is studied using path
2 integrals. The characterization of the initial state, which is assumed to be restricted to one side of
thebarrier, and thetime evolution of thedensity matrix lead to a three-fold path integral which is
n
evaluatedinthesemiclassical limit. Thesemiclassical trajectoriesarefoundtomoveinthecomplex
a
J coordinate plane and barrier penetration only arises due to fluctuations. Both the form of the
semiclassical paths and the relevant fluctuations change significantly as a function of temperature.
6
2 Thesemiclassical analysis leads toadetailed pictureof barrierpenetration in therealtimedomain
and the changeover from thermal activation to quantum tunneling. Deep tunneling is associated
] with quasi-zero modes in the fluctuation spectrum about the semiclassical orbits in the long time
h
limit. The connection between this real time description of tunneling and the standard imaginary
c
timeinstanton approach isestablished. Specificresults aregivenfor adoublewell potentialand an
e
m Eckart barrier.
- PACS numbers:03.65.Sq,82.20.Db,05.40.+j
t
a
t
s
.
t
a I. INTRODUCTION
m
- Semiclassical theories have been found extremely powerful in understanding the dynamics of complex quantum
d
mechanical systems. Special attention has been paid to theories of tunneling processes as they occur in physics,
n
o chemistry, and biology. Currently, a variety of quantum rate theories are in use explaining experimental findings
c for several situations of interest [1]. Thereby, roughly speaking, two different strategies can be distinguished. The
[ first class of approaches constructs the rate from purely thermodynamic considerations. An example is the bounce
1 or instanton method (also called ImF method), originated by Langer [2] and extended by several authors [3]. In
v essence, tunneling rates are derived from the imaginary time dynamics in the inverted potential. Other approaches
0 of this type [4,5] start from periodic orbit theory [6] in imaginary time. Tractable rate formulae are obtained with
8
the centroid method [7] leading to a semi-empirical separationof dynamical and thermal factors. These methods are
3
computationallyveryefficientandhavebeenappliedsuccessfullytosystemsasdiverseastunneling centersinmetals,
1
0 Josephson junctions, or hydrogen bonds, to name but a few. However,closer examination reveals that these theories
0 are based in one way or another on ad hoc assumptions that are not derived from first principles. For instance, the
0 ImF method postulates a relation between the decay rate and the imaginarypart of the free energy. In fact, in some
/
t casesthermodynamic methods failto predictthe correctrate,e.g. they do notreproduce the energydiffusionlimited
a
decay for very weakly damped systems at finite temperature. Moreover, these methods are designed to describe
m
incoherent decay only.
-
The second class of theories describes barrier crossing in terms of dynamical quantities. Perhaps most familiar is
d
n Yamamoto’s rate formula [8], in essence a Kubo type formula relating the rate with a flux-flux correlation function.
o As shown in [9] it is exact only for scattering problems, while in multi- or metastable systems one has to assume
c the existence of a plateau region for times long compared to typical relaxation times but small compared to decay
:
v times. Thisrestrictstheapproachtoincoherentdecay. Recentadvancesintherealtimedescriptionofbarriercrossing
i based on flux-flux correlations were made e.g. by Voth et al. [10], who incorporated analytically known dynamical
X
factors into the rate expressions. Quite recently, Pollakand coworkers[11,12]were able to improvethis idea by using
r
a a thermally symmetrized flux operator. A different way of including dynamical information in an approximate way,
favoredbyMillerandcoworkers,employssemiclassicalInitialValueRepresentation(IVR)forthequantumpropagator
[13]. Althoughquitesuccessfulathightomoderatetemperatures,wherethequantumdynamicsisgovernedbyquasi-
classical above-barrier processes, these approaches usually fail at low temperatures where deep tunneling prevails.
Finally, we mention a kind of hybrid approach, the “real time” instanton theory [14] which includes tunneling in
the real time propagator by means of instantaneous tunneling transitions. This method is restricted to multi-stable
systems in the low temperature limit [15]. Hence, although the semiclassical theory of quantum tunneling is often
regarded to be well settled, this turns out to be true only for some limiting cases. What would be desirable is a
1
semiclassicaltheory starting from firstprinciples that coversthe entire range oftemperatures as wellas coherentand
incoherent tunneling processes.
In the realm of classical physics the theory of thermally activated rates is rather firmly based. In a seminal paper
[16] Kramers determined thermal decay rates from the equation of motion for the phase space distribution function,
i.e., from the real time dynamics of the system. A corresponding treatment of tunneling in the semiclassical limit
seems not to be possible, since all real time minimal action paths connecting two sides of the barrier have energies
largerthanthebarrierenergy. Inthedynamicalapproachesdiscussedabovethesetrajectoriesdoaccountfortunneling
corrections to classicalrates [13,17],but it is usually argued[18,19]that within a semiclassicaltheory deep tunneling
can only be described by incorporating in addition imaginary time trajectories as they are used in thermodynamic
methods. The lack of a first principles semiclassical theory of tunneling is intimately connected with a fundamental
shortcoming of semiclassical real time propagators derived from a single dominant path. Since tunneling arises from
coherentinterferenceofwaves,asatisfactorysemiclassicaltheoryneedstocapturethisinterferencepatternintermsof
anappropriatefamilyofrealtime pathsratherthanmakingad-hoc modificationsofsimplesemiclassicalpropagators.
Very recently, we have proposed a theory for transport across a barrier based solely on the real time dynamics
of the density matrix and the semiclassical approximation [20]. As a notable feature, the method applies equally
well to damped and undamped systems and comprises in a unified way the entire range from thermally activated
decay to low temperature tunneling. In appropriate limits the results of other methods are recovered. Here, we
explain the technicalities of the approach and evaluate the semiclassical propagator in detail for two paradigmatic
systems, namely, a double well potential and an Eckart barrier. Thereby, we mainly concentrate on one-dimensional
models. It turns out that the corresponding theory already reveals the basic structure and that the generalizationto
multidimensional systems, though tedious in detail, is straightforwardwithin the path integral formalism [21].
The article is organized as follows. Next, in Sec. II, we outline the general semiclassical theory which is then used
in Sec. III to derive as a simple example the stationary flux across a parabolic barrier. The main part of the paper
studies the real time dynamics from high down to vanishing temperature for the cases of a double well potential
(Sec. IV) and an Eckart barrier (Sec. V). Finally, in Sec. VI we summarize the main features of the approach and
present our conclusions.
II. GENERAL THEORY
We consider a statistical ensemble of quantum mechanical particles of mass M moving in a barrier potential V(q)
atinversetemperatureβ =1/k T. We chosethe coordinateq sothatthe barriertopislocatedatq =0 andmeasure
B
energies relative to the barrier energy by putting V(0) = 0. The initial nonequilibrium state is assumed to be of
the form of an equilibrium state restricted to the left side of the barrier. Below we will invoke the semiclassical
approximation which is appropriate provided the barrier height V is by far the largest energy scale in the system.
b
The time evolution of the density matrix ρ(t)=exp( iHt/¯h)ρ(0)exp(iHt/¯h) reads in coordinate representation
−
ρ(q ,q′,t)= dq dq′G (q ,q )ρ(q ,q′,0)G (q′,q′)∗ (1)
f f i i t f i i i t f i
Z
where the real-time propagator is given by
G (q,q′)= q exp( iHt/¯h)q′ (2)
t
h | − | i
and ρ(q ,q′,0) describes the initial state. In principle, for our purpose any initial distribution that matches onto
i i
equilibrium on the left side and vanishes on the right side of the barrier top is appropriate. As long as the restricted
equilibrium state gives vanishing probability to find the particle on the right side of the barrier top, different initial
preparations lead to the same long time behavior of the density matrix. Here we put explicitly
ρ(q,q′,0)=Z−1ρ (q,q′)θ( q)θ( q′) (3)
β
− −
for convenience with the proper normalization factor Z and the equilibrium density matrix
ρ (q,q′)= q exp( βH)q′ . (4)
β
h | − | i
Now,employingthepathintegralrepresentationforexp( itH/¯h)andexp( βH),respectively,theaboveintegrand
± −
in Eq. (1) can be written as a three-fold path integral where two real time paths q(u) and q′(u) run in the interval
0 u t from q and q′ to fixed endpoints q and q′, respectively, while those former coordinates are connected by
≤ ≤ i i f f
2
animaginarytime pathq¯(σ) in the interval0 σ ¯hβ, see fig.1. The realtime paths describe the time evolutionof
≤ ≤
the system and the imaginary time path the initial state. Of course, a more complete theory wouldexplicitly include
the coupling to a heat bath environment. In fact, the generalscheme of this approachin the case of damped systems
has already been given elsewhere [21]. Much of the analysis presented below can in principle be extended to this
situation, however, only a limited number of steps can be carried out analytically due to the more complicated form
ofthe effective actionfunctionals. Here,we limit ourselvesto undamped motionwhich allowsto treatdeep tunneling
without resorting to numerical methods. This way, the guiding concepts will become more transparent.
The density matrix ρ(q,q′,t) contains all information about the nonequilibrium quantum process, in particular,
the average of the operator F = [pδ(q)+δ(q)p]/2M gives the flux out of the metastable state, i.e., in coordinate
representation
J(t)=(h¯/2iM)[∂ρ(q , q ,t)/∂q ] . (5)
f − f f qf=0
If the flux becomes quasi-stationary,J(t)=J within a certain“plateauregion”oftime, the escaperatefollows from
fl
Γ=J .
fl
WhileanexactsolutionofEq.(1)foranharmonicbarrierpotentialsisnotpossible,ahighbarriernaturallysuggests
a semiclassical approximation. In the semiclassical limit the above path integrals are dominated by minimal action
paths determined by Hamilton’s equation of motion either in the potential V(q) (for the real time propagators) or
V(q) (for the equilibrium density matrix). Eachpath contributes with anexponential factorcontaining its minimal
−
actionanda prefactorarisingfromthe Gaussianfluctuationsaboutthe minimal actionpaths. Specifically,the action
in real time reads
t
S(q,q′)= du[Mq˙2/2 V(q)] (6)
−
Z0
while its imaginary time version, the so-called Euclidian action, is given by
h¯β
S¯(q,q′)= dσ[Mq¯˙2/2+V(q¯)]. (7)
Z0
Thus, in the Gaussian semiclassics the propagator (2) is approximated as
i π
G (q,q′)= A(q,q′) exp S(q,q′) i ν (8)
t
¯h − 2
cl.paths (cid:20) (cid:21)
X p
whereA(q,q′)=[ ∂2S(q,q′)/∂q∂q′]/2πi¯handν istheMaslovindex. Throughoutthispaperwealsouseanequivalent
−
representation of the prefactor, namely,
iM ∂2W(q,q′) −1
A(q,q′)= q˙(0)q˙(t) (9)
2π¯h ∂E2
(cid:20) (cid:21)
where W(q,q′) = q′dq′′p = S(q,q′)+Et is the short action. The corresponding approximation to the equilibrium
q
density matrix (4) follows by formal analytic continuation t i¯hβ, i.e., by replacing S(q,q′) by iS¯(q,q′) in Eq. (8)
R →−
withν =0. Asaresult,theintegrandinEq.(1)iscompletelydeterminedbyclassicalmechanicsinrealandimaginary
time, respectively, and dominated by an action factor exp[ Σ(q ,q′ q ,q′)/¯h iπ(ν ν′)/2] with
− f f| i i − −
Σ(q ,q′ q ,q′)= iS(q ,q )+S¯(q ,q′)+iS(q′,q′). (10)
f f| i i − f i i i f i
Withtheapproximateintegrandathand,itisconsistenttoevaluatetheordinaryintegrationsinEq.(1)instationary
phase. The stationaryphase points are determined by minimizing Σ with respect to the initial coordinatesq ,q′, i.e.,
i i
∂Σ ∂Σ
=0, =0. (11)
∂q ∂q′
i(cid:12)(cid:12)(qf,qf′) i(cid:12)(cid:12)(qf,qf′)
(cid:12) (cid:12)
Since the endpoints q ,q′ are fixed, the res(cid:12)ulting stationary p(cid:12)hase points q (t) and q′(t) are functions of time with
f f s s
q (0) = q , q′(0) = q′. For finite t these roots are in general complex. The dominant imaginary time path q¯(σ)
s f s f s
connectsq′(t)withq (t),andthetworealtimepathsq(u)andq′(u)connectq (t)andq′(t)withq andq′,respectively.
s s s s f f
3
Hence,thesteepest–descentapproximationnaturallyprovidesamappingfromtheintegrationcontourinthecomplex
timeplaneontoaloopinthe complexcoordinatespaceconnectingtheendpoints[fig.1]. Toavoidpotentialconfusion
withothermethods,weemphasizethatthe appearanceofcomplexpaths hasnothingtodowithtunneling butrather
is merely a consequence of the stationary phase approximation and holds also for systems with no barrier at all. In
fact,itturnsoutthatthecomplexsemiclassicalreal-timetrajectoriesusedherenevercrossthebarriertop,incontrast
to paths emerging from ad hoc complexification procedures occasionally adopted to describe barrier penetration [22].
Starting from the steepest descent conditions (11) and exploiting Hamilton-Jacobi mechanics, one immediately
derives
p (0)=ip¯ (h¯β) , p′(0)=ip¯ (0) , E =E′ =E¯. (12)
s s s s
Here,p (u) [p′(u)]isthemomentumoftherealtimepathq(u) [q′(u)]withenergyE [E′]connectingq [q′]andq [q′];
s s s s f f
accordingly, p¯ (σ) denotes the momentum of the imaginary time path q¯(σ) running from q′ to q with Euclidian
s s s s
energy E¯ = p¯2/2M +V(q¯), see fig. 1. Eq. (12) can also be expressed as dΣ/dt=0 with the solution
− s
Σ(q ,q′ q ,q′)=S¯(q ,q′). (13)
f f| s s f f
Hence along the loop of steepest descent paths the full action is just given by the equilibrium action and thus
independent of time. Differentiating Eq. (13) with respect to q ,q′, one finds
f f
p (t)=ip¯ (h¯β), p′(t)=ip¯ (0) (14)
s 0 s 0
where p¯ (σ) is now the momentum of the imaginary time path q¯ (σ) connecting q′ with q in imaginary time h¯β.
0 0 f f
This path has Euclidean energy E¯ which depends on q ,q′ and h¯β but not on t. Hence, we first deduce that the
f f f
energiesinEq.(12)aregivenbyE¯ ,whichimplies energyandmomentum conservationthroughoutthe loopin fig.1.
f
Secondly, we arriveat the remarkableresult that the sequence oftime-dependent stationaryphase points q (t) [q′(t)]
s s
is itself a minimal action path starting at q (0)=q [q′(0)=q′] with energy E¯ .
s f s f f
To complete the ordinary integrations in Eq. (1) over the initial coordinates q ,q′ we transform to fluctuations
i i
y =q q andy′ =q′ q′ aboutthe stationaryphasepoints. Anexpansionofthefullaction(10)forfixedendpoints
q ,q′ia−rousnd the statiio−nasry phase points up to second order leads to Σ(q ,q′ q ,q′)=S¯(q ,q′)+δΣ(2)(y,y′) with
f f f f| i i f f
1 y
δ(2)Σ(y,y′)= (y,y′) Σ(2) (15)
2 y′
(cid:18) (cid:19)
where
Σ(2) = ΣssΣss′ (16)
Σss′Σs′s′
(cid:18) (cid:19)
is the matrix of second order derivatives, Σ =∂2Σ(q ,q′)/∂q2 etc., to be taken at q =q ,q′ =q′.
ss i i i i s i s
Inserting Eq. (3) into Eq. (1), the integrand now reduces to a product of Gaussian weight factors for deviations
from the stationary phase points and an initial state factor θ( q y)θ( q′ y′) describing deviations from thermal
− s− − s−
equilibriumatt=0. Providedthereisonlyonesemiclassicalpathforeachofthepropagators,weobtainfromEq.(1)
by virtue of Eqs. (13) and (15) the semiclassical time dependent density matrix in the form
1
ρ(q ,q′,t)= ρ (q ,q′)g(q ,q′,t). (17)
f f Z β f f f f
Here, deviations from equilibrium are described by a “form factor”
g(q ,q′,t)= 1 u(qs)dz u′(z,qs′)dz′e−(z2+z′2), (18)
f f π
Z−∞ Z−∞
where
u(q )= q Det[Σ(2)], u′(q′,z)= q′ Σs′s′ +z Σss′ (19)
s − ss 2h¯Σs′s′ s − sr 2h¯ Det[Σ(2)]
with Det[Σ(2)]=ΣssΣs′s′ (Σss′)2. In deriving Eq. (17) we invoked that Hampilton Jacobi mechanics implies [23]
−
4
A(q ,q )A¯(q ,q′)A(q′,q′) 1/2
f s s s f s =A¯(q ,q′). (20)
" Det[Σ(2)] # f f
Note that for an initial equilibrium state, formally θ() 1 in Eq. (3) so that u,u′ in Eq. (18), the form factor
· → →∞
becomes 1 and the semiclassical density matrix is in fact stationary. If there is more than one classical path one
has to sum in Eq. (17) over the contributions of all of them. Certainly, the above formulae (18) and (19) are only
applicableaslongastheGaussiansemiclassicalandstationaryphaseapproximationsarevalid,i.e.asfluctuationsare
sufficiently small. This will be seen to be no longer the case for low temperatures and/or very long times. How the
classical paths in the complex plane can then be used as a skeleton for an extended semiclassical/stationary phase
calculation will be shown below.
Intheremainingpartsofthearticleweapplythegeneralformalismtospecificbarrierpotentials. Thereby,sincewe
areparticularlyinterestedin the flux acrossthe barrier,we restrictourinvestigationto non-diagonalend-coordinates
q and q′ = q close to the barrier top. This does not mean that we may constrain ourselves to study only local
f f − f
dynamics near the barrier top. Especially for lower temperatures, the nonequilibrium state in the barrier region is
predominantly governedby global properties of the potential.
III. PARABOLIC BARRIER
The semiclassical and the stationary phase approximations are always exact for quadratic potentials. Hence, as a
simple test case we consider first a parabolic barrier
1
V(q)= Mω2q2. (21)
−2 b
Accordingly,the imaginarytime dynamics runs in aharmonicoscillatorpotential. For the minimal actionpath q¯ (σ)
0
connecting q with q in time h¯β one obtains
f f
−
q
f
q¯ (σ)= sin[ω (σ ¯hβ/2)]. (22)
0 b
sin(ω ¯hβ/2) −
b
This leads to the well-known equilibrium density matrix
1 q2
f
ρ (q , q )= exp cot(ω ¯hβ/2) (23)
β f − f 4πδb2sin(ωb¯hβ) "− b 2δb2#
p
with the relevant length scale δ = ¯h/2Mω .
b b
The real-time dynamics simply follows. The classical real-time paths q(u) and q′(u) lead to the endpoints q and
p f
q , respectively, and hence obey q(t) = q , q′(t) = q . On the other hand, the stationary phase condition (12)
f f f
− −
implies q˙(t)=iq¯˙(h¯β), q˙′(t)=iq¯˙(0) and we readily find q(u)=q¯(h¯β it+iu), i.e.,
−
q
f
q(u)= sin[ω (h¯β/2 it+iu)],
b
sin(ω ¯hβ/2) −
b
q′(u)=q(u+i¯hβ), 0 u t. (24)
≤ ≤
At time t the imaginary time path q¯ (σ) from q to q is mapped onto the path q¯ (σ) = q(i¯hβ iσ),0 σ ¯hβ
0 f f s
− − ≤ ≤
connecting q′(t)=q′(0) with q (t)=q(0) (fig. 2). The stationary phase points
s s
q
f
q (t)= sin[ω (h¯β/2 it)]
s b
sin(ω ¯hβ/2) −
b
q′(t)= q (t)∗ (25)
s − s
are as functions of t also classical paths moving away from the barrier top as t increases – q (t) to the right and
s
q′(t) to the left for q > 0. For longer times ω t 1, the stationary phase points asymptotically tend towards the
s f b ≫
limiting trajectories starting from q = 0 referred to as asymptotes henceforth. Similar to separatrices in classical
f
phase space, these asymptotes divide the complex plane in regions of negative and positive Euclidian energy: in the
sectorsincluding the realaxisclassicalrealtime motionhasE¯ <0(but E 0orE >0)while inthe remainingparts
≤
5
E¯ >0. However,in contrastto simple classicalseparatricesthe asymptotes are temperature dependent. The angel α
of the q and q′-asymptotes with the positive and negative real axis, respectively, is found as
s s
π ω ¯hβ
b
α= − . (26)
2
Now, for the nonequilibrium preparation (3) the initial coordinates q ,q′ are constrained to Re q , Re q′ 0.
i i { i} { i} ≤
Sinceq (t)andq′(t)areondifferentsidesofthebarrier,ρ(q , q ,t)gainsnonvanishingvaluesonlyduetofluctuations
s s f − f
that effectively shift q away from q and across the barrier top [see Eq. (19)]. For the parabolic barrier potential the
i s
matrix elements in (16) take the simple form
Σss =Mωb[cot(ωb¯hβ)−icoth(ωbt)], Σs′s′ =Σ∗ss,
Mω
b
Σss′ = (27)
−sin(ω ¯hβ)
b
so that the matrix Σ(2) can easily be diagonalized. One finds for the eigenvalues
1/2
λ 1
± =cot(ω ¯hβ) cot(ω ¯hβ)2 . (28)
Mω b ± b − sinh(ω t)2
b (cid:20) b (cid:21)
While,inprinciple,withEq.(18)wecannowevaluatethecompletedynamicsofthedensitymatrix,wewillconcentrate
hereonthelongtimeasymptoticsofthenonequilibriumstate. Then,intheasymptoticregionω t 1theeigenvalue
b
≫
λ tendstozeroasλ Mω exp( ω t),reflectingtheinstabilityoftheparabolicbarrier. Hence,fluctuationsaround
− − b b
∝ −
the stationary phase points with the least action increase occur in the direction of the eigenvector with eigenvalue
λ . These fluctuations are of the form y = y exp[i(α + ω ¯hβ)] and y′ = y exp(iα) so that q and q′ move
− i | i| b i | i| i i
simultaneouslyalongtheirasymptotesmeetingatthebarriertop(seefig.2). Now,insertingthematrixelements(27)
and the stationary phase points (25) into Eq. (19) and considering the limit ω t 1, the relevant form factor turns
b
≫
out to be stationary
g (q, q)= 1 i2qΩdx e−x2 (29)
fl
− √π
Z−∞
withΩ= cot(ω ¯hβ/2)/(8δ2). Thecorrespondingconstantfluxacrossthe barrierisobtainedfromEqs.(5)and(17)
b b
as
p
¯h ∂g (q , q )
fl f f
J = ρ (0,0) − (30)
fl β
2ZM ∂q
f (cid:12)qf=0
(cid:12)
(cid:12)
which leads to the well-known result (cid:12)
ω 1
b
Γ=J = . (31)
fl
4π Zsin(ω ¯hβ/2)
b
Here, Z denotes an appropriate normalization constant which cannot be derived from the pure parabolic potential.
This is not really a problem since realistic potentials always exhibit a well-behaved potential minimum and then
Z follows e.g. as the relative normalization with respect to this minimum. Note that for the quadratic potential
the results (29) and (31) are formally valid for all times ω t 1. However, it was shown in [24] that due to the
b
≫
lack of a well-behaved ground state it makes physically only sense to use the parabolic barrier in T > T where
c
ω ¯h/k T = π. For lower temperatures ω ¯hβ π large quantum fluctuations render the Gaussian approximation
b B c b
→
insufficient. Interestingly, the rate expression (31) diverges at the lower temperature T = T /2 only where the
0 c
parabolic density matrix (23) ceases to exist.
IV. DOUBLE WELL POTENTIAL
A model well-behaved for the entire range of temperatures with many applications is the bistable dynamics of a
particle moving in a double well potential
6
Mω2 q2
V(q)= b q2 1 (32)
− 2 − 2q2
(cid:20) a(cid:21)
Here, the barrier is located at q = 0, the wells at q = q , and the barrier height is V = V(q ) = (Mω2/4)q2.
± a b − a b a
This potential exhibits rich quantum dynamics, namely, incoherent hopping between the wells over a broad range of
temperatures that changes to coherentoscillations for T 0. Due to the complexity of the dynamics this is a highly
→
nontrivial problem for the semiclassical approachwhere the ratio δ /q serves as the small parameter.
b a
A. Thermal equilibrium
The Euclidian mechanics in the inverted potential V(q) can be solved exactly using Jacobian elliptic functions
−
[25]. For the general solution one obtains
q¯ (q ,σ)=Bsn[ω(B)σ φ m], 0 σ ¯hβ (33)
0 f f
− | ≤ ≤
where the boundary conditions q¯ (q ,0) = q and q¯ (q ,¯hβ) = q fix the amplitude B and phase φ . Since the
0 f f 0 f f f
−
potential is no longer purely quadratic – depending on temperature – there may be several solutions each of them
with another amplitude. In Eq. (33) the frequency is given by
B2
ω(B)=ω 1 η2, η2 = , (34)
b − 2q2
a
p
and the phase can be represented as an incomplete elliptic integral
qf/B 1
φ =F(q /B m)= dx (35)
f f
| (1 x2)(1 mx2)
Z0 − −
with the so-called modul m = η2/(1 η2). From the bounpdary condition q¯ (h¯β) = q¯ (0) and the periodicity of
0 0
− −
the Jacobian function, sn[z+2rK(m)m] = ( 1)rsn[z m],r = 1,2,3,... with K(m) = F(1,m), the amplitude B is
| − |
determined by
ω(B)h¯β =2rK(m)+[1+( 1)r]φ . (36)
f
−
Since K(m),φ >0, for fixed ω ¯hβ real solutions to this equation exist only for a finite number of integers r 0.
f b
≥
Letusbrieflydiscussthetrajectoriesq¯ (q ,σ)asthetemperatureislowered. Forhightemperaturesonlysolutionsof
0 f
Eq.(36)withr =0existcorrespondingtodirectpaths from q to q ;particularly,q¯ (0,σ)=0. Asthe temperature
f f 0
−
drops below the critical temperature
T =h¯ω /πk , (37)
c b B
i.e. ω ¯hβ > π, solutions of Eq. (36) with r = 1 arise. Then, for q = 0 the barrier top can be joined with itself also
b f
by two nonlocal paths denoted by q¯ (0,σ) oscillating in V(q) to the right and to the left with amplitudes q ,
± 1
respectively, and energy E¯ = V(q ). With further decrea−sing temperature q grows and eventually saturates a±t q
1 1 1 a
for T 0. For finite q the situation is rather similar: oscillating paths q¯ (q ,σ) exist for all q < q . These paths
f ± f f 1
→
connect q with q via a turning point at q , thus, differing from q¯ (0,σ) only by a phase shift. The described
f f 1 ±
− ±
scenariorepeatsinananalogwayatallT =T /r,r =2,3,4,...,wherercountsthenumberofturningpoints. Atzero
c
temperature all these oscillating paths reach q with the same energy E¯ =V(q ) and are then called instantons.
a a a
±
NowthatallproperEuclidiantrajectoriesareidentified,thesemiclassicalequilibriumstatefollowsreadily. Forhigh
temperatures T >T and end-coordinates q near the barrier top, ρ (q , q ) basically coincides with the parabolic
c f β f f
−
result (23) and anharmonic corrections are negligible. This situation changes drastically for temperatures near T .
c
Then,thebifurcationofnewclassicalpathsleadstolargequantumfluctuationsandonehastogobeyondtheGaussian
approximation of the fluctuation integral. Slightly below T a caustic appears for q = q . Since the region around
c f 1
T was already studied in detail elsewhere [24] we omit here this crossover region and proceed with temperatures
c
sufficiently below T where near the barrier top Gaussiansemiclassics is againapplicable. It turns out that the paths
c
newly emerging near T are stable and dominate ρ (q,q′) for all T > T > 0 while the unstable “high temperature”
c β c
path and those springing up at lower T give negligible contributions. Since for q <q all paths q¯ (q ,σ) differ only
f 1 ± f
by a phase shift, one has for the corresponding actions
7
S¯ (q , q )=S¯ (0,0)=S¯ (0,0) (38)
± f f + −
−
so that
ρ (q , q )=2[A¯(q , q )]1/2exp[ S¯ (0,0)/¯h], q <q . (39)
β f f f f + f 1
− − −
Accordingly, the matrix element ρ (q , q ) changes to a non-Gaussiandistribution with a local minimum at q =0
β f f f
and two maxima at q = q . Thereby−S¯ (0,0)< 0, so that the probability ρ (0,0) to find the particle at q = 0 is
f 1 + β
±
substantially enhanced compared to its classical value.
For T 0 it is no longer sufficient to include only the trajectories with r = 1 in the semiclassical analysis but
→
rather all other paths with r > 1 must also be taken into account. This is due to the fact that the smaller action
factors ofthese latter paths arecompensated for by zero mode phase factorsfrom the correspondingfluctuation path
integrals. Accordingly, all instanton contributions are summed up to yield e.g. for coordinates near the barrier top
8 ¯hβ∆ q q ¯hβ∆
ρ (q , q )= e−β[V(qa)+h¯ωa/2]−W¯a/h¯ cosh +cosh f a sinh . (40)
β f − f δ √2π 2 2δ2 2
a (cid:20) (cid:18) (cid:19) (cid:18) a (cid:19) (cid:18) (cid:19)(cid:21)
Here, W¯ W¯( q ,q )= ¯hβE¯+S¯( q ,q ) is the short action for an instanton from q to q . Further,
a a a a a a a
≡ − − − −
4q
∆=ω a exp( W¯ /¯h) (41)
a a
√2πδ −
a
denotes the WKB tunnel splitting with the well frequency ω =ω √2 and δ2 =h¯/2Mω .
a b a a
B. Dynamics of stationary phase points
Asincaseoftheparabolicbarrier,thestationaryreal-timepathscanbedirectlyinferedfromtheEuclidiandynamics
at t=0. From Eq. (33) and the stationary phase condition we have
q (t)=Bsn[φ iω(B)tm], q′(t)= q [( 1)r+1t] (42)
s f − | s − s −
and q¯ (σ) is mapped at time t onto q¯(σ) = q¯ [σ +i( 1)r+1t] where r follows from Eq. (36). In the sequel we
0 s 0
−
always formulate the semiclassical theory in terms of the real time paths q ,q′ that “start” at the endpoints q , q ,
s s f − f
respectively, and reach the initial points q ,q′ after time t. Since the endpoints q , q are fixed, while the most
i i f − f
relevantinitial coordinatesdepend ontime, this backwardview of the dynamics is in factmoretransparent. The real
time trajectories now start from the end coordinates we are interested in and lead to the relevant initial coordinates
that need to be integrated over with an integrand weighted according to the initial deviations from equilibrium. The
path q (t) runs in the complex coordinate plane as a periodic orbit with period t (q ) = 2K(1 m)/ω(B) (fig. 3).
s p f
−
Within one period it connects q with q via a loop crossing the real axis also after time t (q )/2 at the point
f f p f
q (q ) = q [q ,t (q )/2] q . Thus q (t) stays always on the same side of the barrier top and likewise q′(t) on the
c f s f p f ≥ a s s
othersidesothatthecomplexdynamicsofthestationaryphasepointsstartingfromq and q ,respectively,reflects
f f
−
a bounded motion in either of the potential wells.
Let us consider the stationary orbits as the temperature decreases. For high temperatures T >T , i.e. r =0, each
c
q -dependent loop carriesits own period t (q ) and energy E(q ). If q =0, t is small for T T and the real time
f p f f f p c
6 ≫
dynamics corresponds to a fast bouncing back and forth in the well. As the temperature is loweredthe period grows
while simultaneously the “width” of the loop q (q ) shrinks. In the special case q =0 the real time path reduces to
c f f
aconstantq (0,t)=0. Fortemperatures T <T the situationchangesaccordingto the appearenceofnew oscillating
s c
Euclidian paths q¯ (q ,σ) for q <q . In contrast to the high temperature case all stationary phase point paths with
± f f 1
q <q have then the same period t (q )=t (q ) and energy E(q )=E¯ , and differ only in their respective phases.
f 1 p f p 1 f 1
Specialcasesareq =0andq =q : Thepathq (q ,t),q 0runsalongtheimaginaryaxis,whiletheorbitq (q ,t)
f f 1 s f f s 1
→
degenerates to a usual well oscillation along the real axis. These properties have a direct effect on the corresponding
actions. One finds by employing Cauchy’s theorem that after each period
S[q (q ,nt ),q ]=S[q (q ,nt ),q ], q q . (43)
s f p f s 1 p 1 f 1
≤
Hence, all q (q ,t) for q < q can be seen as phase shifted copies of the specific real path q (q ,t) having the same
s f f 1 s 1
energy, period, and action increase during one period. In particular, the period t (q ) is large for T<T when q is
p 1 c 1
still small, t (q ) ln(q /q )/ω , and drops down to t (q )=2π/ω in the limit T 0. ∼
p 1 a 1 b p a a
≈ →
8
C. Nonequilibrium dynamics for high and moderately low temperatures
For t=0 the density matrix is given by the initial state (3). The semiclassical time evolution of this state follows
by inserting the proper classical paths into Eq.(17). In the sequel, we mainly focus on the long time dynamics and
are especially interested in a plateau region where the time evolution becomes quasi-stationary.
We start by addressing the question when a plateau region does exist at all. Inserting the classical paths into
g(q , q ,t)inEq.(19),the detailedanalysisrevealsthatthis functionbecomes stationarywhenthe ratioq (t)/p (t)
f f s s
−
reducestoaconstant. Sincethis willonlyoccurwithintheparabolicbarrierregion,aleastupperboundforaplateau
regionfollowsfromthetimeintervalwithinwhichaparticlestartingatatypicalpointq nearthebarriertopcontinues
b
to experience a nearly parabolic potential. This leads to t t (q ). The lower bound is obvious: it is given by the
p b
≪
transient time near the top, i.e. by 1/ω . Hence, a plateau region can be estimated to occur as long as a there is
b
clear separation of time scales between local barrier motion and global well oscillations, i.e., 1 ω t ω t (q ).
b b p b
≪ ≪
Particularly,inthetemperaturedomainwhereaplateauregionexists,theGaussian/stationaryphaseapproximations
are valid and we can actually calculate a rate. These approximations break down when the stationary phase points
move far from the barrier top and times of order t (q ) become relevant for the barrier crossing. This range will be
p b
addressed in the next section.
For high temperatures T > T the typical length scale q can be identified with δ = ¯h/2Mω . Then, the
c b b b
separationof time scales fails for very high temperatures where k T>8V meaning that the thermal energy is of the
B b p
sameorderorlargerthanthebarrierheight. Withdecreasingtemper∼aturet growssothatforω ¯hβ oforder1awide
p b
plateaurangeappearswitht ln(q /δ )/ω . Inthe correspondingdensity matrixanharmoniccorrectionsaresmall,
p a b b
≈
and we obtain approximately the parabolic result (29). To get the rate, here, the proper normalizationconstant Z is
taken as the partition function of the harmonic well oscillator
1
Z = eβVb. (44)
2sinh(ω ¯hβ/2)
a
Hence, from Eq. (30) one regains the known result
ω sinh(ω ¯hβ/2)
Γ= b a e−βVb (45)
2π sin(ω ¯hβ/2)
b
with the exponentialArrhenius factor and a characteristic¯h-dependent prefactorthat formally tends to ω /ω in the
a b
classical limit and describes the quantum enhancement of the rate as T is approached.
c
At this point we have to be very careful: a detailed analysis [21] of the full density matrix ρ(q ,q′,t), not only of
f f
its nondiagonal part, reveals that the nonequilibrium effects described by the flux state are restricted to the barrier
region only in the presence of damping, consistent with the fact that finite temperature decay rates require coupling
to a heat bath. In the absence of damping the full density matrix does not become quasistationaryand the real time
trajectories explore the strongly anharmonic range of the potential. Hence, an evaluation of the rate based upon a
supposedly quasi-stationary flux state ρ (q , q ) for the undamped case corresponds to the transition state theory
fl f f
−
result. We refer to [21] for a detailed discussion of this point.
As the temperature reaches T large quantum fluctuations occur and the impact of anharmonicities becomes sub-
c
stantial. Adetailedstudy ofthebifurcationrangearoundT isquite tediousandwasalreadypresentedin[26]. Thus,
c
we omit explicit results here and proceed with temperatures T<T where for coordinates close to the barrier top a
c
Gaussian approximation – then around the newly emerging pa∼ths with amplitudes q – is again appropriate. As
1
±
discussed above all real time paths with q < q have now the same oscillation period t (q ) ln(q /q )/ω . One
f 1 p 1 a 1 b
≈
observes that even though they are influenced by the anharmonicity of the potential via the Euclidian amplitude q ,
1
their time evolutionfor T<T is still dominated by parabolic properties. Then, a somewhat lengthy algebra leads to
c
the quasi-stationary densit∼y matrix
1
(+) (−)
ρ (q , q )= ρ (q , q ) g (q , q )+g (q , q ) (46)
fl f − f 2 β f − f fl f − f fl f − f
h i
whereg(±) describethe contributionsfromeachofthe twooscillatingEuclidianpaths. Note thatdue to symmetryin
fl
ρ (q , q ), these contributions are identical and just lead to a factor of 2. In the temperature domain studied here,
β f f
−
q can be gained analytically from Eq. (36) as
1
1/2
2q π2
a
q = 1 . (47)
1 √3 − ω2¯h2β2!
b
9
Accordingly, one finds with some algebra for the thermal distribution
2
1 ω ¯hβq2 π2
ρ (0,0)= exp b a 1 . (48)
β 2πδb2|sin(ωb¯hβ)| 12δb2 − ωb2¯h2β2!
p
The form factor now has two contributions of the form (18) and for κ = ∂ln(q )/∂(ω ¯hβ) 1 the corresponding
1 1 b
≫
integration boundaries read
q +iq
u(±)(q )= ± 1 f, u′(±)(q ,z)=κ u(±)(q ) z . (49)
f f 1 f
4δ √κ −
b 1
h i
Here,weemployedthatnearT the derivativeΣ is dominatedby∂2S¯/∂q2 whichis proportionalto∂E /∂(ω ¯hβ)
c ss i 1 b ∝
κ . This way, using the normalization (44), the result for the rate is
1
ω sinh(ω ¯hβ/2) π2
Γ= b a ω ¯hβ e−βVb. (50)
b
2π 2 sin(ωb¯hβ) s − ωb¯hβ
| |
This expressionis valid for temperaturespT <T where still κ 1, a region which can be estimated as T somewhat
c 1
≫
larger than T /2. There are two interesting observations to mention: first, the exponentially large term in the
c
thermal distribution (48) – a consequence of the new Euclidian paths – is exactly canceled by a corresponding
term which arises from the derivative of the form factor. This way, the rate is still dominated by the characteristic
“Arrheniusfactor”. Second,inthelimitT T theabove-T formula(45)andthebelow-T result(50)bothapproach
c c c
→
Γ =(ω /2π)sinh(ω ¯hβ/2)exp( βV ), however, the derivatives ∂Γ/∂T are different. This discontinuity in the slope
c b a b
−
of the temperature dependent rate is removed by the full semiclassical theory [26] which takes the non-Gaussian
fluctuations near T into account and leads to a smooth changeoverbetween the rate formulas (45) and (50).
c
With further decreasing temperature the amplitude q tends to saturate at q so that κ 0 and the above
1 a 1
→
rate expression is no longer applicable. Furthermore, the plateau region shrinks and eventually vanishes so that the
assumption of a quasistationary flux state becomes inadequate even in the sense of a transition state theory limit
of a weak damping theory. We note that in case of finite damping a meaningful rate can be found for much lower
temperatures,thendescribingincoherentquantumtunneling. Toinvestigatethetimedependenceofρ(q , q ,t)with
f f
−
no damping in the limit of deep tunneling, we consider the case T =0 in the next section.
D. Nonequilibrium dynamics for zero temperature
Any Gaussian semiclassics to the real time propagator is expected to break down for very low temperatures and
very long times where quantum tunneling comes into play. In fact, as yet no satisfactory semiclassicalprocedure was
found to describe deep tunneling in the time domain. The crucial question thereby is: how can classical trajectories
that either oscillate in one of the potential wells, here with energy E < 0, or move over the barrier, here with
E > 0, produce exponentially small contributions to the semiclassical propagator which originate from quantum
states connectingthe twowells underthe barrier. Here,we presenta mechanismthatis basedonly uponthe complex
plane mechanics discussed above and avoids any additional ad hoc insertion of “barrier paths”. Since the complex
planedynamicsbehavesastheusualclassicalrealtimemechanics,pathswithE <0nevercrossthebarrier. However,
a full semiclassicaltreatment needs to account for the dominant fluctuations about the semiclassicalpaths. Now, for
T <T there is a whole family of loop-like orbits in the complex plane; all with the same energy, period, and action
c
increaseafteroneperioddifferingfromeachotheronlybytheir respectivephases,i.e.bytheircrossingpointsq q
f 1
≤
with the realaxis. It turns outthat eachtime these orbitspass their end-coordinateq there areother trajectoriesof
f
this family arbitrarily close in phase space (see fig. 4). The role of quantum mechanics then is to induce transitions
betweenthese orbits viasmallfluctuations. For sufficiently longtimes a path startingata certainq nearthe barrier
f
top may successively slip down to an orbit with another phase q′, eventually reach the stable regions around q ,
f ± a
and fluctuate in the long time limit between these regions. That this scenario actually describes the low temperature
coherent tunneling dynamics has been discussed briefly in [20] and will be described in some detail in the sequel.
For T =0 the amplitude of the Euclidian time paths is q =q . Thus, all stationary paths q (q ,t), q′(q ,t) have
1 a s f s f
energyE =V(q )andperiodt t (q )=2π/ω . Then,theEuclidianactionS¯(q ,q′)suppressesenergyfluctuations
a a ≡ p a a i i
around E = V(q ) exponentially, so that classical paths running in time t from q = q and q′ = q′ to q and q ,
respectively,i.e.waithE =V(q ),arenegligible. Further,studyingtheshortaction,iW6 (q,sq′)= iq′6 dqsp Sf(q,q′)+−Eft,
6 a q ≡
R
10
