Higher-order solutions to non-Markovian quantum dynamics via hierarchical functional derivative

0.36 MB·English

by  Da-Wei Luo

Higher-order solutions to non-Markovian quantum dynamics via hierarchical
functional derivative
Da-Wei Luo,1 Chi-Hang Lam,2 Lian-Ao Wu,3,4 Ting Yu,1,5 Hai-Qing Lin,1 and J. Q. You1
1Beijing Computational Science Research Center, Beijing 100094, China
2Department of Applied Physics, Hong Kong Polytechnic University, Hung Hom, Hong Kong, China
3Department of Theoretical Physics and History of Science,
The Basque Country University (EHU/UPV), PO Box 644, 48080 Bilbao, Spain
4Ikerbasque, Basque Foundation for Science, 48011 Bilbao, Spain
5Center for Controlled Quantum Systems and Department of Physics and Engineering Physics,
Stevens Institute of Technology, Hoboken, New Jersey 07030, USA
(Dated: September9, 2015)
5
1 Solving realistic quantum systems coupled to an environment is a challenging task. Here we
0 developahierarchicalfunctionalderivative(HFD)approachforefficientlysolvingthenon-Markovian
2 quantumtrajectoriesofanopenquantumsystemembeddedinabosonicbath. Anexplicitexpression
for arbitrary order HFD equation is derived systematically. Moreover, it is found that for an
p
analytically solvable model, this hierarchical equation naturally terminates at a given order and
e
S thus becomes exactly solvable. This HFD approach provides a systematic method to study the
non-Markovian quantumdynamics of an open system coupled to a bosonic environment.

