Environmental engineering for quantum energy transport

Transport phenomena are ubiquitous throughout the science, engineering, and technology disciplines as it concerns energy, mass, charge, and information exchange between systems. In particular, energy transport in the nanoscale regime has attracted significant attention within the physical science community due to its potential to explain complex phenomena like the electronic energy transfer in molecular crystals or the Fenna-Matthews-Olson (FMO)/light harvesting complexes in photosynthetic bacteria with long time coherences. Energy transport in these systems is highly affected by environmental noise but surprisingly not always in a detrimental way. It was recently found that situations exist where noise actually enhances the transport phenomena. Such noise can take many forms, but can be characterised in three basic behaviors: quantum, correlation in time, or space. All have been shown potential to offer an energy transport enhancement. The focus of this work is on quantum transport caused by stochastic environment with spatio-temporal correlation. We consider a multi-site nearest neighbor interaction model with pure dephasing environmental noise with spatio-temporal correlation and show how an accelerated rate for the energy transfer results especially under negative spatial correlation (anti-correlation). Spatial anti-correlation provides another control parameter to help one establish the most efficient transfer of energy and may provide new insights into the working of exciton transport in photosynthetic complexes. Further the usage of spatio-temporal correlated noise may be a beneficial resource for efficient transport in large scale quantum networks. Propagation of energy through microscopic noisy systems is more efficient when the noise shows site-to-site anticorrelation. Chikako Uchiyama and collaborators from Japan modeled the transport of energy in two-sites and three-sites chains when the environmental noise is not completely random, but rather can show some level of spatial correlation (meaning that the fluctuation is the same for different sites) or anticorrelation (the fluctuation is opposite in two neighbor sites). What they found is that spatial anticorrelation gives a clear advantage for energy transport. They also explored the role played by the time scale at which these correlations are present, finding an optimal coherence time which maximizes energy transport. These studies could help understanding processes such as photosynthesis, both at a fundamental level and for applications.


INTRODUCTION
Energy transport is a fundamental primitive in our world operating on length scales ranging from the atomic scale to cosmological ones. It is at the core of natural life as well as our current technology. Any, even slight, improvements in transport efficiency can bring profound effects. In the natural world this can lead to a species dominating another, while in the technological world it can lead to lower energy consumption devices. Many of these improvements in the technological arena have been achieved by better device engineering reducing noise and imperfections. While this seems a perfectly logical approach, recently however it has been found, though nature already knew, that energy transport can be enhanced by adding environmental noise. This counter intuitive behavior has shed new light on the conventional thinking that noises in transport phenomena should be removed, and it was indeed triggered by energy transport discussions of light harvesting complex of photosynthetic bacteria. 1 The intensive research on four wave mixing in the light harvesting complex, such as the Fenna-Matthews-Olson (FMO) complex in green sulfur bacteria, 2-4 light harvesting complex in marine algae 5 or molecular crystals, 6 has given us clues to understand the mechanisms for efficient quantum energy transport. The striking long lived coherence observed in these experiments strongly suggests the positive effect of environmental noise. Extending this concept to the dynamics of exciton motion with delta-correlated stochastic noise (white noise) proposed by Haken and Strobl, 7 a simple theoretical model using pure dephasing was introduced to assist the excitation transport, referred as the environment-assisted quantum transport (ENAQT) 8 or dephasing-assisted transport. 9 To describe the long-lived coherence observations, the model was further extended to include finite noise correlation times and lengths, including a colored noise approach using dichotomic stochastic process, 10 a reservoir modeling with an infinite number of quantum harmonic oscillators, [11][12][13][14][15] and an analysis with finite correlation lengths. [16][17][18][19] Positive spatial correlations were shown to reduce the environmental effects, assisting the long-lived coherence. Both effects of temporal and spatial correlation together 20,21 were shown to be beneficial, and recently an artificial realization of the ENAQT model has been proposed. 22,23 The long-lived coherence implies the lasting energy exchange between sites, which however may affect the efficiency on the energy transfer rate. If the long-lived coherence can contribute to a higher transport efficiency, these contradicting evidences need to be somehow resolved. Spatial correlation might give us a solution to this. [17][18][19]24,25 Cao and Silbey 17 extended the interaction between the sites to include a phase relation and found the dependency of transfer efficiency on the phase, while the spatial correlation by propagation of environmental phonon modes 24 and the application of the extended ENAQT model with the spatial correlation to the photosynthetic bacteria 18,19 have been considered. Anticorrelation in noise has also been taken into account as the effects of anti-correlated (anti-phase) motion of two harmonic oscillators were investigated, 25 and the anti-spatial correlations between the bacteriochlorophylls in the FMO complex could both positively and negatively influence the energy transport. 19 These results imply that the phase relation can affect the quantum energy transport, which leads us to the natural question of what is the best spatio-temporal correlation function.
In this manuscript, we consider a simple multi-site model for the ENAQT using spatio-temporal correlated stochastic noise processes and show how we can engineer the environment to enhance quantum energy transport. The fundamental difference with the conventional treatment is to exploit the negative parameter regime (anti-correlation) of the spatial correlation as well as including finite temporal correlations with exponential decay. For the two site model, using the 2nd order of timeconvolutionless (TCL) master equation, we analytically find that the time evolution of the population under negative spatial correlation is significantly faster than uncorrelated and positively correlated spatial noise. Extending the spatial correlations to include anti-correlations provides another control parameter to find efficient energy transfer regimes. Using two independent measures for transfer efficiency, we find a significant improvement in transported quantity and elapsed time for negative spatialcorrelation compared to noise with delta-correlation both in space and time as in the original ENAQT model. Extending to the three site model, we find that the negative spatial correlation between every other site shows the best transfer efficiency, while the negative correlation between the next nearest neighbors is worse.

