Subscriber access provided by UNIV OF PITTSBURGH

Article

Chebyshev Expansion Applied to Dissipative Quantum Systems Bogdan Popescu, Hasan Rahman, and Ulrich Kleinekathoefer J. Phys. Chem. A, Just Accepted Manuscript • DOI: 10.1021/acs.jpca.5b12237 • Publication Date (Web): 04 Feb 2016 Downloaded from http://pubs.acs.org on February 7, 2016

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

The Journal of Physical Chemistry A is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Chebyshev Expansion Applied to Dissipative Quantum Systems Bogdan Popescu,†,‡ Hasan Rahman,† and Ulrich Kleinekath¨ofer∗,† †Department of Physics and Earth Sciences, Jacobs University Bremen, Campus Ring 1, 28759 Bremen, Germany ‡Current address: Max-Planck-Institute for the Physics of Complex Systems, N¨othnitzer Str. 38, 01187, Dresden, Germany E-mail: [email protected]

1

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Abstract To determine the dynamics of a molecular aggregate under the influence of a strongly time-dependent perturbation within a dissipative environment is still, in general, a challenge. The time-dependent perturbation might be, for example, due to external fields or explicitly treated fluctuations within the environment. Methods to calculate the dynamics in these cases do exist though some of these approaches assume that the corresponding correlation functions can be written as a weighted sum of exponentials. One such theory is the Hierarchical Equations of Motion approach. If the environment, however, is described by a complex spectral density or if its temperature is low, these approaches become very inefficient. Therefore we propose a scheme based on a Chebyshev decomposition of the bath correlation functions and detail the respective quantum master equations within second-order perturbation theory in the environmental coupling. Similar approaches have recently been proposed for systems coupled to fermionic reservoirs. The proposed scheme is tested for a simple two-level system and compared to existing results. Furthermore, the advantages and disadvantages of the present Chebyshev approach are discussed.

Introduction To tackle the external influence of, e.g., a protein environment on the excitation energy transfer in systems, such as pigment-protein complexes, various methods based on the reduced density operator (RDO) are available. 1–3 The time evolution of the latter can, for instance, be determined using quantum master equations (QMEs). 1–3 However, these have severe limitations as they are commonly derived in the lowest-order perturbation theory. To overcome this bottleneck the hierarchical equations of motion (HEOM) approach is one possible and popular possibility. 4–8 The HEOM theory was initially formulated for dissipative systems by Tanimura and Kubo 4 and extended to fermionic systems by Yan and coworkers. 5 By now it has already been applied to a wide range of bosonic 6,9,10 and fermionic scenarios. 5,11–13 2

ACS Paragon Plus Environment

Page 2 of 31

Page 3 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Most notably, for the fermionic case its second-tier truncation is exact for non-interacting particles and was connected to Keldysh non-equilibrium Green’s functions, whereas the firsttier corresponds to time-non-local master equations. 13 However, HEOM calculations for large systems are still challenging, merely due to their large computational demands. 7,8 Some RDO theories such as the HEOM approach involve a so-called multi-Lorentzian parametrization of the reservoir spectral densities 14 as well as a Matsubara or equivalent expansion of the Fermi function. 15,16 These expansions turn the bath correlation functions (BCFs) or, equivalently, the self-energies of the leads into a sum of weighted exponentials. Although feasible in principle, the reconstruction of environmental spectral densities in terms of Lorentzian functions can be rather cumbersome in practice. 17 For example, as discussed in Ref. 18 for the case of graphene nanoribbons, the Lorentzian fitting of a given lead self-energy can be pursued by employing a least-square solver with the goal of minimizing deviation functions between the fitted object and its Lorentzian reconstruction. Moreover, in particular cases such as light-absorbing aggregates, where the spectral densities show a non-linear behaviour for small frequencies, multi-Lorentzian fitting fails in that range. Therefore new types of parametrizations were proposed 19 using other classes of fit functions better suited to mimic super-ohmic behaviour. The fitting procedure becomes even more complicated when trying to embed experimental or numerically obtained spectral densities, which, at least for certain systems, are known to contain a strongly varying peak structure. 20,21 For fermionic reservoirs, Tian and Chen 11,12 employed the Jacobi-Anger identity to expand the HEOM (truncated at the second-tier) in terms of Chebyshev polynomials. This scheme has the advantage that the correlation function does not have to be a weighted sum of exponentials, which also implies that very low temperatures can be used without any additional problems. At the same time, Chebyshev expansions do have the drawback that the simulation time is limited by the truncation order of the series that is infinite in principle. In Ref. 11 the auxiliary density operators (ADOs) and the lead spectral density are decomposed in terms of Chebyshev polynomials. Recently, we took up the idea of using a Chebyshev

3

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

expansion but used it directly to decompose the BCF which is much more in the spirit of the original perturbative or HEOM schemes. 22–24 To this end, a scheme for molecular wire scenarios, i. e. for fermionic leads, has been detailed to second-order in the wire-lead coupling. 25 Compared to the scheme by Tian and Chen 11,12 the scheme by the present authors only includes a single expansion instead of a twofold one. So there is only one truncation step which might turn out to be advantageous. The goal of the present contribution is to adapt the QME derived in Ref. 25 to the case of a system in contact with bosonic reservoirs. Note that polynomial expansions and especially Chebyshev expansions have been applied extensively also in the context of wave packet propagation 26–29 also for inhomogeneous 30 or explicitly time-dependent cases. 31 For open quantum systems Chebyshev polynomials 32 as well as their generalizations, the Faber polynomials have been employed, as have been Newton polynomials. 33 Recently, a decomposition scheme termed time-evolving density with orthogonal polynomials algorithm (TEDOPA) has been reported which is able to treat dissipative systems with arbitrary spectral densities. 34 The article is organized as follows. In the subsequent section the tight-binding modeling, typically used for open quantum systems is discussed whereas in Sec. III the Chebyshev expansion is detailed. In Sec. IV we employ the Chebyshev expansion to derive QMEs in second-order perturbation theory. Finally in Sec. V, numerical results are shown and discussed. In particular, we apply the new equations to a typical benchmarking two-level model and compare the results with other standard methods such as the HEOM.

Model Hamiltonian In theories for open quantum systems, the Hamiltonian of a composite system is typically decomposed into three terms according to H = HS + HR + HSR . Herein, the Hamiltonian of the quantum system is denoted by HS (t), the one of the bosonic reservoir, playing the role of the environment, by HR , whereas HSR stands for the system-reservoir coupling Hamiltonian.

4

ACS Paragon Plus Environment

Page 4 of 31

Page 5 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Since a particular focus is placed on the transfer of a single electronic excitation, we restrict ourselves to the single exciton manifold. Thus the system basis states |ii assume one exciton at site i, corresponding to a molecular orbital, whereas the other sites are supposed to be in their ground states. The ground state without any excitations |0i is not explicitly considered here. Moreover, the electronic system is modeled by N sites with site energies Ei at site i and electronic couplings Vij between any site i and site j, leading to the following tight-binding description for HS HS =

