Energy transfer in quantum molecular chain -- two models of inhomogeneity

We study a linear chain of oscillators with inhomogeneity in their interactions with phonon bath. In a previous work on the Markovian master equation of the system, we investigated a model in which the difference in the site-phonon coupling between adjacent oscillators is the same throughout the chain. Here we look into another model in which the oscillators are coupled to the phonon bath with alternating strength at successive sites. Whereas in the first model all exciton modes are connected, in the second model they are coupled in pairs that are not connected to each other. Owing to this special structure in the coupling, the excitation numbers of different modes can be solved exactly in the steady state. In the first model, the minima of the excitation profile in the site basis occur at the edges of the chain, whereas in the second model the maxima occur at the edges. The energy transfer efficiency in the first model is affected by the source power whereas in the second model the efficiency is independent of it. A distinct feature in the second model is that a sink placed at the middle of the chain is able to distinguish between chains with even and odd number of sites. The energy transfer efficiency in a chain with even number of sites is higher than a chain with odd number of sites. Therefore, it reveals the discrete nature of the chain. In the limit of very long chain when the discreteness of the chain is less evident, the efficiencies approach each other.


MOTIVATIONS
In this paper, we summarise our work on the energy transfer in quantum molecular chain [1] under the framework of open quantum systems [2,3], and introduce a new model of site-phonon coupling as a further application to the results. The Hamiltonian we study describes the transfer of energy in photosynthetic systems [3,4]. It also contained solitonic solutions when applied to the energy transfer in α-helix protein chain [5]. Owing to the complexity in real systems, studies are mostly carried out numerically with the approximation that each site contains a single excited level [4,6,7,8,9,10,11,12]. The time evolution of the excitation is followed and the energy yield is obtained.
To carry out an analytic study on the system, we study a simplified version of a linear chain of oscillators interacting with a phonon bath. In our approach, we first simplify the Hamiltonian through a unitary transformation into a dressed basis. A further transformation into an exciton basis then diagonalises the Hamiltonian of the chain. Supplemented by a few assumptions on the system, we obtain a reduced dynamics of the chain in a phonon bath under thermal equilibrium.
In our previous studies [1], we assumed a model in which a difference in the site-phonon coupling between nearest neighbours is the same. The relaxation rate could be calculated analytically and all the exciton modes are coupled. We also studied the excitation profile of the chain and the energy transfer efficiency in the steady state in detail. Here, we explore a different model in which the oscillators are coupled to the phonon bath at alternating strength in the subsequent sites. This implies that the difference in the site-phonon coupling changes sign between pairs of subsequent nearest neighbours. In this model, pairs of exciton modes are coupled separately. Each pair does not interact with the rest. The excitonic structure is very similar to a system that consists of many independently coupled two-level systems with different excitation energy. This model can be solved analytically in the steady state. It leads to very different behaviours compared to the previous model.
In the following, we will first summarise the main steps taken and results obtained in our previous work. We then calculate the relevant quantities of the second model and compare the results with the first model. We hope that our analytic work on such simplified systems can provide insights into the energy transfer mechanism in quantum molecular chains.