Model
Let us begin by describing a simple multi-site model for the ENAQT whose energy level diagram we schematically depict in Fig. 1. For N sites, the total Hamiltonian can be written as H t where ε is the cumulant expansion parameter of the timeconvolutionless master equation with |n〉 being the n-th excitation site, ω n the n-th site Larmor frequency, V nm the transition frequency between the n-th and m-th site, and f n (t) the fluctuating frequency on the n-th site which we consider as a stochastic process. We assume that the average of the frequency fluctuation for each site n is zero as 〈f n (t)〉 = 0, and the correlation function 〈f n (0)f m (t)〉 of these fluctuations can for convenience be described by a simple exponential decay as i¼ c n;m Δ 2 n;m exp À t j j=τ c;fn;mg Â Ã ; where to allow for both positive and negative (anti) spatial correlations we introduce the quantity c n,m with n; m ¼ 1; 2 N f g . It is defined over the range by −1 ≤ c n,m ≤ 1 with the extremal values −1 (+1) corresponding to perfectly anticorrelated (correlated) noise respectively. Δ n,m is the amplitude of the fluctuation, whereas τ c;fn;mg is the correlation time of the fluctuation. We also include damping with the rate κ at the end of the linear chain to trap the system in its ground state. 17 This corresponds to the excitation leaving the chain. Now averaging the density operator over the fluctuation using a timeconvolutionless decomposition approach, we obtain the master equation where for an arbitrary operator X withL 1 t ð Þ ¼ e iL0t L 1 t ð Þe ÀiL0t . The master equation given by (4) is inherently non-Markovian and incorporates both spatial and temporal correlations. In the limit that such correlations vanish, it reduces to a Markovian master equation (see "Methods"). Let us now consider a simple example. Two-site model As the description of our multi-site model has been completed, we will now consider the simplest situation, N = 2. In the two-site model, when initially only the site 1 is fully excited, the time evolution of the probability of finding the second site excited is shown in Fig. 2. Figure 2a represents the case where the fluctuation has no cross correlation between sites c n;m ¼ δ n;m À Á with the inset showing the short-time behavior. The uncorrelated situation (Fig. 2a) shows that there is an optimum time to transport the excitation to the second site depending on τ c , where we set the amplitude and correlation time of fluctuation at each site to be the same as Δ n,n = Δ and τ c; n;n f g ¼ τ c for n = {1,2}. As the correlation time τ c increases, the change of fluctuation carries a longer time memory effect and the optimum time decreases. Next in Fig. 2(b) we consider the situation where the fluctuation is anticorrelated between sites 1 and 2 by setting c 1;2 ¼ c 2;1 ¼ À1.
Assuming that the amplitude and correlation time of the fluctuation of each and between sites are the same, Δ n,m = Δ and τ c; n;m f g ¼ τ c for n,m = {1, 2}, we show in Fig. 2b that the transport finishes faster as the correlation time becomes longer. Comparing Fig. 2a with b, we find that the anti-correlation between the sites makes transport to finish faster than the uncorrelated case. This is not unexpected when one examines the two time noise correlation function given by iif the noise is the same for each site, while for the anti-correlated noise, the quantity becomes, ϕ t ð Þ ¼ 4 f 1 0 ð Þf 1 Àτ ð Þ h i , to be twice the size of the uncorrelated case.
To quantify this improvement in more details, let us examine two independent measures: (1) average trapping time and (2) ratio of transported quantity. The average trapping time is defined by 17 The average trapping time indicates how long the population of the initial excitation remains in the system. Since the excitation Fig. 1 Schematic energy level diagram of our multi-site system illustrates the case of three sites. Each site n has an associated excited state |n〉 with Larmor frequency ω n , and the excited state at each site fluctuates with the frequency fluctuation f n (t) around a fixed level. Each site has a transition frequency V ij between the adjacent sites i = n and j = n + 1. The final site in the chain, i.e., site n = 3 in this illustration, includes a decay channel with rate κ to trap the system in its ground state. This corresponds to the excitation leaving the chain Environmental engineering for quantum energy transport C Uchiyama et al. trapping only at the second site with rate κ leads to R 1 0 dtρ 22 t ð Þ ¼ κ À1 , we focus on the average trapping time 〈t〉 minus an offset as κ −1 . This is drawn in Fig. 3a for varying degree of spatial correlation c 1;2 ¼ c 2;1 ¼ c. This indicates that the transport feature depends on the correlation between sites. Especially, as c approaches to −1, the average trapping time decreases, i.e., the transport finishes faster. Moreover, the quantum yield becomes larger as c approaches to~1, since the average trapping time is found to be inversely related to the quantum yield 17 Þwhere k d is the decay rate of recombination (see "Methods").
Next, we introduce a quantity η to define an integrated probability that the excitation has been transported to the second site upto a time t u as Since the trapping rate κ is given by η indicates the ratio of transported quantity between upto a finite time t u and completion of the transport. In Fig. 3b, we chose t u to be 2000Δ when the transport to the second site is 98% completed for α = 1 and c = −1. We used it for all other evaluation in the figure. It is straightforward to observe from Figs. 2 and 3 the dependence of the energy transfer on both the degree of spatial correlation c and temporal correlation time α ¼ Δ Á τ c . In both measure, there is a clear monotonic dependence on the degree of spatial correlation, with the best results occurring as we approach perfect anti-correlation (ϕ(t) is larger in this case). We also observe a non-monotonic dependence of the populations time evolution and the measures of transfer efficiency for temporal correlations. This indicates that there will be a condition for optimal energy transfer on the correlation time. For small and decreasing α, the energy transfer takes longer time due to the fact that fluctuation becomes too fast making it difficult for the energy gaps to be close and the transition probability between sites becomes smaller (this moves us towards the Markovian limit). For α larger than the optimal value, the correlation time becomes larger making the transition probability smaller again. Thus for the best performance one should operate near this optimal value. Now what happens when we have more sites? Three-site model The two-site model has shown how spatial correlations are potentially an important resource for efficient energy transfer. The natural question which follows to this would be whether this is true for when the system has more than two sites. We extend the model to three-site linear chain under the nearest neighbor interaction with the interaction strengths given by V 12 = V 21 = V 23 = V 32 = V and V 13 = V 31 = 0. In Fig. 4 ig. 3 Two-site transport properties and their dependence on the degree of spatial correlation. In a, we show the average trapping timet minusκ À1 indicating the average time for the excitation trapped in the system, i.e., in the site 1 for the two-site model. In b, we show the ratio of transported probability. In both situations, we find a monotonic dependence on degree of spatial correlation parameter c, which indicates that transport is mostly accelerated for the anti-correlated energy fluctuations. The dependence of both measures on the correlation time is shown with α = Δ·τ c as 0.1,0.3,1 and 10. The dashed lines corresponds to α ¼ 0:1, dotdashed lines to α ¼ 0:3, solid lines to α ¼ 1 and dotted lines to α ¼ 10. Focusing on the uncorrelated case, as α increases up to 1, the average trapping time (ratio of transported quantity) decreases (increases). However, for α ¼ 10, the former (the latter) increases (decreases), which shows that a suitable correlation time of fluctuation can accelerate the transport. While the average trapping time (ratio of transported quantity) increases (decreases), for the positive degree of correlation, the acceleration of transport becomes negative. Other parameters are the same as in Fig. (2) Environmental engineering for quantum energy transport C Uchiyama et al.
spatially-uncorrelated noise and two spatially-anti correlated noise cases. Figure 4a corresponds to the uncorrelated case with c n;m ¼ δ n;m for n,m = {1, 2, 3}, while Fig. 4b represents the case for the anti-correlated noise between the nearest neighbor sites with c 1,2 = c 2,1 = c 2,3 = c 3,2 = −1 and correlated between the end sites c 1,3 = c 3,1 = 1. Finally, Fig. 4c represents the case for the anticorrelated noise between the sites 1 and 3, and the sites 2 and 3, with c 1,3 = c 3,1 = c 2,3 = c 3,2 = −1 and correlated noise between the sites 1 and 2, c 1,2 = c 2,1 = 1. In each sub-figure, we find that the dependence of the dynamics on the correlation time for the twosite model remains for the three-site case: transport finishes faster as the correlation time becomes longer particularly in the short correlation time regime, αt1. Comparing Fig. 4a-c, we find that the anti-correlated between the nearest neighbor sites typically shows the most efficient transport.
To explore these features systematically, we evaluate the average trapping time 〈t〉. In Fig. 5, we show the contour plot of the dependence of t h i À κ À1 on the degree of correlation c n,m . Comparing Fig. 5a, b, we find that "anti-ferromagnetic" correlation where the fluctuation is anti-correlated between the adjoint site shows the most efficient transport. The origin of the efficient transport lies at the higher transition probability between the adjoint cites with anti-correlation, which in the three site model corresponds to the"anti-ferromagnetic" correlation. There the average trapping time t h i À κ À1 for the uncorrelated case is 12% longer than the "anti-ferromagnetic" case. This clearly shows the acceleration in the energy transport occurs with "antiferromagnetic" correlations, while the difference is smaller as α increases (5% for α = 0.5). The question follows is how robust this efficiency in energy transport.
The most important feature is the stability against inhomogeneity in coupling between adjoint cites. To see this, we deviate V 23 /Δ from 0.1 to be 0.15 and 0.3, keeping V 12 =Δ ¼ 0:1 in the three site model. The average trapping time decreases for any spatial correlation as increasing V 23 /Δ, due to the stronger coupling strength accelerating the energy transfer. The difference of the average trapping time 〈t〉 between uncorrelated and "antiferromagnetic" case decreases from~12% to~8.5% for V 23 /Δ = 0.15 and~4% for V 23 /Δ = 0.3. This indicates a significant robustness in the acceleration mechanism, yielding the robustness around~10% with the inhomogeneous coupling upto 150%. The inset in each figure shows the short-time behavior. Comparing the behaviour between the different case, we find that the "anti-ferromagnetic" correlation (case b) accelerates the transport most. We have set the paramters as, ω 1 =Δ ¼ 1:5; ω 2 =Δ ¼ 1:2; ω 3 =Δ ¼ 1:0; V=Δ ¼ 0:1 and κ=Δ ¼κ ð Þ¼0:005 while other parameters are the same as in Fig. (2)