X

Ei |ii hi| +

X

Vij |ii hj| .

(1)

i6=j

i

Note that in the case of exciton transfer, which is often described by the above Hamiltonian, the exciton energies can be determined from electronic structure calculations together with the electronic couplings and are mapped onto the above tight-binding model. 20 The bosonic reservoir is treated as an infinite collection of displaced harmonic oscillators, linearly coupled to the system degrees of freedom and acting as a heat bath

HR =

X p2ξ mξ ωξ2 x2ξ  . + 2m 2 ξ ξ

(2)

In this equation pξ , xξ , ωξ and mξ denote the momenta, displacements around the equilibrium, frequencies, and masses of the bath oscillators, respectively. Finally, the system-reservoir coupling Hamiltonian HSR can be written as a sum of products of system and bath operators. The latter are denoted by Φj , while the former are supposed to be site diagonal for excitonic systems 3 HSR =

X j

Φ j Kj =

XX j

cjξ xξ |ji hj| .

(3)

ξ

Here the cjξ define the coupling strengths of the bath oscillators to the system states and the system part of the system-reservoir coupling is given by a projector onto site j, i.e., Kj = |ji hj|. In this form an interaction between the bath and the system at site j exists only if an excitation is present. Also note that we restricted ourselves to situations in which

5

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 31

each site has its own thermal bath which is assumed to be uncorrelated to the bath modes at other sites. The effect of the environment on the system is typically given in form of the spectral density, which contains information concerning the spectrum of the environment and the system-bath coupling Jj (ω) =

X ξ

c2jξ δ(ω − ωξ ) . 2mξ ωξ

(4)

For a set of bath oscillators the spectral density consists of a sequence of δ-functions, assumed sufficiently dense in order to form a continuous spectrum. Quantitatively, Jj (ω) is related to the bath correlation functions by 3 Z Z ~ ∞ e−iωt Jj (ω) ~ ∞ Cj (t) = dω dω e−iωt Jj (ω)[1 + n(~ω)] = π −∞ 1 − e−β~ω π −∞ Z ~ ∞ dω Jj (ω)[cos(ωt) coth(β~ω/2) − i sin(ωt)] , = π 0

(5)

where β = 1/(kB T ) denotes the inverse temperature and n(~ω) the Bose-Einstein distribution.

Chebyshev decomposition of correlation functions In this section, we sketch the line of reasoning that we followed in Ref. 25 to expand the BCF in terms of Chebyshev polynomials, however, adapted here to bosonic environments. To make use of this type of orthogonal polynomials, a mapping of the spectral range belonging to the Hamiltonian to their interval of definition is in order, i.e. to the interval [−1, 1]. To this end, one first needs to restrict the infinite integration range in Eq. (5) to a finite range [ωmin , ωmax ]. Then, by introducing the dimensionless variable x = (ω − ω ¯ )/Ω, where ω ¯ = (ωmax + ωmin )/2 and Ω = (ωmax − ωmin )/2, the integral expression for the correlation function turns into ~Ω Cj (t) = π

Z

1

dx −1

e−i(Ωx+¯ω)t Jj (Ωx + ω ¯) . −β~(Ωx+¯ ω ) 1−e 6

ACS Paragon Plus Environment

(6)

Page 7 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Next, we expand the Fourier kernel in Eq. (6) by means of the Jacobi-Anger identity 35 given by e−iΩxt = J0 (Ωt) +

∞ X

2(−i)n Jn (Ωt)Tn (x) , ∀t ∈ R , ∀x ∈ [−1, 1] ,

(7)

n=1

with the Chebyshev polynomials Tn and the Bessel functions of the first kind Jn . Subsequently, the BCF can be written as

Cj (t) =

∞ X n=0

J (Ωt)e | n {z

~Ω } π |

−i¯ ωt

Cn (t)

Z



1

X Jj (Ωx + ω ¯) = Cn (t)Ij,n , (8) dx (2 − δ0,n )(−i) Tn (x) 1 − e−β~(Ωx+¯ω) −1 n=0 {z } n

Ij,n

where δ0,n denotes the Kronecker delta. In Eq. (8) the time-dependent quantities Cn (t) as well as the time-independent integrals Ij,n were defined. It is crucial to note that the quantities Ij,n need to be computed only once at runtime, whereas the time-dependence of the BCF is encapsulated in the simpler expressions Cn (t). Moreover, up to this point there were no assumptions on the form of spectral densities nor on the Bose-Einstein distribution, hence also on the involved temperature. This suggests that such an expansion might be advantageous for embedding complex spectral densities and low temperatures. This fact is substantially in contrast with the exponential decomposition of the BCF, which becomes very demanding for complex spectral densities and/or low temperatures. We would like to note that the numerical evaluation of the integrals Ij,n by standard techniques turns out to be non-trivial due to the strongly oscillating Chebyshev polynomials for larger values of n. However, this can be circumvented by the substitution x = cos(Θ) and the subsequent use of the property Tn (cos(Θ)) = cos(Θ), which renders the integrands to the general form f (x) · cos(ωx), thus amenable to specifically designed integration subroutines, e.g., from the QUADPACK library. 36 The key point in deriving a QME in the exponential decomposition of the correlation function 22 combined with a complex-pole expansion of the Bose-Einstein distribution is to be able to write the BCF as a weighted sum over exponentials. This decomposition guarantees

7

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 31

that time derivatives of partial BCFs can be expressed in terms of themselves. However, to achieve the latter goal one can instead benefit from the Bessel functions of the first kind entering the time-dependent part of the BCF in Eq. (8) since these have a similar property 35  Ω d Jn−1 (Ωt) − Jn+1 (Ωt) . Jn (Ωt) = dt 2

(9)

Consequently, the time derivative of Cn (t) = Jn (Ωt)e−i¯ωt can be given by  d Ω d Cn−1 (t) − Cn+1 (t) , Cn (t) = −i¯ ω Jn (Ωt)e−i¯ωt + e−i¯ωt Jn (Ωt) = −i¯ ω Cn (t) + dt dt 2

(10)

which connects the derivative of Cn (t) to itself and its nearest-rank neighbours. This expression is the key point in deriving a QME in the next section.

Quantum master equations As discussed in previous studies 2,3,22 the second-order time-local quantum master equation can be written as d 1 ρ˜S (t) = − 2 dt ~

Z