EQUATIONS OF MOTION IN EXCITATION NUMBER
In photosynthetic systems, light is absorbed by protein pigments held in a protein scaffold [13]. Light creates excitations that are transferred through the chlorophyll molecules in the protein pigments to reaction centers, where energy is trapped in a useful form of chemical compound. During the process, the excitations are influenced by the vibrations of the molecules in the protein scaffold, which is described by phonon fields. In the well-studied Fenna-Matthews-Olson (FMO) pigment protein complex in green sulphur bacteria, each of the complex consists of three monomers. Each monomer contains seven bacteriochlorophyll-a (BChla) molecules that are connected to each other in a complicated way.
To facilitate our analytical works, we consider a simplified version of a linear chain of coupled oscillators that interact with their nearest neighbours where we use the unitsh = c = 1. a † x , a x are the creation and annihilation operators at site-x. The oscillators have a same natural frequency ω 0 and an intersite coupling J. The chain lies in a protein scaffold whose vibrations are effectively described by a field of phonons H q = ∑ q ω q b † q b q . The phonons interact with the oscillators through the potential where q is a real dimensionless coupling constant. In most works exciton basis is introduced at this stage. However, in our approach we first introduce a new basis through the unitary transformation [3] In the new basis, the number operator of the dressed oscillator remains the same as the number operator a † x a x of the bare oscillator [1]. This simplifies our discussion when we consider the excitation numbers of the oscillator since they are equivalent in both basis. The natural frequency is modified by the interaction between the oscillators and phonons. We assume that this modified term is negligible in order to diagonalize H 0 .
The transformation introduces into the interaction a difference in the interaction χ (x) q between neighbouring sites with phonons This fact enables us to consider different models of inhomogeneity in χ (x) q . After this transformation, we introduce the exciton basis through a discrete sine transform [3] A µ ≡ 2 where µ = 1, 2, · · · , . The exciton basis diagonalizes the Hamiltonian of the chain into ⁞ ⁞ with the exciton energy assuming that the correction to ω 0 is negligible. The spectrum is modulated by the cosine function as shown in Fig. 1. The advantage of the basis (5) is that the exciton energy is a monotonically decreasing function of µ. Owing to this fact, we can interpret the relaxation process among the excitons in a straight-forward way in the next section. As already explained in Ref. [3], there are two modes of energy transfer in quantum molecular chains. The first mode occurs through excitation hopping from one site to the next. The second mode involves transfer of energy between delocalised states that extend over several sites in the form of wave packets. The second mode is the dominant mode of transfer when intramolecular relaxations occur at a slower rate than intermolecular transitions. It is the second mode of transfer that we try to describe through the exciton basis.
In the exciton basis, the interaction Hamiltonian takes the form where L † µν ≡ A † µ A ν is an exciton raising operator. c µν (x) is a coefficient whose sum over x can be calculated explicitly. We will give its expression in the next section.
We assume that the difference in the inhomogeneity is sufficiently small that it can be used as a perturbation parameter to obtain the Markovian master equation of the chain using standard method. The phonon field is assumed to be in thermal equilibrium. The average number of excitation in the σ -mode is given by n σ ≡ A † σ A σ , where · denotes a trace over the density matrix of the exciton. The equation of motion of n σ can be worked out to be We will present the expression of the relaxation rate of exciton Γ µν in the next section. The phonon bath has the canonical distribution with inverse temperature β = 1/(k B T ). The energy gap between two exciton modes is ω µν ≡ ω µ −ω ν . Eq.(10) is a set of coupled nonlinear equations in n σ . We note that in the steady state, the off-diagonal components of the excitations vanish [1].
Assuming that there is a constant source supplying a rate of energy s at site-1, the σ -mode exciton gains a portion of power from the source, where (with z = 1 in Eq. (12)) is a weight factor ∑ σ α (z) σ = 1. Later on in our discussion, we denote the site where the sink is attached to by z. We consider two types of loss mechanisms. One of them is the loss to radiation, where an exciton relaxes to the ground state of the system and radiates off a photon. We assume that the rate of radiation γ r at each site is the same. The other loss mechanism is absorption of excitons by the sink. The energy absorbed turns into useful energy. If the sink is prepared at site-z, then each mode will experience a loss rate of α (z) σ γ s to the sink. Effects of the interactions between the chain and phonon bath are already included in Eq. (10).
After including all the effects, we obtain the following set of rate equations for the number of excitations where the sink is prepared at site-z. We denote the excitation number in the steady state byn σ . When we sum Eq. (14) over all the modes in the steady state, we obtain a consistency equation There are two clearly separated time scales in the reduced dynamics. During the shorter time scale of the order of 1/Γ µν , the excitations will achieve approximate equilibrium before the longer one of the order of 1/ε (z) σ sets in as the system approaches the steady state.
This set of equation is similar to those obtained in in biological systems [14,15] in an attempt to explain energy storage and energy transfer. It was found that when the source power exceeds certain limit, the excitations accumulate to the lowest exciton mode. An example of this collective oscillation to the lowest mode was reported recently in protein macromolecules [16].