Fig. 5 Average trapping timet
Àκ À1 for three-site transport depending on the degree of spatial correlation for α = 0.3. a Shows this quantities dependency onc 2;3 ¼ c 3;2 À Á and c 1;2 ¼ c 2;1 À Á while setting c 1;3 ¼ c 3;1 ¼ c 2;3 c 1;2 . Similarly b shows our quantities dependency on c 1;3 ¼ c 3;1 À Á and c 1;2 ¼ c 2;1 À Á while setting c 2;3 ¼ c 3;2 ¼ c 1;3 c 1;2 . Other parameters are the same as in Fig. 4. The shorter values of average trapping time occurs at the left hand bottom side corner in (a) and right hand bottom side corner in (b). This means the "anti-ferromagnetic" correlation case is the best choice for acceleration of transport in this case Environmental engineering for quantum energy transport C Uchiyama et al.

DISCUSSION AND CONCLUSION
Energy transport in molecular complexes/light harvesting complexes in photosynthetic bacteria has raised many fundamental questions on how our nature operates and our understanding of noise effects in such processes has been challenged. The recent research indicates that noise can enhance transport rates through temporal correlations [11][12][13][14][15] or spatial correlations. 18 In this work we have considered the transport of excitation for a multi-site linear chain model whose energy levels are affected by spatio-temporal correlated stochastic noise processes. 16,17 We find that the energy transport can be accelerated by extending the spatial correlation into the negative region (anti-correlation). In the two-site model, our numerical analysis showed the significant acceleration in the energy transfer with the negative spatial correlation. By extending the model to three sites, we numerically demonstrated that the optimal efficiencies can be obtained with the "anti-ferromagnetic" correlation. The difference in transport with anti-correlated and uncorrelated noises in the four site model becomes even larger. For the temporal correlation dependencies, we explored the energy transfer dependence by changing the parameter α ¼ Δ Á τ c . A non-monotonic dependence on both the population time evolution and the transfer efficiency measures were observed, and hence the correlation time needs to be chosen to achieve the optimum energy transport. These results show new possibilities to understand efficient energy transport in nature and engineer it to our technologies.