t 0

 ˜ SR (t), [H ˜ SR (τ ), ρ˜S (t)ρeq ]] , dτ TrR [H R

(11)

where ρeq R denotes the equilibrium density operator of the reservoir, ρS (t) the system reduced density operator and the tilde on top of the operator symbols marks the interaction picture. Eq. (11) can be expressed in a more convenient form by inserting the expression of the coupling Hamiltonian, Eq. (3), and subsequently changing to the Schr¨odinger picture d i 1X ρS (t) = − LS (t)ρS (t) − [Kj , Λj (t)ρS (t) − ρS (t)Λ†j (t)] . dt ~ ~ j

8

ACS Paragon Plus Environment

(12)

Page 9 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Here LS = [HS (t), •] denotes the Liouville operator and the Λj (t) are auxiliary density operators (ADOs) defined as 22

Λj (t) =

Z

t 0

dτ Cj (t − τ )US (t, τ )Kj US† (t, τ ) ,

(13)

where US (t, τ ) stands for the time-evolution operator of the system. The QME in Eq. (11) is given in the so-called time-local (TL) form 22 rather than in its time-nonlocal form. 14 In Ref. 22 this TL form has been combined with the exponential decomposition of the BCF leading to the exponential decomposition time-local (EDTL) formalism. Alternatively we now make use of the Chebyshev analysis of the BCF detailed in the previous section. This then leads to the Chebyshev decomposition time-local (CDTL) formalism, a notation introduced previously. 25 When the respective decomposition is converged, the EDTL and CDTL approaches yield the same results as expected, since they are based on the same time-local QME. Using the Chebyshev expansion turns Eq. (13) into

Λj (t) =

∞ X n=0

Ij,n

Z |

t

dτ Cn (t − 0

τ )Us (t, τ )Kj Us† (t, τ ) {z

Λj,n (t)

}

=

∞ X

Ij,n Λj,n (t) .

(14)

n=0

This result consists of a time-dependent part of the ADO encapsulated in the new quantities Λj,n (t), which we term partial ADOs. For later usage, one can derive simple master equations for these partial ADOs by directly differentiating their definition expression, Eq. (14), and by using the relation, Eq. (10). In addition, the infinite sum over the partial ADOs is truncated at a finite order NC . As detailed in Ref. 25 and similar to the method of Ref. 11, truncating at a certain NC limits the approach to a NC -dependent simulation time. An estimate for the number of terms needed to obtain converged results can be given by 25,37,38

N = ΩT + 10 ln(ΩT ) .

9

ACS Paragon Plus Environment

(15)

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 31

If this relation is fulfilled, the uniform Chebyshev approximation shows an exponential convergence leading to very accurate results. 27 Summarizing the final set of QMEs for the quantities Λj,n (t) together with the modified forms of the QMEs for the boundary cases, the zeroth and NC -th order, one obtains  i   Kj − Ls (t)Λj,0 (t) − i¯ ω Λj,0 (t) − ΩΛj,1 (t) ,   ~   d Ω i Λj,n (t) = ω Λj,n (t) + [Λj,n−1 (t) − Λj,n+1 (t)] , − Ls (t)Λj,n (t) − i¯  dt ~ 2    i Ω   − L (t)Λ ω Λj,NC (t) + Λj,NC −1 (t) , s j,NC (t) − i¯ ~ 2

if n = 0 , if 0 > n > NC , if n = NC . (16)

The above system along with the summation rule (14) and the equation of motion for the RDO, Eq. (12), form a complete set of coupled differential equations, amenable to numerical treatment. A note is in order concerning the Markovian limit of the obtained set of QMEs. In the present time-local approach, the Markovian limit is identical to moving the upper integration limit in Eqs. (13) or (14) to infinity, i.e., using Λj (t = ∞) while propagating the density matrix ρS . Since the time-evolution of the Λj (t) is independent of ρS , Eq. (16) can be used to first determine Λj (t = ∞) and then be employed in Eq. (12) to propagate ρS . The nonMarkovian results together with this Markovian limit are tested below for a simple model, i.e. a two-level dissipative system, showing its performance, in comparison to the well-established HEOM theory. It is worth stressing again the advantage of the present method over previous, similar approaches, 23 i.e. the fact that the Bose-Einstein distribution and the bath spectral densities only enter the time-independent quantities Ij,n , which need to be computed only once, namely at the beginning of the simulation. Therefore the present scheme has no restriction, neither on the form of the respective spectral densitiy nor on the system temperature, thus being able to cope with numerically or experimentally obtained spectral densities as well as low temperatures.

10

ACS Paragon Plus Environment

Page 11 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Hierarchical equations of motion and ensemble-averaged wave packet dynamics The HEOM scheme can only be employed in cases where the correlation function can be written as a weighted sum of exponentials, implying that it uses the exponential decomposition. 7,10,24 The EDTL scheme 22 is actually equivalent to the first tier HEOM formalism if properly terminated. This restriction implies that the spectral density can be written as a sum of Lorentzian functions and that the Fermi function is written as a series expansion such as the Matsubara series. To this end, auxiliary density operators take into account the non-Markovian dynamics induced by the environment. The exponential decomposition together with the hierarchical scheme leads, in principle, to an infinite number of terms which can be truncated at a finite level, however, depending on the required accuracy. Appropriate approximations can be employed in the case of low temperatures. 39,40 In the present case a code described earlier 7 has been used to obtain results converged with respect to the number of auxiliary terms. The numerical findings coincide with results already detailed in Refs. 41,42 and 43 where also some more details on the convergence have been reported. The other method to which we compare the present CDTL results is the ensembleaveraged wave packet dynamics in its classical-path based version. In this approach, other than in the standard Ehrenfest dynamics as, for example, employed in Ref. 44, the effect of the system on the bath is neglected. In other words, the bath always stays in equilibrium if it is initially in equilibrium. This assumption is common in perturbative density matrix schemes. 3 At the same time, this assumption clearly indicates that the approach will not be able to handle strongly coupled system-bath configurations. In the ensemble-averaged wave packet dynamics, the effect of the environment is given in terms of fluctuations in the site energies ∆Em (t) which are determined as specific matrix elements of bath operators. As a starting point one determines time-dependent trajectories of site energies including the ∆Em (t) and couplings. Subsequently, the time-dependent Schr¨odinger equation is solved for

11

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 31

the state vector |ΨSα (t)i of ensemble member α along N pieces of the trajectory or an ensemble of trajectories. Finally, an ensemble average over the different pieces of the trajectory has to be performed to obtain the expectation value of an observable A leading to Nα 1 X hΨSα (t)|A|ΨSα (t)i = hΨS (t)|A|ΨS (t)i . hA(t)i = Nα α=1

(17)

A comparison of the HEOM and ensemble-averaged wave packet has been reported earlier 43 though not with a second-order perturbative treatment.

Results and Discussion In this section we employ a simple model system, i.e., a dissipative two-site system, to validate the Chebyshev QME and to compare the results to the aforementioned established methods. To this end, we consider the system Hamiltonian of form discussed before, Eq. (1), with site energies E1 and E2 and coupling elements V12 = V21 = V . Each site is coupled independently to its own thermal bath, where the coupling is governed by Drude-Lorentz spectral densities defined as Ji (ω) = 2πλi γi

ω2

ω , + γi2

(18)

with the real parameters γi and λi . Specifically, γi affects the width of the function and represents the inverse correlation time, whereas the reorganization energy λi influences its amplitude. For the sake of simplicity, however, we chose the same spectral density for both sites. This system has been used in several benchmark studies in literature 41–43 in which the simulation parameters were chosen to be V = 100 cm−1 for the electronic coupling, T = 300 K for the temperature, and 1/γi = 100 fs for the correlation time. Moreover, we distinguish two cases, i.e. in Fig. 1 equal site energies E1 = E2 are considered , whereas Fig. 2 shows the case of unequal site energies, i.e. an energy splitting of ∆E = E2 −E1 = 100 cm−1 . For both cases, i.e., Figs. 1 and 2, the reorganization energy, λi assumes four different

12

ACS Paragon Plus Environment

Page 13 of 31

-1

-1

a) λ = 2 cm