8
PACSnumbers: 03.65.Yz,42.50.Lc

]
h
p
I. INTRODUCTION is almost impossible to find any useful analytical solu-
-
t tions. Therefore, one has to develop numerical methods
n
to study the quantum dynamics of the open systems.
a The theory of open quantum systems [1] has received
u great interest because environment-induced effects are However, the application of the non-Markovian QSD is
q greatly hindered unless one can cast it to a numerically
importantinawiderangeofresearchtopicssuchasquan-
[ implementable time-local form. Recently, two hierarchi-
tum information [2] and quantum optics [3]. There were
calapproacheshave been proposed;one is the stochastic
2 considerable studies involving environments modeled by
v either bosonic or fermionic baths (see, e.g., Refs. [1, 4– differential equation (SDE) method [9] based on a func-
9 13]), as well as structured environments such as spin- tional expansion of a system operator, and the other is
8 anapproachusingahierarchyofpurestates(HOPS)[10]
chain baths [14]. Conventionally,Markovapproximation
2 to calculate a related functional derivative. Apart from
has been extensively used because of its simplicity. In-
6
these two approaches, an earlier attempt at a perturba-
deed, this is valid only when memory effects of the en-
0
tive solution of the non-Markovian QSD equation was
. vironment are negligible. However, this approximation
1 based on hierarchical functional derivative (HFD) [7].
becomes invalid when the system-environment coupling
0 However,because of its apparent complexity, the hierar-
is strong or when the environment is structured [1, 15].
5
chical equations presented there were implemented only
1 Therefore, non-Markovianenvironments have to be con-
uptothesecondorderwherethehigher-ordertermswere
: sideredinexplainingnewexperimentaladvancesinquan-
v approximated by a simplified operator. Here we develop
tum optics. Also, they must be considered in quantum
Xi information manipulations in which the environmental a systematic andefficient higher-orderHFD approachto
solvingthenon-Markovianquantumdynamicsofanopen
r memory is utilized to control the entanglement dynam-
a ics [16]. Therefore, it is vital to have a non-Markovian system coupled to a bosonic environment. Remarkably,
a compact explicit expression for an arbitrary order hi-
descriptionofthe system’s dynamics under the influence
erarchicalequationis derived. Moreover,itis foundthat
of the memory effects and the backactionof the environ-
for ananalytically solvablemodel, the hierarchicalequa-
ment. Actually,thishaslongbeenachallengingtaskand
tion naturally terminates at a givenorder,so it becomes
many theoretical approaches have been developed (see,
exactlysolvable. Thus,our method providesnot onlyan
e.g., Refs. [1, 4–10, 17–21]). Among these approaches,
approachforefficiently solvingthe non-Markovianquan-
the non-Markovian quantum state diffusion (QSD) [4–
tum dynamics of an open system, but also a systematic
6] has been proven to be a powerful tool to study the
method for the exact solution of analytically tractable
quantum dynamics of the system and exact analytical
open systems.
resultswerederivedformanyinterestingsystemssuchas
dissipative multi-level atoms [8] and quantum Brownian Thepaperisorganizedasfollows. InSec.II,wepresent
motion [4, 22] which was also analytically solved via a our HFD method for studying the non-MarkovianQSD.
path-integral approach [17]. Quantum continuous mea- Then, we show the relationship of the HFD method to
surement employing the QSD technique was also stud- therecentlydevelopedSDEandHOPSmethodinSec.III
ied [23, 24]. andSec. IV,respectively. Finally,discussionandconclu-
For most realistic open quantum system problems, it sion are given in Sec. V. Moreover, in Appendix A, we
2
present a detailed derivation of the HFD equation for the O operator. It is difficult to analytically obtain the
a general bath correlation function. The mathematical explicitexpressionoftheOoperatorexceptforafewsim-
relationship between the SDE approach and our HFD ple models such as dissipative multi-level atoms [8], the
method is explicitly demostrated in Appendix B. quantum Brownian motion [4] and dissipative multiple
qubits [25]. Owing to the importance of open systems,
efforts [7, 9, 10] have been devoted to develop numeri-
II. THE HFD METHOD FOR cal methods for solving Eq. (3). In Ref. [7], an approach
NON-MARKOVIAN QSD was proposed to calculate the functional derivative by
defining a set of Q operators as
k
We study a generic open system with the Hamilto- t δ
nian [4, 5] Q (t,z˜∗)= α(t,s) Q (t,z˜∗)ds, (4)
k δz˜∗ k−1
Z0 s
H =H + g Lb† +g∗L†b + ω b†b , (1) where Q (t,z˜∗) = O¯(t,z˜∗). Because of the complexity
s k k k k k k k 0
Xk (cid:16) (cid:17) Xk involved, an explicit hierarchical equation of motion for
this set of operators was obtained in Ref. [7] only up
where Hs is the Hamiltonianof the system ofinterest, L to the second-order Q operator while the higher-order
2
is the coupling operator called Lindblad operator, bk is termswereapproximatedbyasimplifiedoperator. Below
the kth mode annihilation operator of the bosonic bath we give an explicit hierarchicalequation of motion up to
withafrequencyωk andgk denotesthecouplingstrength any high orders and also applicable to an arbitrary bath
between the system and the bosonic bath. The bath correlationfunction α(t,s). For this purpose, we include
state can be specified by a set of complex numbers {zk} the time derivatives of α(t,s) and define
labeling the coherent state of all bath modes and the ef-
fect of the bathis characterizedby the zero-temperature Q(j)(t,z∗)=P(jk)P(jk−1)...P(j1)
k
bath correlation function α(t,s) = k|gk|2e−iωk(t−s). t
Defining zt∗ ≡ −i kgk∗zk∗eiωkt, one can interpret zk as × ds0α(j0)(t,s0)O(t,s0,z∗), (5)
a Gaussian random variable and zt∗ bPecomes a Gaussian Z0
P
process with its statistical mean given by the bath cor- where j= j e ≡(j ,j ,...,j ) with e being the ith
i i i 0 1 k i
r|tψuelzma∗t(istot)nait=efu|nΨhzct∗toi|toΨ(ntt)oiαt(o(ttf),its,h)eob=totathaihnlzestdyzss∗btiyei.mprToonhjteeoctwtihnaegvebtahftuehnqcsuttaiaontne- u∂∂ntjjiαt(vte,cst)o.rPI,tPis(jd)er=iveRd0t(dsseeα(Aj)p(pt,esn)dδizδxs∗,Aa)ntdhaαt(tjh)(ist,sse)to=f
operators satisfy
|zi, is called a quantum trajectory and obeys a linear
QSD equation k
∂
Q(j)(t,z∗)= α(ji)(0) L,Q(D(j,i))(t,z∗)
∂ ∂t k k−1
|ψz∗(t)i= −iHs+Lz∗−L†O¯(t,z) |ψz∗(t)i, (2) Xi=1 h i
∂t
k
(cid:2) (cid:3) + Q(j+ei)(t,z∗)−L†Q(0,j)(t,z∗)
where O is an operator defined by the functional deriva- k k+1
tive δzδ∗|ψz∗(t)i = O(t,s,z∗)|ψz∗(t)i and O¯(t,z∗) = Xi=0
tα(t,ss)O(t,s,z∗)ds [4]. Non-Markovian master equa- + −iHsys+Lzt∗,Qk(j)(t,z∗) +δ0,kα(j)(0)L
0
tions can also be obtained using this approach [6, 22]. hk i
FRor this linear QSD, |ψz∗(t)i is an unnormalized wave − L†Q(0,ci)(t,z∗),Q(j0,¯ci)(t,z∗) , (6)
i k−i
function, so it can be of different orders of magnitude at
Xi=0Xci h i
different evolution times. This gives rise to an inefficient
Monte Carlo sampling in numerical calculations. Thus, where Q(k0)(t,z˜∗) = Qk(t,z˜∗), α(j)(0) = α(j)(t,t),
oneutilizesanon-linearQSD[4]forthenormalizedstate D(j,i)≡(j0,...ji−1,ji+1,...,jk) excludes ji in the vec-
|ψ˜z˜∗i=|ψz∗(t)i/|||ψz∗(t)i||: tor j, and ci is the sum of k!/[i!(k−i)!] terms which
includeallpossiblecasesofchoosingielementsfromak-
P
d|ψd˜tz˜∗i = −iHs+∆t(L)z˜t∗−∆t(L†)O¯(t,z˜∗) dsiiomnesn[s6i]o,nwael vceacntoarls(oj1i,m..m.e,djkia).teUlysipnrgovfuentchtaiotnQal(kje)(xtp,azn∗)-
(cid:2)+h∆t(L†)O¯(t,z˜∗)it |ψ˜z˜∗(t)i, (3) operators have such a symmetry that Q(kj0,j1,...,jk)(t,z∗)
coincide for all permutations of (j ,j ,...,j ). This
1 2 k
which is derived via the Girsanov(cid:3)transformation z˜∗ = enables us to greatly reduce the number of equations
t
z∗+ tα∗(t,s)hL†i ds. Here ∆ (L)≡L−hLi , and the that should be solved, because we only need to calcu-
retduce0d density opserator ρ (t)t≡ Tr |Ψ ihtΨ | can late those Q(j)(t,z∗) operators with normal-ordered j
R s env tot tot k i
beobtainedfromtheensembleaverageofthenormalized (j ≤j ≤...≤j ).
1 2 k
quantum trajectories as ρs(t)=hh|ψ˜z˜∗(t)ihψ˜z˜(t)|ii. It is alsoworthpointing out that since the QSD equa-
Thenon-MarkovianQSDequation(3)isexact,butthe tionfor a finite-temperature bathwith self-adjointLind-
key challenge in solving it relies on the determination of bladoperatorL(i.e.,L=L†)isformallythesameasthe
3
zero-temperature one with a different finite-temperature
bath correlation function α(t,s) [4], so Eq. (6) also ap-
pliestothatfinite-temperaturecase. Notethatthefinite-
temperaturecasewithaself-adjointLcoversalargenum-
ber of open-system problems, including the important
spin-boson model without the rotating-waveapproxima-
tion (see, e.g., [9]).
When modeling the structure of the bosonic bath,
an important and widely used choice is the Lorentzian
spectrum, which corresponds to an environmental noise
z of the Ornstein-Uhlenbeck type with autocorrelation
t
α(t,s) = Γγexp(−γ|t − s|)/2. This kind of bosonic
bath has been used to describe many interesting prob-
lems [1, 3], and it is easy to observe a non-Markovianto
Markovian crossover by just increasing γ. In addition,
this can greatly simplify the hierarchicalequations since
α(j)(t,s) = (−γ)jα(t,s) and Q(kj) = (−γ)JsQk, where
J = j . In this case, by utilizing our general re-
s i i
sult (6), a compact hierarchical equation for the Q op-
k
P
erators given by Eq. (4) can be explicitly written as
∂
Q (t,z˜∗)=kα(0)[L,Q (t,z˜∗)]+δ α(0)L
k k−1 0,k
∂t
−(k+1)γQ (t,z˜∗)+[−iH +Lz˜∗,Q (t,z˜∗)] (7)
k s t k
k
−L†Q (t,z˜∗)− Ck L†Q (t,z˜∗),Q (t,z˜∗) ,
k+1 i i k−i
i=0
X (cid:2) (cid:3)
FIG.1. (Coloronline)(a)Ensembleaveragevaluesofthean-
where Cik ≡ k!/[i!(k − i)!] is the binomial coefficient gular momentum hJxi, hJyi, and hJzi versus time t. (b) En-
and the initial conditions are Qk(0,z˜∗) = 0 for k ≥ sembleaveragevaluesofthetracenormofQ0(t,z˜∗),Q1(t,z˜∗),
0. This set of hierarchical equations can be numeri- andQk(t,z˜∗),with k≥2, asafunction oft. Inboth (a)and
cally solved perturbatively if terminated at order N +1 (b), we use 1000 noise realizations and γ =Γ = ω =1. The
by putting either Q (t,z˜∗) = 0 or Q (t,z˜∗) = solid curves correspond to the analytic solutions; the solid,
N+1 N+1
tdsα(t,s)[L,Q (t,z˜∗)] [7]. open circles and triangles correspond to the results obtained
0 N using ourHFD approach.
Now we also explore the use of this approachas a sys-
R
tematic method for models with exact solutions. For
open systems that are known to be exactly solvable by
methodtodealwithanopensystemwithunknownprop-
QSD, the number of expansion terms in
erties by just implementing the hierarchical equation; if
t the hierarchical equation has a natural termination at a
O¯(t,z˜∗)=O¯(0)(t)+ O¯(1)(t,v1)z˜v∗1dv1 givenorder,theconsideredmodelisanalyticallysolvable
Z0 and our results will be essentially accurate.
t t
+ O¯(2)(t,v ,v )z˜∗ z˜∗ dv dv +.... (8) Asanillustrativeexample,weconsiderananalytically
Z0 Z0 1 2 v1 v2 1 2 solvable three-level system [8], with Hsys = ωJz, and
L = J . The functional expansion of its O operator
−
must be finite [6]. Thus, O¯(n) ≡0 for all n larger than a is only up to the first-order term, i.e., N = 1. Thus,
c
givenfinite integer Nc. In this case, we can readily show only Q0(t,z˜∗) and Q1(t) are involved in the hierarchical
that equation, and Q ≡ 0 for k ≥ 2. We numerically solve
k
the hierarchical equation (7) up to order N = 10. From
t t
Q (t)=N ! ... α(t,v )...α(t,v ) Fig. 1(a), it can be seen that the numerically calculated
Nc c 1 Nc
Z0 Z0 ensemble averages hJii, i = x, y and z, agree well with
×O¯(Nc)(t,v1,...vNc)dv1...dvNc, (9) theexactlysolvedresults. Also,thetracenorms||Qk||≡
Tr Q†Q are calculated in Fig. 1(b), which shows
whichisindependentofthenoisez˜∗. Therefore,itsfunc- k k
tionalderivativewithrespecttoz˜∗tiszero. FromEq.(4), tha(cid:18)tqfor k ≥(cid:19)2, the Q operators remain zero and the
s k
itfollowsthatQ =0forallk ≥N +1,sothatthehier- hierarchical equation naturally terminates at order k =
k c
archicalequation naturally terminates and the approach 2 (i.e., N = 1), in full consistency with the analytical
c
becomes exact. This provides a useful and systematic derivations.
4
proximation (RWA) was studied using this approach by
consideringanOrnstein-Uhlenbecknoise(see[9]). Thisis
a typicalexample where the explicit form of the O oper-
ator is unknownand the hierarchicalapproachcan serve
as a powerful numerical tool. In Fig. 2 (a), we display
the expectation value of the angular momentum hσ i as
z
a function of time using both the SDE and our HFD ap-
proaches. HerewealsouseanOrnstein-Uhlenbecknoise,
so as to directly compare with the SDE results. An ex-
cellentagreementisreachedforthenumericalresultsob-
tainedbythem. Asamatteroffact,these twohierarchi-
calmethodsareverycloselyrelatedmathematically,and
we provein Appendix B that the operatorsQ andQ(n)
k m
are related by
N
n!
Q (t,z˜∗)= Q(n)(t,z˜∗). (11)
k (n−k)! k
n=k
X
The key advantage of the HFD method over the SDE
approachisthattheHFDequation(7)effectivelygroups
theQ(n)operatoraccordingtoEq.(11)andsumsuptheir
m
contributions. This makes the HFD method much more
efficient thanthe SDE approach. In fact, for a giventer-
minationorderN,theSDEapproachneedstosimultane-
ouslysolvecoupledequationsforQ(n) wheren+m≤N,
m
resultinginatotalnumberof 1(N+2)(N+1)equations.
FIG.2. (Coloronline)(a)Ensembleaveragehσziversustime 2
t for γ =0.2, 0.4 and 0.8, obtained using the SDE approach On the other hand, the HFD approach greatly reduces
(solid curves) and the HFD method (solid, open circles and the total number of equations to N + 1. Figure 2(b)
triangles). Here N = 30, Γγ = 0.2, and 1000 noise realiza- shows the total number of coupled differential equations
tions are used. (b) Total number of equations versus termi- that shouldbe solvedin both the HFD and the SDE ap-
nation order N, where the curve marked by triangles (solid proaches. The very high efficiency of the HFD method
circles)correspondstothenumberofcoupledequationsofmo- over the SDE approach is clearly seen when increasing
tion thatshould besimultaneously solvedin theSDE(HFD) N.
approach.
III. THE RELATIONSHIP TO THE SDE
METHOD IV. THE RELATIONSHIP TO THE HOPS
METHOD
Recently, a numerical approach [9] was developed to
solve the non-Markovian QSD equation via a set of Previously, almost all QSD approaches are focused on
stochasticdifferentialequations(SDE).Akeypartofthe howtocalculatethefunctionalderivativeassociatedwith
SDE formulation is the introduction of a Q operator the O operator. Recently, a hierarchy of pure states
(HOPS)approach[10]wasdevelopedasanumericaltool
t t
which, instead of using the O operator, introduces a set
Q(n)(t,z˜∗)= ... α(t−v )...α(t−v )z˜∗
m 1 m vm+1 of pure states
Z0 Z0
...z˜∗ O¯(n)(t,v ,...,v )dv ...dv , (10)
vn 1 n 1 n t δ
|ψ (t)i= α(t,s) ds|ψ (t)i, (12)
fwuhnecrteionO¯a(ln)e(xt,pva1n,s.i.o.n,vtne)rmcorinresEpqo.nd(s8)toantdheO¯n(tt,hz˜-o∗)rde=r k Z0 δz˜s∗ k−1
tribn=ut0eQt(0onO)¯(t(,t,z˜z˜∗∗).)bFuotrfmorm6=a0s,eQto(mnf)hideoranrcohticdairleecqtulyatcioonns- wchhicearel e|ψqu0(att)iion≡s|oψf(mt)io.tiIonntwhaissafopupnrodafcohr,|aψse(tt)oi.f hItiseraadr--
wPith Q(n) and need to be solved simultaneously. Within vantageisthatthe hierarchicalequationsdekalwithstate
0
thisframework,onecancalculatethequantumtrajectory vectors of size dim(H )×1 rather than the operators of
s
up to an arbitrary order of the environmental noise. sizedim(H )×dim(H ),wheredim(H )isthedimension
s s s
The spin-boson model without the rotating-wave ap- of the system’s Hilbert space. From the definition, it is
5
easy to see that naturally terminates at a given order and becomes ex-
actly solvable for an analytically solvable model, it pro-
|ψ1(t)i=Q0(t,z˜∗)|ψ0(t)i, videsasystematicperturbationforagenericopensystem
|ψ (t)i=Q (t,z˜∗)|ψ (t)i+Q (t,z˜∗)|ψ (t)i, irrespectiveoftheexistenceofthetime-localO operator.
2 1 0 0 1
OurHFDmethodappliestoanarbitrarybathcorrelation
|ψ (t)i=Q (t,z˜∗)|ψ (t)i+2Q (t,z˜∗)|ψ (t)i
3 2 0 1 1 function. Also,itcanbenaturallyextendedtothecaseof
+Q (t,z˜∗)|ψ (t)i, ......,
0 2 afinite-temperaturebathwhenthe Lindbladoperatorin
the interaction Hamiltonian is self-adjoint, including the
and in general,
important spin-boson model without the rotating-wave
approximation.
k−1
|ψ i= Ck−1Q (t,z˜∗)|ψ i. (13)
k i i k−i−1
ACKNOWLEDGMENTS
i=0
X
Therefore, it is also possible to formulate the HOPS ap- This work is supported by the National Natural Sci-
proach using the Qk operators in Eq. (4). As shown enceFoundationofChinaNo.91421102andtheNational
above, for an exactly solvable model, Qi(t,z˜∗) = 0 for Basic Research Program of China No. 2014CB921401.
all i ≥ Nc + 1, so that the hierarchical equations for C.H.L. is supported by HK GRF (Project No. 501213).
Qi(t,z˜∗)terminateatthe Ncthorder. Itisinterestingto L.W. is supported by the Spanish MICINN (Project
note that there are infinite number of hierarchical equa- No. FIS2012-36673-C03-03). T.Y. is supported by the
tions for |ψki in the HOPS approach, but only Nc +1 NSF PHY-0925174and DOD/AF/AFOSR No. FA9550-
nonzero operators Q0(t,z˜∗),Q1(t,z˜∗),...,QNc(t,z˜∗) are 12-1-0001. L.W. also acknowledges CSRC, Beijing for
neededtoexpressallofthese|ψki. Thisreflectstheexact the hospitality during his visit to it.
solvability of the considered model.
Appendix A: Derivation for the HFD equation with
V. DISCUSSION AND CONCLUSION a general bath correlation function
Solving non-Markoviandynamics of an open quantum In this appendix, we derive the HFD equation [i.e.,
system has long been a challenge. The conventional Eq. (6)] for an open system with a general bath correla-
master equation for the system’s reduced density ma- tion function α(t,s). Define
trix is driven by a super-operator K, which is given by
the perturbation expansion with respect to the strength Q(j)(t,z∗)=P(jk)P(jk−1)...P(j1)
k
of system-bath interaction [1]. Whereas the expansion t
could be done manually up to the orders beyond 2, it × ds α(j0)(t,s )O(t,s ,z∗), (A1)
0 0 0
is difficult to find an automatic numerical way to obtain Z0
the higher-order K operator. As such, the conventional
where j= j e ≡(j ,j ,...,j ) with e being the ith
master equation is used in the weakly-coupled scenario. i i i 0 1 k i
unit vector, P(j) = tdsα(j)(t,s) δ , and α(j)(t,s) =
Likewise,thequantumstatediffusionequationdrivenby P 0 δz∗
s
an O operator (which plays a role similar to K) encoun- ∂∂tjjα(t,s). At the zeroRth order, using the quantum state
ters the same problem. As a breakthrough, this work diffusion (QSD) equation
developsasystematicandefficienthigher-orderHFDap-
proachfor solvingthe non-Markovianquantumtrajecto- O˙(t,s,z∗)]= −iH +Lz∗−L†Q(0)(t,z∗),O(t,s,z∗)
sys t 0
riesofanopensystemcoupledtoabosonicenvironment.
A compact explicit expression for an arbitrary order hi- h δQ(0)(t,z∗) i
−L† 0 , (A2)
erarchical equation of motion is derived and it can be δz∗
s
efficiently implemented numerically. As a distinctive ad-
vantage of this method, while this hierarchical equation we have
∂ ∂ t
Q(j0)(t,z∗)= ds α(j0)(t,s )O(t,s ,z∗)
∂t 0 ∂t 0 0 0
(cid:20)Z0 (cid:21)
t ∂ t
=α(j0)(0)O(t,t,z∗)+ ds α(j0)(t,s ) O(t,s ,z∗)+ ds α(j0)(t,s)O˙(t,s ,z∗)
0 0 0 0 0
∂t
Z0 (cid:16) (cid:17) Z0
=α(j0)(0)L+Q(j0+1)(t,z∗)−L†Q(0,j0)(t,z∗)+ −iH +Lz∗−L†Q(0)(t,z∗),Q(j0)(t,z∗) . (A3)
0 1 sys t 0 0
h i
6
For k ≥1, it can be seen that i elements from (j ,...,j ) and its complement ¯c that
1 k i
contains the remaining elements. For example, for
(1,2,3,4,5), c can be (1,4) and then¯c =(2,3,5).
∂ ∂ t δ 2 2
P(j) = dsα(j)(t,s) We thus have
∂t ∂t δz∗
(cid:20)Z0 s(cid:21)
=α(j)(0) δ + tdsα(j+1)(t,s) δ ∂∂tQk(j)(t,z∗)= ∂∂t P(jk)P(jk−1)...P(j1)Q0(j0)(t,z∗)
δz∗ δz∗
=α(j)(0)δδzt∗ +ZP0(j+1). s (A4) = k hP(jk)P(jk−1)...∂P∂(tji) ...P(j1)Qi(0j0)
t i=1
X
∂Q(j0)(t,z∗)
Here we define D(j,i) ≡ (j ,...j ,j ,...,j ), which +P(jk)P(jk−1)...P(j1) 0 .
0 i−1 i+1 k ∂t
excludes ji in the vector j. For any n-dimensional (A5)
vector v = (v ,...,v ) and integer x, the notation
1 n
(x,v) represents a new (n+1)-dimensional vector v′ = UsingEq.(A4)forthe firsttermontheRHSofEq.(A5)
(x,v ,...,v ). Wealsointroduceavectorc thatchooses and Eq. (A3) for the second term, we have
1 n i
k
∂ δ
Q(j)(t,z∗)= P(jk)P(jk−1)...P(ji+1)α(ji)(0) P(ji−1)...P(j1)Q(j0)(t,z∗)
∂t k δz∗ 0
i=1 t
X
k
+ P(jk)P(jk−1)...P(ji+1)...P(j1)Q(j0)(t,z∗)
0
i=1
X
+Q(j0+1,j1,...,jk)(t,z∗)+ −iH +Lz∗,Q(j0,j1,...,jk)(t,z∗) −L†Q(0,j0,j1,...,jk)(t,z∗)
k s t k k+1
−P(jk)P(jk−1)...P(j1) Lh†Q(0)(t,z∗),Q(j0)(t,z∗) . i (A6)
0 0
h i
Note that
P(jk)P(jk−1)...P(ji+1)...P(j1)Q(j0)(t,z∗)=Q(j+ei)(t,z∗),
0 k
δ
P(jk)P(jk−1)...P(ji+1)α(ji)(0) P(ji−1)...P(j1)Q(j0)(t,z∗)=α(ji)(0) L,Q(D(j,i))(t,z∗) . (A7)
δz∗ 0 k−1
t
h i
Then, we finally have
k k
∂
Q(j)(t,z∗)= α(ji)(0) L,Q(D(j,i))(t,z∗) + Q(j+ei)(t,z∗)+Q(j+e0)(t,z∗)−L†Q(0,j)(t,z∗)
∂t k k−1 k k k+1
Xi=1 h i Xi=1
+ −iH +Lz∗,Q(j)(t,z∗) −P(jk)P(jk−1)...P(j1) L†Q(0)(t,z∗),Q(j0)(t,z∗)
sys t k 0 0
kh i k h i
= α(ji)(0) L,Q(D(j,i))(t,z∗) + Q(j+ei)(t,z∗)−L†Q(0,j)(t,z∗)
k−1 k k+1
Xi=1 h i Xi=0
k
+ −iH +Lz∗,Q(j)(t,z∗) − L†Q(0,ci)(t,z∗),Q(j0,¯ci)(t,z∗) , (A8)
sys t k i k−i
h i Xi=0Xci h i
where is the sum of Ck ≡ k!/[i!(k−i)!] terms that wherewehavechosenioperatorsinP(jk)P(jk−1)...P(j1)
ci i
include all possible cases of choosing i elements from a to act on Q(0)(t,z∗), and the remaining (k − i)
P 0
k-dimensional vector (j1,...,jk). This comes from the P(jx) operators act on Q(j0)(t,z∗). For example,
term 0
for P(j5)P(j4)P(j3)P(j2)P(j1), one possible case with
i = 2 is that P(j3)P(j1) acts on Q(0)(t,z∗) to give
0
P(jk)P(jk−1)...P(j1) L†Q(00)(t,z∗),Q(0j0)(t,z∗) , (A9) Q2(0,j1,j3)(t,z∗),andP(j5)P(j4)P(j2)actsonQ0(j0)(t,z∗)to
h i
7
give Q(j0,j2,j4,j5)(t,z∗). For the Q(j) operator, termina- where N is a suitably chosen integer.
3 k
tionoftheindexkisdeterminedbythefunctionalderiva- As an example, we use an exactly solvable model to
tive of the O operator,whereas termination of the index show what the HFD equations look like in the pres-
j is determined by the bath correlation function α(t,s). ence of a general bath correlation function. Consider
In numerical calculations, this set of equations can be a three-level atom in a multi-mode bosonic bath [8],
terminated by imposing Q(kj) ≡ 0 for all k+ iji > N, where Hsys = ωJz, and L = J−. The HFD equation
for Q(j)(t,z∗), with j≡(j ), is given by
P 0 0
∂
Q(j0)(t,z∗)=α(j0)(0)L+Q(j0+1)(t,z∗)−L†Q(0,j0)(t,z∗)+ −iH +Lz∗−L†Q(0)(t,z∗),Q(j0)(t,z∗) , (A10)
∂t 0 0 1 sys t 0 0
h i
where j =0,1,2,.... The HFD equation for Q(j)(t,z∗), with j≡(j ,j ), is given by
0 1 0 1
∂
Q(j0,j1)(t,z∗)=α(j1)(0) L,Q(j0)(t,z∗) +Q(j0+1,j1)(t,z∗)+Q(j0,j1+1)(t,z∗)+ −iH +Lz∗,Q(j0,j1)(t,z∗)
∂t 1 0 1 1 sys t 1
h i h i
− L†Q(0)(t,z∗),Q(j0,j1)(t,z∗) − L†Q(0,j0)(t,z∗),Q(j0)(t,z∗) , (A11)
0 1 1 0
h i h i
where j ,j = 0,1,2,... and j ≤ j . Because the func- Appendix B: The relationship between SDE and
0 1 0 1
tionalexpansionoftheOoperatorterminatesatthefirst HFD methods
orderforthisexactlysolvablethree-levelmodel,itfollows
from Eqs. (A1) and (8) that In this appendix, we explicitly show the relationship
between the stochastic differential equation (SDE) ap-
Q(j)(t,z∗)=0, j≡(j ,j ,...,j ), k ≥2. (A12) proach [9] and our generalized hierarchical functional
k 0 1 k
derivative (HFD) approach.
Solving HFD equations in Eqs. (A10) and (A11), we From the definition of the Qm(n) operator in Eq. (10),
can obtain O¯ for an arbitrary bath correlation function we have
α(t,s) and then determine the quantum-dynamical be- t δ
dsα(t,s) Q(n)(t,z˜∗)=(n−m)Q(n) . (B1)
havior of the system using the QSD equation. To ter- δz˜∗ m m+1
minate Eqs. (A10) and (A11) in numerical calculations, Z0 s
wecanimposeQ(jmax+1) =aQ(jmax)andQ(jmax,jmax+1) = Thus, from Eqs. (A1) and (4), it follows that
0 0 1
aQ(jmax,jmax),wherej isasuitablychosenintegerand N
1 max Q (t,z˜∗)≡ O¯(t,z˜∗)= Q(n)(t,z˜∗),
aisaparameterdeterminedbythebathcorrelationfunc- 0 0
tion. nX=0
δ
Inparticular,foranOrnstein-Uhlenbecknoise,wehave Q (t,z˜∗)= dsα(t,s) Q (t,z˜∗)
Q(01) =−γQ(00) =−γQ0 andQ(10,1) =−γQ(10,0) =−γQ1. 1 Z δz˜s∗ 0
Thus,j =0anda=−γ. Then,Eqs.(A10)and(A11) N
simplymreadxuce to = nQ1(n), (B2)
n=1
X
δ
Q (t,z˜∗)= dsα(t,s) Q (t,z˜∗)
∂ 2 δz˜∗ 1
Q (t,z∗)=α(0)L−γQ (t,z∗)−L†Q (t,z∗) Z s
0 0 1
∂t N
+ −iHsys+Lzt∗−L†Q0(t,z∗),Q0(t,z∗) , (A13) = n(n−1)Q2(n),
n=2
X
(cid:2) (cid:3)
and ......
By mathematical induction, it can be obtained that if
∂
Q1(t,z∗)=α(0)[L,Q0(t,z∗)]−2γQ1(t,z∗) N
∂t Q (t,z˜∗)= n(n−1)...(n−k+1)Q(n)
+[−iH +Lz∗,Q (t,z∗)]− L†Q (t,z∗),Q (t,z∗) k k
sys t 1 0 1 n=k
X
− L†Q1(t,z∗),Q0(t,z∗) , (cid:2) (A1(cid:3)4) N n!
= Q(n)(t,z˜∗), (B3)
(cid:2) (cid:3) (n−k)! k
which give the results in Ref. [8]. n=k
X
8
then for Q , we have This proves the result given in Eq. (11).
k+1
N n! t δ
Q (t,z˜∗)= dsα(t,s) Q(n)(t,z˜∗)
k+1 (n−k)! δz˜∗ k
n=k Z0 s
X
N When an Ornstein-Uhlenbeck noise is considered, the
n!
= (n−k−1)!Q(kn+)1(t,z˜∗). (B4) hierarchical equation of motion for the Qk(n) operator is
n=k+1 derived in Ref. [9] as
X
∂ k (n−k)
Q(n)(t,z˜∗)=δ α(0)L+ α(0) L,Q(n−1)(t,z˜∗) + z˜∗ L,Q(n−1)(t,z˜∗)
∂t k n,0 max{1,n} k−1 max{1,n} t k
h i h i
−(k+1)γQ(n)(t,z˜∗)+ −iH ,Q(n)(t,z˜∗) −(n+1)L†Q(n+1)(t,z˜∗)
k sys k k+1
n CpCn−p h i
− l n−k−l L†Q(p) (t,z˜∗),Q(n−p)(t,z˜∗) . (B5)
Cn p−l k−p+l
Xp=0Xl k h i
Substituting Eq. (B5) into ∂ Q (t,z˜∗)= N n! ∂ Q(n)(t,z˜∗), it follows that for k =0, one has
∂t k n=k (n−k)!∂t k
N P
∂ ∂
Q (t,z˜∗)= Q(n)(t,z˜∗)
∂t 0 ∂t 0
n=0
X
N
=α(0)L+ z˜∗ L,Q(n−1)(t,z˜∗) −γQ(n)(t,z˜∗)+ −iH ,Q(n)(t,z˜∗)
t 0 0 sys 0
nX=0(cid:26) h i h i
n p CpCn−p
− l n−l L†Q(p) (t,z˜∗),Q(n−p)(t,z˜∗) −(n+1)L†Q(n+1)(t,z˜∗)
Cn p−l l−p 1
Xp=0Xl=p 0 h i (cid:27)
=α(0)L+z˜∗[L,Q (t,z˜∗)]−γQ (t,z˜∗)+[−iH ,Q (t,z˜∗)]− L†Q (t,z˜∗),Q (t,z˜∗) −L†Q (t,z˜∗),
t 0 0 sys 0 0 0 1
(B6)
(cid:2) (cid:3)
and for k≥1, one has
N
∂ n! ∂
Q (t,z˜∗)= Q(n)(t,z˜∗)
∂t k (n−k)!∂t k
n=k
X
N
n! k (n−k)
= α(0) L,Q(n−1)(t,z˜∗) + z˜∗ L,Q(n−1)(t,z˜∗) −(k+1)γQ(n)(t,z˜∗)
(n−k)! n k−1 n t k k
nX=k (cid:26) h i h i
n CpCn−p
+ −iH ,Q(n)(t,z˜∗) −(n+1)L†Q(n+1)(t,z˜∗)− l n−k−l L†Q(p) (t,z˜∗),Q(n−p)(t,z˜∗)
sys k k+1 Cn p−l k−p+l
h i Xp=0Xl k h i)
=kα(0)[L,Q (t,z˜∗)]+z˜∗[L,Q (t,z˜∗)]−(k+1)γQ (t,z˜∗)−L†Q (t,z˜∗)
k−1 t k k k+1
N n! n CpCn−p
+[−iH ,Q (t,z˜∗)]− l n−k−l L†Q(p) (t,z˜∗),Q(n−p)(t,z˜∗) . (B7)
sys k (n−k)! Cn p−l k−p+l
nX=k Xp=0Xl k h i
The last term on the RHS of Eq. (B7) can be rewritten as
N n! CpCn−p
l n−k−l L†Q(p) (t,z˜∗),Q(n−p)(t,z˜∗)
(n−k)! Cn p−l k−p+l
nX=k Xp,l k h i
k! p! (n−p)!
= L† Q(p) (t,z˜∗), Q(n−p)(t,z˜∗)
(p−l)!(k+l−p)! l! p−l (n−k−l)! k−p+l
n=k p,l (cid:20) (cid:21)
XX
= Ck L†Q (t,z˜∗),Q (t,z˜∗) , (B8)
i i k−i
i
X (cid:2) (cid:3)
9
where the substitution p−l →i is taken in the last line. Eq. (7), i.e., the hierarchical equation of motion for the
Finally, substituting Eq. (B8) into Eq. (B7), we obtain Q operator. This proves that Eq. (B5) is mathemati-
k
cally equivalent to Eq. (7).
[1] H.P.Breuerand F.Petruccione, Theory of Open Quan- [14] H. T. Quan et al., Phys. Rev. Lett. 96, 140604 (2006);
tum Systems (Oxford,New York,2002). D.-W. Luo, H.-Q. Lin, J.-B. Xu, and D.-X. Yao, Phys.
[2] M. A. Nielsen and I. L. Chuang, Quantum Computation Rev. A 84, 062112 (2011).
andQuantumInformation (CambridgeUniversityPress, [15] T. Ma et al., Phys.Rev.A 90, 042108 (2014).
2000). [16] B.Bellomo,R.LoFranco,andG.Compagno,Phys.Rev.
[3] M.O.ScullyandM.S.Zubairy,Quantum Optics (Cam- Lett. 99, 160502 (2007).
bridge UniversityPress, Cambridge, England, 1997). [17] B. L. Hu,J. P. Paz, and Y. H.Zhang, Phys. Rev.D 45,
[4] L. Di´osi, N. Gisin, and W. T. Strunz, Phys. Rev. A 58, 2843 (1992)
1699 (1998). [18] H. Li, J. Shao, and S. Wang, Phys. Rev. E 84, 051112
[5] W. T. Strunz, L. Di´osi, and N. Gisin, Phys. Rev. Lett. (2011).
82, 1801 (1999). [19] A. Bassi, Phys.Rev.A 67, 062101 (2003).
[6] T. Yu,L. Di´osi, N.Gisin, and W. T. Strunz,Phys. Rev. [20] H.-P. Breuer, Phys. Rev.A 70, 012106 (2004).
A 60, 91 (1999). [21] J. Piilo, K. H¨ark¨onen, S. Maniscalco, and K.-A. Suomi-
[7] J. Gambetta and H. M. Wiseman, Phys. Rev. A 66, nen, Phys.Rev. A 79, 062112 (2009).
052105 (2002). [22] W.T.StrunzandT.Yu,Phys.Rev.A69,052115(2004).
[8] J. Jing andT. Yu,Phys.Rev.Lett.105, 240403 (2010). [23] L. Diosi, Phys. Rev. Lett. 100, 080401 (2008); 101,
[9] Z.-Z.Li et al.,Phys.Rev.A 90, 022122 (2014). 149902 (2008).
[10] D.Suess,A.Eisfeld, andW.T.Strunz,Phys.Rev.Lett. [24] H.M.WisemanandJ.Gambetta,Phys.Rev.Lett.101,
113, 150403 (2014). 140401 (2008).
[11] L.Diosi,L.Ferialdi, Phys.Rev.Lett.113,200403 (2014). [25] J.Jing,X.Zhao,J.Q.You,andT.Yu,Phys.Rev.A85,
[12] M.ChenandJ.Q.You,Phys.Rev.A87,052108(2013). 042106 (2012).
[13] X.Zhao,W.Shi,L.-A.Wu,andT.Yu,Phys.Rev.A86,
032116 (2012).