An open systems approach
As shown in Fig. 1 our multi-site model contains fluctuating excited state energy levels as well as an energy trap on the last site (energy decay). In principle this means an open systems approach must be used-especially as our energy level fluctuations are stochastic in nature. It would thus be natural to write down a master equation in Lindblad form 26,27 (which can easily handle the last site energy decay), which would force us to use a white noise model where the fluctuations are not correlated in time. However in this case we want to examine temporal and spatial correlation effects. This means the typical master equation is not appropriate, however a master equation can be derived using time-convolutionless decomposition techniques. [28][29][30][31][32][33][34][35][36] Derivation of the time convolution type of master equation The master equation given by Eq. (4) is obtained by extracting the time evolution of the excitations in each site from the one of the total system which is written by the Liouville-von Neuman equation as where W(t) is the density operator for the total system and L t ð Þ is the Liouville operator defined as L t for an arbitrary operator X. In such a case, the extraction averages W(t) over the stochastic process of the fluctuation. Our purpose is to obtain the time evolution of the reduced density operator 〈W(t)〉(≡ρ(t)) where we denote the average operation as hÁ Á Ái. For this purpose, it is convenient to use the projection operator method. [28][29][30][31][32][33][34][35][36] Introducing a projection operator, P, which is an idempotent operator with a property as P 2 ¼ P, we describe the reduced density Þ . To obtain the time evolution of PWðtÞ, we use the formal solution of Eq.
where T + is an increasing time ordering operator from the right to the left. The procedure to obtain the master equation, Eq. (4), is roughly divided into the following six steps: 1. First we define the relevant P part and the irrelevant Qð 1 À PÞ part of the time evolution operatorÛ þ ðt; t 0 Þ as where we use interaction picture with the definition asÛ þ t; t 0 ð Þ¼ 2. Then we derive simultaneous differential equations for x(t) and y(t) as d dt 3. Next the formal solution of the irrelevant Q part can be written as, 4. We then rewrite x(τ) in Eq. (11) with x(t) and y(t) using the relation where we define þ ðt; τÞQðÀiL 1 ðτÞÞPÛ À ðt; τÞdτ 1 À σðtÞ: 5. Next we substitute the formal solution of y(t) into Eq. (9) and multiply W(t 0 ) from the right on the both hand sides of Eq. (9), withρðtÞ ¼ e ÀiL0ðtÀt0Þ ρðtÞ, Ξðt; t 0 Þ ¼ PðÀiL 1 ðtÞÞð1 À θðtÞ À1 Þ and Iðt; t 0 Þ ¼ PðÀiL 1 ðtÞÞθðtÞ À1V ðt; t 0 ÞQ. 6. Eq. (15) is then expanded with using the relation as θ t ð Þ À1 ¼ P 1 n¼0 σ t ð Þ n . This gives 7. Expansion of Eq. (16) up to the second order of cumulant for L 1 t ð Þ with using an assumption as f m t ð Þ h i¼0 and QW t 0 ð Þ ¼ 0 and transformation into the original picture from the interaction picture gives Eq. (4). Now let us consider the specific example, two site model.