b) λ = 20 cm

1

CDTL scheme HEOM Ensemble Average CDTL Markov Limit

0.8 0.6

Population of Excited Site (ρ11)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

1 0.8 0.6

0.4

0.4

0.2

0.2

0

0

200

400

600

800

1000 0

-1

1

200

400

600

800

0 1000

-1

c) λ = 100 cm

d) λ = 500 cm

1

0.8

0.8

0.6

0.6

0.4

0.4

0.2

0.2

0

0

200

400

600

800

1000 0

200

400

600

800

0 1000

Time [fs] Figure 1: Time evolution of the population on the first site which is intially excited for different reorganization energies. The site energies are assumed the same, i.e., ∆E = E2 − E1 = 0. values, indicated in each panel spanning the system-bath coupling regime from low to high values. The plots show the population dynamics of the first site, i.e., element (1, 1) of the density operator, given by the present Chebyshev QME as compared to its counterparts computed by the HEOM scheme. The HEOM results including details of their calculation have been reported in Ref. 43 already. The important point for the present study is that these findings can serve as benchmark results since they are converged in the system-bath coupling. In Fig. 1 the results for the case of equal site energies are shown. For the lowest reorganization energy (i.e. λ = 2 cm−1 ) the agreement with the HEOM results is good which points out that the present scheme correctly describes the underlying dynamics. Moving on

13

ACS Paragon Plus Environment

The Journal of Physical Chemistry

-1

-1

a) λ = 2 cm

b) λ = 20 cm

1

CDTL Scheme Ensemble Average HEOM CDTL Markov Limit

0.8 0.6

Population of Excited Site (ρ11)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 31

1 0.8 0.6

0.4

0.4

0.2

0.2

0

0

200

400

600

800

1000 0

-1

1

200

400

600

800

0 1000

-1

c) λ = 100 cm

d) λ = 500 cm

1

0.8

0.8

0.6

0.6

0.4

0.4

0.2

0.2

0

0

200

400

600

800

1000 0

200

400

600

800

0 1000

Time [fs] Figure 2: Same as in Fig. 1 but for unequal site energies ∆E = E2 − E1 = 100 cm−1 . to intermediate values of λ, i.e. Figs. 1b and 1c the damped oscillations in the Chebyshev dynamics start to differ from those of the HEOM results mainly in amplitude. Moreover, a slight shift in frequency can be observed as well. This finding indicates that this medium magnitude range of the system-bath coupling already marks the start of the breakdown of second-order perturbation theory which is the basis of the CDTL QME. However, it is interesting to note that the scheme still correctly describes the long-time limit which corresponds to equal population of both sites for equal site energies. Using the CDTL approach, the correct thermal state is also obtained for the largest reorganization energy of λ = 500 cm−1 in Fig. 1. This second-order scheme reaches this limit, however, much too fast compared to the converged HEOM results. The case of unequal site energies is shown in Fig. 2. The agreement between the two approaches is quite similar to that for the case of equal site ener-

14

ACS Paragon Plus Environment

Page 15 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

gies. Again, the CDTL scheme almost reaches the same equilibrium value as the converged results. Small deviations in the thermal state can be explained by the second-order nature of the present approach. 45 In addition, the Markovian limit of the present scheme is shown in Figs. 1 and 2. The procedure how this limit has been obtained, was detailed in the methods section. For the case of equal site energies, the Markovian limit actually shows larger oscillations than the non-Markovian results and is less accurate at the same time. This behavior is different for the case of unequal site energies. Not many differences are visible for the smallest reorganisation energy. For the reorganisation energy of λ = 20 cm−1 the agreement with the HEOM result is better for the Markovian than for the non-Markovian dynamics. For even larger reorganisation energies, the reminiscences of oscillations, present in the non-Markovian case, vanish completely in the Markovian limit. This behavior of the Markovian case agrees well with the Markovian Redfield dynamics as reported earlier. 41,46 Moreover, shown in Figs. 1 and 2 are the results for the classical path-based ensembleaverage wave packet dynamics outlined in the previous section. It is interesting to see how different this approximation behaves compared to a perturbative treatment. For the case of equal site energies in Fig. 1, the wave packet dynamics adequately agrees with the HEOM results even for reorganisation energies up to λ = 100 cm−1 . Only for the largest reorganisation energy quantitative differences are visible, though qualitatively the behavior is similar again. So for equal site energies, the ensemble-average wave packet dynamics approach outperforms the present Chebyshev scheme which is based on a second-order treatment. In the wave-packet dynamics approach the dephasing is actually included properly in the theory by the ensemble average. This theory is missing a relaxation term, however. A fact which is invisible for equal site energies. Moving on to a bias between the site energies this deficiency of the wave packet dynamics approach becomes obvious since the long term limit is not properly reproduced, a fact which is well known. 43 Nevertheless, it is clearly visible that the dephasing part of the dynamics is still quite accurately reproduced, although the dynamics

15

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1.2 Jmax

1

J(w)

Population of Excited Site ( r11)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 31

0.8 0

0.6

0

w max

w max/2

w

0.4

0.2 0

Unequal Site Energies Equal Site Energies

0

500

1000

1500

2000

2500

3000

3500

Time [fs] Figure 3: Population dynamics of the first site determined by the CDTL scheme, for both configurations, i.e., ∆E = 0 and ∆E = 100 cm−1 . The employed spectral density is shown in the inset. It consists of a linearly increasing function which drops to zero at half the maximum frequency. The maximum of the spectral density Jmax is chosen such that it should lie within the validity range of perturbation theory. lead to the wrong thermal state. Finally, in Fig. 3 we deal with a non-smooth spectral density as depicted in the inset. The spectral density linearly increases with the frequency and then suddenly drops to zero. Such a spectral density is practically not reproducible by a sum of Lorentzian functions, i.e., it can neither be handled properly by the EDTL nor by standard HEOM approaches. Note that the spectral density is extended to negative frequencies by Ji (−ω) = −Ji (ω) before being inserted into Eq. (5). The thermal equilibria are the same as in Fig. 1 and 2, respectively, since these do not depend on the spectral density when neglecting higher-order corrections. 45 This example clearly shows the advantage of the present scheme in being able to use any