TWO MODELS OF INHOMOGENEITY IN SITE-PHONON COUPLING
(1) First model In the first model of inhomogeneity in site-phonon coupling, we assume that the difference between consecutive sites is constant along the chain where η is a parameter that estimates the magnitude of the difference. We assume that it is a small parameter and use it as a perturbation parameter to obtain the Markovian master equation. The sum of the coefficient c µν (x) over x in Eq.(4) can then be calculated analytically to yield and 0 otherwise. An interesting feature of this model is that only those modes satisfying the condition ν − µ = odd are coupled, see Fig. 1(a) for illustrations. It follows that the relaxation rate of excitons between two modes is [1] and 0 otherwise, where γ d is the dephasing rate of a single site under the influence of a phonon bath. Notice that the denominator is proportional to the fourth-order of the energy gap between two excitons. Therefore, the larger the difference in the index m = ν − µ, the wider the energy gap (cf. Fig. 1) and the smaller the relaxation rate. We can then group the transitions into sets labelled by m, as illustrated in Fig. 1(a). The transitions between adjacent levels labeled by m = 1 are fastest. The square of sine functions in the numerator in Eq. (19) indicates that the rate is fastest for excitons that lie at the center of the spectrum, highlighted by the circle in Fig. 1(a). Based on this information, we infer that the excitons relax by cascading down the energy spectrum. The lowest level of the excitons then acts as a metastable state before the excitation finally relaxes to the true ground state by losing energy to the sink or radiation. An analysis shows that the relaxation rate increases as 2 in the long chain limit [1].
(2) Second model In the second model, we assume that the sites interact with the phonon bath in alternating strength, depending on whether they are odd or even sites in the chain. Consequently, we have the following difference in the site-phonon coupling It leads to the coefficient and 0 otherwise. Notice that in this model only two modes satisfying the condition µ + ν = + 1 are coupled. Examples of energy levels coupled in pairs are illustrated in Fig. 1(b). The relaxation rate turns out to be constant for all pairs, and 0 otherwise. It suggests that the system behaves like a set of independent two-level systems (with different energy gaps) separately connected to the source, phonon bath and sink. In this case, the set of equations (14) together with Eq.(10) reduces to a set of equations coupled in pair, where µ + ν = + 1 only. For odd , the unpaired mode ( + 1)/2 evolves alone, i.e., it is not coupled to other modes, Denotingn as the excitation number at the steady state, the total excitation number of two coupled modes is   Fig. 2(a) shows the whole range of excitation profile. Fig. 2(b) is a magnified view of Fig. 2(a) in the range 0.036 ≤n To obtain Eq.(26), we use s µ = s ν and ε ν for µ + ν = + 1. Owing to the simplicity in the coupling, the model has exact solution in the steady state. Using Eq.(26), we solve the simultaneous equations (23)-(24) at the steady state to get In Eq.(28), we have chosen the negative branch of the square root, otherwisen ν can be negative, which is not permitted.

EXCITATION PROFILES IN THE SITE BASIS
The excitation number at site-x,n x a x , can be obtained by summing over all the exciton modes weighted by α Let us compare the excitation profiles of the two models using typical parameters from the well-studied Fenna-Matthews-Olson (FMO) pigment protein complex in green sulphur bacteria [1,6,13]. We have approximately J = 100 cm −1 and ω 0 = 12, 500 cm −1 . For decay rate, we use the dephasing rate γ d = 20 ps −1 , radiation rate γ r = 0.001 ps −1 , and a trapping power of the sink γ s = 1 ps −1 . We choose a source power of one excitation per nanosecond s = 0.001 ps −1 as a reference of low power source, and an inhomogeneity parameter of η = 0.1. We start with a chain at the ground state as an initial condition. (1) First model The excitation profile in the site basis (29) and the efficiency of the first model had been studied in detail in Ref. [1]. In this model, the position of the sink does not affect the overall feature of the profile much, see Fig. 2(a) for examples of low source power s = 0.001 ps −1 , where the excitation profiles are labelled by the -and -curves. We notice that the excitations distribute almost uniformly across the chain and the two curves almost coincide. Fig. 2(b) gives a magnified view of the profiles. The minima of the excitations occur at the edges of the chain.
When the power of the source exceeds certain limit, the excitations will "condense" to the lowest mode [14,15]. Now the excitation profile in the site basis has a prominent maximum at the middle of the chain, illustrated in Fig. 2(a) by the •and •-curves for a source power of s = 10 ps −1 . Intuitively, we expect that positioning the sink at the middle of the chain will lead to a higher energy transfer efficiency. This is confirmed by the results in the next section.
where we use α µ for µ + ν = + 1. The symbol "mid" equals /2 or ( + 1)/2 depending on whether is even or odd, respectively. The profile function is independent of s. It is determined by the position of the sink z. Since the total excitation numberN is also proportional to s, the normalize excitation populationn   Fig. 2(b). This is in sharp contrast to the first model where minima occur at the edges of the chain.
When the sink is shifted to the middle of the chain, the sink power α (z) σ γ s vanishes at alternate exciton indices for a chain with odd number of sites. This means that half of the modes lose energy through radiation (and the phonon bath) only. As a result, the efficiency is expected to be lower, see the discussion in the next section. It turns out that the profile in the site basis exhibits maxima at the edges too, see the -curve in Fig. 2(a). Most interestingly, the middle site-13 has a sharp minimum with almost zero excitation. On either side of the sharp minimum, we find two hollow excitation profiles. In the case of even , the profile is similar to odd except now there are two minima at the middle of the chain (not shown in the figure). By carefully arranging the sink at a distance of three quarter z = 17 from the edge of the chain, it is possible to create an excitation profile with three hollow structures, as shown in Fig. 3.