TCL type master equation for the 2-site model
The concrete form of the master equation for the two-site model is written as where ρ nm t ð Þ is the (n, m) element of the reduced density operator ρ t ð Þ, κ is the trap frequency at the 2nd site and we F n ω 1 ; ω 2 ; V 12 ; t ð Þ for n = 1,2 are defined as Environmental engineering for quantum energy transport C Uchiyama et al.
In Eq. (20), f n 0 ð Þf m Àτ ð Þ h i with n; m ¼ 1; 2 f gis the correlation function of the fluctuation of the frequency. The obtained master equation shows that the transport is controlled by adjusting correlation function of the fluctuation on each site as well as the one between the fluctuation on the different cite such as f 1 0 ð Þf 2 Àτ ð Þ h i . Assuming that the amplitude and correlation time of the fluctuation of each and between site are the same, Δ n;m ¼ Δ and τ c; n;m f g ¼ τ c for n; m ¼ 1; 2 f g, ϕ(t) in Eq. (20) is given by ifor the uncorrelated case (c = 0). The dynamics described by time convolutionless(TCL) type of master equation is compared to the one by hierarchical equations of motion (HEOM) for the spin-boson system in 37 where they find that the second order and fourth order of TCL equation and HEOM shows almost the same dynamics for weak coupling case.