16

ACS Paragon Plus Environment

Page 17 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

form of spectral density. Similar examples can be tested for low temperatures though one runs the risk of leaving the validity range of the second-order perturbation theory.

Concluding remarks In this contribution we have introduced a new quantum master equation based on a secondorder treatment in the environmental coupling employing a Chebyshev decomposition of the bath correlation function. The CDTL scheme is an extension of a recently reported approach for coupling to fermionic reservoirs. 25 As in earlier approaches 14,23,24 including the EDTL one, 22 the present scheme does not make any additional assumptions and approximations in the time-dependence of the system Hamiltonian. So the approach is able to treat systems with rapidly fluctuating and/or strong external or internal perturbations. In the fermionic case this has been employed, for example, in studying the coherent destruction of tunneling in molecular wire systems. 47 In these investigations we used the EDTL scheme which is only efficiently usable when the spectral density can be decomposed in terms of a few Lorentzian functions and if the temperature is not too low since otherwise too many terms are needed in the decomposition of the Fermi function. The CDTL approach lifts the restriction of having to use simple forms of spectral densities and not too low temperatures. The drawback is that the simulation time is limited by the number of terms used in the Chebyshev expansion. As detailed in Eq. 15 the present scheme depends quasi-linearly on the simulation time, i.e., when the damping is small resulting in a longer simulation time needed until the steady state is reached, the scheme becomes less efficient with regard to computational time. At the same time, the proposed scheme depends quasi-linearly on the spectral range Ω of the underlying system as well. So it becomes more efficient for systems with a narrow spectral range. Moreover, a slight dependence on the reorganisation energy of the environmental coupling was observed in the numerical studies which needs to be analyzed in more detail in the future.

17

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 31

In recent years complex spectral densities have been extracted from experiments 48,49 or from a combination of classical molecular dynamics simulations and quantum chemical calculations. 20,21,50–53 These spectral densities show a very complex structure which can only approximately be reproduced by a sum of Lorentzian functions. 8 To be able to treat systems with these spectral densities one probably first has to extend the present scheme to higher orders or to the HEOM scheme, but we have clearly shown that such an approach based on a Chebyshev decomposition is feasible and a practical alternative to the exponential decomposition schemes.

Acknowledgements Financial support by the Deutsche Forschungsgemeinschaft (DFG) through project KL1299/141 is gratefully acknowledged.

References (1) Weiss, U. Quantum Dissipative Systems, 2nd ed.; World Scientific: Singapore, 1999. (2) Breuer, H. P.; Petruccione, F. The Theory of Open Quantum Systems; Oxford University Press: Oxford, 2002. (3) May, V.; K¨ uhn, O. Charge and Energy Transfer in Molecular Systems, 3rd ed.; Wiley– VCH: Berlin, 2011. (4) Tanimura, Y.; Kubo, R. Time Evolution of a Quantum System in Contact with a Nearly Gaussian-Markoffian Noise Bath. J. Phys. Soc. Jpn. 1989, 58, 101. (5) Jin, J. S.; Zheng, X.; Yan, Y. J. Exact Dynamics of Dissipative Electronic Systems and Quantum Transport: Hierarchical Equations of Motion Approach. J. Chem. Phys. 2008, 128, 234703. 18

ACS Paragon Plus Environment

Page 19 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(6) Ishizaki, A.; Fleming, G. R. Theoretical Examination of Quantum Coherence in a Photosynthetic System at Physiological Temperature. Proc. Natl. Acad. Sci. USA 2009, 106, 17255–17260. (7) Str¨ umpfer, J.; Schulten, K. Open Quantum Dynamics Calculations with the Hierarchy Equations of Motion on Parallel Computers. J. Chem. Theory Comput. 2012, 8, 2808– 2816. (8) Kreisbeck, C.; Kramer, T.; Aspuru-Guzik, A. Scalable High-Performance Algorithm for the Simulation of Exciton Dynamics. Application to the Light-Harvesting Complex II in the Presence of Resonant Vibrational Modes. J. Chem. Theory Comput. 2014, 10, 4045–4054. (9) Schr¨oder, M.; Schreiber, M.; Kleinekath¨ofer, U. Reduced Dynamics of Coupled Harmonic and Anharmonic Oscillators Using Higher-order Perturbation Theory. J. Chem. Phys. 2007, 126, 114102–1–10. (10) Str¨ umpfer, J.; Schulten, K. Excited State Dynamics in Photosynthetic Reaction Center and Light Harvesting Complex 1. J. Chem. Phys. 2012, 137, 065101–8. (11) Tian, H.; Chen, G. An Efficient Solution of Liouville-von Neumann Equation That Is Applicable to Zero and Finite Temperatures. J. Chem. Phys. 2012, 137, 204114. (12) Tian, H.; Chen, G. Application of Hierarchical Equations of Motion (HEOM) to Time Dependent Quantum Transport at Zero and Finite Temperatures. Eur. Phys. J. B 2013, 86, 411. (13) Popescu, B.; Kleinekath¨ofer, U. Treatment of Time-dependent Effects in Molecular Junctions. phys. stat. sol. (b) 2013, 250, 2288–2297, Special issue ’Quantum Transport at the Molecular Scale’.

19

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(14) Meier, C.; Tannor, D. J. Non-Markovian Evolution of the Density Operator in the Presence of Strong Laser Fields. J. Chem. Phys. 1999, 111, 3365. (15) Hu, J.; Xu, R. X.; Yan, Y. J. Pade Spectrum Decomposition of Fermi Function and Bose Function. J. Chem. Phys. 2010, 133, 101106. (16) Croy, A.; Saalmann, U. Partial Fraction Decomposition of the Fermi Function. Phys. Rev. B 2010, 82, 159904. (17) Liu, H.; Zhu, L.; Bai, S.; Shi, Q. Reduced Quantum Dynamics with Arbitrary Bath Spectral Densities: Hierarchical Equations of Motion Based on Several Different Bath Decomposition Schemes. J. Chem. Phys. 2014, 140, 134106. (18) Xie, H.; Kwok, Y.; Zhang, Y.; Jiang, F.; Zheng, X.; Yan, Y.; Chen, G. Time-dependent Quantum Transport Theory and Its Applications to Graphene Nanoribbons. Phys. Status Solidi B 2013, B, 250. (19) Ritschel, G.; Eisfeld, A. Analytic Representation of Bath Correlation Functions for Ohmic and Superohmic Spectral Densities Using Simple Poles. J. Chem. Phys. 2014, 141, 094101. (20) Olbrich, C.; Kleinekath¨ofer, U. Time-dependent Atomistic View on the Electronic Relaxation in Light-harvesting System II. J. Phys. Chem. B 2010, 114, 12427–12437. (21) Jurinovich, S.; Viani, L.; Curutchet, C.; Mennucci, B. Limits and Potentials of Quantum Chemical Methods in Modelling Photosynthetic Antennae. Phys. Chem. Chem. Phys. 2015, 17, 30783–30792. (22) Kleinekath¨ofer, U. Non-Markovian Theories Based on the Decomposition of the Spectral Density. J. Chem. Phys. 2004, 121, 2505.