ENERGY TRANSFER EFFICIENCY
The efficiency of energy transfer along the chain is defined to be the ratio of the rate of useful energy absorbed by the sink to the rate of total energy supplied by the source We use the fact that a σ -mode exciton loses ω σ of energy to the sink at a rate of γ s α (z) σnσ when it relaxes to the ground state.
We now summarise the results of the energy transfer efficiency of the first model already studied extensively in Ref. [1]. In general, the efficiency reduces as the length of the chain increases. In Ref. [1], it was shown that the external parameters, such as γ s and γ r , play a dominant role in determining the energy transfer efficiency. Internal parameters related to the chain, such as J, χ (x) q , ω 0 and model related parameters such as η, have mild influence on the efficiency.
The analytical expressions of energy transfer efficiency for low power and high power sources were already obtained in Ref. [1]. In general, the efficiency depends on the position of the sink. Below we study the efficiency for two positions of the sink, either it lies at the end or at the middle of the chain.
(1) First model Fig. 4 shows that the energy transfer efficiency using the parameters listed in the first paragraph of the last section. In the first model at low source power s = 0.001 ps −1 , -and -curve label a sink prepared at the end and at the middle of the chain, respectively. Though the efficiency of the former configuration is slightly greater than the latter, they almost coincide, see the inset of Fig. 4. At high source power, we already learned from Ref. [1] that positioning a sink at the middle yields better efficiency (not shown in the figure). This is because at high source power, excitations will accumulate to the lowest exciton level [14,15]. The excitation profile of the lowest exciton level in the site basis turns out to have a maximum at the middle of the chain, as illustrated by the •-curve in Fig. 2(a). As a result, a sink placed at the middle is more efficient in trapping energy from the chain.
(2) Second model In the second model, a sink prepared at the end of the chain ( -curve) gives very similar efficiency as the first model, see the inset in Fig. 4. However, when the position of the sink shifts to the middle of the chain, the efficiency becomes smaller and depends on whether the chain has odd or even number of sites. The efficiency is almost halved in odd-chain ( -curve) compared to even-chain ( -curve). This is dictated by the coupled-in-pair excitation structure of the second model, as well as by the distribution in the power of the sink across the exciton modes. In fact, half of the sink power α (z) σ γ s vanishes in odd-chain for a sink at the middle z = ( + 1)/2. It is interesting to note that this reveals the discrete nature of the chain in that it distinguishes between odd-chain that is less efficient in transferring energy to the sink compared to even-chain.
We expect that the discrete nature becomes less prominent in very long chain. In fact, when the site number increases, the efficiency in even-chain approaches the one in odd-chain. This is illustrated in Fig. 5, which shows the efficiencies for two values of sink power γ s when energy is supplied at two different source power. The rest of the parameters used in the plots are listed in the first paragraph of the previous section. In the figure, filled shapes (•, , , ) denote chains with odd = 3, 49, 99, 149, · · · , and empty shapes (•, , , ) denote chain with even = 2, 50, 100, 150, · · · . The sinks are located at the middle of the chain in all the plots. At a rate of sink γ s = 1 ps −1 , (•, •) and ( , ) label s = 0.001 ps −1 and 10 ps −1 , respectively. Whereas when γ s is reduced to 0.1 ps −1 , ( , ) and ( , ) label s = 0.001 ps −1 and 10 ps −1 , respectively.
The plots verify the fact that in the second model, the energy transfer efficiency is independent of the source power, a result that follows from Eq.(30). We also notice that in all cases the efficiency of even-chain approaches the one of odd-chain as increases. The approach is faster for smaller sink power γ s . In very long chain, the efficiency should be insensitive to whether the chain has odd or even number of sites.

CONCLUSION
We show that by designing different models of coupling between the oscillators and the phonon bath in a quantum molecular chain, it is possible to produce different excitation profiles along the chain with different efficiencies. When the difference in the site-phonon coupling alternates in sign between successive sites along a chain, the structure in the exciton couplings is very different compared to a model introduced in our previous work, where the difference in the site-phonon coupling is the same in successive sites. In the former model, only pairs of excitons are coupled, and all transitions have the same relaxation rate. Furthermore, the excitation profiles in the site basis are independent of the source power. It is also interesting to learn that a sink prepared at the middle of the chain is able to tell apart whether the chain has odd or even number of sites in that they yield very different efficiency. As the number of sites increases, the evenness or oddness of the chain becomes less obvious and the efficiency approaches each other.