Relation with the traditional Lindblad master equation
The TCL master equation given by Eq. (4) in a specific limit reduces to the typical Lindblad master equation. In such a case we set the upper integral limit t ! 1 giving where all of the coefficients are time-independent. This is the Born-Markov approximation. 30,31,36 For example, the time-dependent coefficients of the two site model ignoring spatial correlation (c = 0) can be approximated as Taking the limit of the correlation time τ c to approach zero with maintaining Δ 2 τ c finite, 31 we obtain lim τ c ! 0 Defining lim τ c ! 0 F 2 ω 1 ; ω 2 ; Γ 12 ; 1 ð Þ γ, our master equation given by Eq. (4) reduces to the Lindblad form: 26,27 where A m ¼ jmihmj. We find that Eq. (24) is the same as the master equation obtained in 8 for the Haken-Strobl model. Moreover, when we substitute the stationary value of ρ 12 t ð Þ into the differential equation, we obtain the same differential equation as in. 17 Besides the above, let us note that, by using the time convolutionless type of master equation, we need only to take long-time limit as the approximation procedure, which is much simpler than using the time-nonlocal (time-convolution) type of master equation. 30 Numerical methods The effect of the fluctuation on the time evolution of the density operator is described with the time dependent coefficient of the second term in the right hand side of the time-convolutionless type of master equation given by Eq. (4). To obtain the time evolution of the density operator, we numerically solved the master equation by evaluating the time dependent coefficient and iterating the equation step by step with the coefficient. In the evaluation, we scaled the time variable and parameters with the strength of fluctuation Δ.

Quantum yield
The quantum yield is defined as q ¼ P n k t;n τ n P n k t;n τ n þ P n k d;n τ n ; (25) where k d;n is the the decay rate at the n-th site by recombination, k t;n is the one by trap and τ n ¼ R 1 0 ρ nn t ð Þdt in. 17 q indicates the ratio of trapped quantity to the total loss by recombination and trap. Cao and Silbey showed that the decay rate by recombination is necessary to be much smaller than the trapping rate, k d;n ( k t;n to obtain a high quantum yield. 17 In such situation, the dependence of τ n on k d;n can be neglected to give P n k t;n τ n k d;n ¼ 0 À Á % 1. Thus q % 1 1 þ P n k d;n τ n k d;n ¼ 0 where the last form is obtained by setting the values of k d;n to be the same as k d for all of the state n and with t h i ¼ P n τ n , which is called as the average trapping time.

Data availability
The data that support the findings of this study are available from the corresponding author upon reasonable request.