20

ACS Paragon Plus Environment

Page 20 of 31

Page 21 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(23) Welack, S.; Schreiber, M.; Kleinekath¨ofer, U. The Influence of Ultra-fast Laser Pulses on Electron Transfer in Molecular Wires Studied by a Non-Markovian Density Matrix Approach. J. Chem. Phys. 2006, 124, 044712–1–9. (24) Tanimura, Y. Stochastic Liouville, Langevin, Fokker-Planck, and Master Equation Approaches to Quantum Dissipative Systems. J. Phys. Soc. Jpn. 2006, 75, 082001. (25) Popescu, B.; Rahman, H.; Kleinekath¨ofer, U. Using the Chebyshev Expansion in Quantum Transport Calculations. J. Chem. Phys. 2015, 142, 154103. (26) Kosloff, R. Time-dependent Quantum-mechanical Methods for Molecular Dynamics. J. Phys. Chem. 1988, 92, 2087–2100. (27) Kosloff, R. Propagation Methods for Quantum Molecular Dynamics. Annu. Rev. Phys. Chem. 1994, 45, 145. (28) Chen, R.; Guo, H. The Chebyshev Propagator for Quantum Systems. Comp. Phys. Commun. 1999, 119, 19. (29) Vijay, A.; Metiu, H. A Polynomial Expansion of the Quantum Propagator, the Green’s Function, and the Spectral Density Operator. J. Chem. Phys. 2002, 116, 60–68. (30) Ndong, M.; Tal Ezer, H.; Kosloff, R.; Koch, C. P. A Chebychev Propagator for Inhomogeneous Schr¨odinger Equations. J. Chem. Phys. 2009, 130, 124108. (31) Ndong, M.; Tal Ezer, H.; Kosloff, R.; Koch, C. P. A Chebychev Propagator with Iterative Time Ordering for Explicitly Time-dependent Hamiltonians. J. Chem. Phys. 2010, 132, 064105. (32) Guo, H.; Chen, R. Short-time Chebyshev Propagator for the Liouville-von Neumann Equation. J. Chem. Phys. 1999, 110, 6626.

21

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(33) Huisinga, W.; Pesce, L.; Kosloff, R.; Saalfrank, P. Faber and Newton Polynomial Integrators for Open-system Density Matrix Propagation. J. Chem. Phys. 1999, 110, 5538. (34) Woods, M. P.; Groux, R.; Chin, A. W.; Huelga, S. F.; Plenio, M. B. Mappings of Open Quantum Systems onto Chain Representations and Markovian Embeddings. J. Math. Phys. 2014, 55, 032101. (35) Arfken, G. Mathematical Methods for Physicists, 3rd ed.; Academic Press: San Diego, 1985. ¨ (36) de Doncker-Kapenga, E.; Kahaner, D. K.; Piessens, R.; Uberhuber, C. W. Quadpack: A Subroutine Package for Automatic Integration; Springer, 1983. (37) Coifman, R.; Rokhlin, V.; Wandzura, S. The Fast Multipole Method for the Wave Equation: A Pedestrian Prescription. Antennas and Propagation Magazine, IEEE 1993, 35, 7–12. (38) Landa, Y.; Tanushev, N.; Tsai, R. Discovery of Point Sources in the Helmholtz Equation Posed in Unknown Domains with Obstacles. Comm. in Math. Sci. 2011, 9, 903–928. (39) Ishizaki, A.; Tanimura, Y. Quantum Dynamics of Systems Strongly Coupled to LowTemperature Colored Noise Bath: Reduced Hierachy Equations Approach. J. Phys. Soc. Jpn. 2005, 74, 3131–3134. (40) Amin, A. F.; Li, G.-Q.; Phillips, A. H.; Kleinekath¨ofer, U. Coherent Control of the Spin Current Through a Quantum Dot. Eur. Phys. J. B 2009, 68, 103–109. (41) Ishizaki, A.; Calhoun, T. R.; Schlau Cohen, G. S.; Fleming, G. R. Quantum Coherence and Its Interplay with Protein Environments in Photosynthetic Electronic Energy Transfer. Phys. Chem. Chem. Phys. 2010, 12, 7319–7337.

22

ACS Paragon Plus Environment

Page 22 of 31

Page 23 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(42) Ishizaki, A.; Fleming, G. R. On the Interpretation of Quantum Coherent Beats Observed in Two-Dimensional Electronic Spectra of Photosynthetic Light Harvesting Complexes. J. Phys. Chem. B 2011, 115, 6227–6233. (43) Aghtar, M.; Liebers, J.; Str¨ umpfer, J.; Schulten, K.; Kleinekath¨ofer, U. Juxtaposing Density Matrix and Classical Path-based Wave Packet Dynamics. J. Chem. Phys. 2012, 136, 214101. (44) van der Vegte, C. P.; Dijkstra, A. G.; Knoester, J.; Jansen, T. L. C. Calculating TwoDimensional Spectra with the Mixed Quantum-Classical Ehrenfest Method. J. Phys. Chem. A 2013, 117, 5970–5980. (45) Geva, E.; Rosenman, E.; Tannor, D. J. On the Second-order Corrections to the Quantum Canonical Equilibrium Density Matrix. J. Chem. Phys. 2000, 113, 1380. (46) Ishizaki, A.; Fleming, G. R. Unified treatment of quantum coherent and incoherent hopping dynamics in electronic energy transfer: Reduced hierarchy equation approach. J. Chem. Phys. 2009, 130, 234111. (47) Li, G.-Q.; Schreiber, M.; Kleinekath¨ofer, U. Coherent Laser Control of the Current Through Molecular Junctions. EPL 2007, 79, 27006–1–6. (48) Adolphs, J.; Renger, T. How Proteins Trigger Excitation Energy Transfer in the Fmo Complex of Green Sulfur Bacteria. Biophys. J. 2006, 91, 2778–2797. (49) Kell, A.; Feng, X.; Reppert, M.; Jankowiak, R. On the Shape of the Phonon Spectral Density in Photosynthetic Complexes. J. Phys. Chem. B 2013, 117, 7317–7323. (50) Valleau, S.; Eisfeld, A.; Aspuru Guzik, A. On the Alternatives for Bath Correlators and Spectral Densities from Mixed Quantum-Classical Simulations. J. Chem. Phys. 2012, 137, 224103–13.

23

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(51) Aghtar, M.; Str¨ umpfer, J.; Olbrich, C.; Schulten, K.; Kleinekath¨ofer, U. Different Types of Vibrations Interacting with Electronic Excitations in Phycoerythrin 545 and FennaMatthews-Olson Antenna Systems. J. Phys. Chem. Lett. 2014, 5, 3131–3137. (52) Viani, L.; Corbella, M.; Curutchet, C.; O’Reilly, E. J.; Olaya Castro, A.; Mennucci, B. Molecular Basis of the Exciton-phonon Interactions in the PE545 Light-harvesting Complex. Phys. Chem. Chem. Phys. 2014, 16, 16302–16311. (53) Jurinovich, S.; Curutchet, C.; Mennucci, B. The Fenna-Matthews-Olson Protein Revisited: A Fully Polarizable (TD)DFT/MM Description. ChemPhysChem 2014, 15, 3194–3204.

References (1) Weiss, U. Quantum Dissipative Systems, 2nd ed.; World Scientific: Singapore, 1999. (2) Breuer, H. P.; Petruccione, F. The Theory of Open Quantum Systems; Oxford University Press: Oxford, 2002. (3) May, V.; K¨ uhn, O. Charge and Energy Transfer in Molecular Systems, 3rd ed.; Wiley– VCH: Berlin, 2011. (4) Tanimura, Y.; Kubo, R. Time Evolution of a Quantum System in Contact with a Nearly Gaussian-Markoffian Noise Bath. J. Phys. Soc. Jpn. 1989, 58, 101. (5) Jin, J. S.; Zheng, X.; Yan, Y. J. Exact Dynamics of Dissipative Electronic Systems and Quantum Transport: Hierarchical Equations of Motion Approach. J. Chem. Phys. 2008, 128, 234703. (6) Ishizaki, A.; Fleming, G. R. Theoretical Examination of Quantum Coherence in a Photosynthetic System at Physiological Temperature. Proc. Natl. Acad. Sci. USA 2009, 106, 17255–17260. 24

ACS Paragon Plus Environment

Page 24 of 31

Page 25 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(7) Str¨ umpfer, J.; Schulten, K. Open Quantum Dynamics Calculations with the Hierarchy Equations of Motion on Parallel Computers. J. Chem. Theory Comput. 2012, 8, 2808– 2816. (8) Kreisbeck, C.; Kramer, T.; Aspuru-Guzik, A. Scalable High-Performance Algorithm for the Simulation of Exciton Dynamics. Application to the Light-Harvesting Complex II in the Presence of Resonant Vibrational Modes. J. Chem. Theory Comput. 2014, 10, 4045–4054. (9) Schr¨oder, M.; Schreiber, M.; Kleinekath¨ofer, U. Reduced Dynamics of Coupled Harmonic and Anharmonic Oscillators Using Higher-order Perturbation Theory. J. Chem. Phys. 2007, 126, 114102–1–10. (10) Str¨ umpfer, J.; Schulten, K. Excited State Dynamics in Photosynthetic Reaction Center and Light Harvesting Complex 1. J. Chem. Phys. 2012, 137, 065101–8. (11) Tian, H.; Chen, G. An Efficient Solution of Liouville-von Neumann Equation That Is Applicable to Zero and Finite Temperatures. J. Chem. Phys. 2012, 137, 204114. (12) Tian, H.; Chen, G. Application of Hierarchical Equations of Motion (HEOM) to Time Dependent Quantum Transport at Zero and Finite Temperatures. Eur. Phys. J. B 2013, 86, 411. (13) Popescu, B.; Kleinekath¨ofer, U. Treatment of Time-dependent Effects in Molecular Junctions. phys. stat. sol. (b) 2013, 250, 2288–2297, Special issue ’Quantum Transport at the Molecular Scale’. (14) Meier, C.; Tannor, D. J. Non-Markovian Evolution of the Density Operator in the Presence of Strong Laser Fields. J. Chem. Phys. 1999, 111, 3365. (15) Hu, J.; Xu, R. X.; Yan, Y. J. Pade Spectrum Decomposition of Fermi Function and Bose Function. J. Chem. Phys. 2010, 133, 101106. 25

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(16) Croy, A.; Saalmann, U. Partial Fraction Decomposition of the Fermi Function. Phys. Rev. B 2010, 82, 159904. (17) Liu, H.; Zhu, L.; Bai, S.; Shi, Q. Reduced Quantum Dynamics with Arbitrary Bath Spectral Densities: Hierarchical Equations of Motion Based on Several Different Bath Decomposition Schemes. J. Chem. Phys. 2014, 140, 134106. (18) Xie, H.; Kwok, Y.; Zhang, Y.; Jiang, F.; Zheng, X.; Yan, Y.; Chen, G. Time-dependent Quantum Transport Theory and Its Applications to Graphene Nanoribbons. Phys. Status Solidi B 2013, B, 250. (19) Ritschel, G.; Eisfeld, A. Analytic Representation of Bath Correlation Functions for Ohmic and Superohmic Spectral Densities Using Simple Poles. J. Chem. Phys. 2014, 141, 094101. (20) Olbrich, C.; Kleinekath¨ofer, U. Time-dependent Atomistic View on the Electronic Relaxation in Light-harvesting System II. J. Phys. Chem. B 2010, 114, 12427–12437. (21) Jurinovich, S.; Viani, L.; Curutchet, C.; Mennucci, B. Limits and Potentials of Quantum Chemical Methods in Modelling Photosynthetic Antennae. Phys. Chem. Chem. Phys. 2015, 17, 30783–30792. (22) Kleinekath¨ofer, U. Non-Markovian Theories Based on the Decomposition of the Spectral Density. J. Chem. Phys. 2004, 121, 2505. (23) Welack, S.; Schreiber, M.; Kleinekath¨ofer, U. The Influence of Ultra-fast Laser Pulses on Electron Transfer in Molecular Wires Studied by a Non-Markovian Density Matrix Approach. J. Chem. Phys. 2006, 124, 044712–1–9. (24) Tanimura, Y. Stochastic Liouville, Langevin, Fokker-Planck, and Master Equation Approaches to Quantum Dissipative Systems. J. Phys. Soc. Jpn. 2006, 75, 082001.

26

ACS Paragon Plus Environment

Page 26 of 31

Page 27 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(25) Popescu, B.; Rahman, H.; Kleinekath¨ofer, U. Using the Chebyshev Expansion in Quantum Transport Calculations. J. Chem. Phys. 2015, 142, 154103. (26) Kosloff, R. Time-dependent Quantum-mechanical Methods for Molecular Dynamics. J. Phys. Chem. 1988, 92, 2087–2100. (27) Kosloff, R. Propagation Methods for Quantum Molecular Dynamics. Annu. Rev. Phys. Chem. 1994, 45, 145. (28) Chen, R.; Guo, H. The Chebyshev Propagator for Quantum Systems. Comp. Phys. Commun. 1999, 119, 19. (29) Vijay, A.; Metiu, H. A Polynomial Expansion of the Quantum Propagator, the Green’s Function, and the Spectral Density Operator. J. Chem. Phys. 2002, 116, 60–68. (30) Ndong, M.; Tal Ezer, H.; Kosloff, R.; Koch, C. P. A Chebychev Propagator for Inhomogeneous Schr¨odinger Equations. J. Chem. Phys. 2009, 130, 124108. (31) Ndong, M.; Tal Ezer, H.; Kosloff, R.; Koch, C. P. A Chebychev Propagator with Iterative Time Ordering for Explicitly Time-dependent Hamiltonians. J. Chem. Phys. 2010, 132, 064105. (32) Guo, H.; Chen, R. Short-time Chebyshev Propagator for the Liouville-von Neumann Equation. J. Chem. Phys. 1999, 110, 6626. (33) Huisinga, W.; Pesce, L.; Kosloff, R.; Saalfrank, P. Faber and Newton Polynomial Integrators for Open-system Density Matrix Propagation. J. Chem. Phys. 1999, 110, 5538. (34) Woods, M. P.; Groux, R.; Chin, A. W.; Huelga, S. F.; Plenio, M. B. Mappings of Open Quantum Systems onto Chain Representations and Markovian Embeddings. J. Math. Phys. 2014, 55, 032101.

27

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(35) Arfken, G. Mathematical Methods for Physicists, 3rd ed.; Academic Press: San Diego, 1985. ¨ (36) de Doncker-Kapenga, E.; Kahaner, D. K.; Piessens, R.; Uberhuber, C. W. Quadpack: A Subroutine Package for Automatic Integration; Springer, 1983. (37) Coifman, R.; Rokhlin, V.; Wandzura, S. The Fast Multipole Method for the Wave Equation: A Pedestrian Prescription. Antennas and Propagation Magazine, IEEE 1993, 35, 7–12. (38) Landa, Y.; Tanushev, N.; Tsai, R. Discovery of Point Sources in the Helmholtz Equation Posed in Unknown Domains with Obstacles. Comm. in Math. Sci. 2011, 9, 903–928. (39) Ishizaki, A.; Tanimura, Y. Quantum Dynamics of Systems Strongly Coupled to LowTemperature Colored Noise Bath: Reduced Hierachy Equations Approach. J. Phys. Soc. Jpn. 2005, 74, 3131–3134. (40) Amin, A. F.; Li, G.-Q.; Phillips, A. H.; Kleinekath¨ofer, U. Coherent Control of the Spin Current Through a Quantum Dot. Eur. Phys. J. B 2009, 68, 103–109. (41) Ishizaki, A.; Calhoun, T. R.; Schlau Cohen, G. S.; Fleming, G. R. Quantum Coherence and Its Interplay with Protein Environments in Photosynthetic Electronic Energy Transfer. Phys. Chem. Chem. Phys. 2010, 12, 7319–7337. (42) Ishizaki, A.; Fleming, G. R. On the Interpretation of Quantum Coherent Beats Observed in Two-Dimensional Electronic Spectra of Photosynthetic Light Harvesting Complexes. J. Phys. Chem. B 2011, 115, 6227–6233. (43) Aghtar, M.; Liebers, J.; Str¨ umpfer, J.; Schulten, K.; Kleinekath¨ofer, U. Juxtaposing Density Matrix and Classical Path-based Wave Packet Dynamics. J. Chem. Phys. 2012, 136, 214101.

28

ACS Paragon Plus Environment

Page 28 of 31

Page 29 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(44) van der Vegte, C. P.; Dijkstra, A. G.; Knoester, J.; Jansen, T. L. C. Calculating TwoDimensional Spectra with the Mixed Quantum-Classical Ehrenfest Method. J. Phys. Chem. A 2013, 117, 5970–5980. (45) Geva, E.; Rosenman, E.; Tannor, D. J. On the Second-order Corrections to the Quantum Canonical Equilibrium Density Matrix. J. Chem. Phys. 2000, 113, 1380. (46) Ishizaki, A.; Fleming, G. R. Unified treatment of quantum coherent and incoherent hopping dynamics in electronic energy transfer: Reduced hierarchy equation approach. J. Chem. Phys. 2009, 130, 234111. (47) Li, G.-Q.; Schreiber, M.; Kleinekath¨ofer, U. Coherent Laser Control of the Current Through Molecular Junctions. EPL 2007, 79, 27006–1–6. (48) Adolphs, J.; Renger, T. How Proteins Trigger Excitation Energy Transfer in the Fmo Complex of Green Sulfur Bacteria. Biophys. J. 2006, 91, 2778–2797. (49) Kell, A.; Feng, X.; Reppert, M.; Jankowiak, R. On the Shape of the Phonon Spectral Density in Photosynthetic Complexes. J. Phys. Chem. B 2013, 117, 7317–7323. (50) Valleau, S.; Eisfeld, A.; Aspuru Guzik, A. On the Alternatives for Bath Correlators and Spectral Densities from Mixed Quantum-Classical Simulations. J. Chem. Phys. 2012, 137, 224103–13. (51) Aghtar, M.; Str¨ umpfer, J.; Olbrich, C.; Schulten, K.; Kleinekath¨ofer, U. Different Types of Vibrations Interacting with Electronic Excitations in Phycoerythrin 545 and FennaMatthews-Olson Antenna Systems. J. Phys. Chem. Lett. 2014, 5, 3131–3137. (52) Viani, L.; Corbella, M.; Curutchet, C.; O’Reilly, E. J.; Olaya Castro, A.; Mennucci, B. Molecular Basis of the Exciton-phonon Interactions in the PE545 Light-harvesting Complex. Phys. Chem. Chem. Phys. 2014, 16, 16302–16311.

29

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(53) Jurinovich, S.; Curutchet, C.; Mennucci, B. The Fenna-Matthews-Olson Protein Revisited: A Fully Polarizable (TD)DFT/MM Description. ChemPhysChem 2014, 15, 3194–3204.

30

ACS Paragon Plus Environment

Page 30 of 31

Page 31 of 31

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 ACS Paragon Plus Environment

Chebyshev Expansion Applied to Dissipative Quantum Systems.

To determine the dynamics of a molecular aggregate under the influence of a strongly time-dependent perturbation within a dissipative environment is s...
1MB Sizes 3 Downloads 18 Views