Simulated two-dimensional electronic spectroscopy of the eight-bacteriochlorophyll FMO complex Shu-Hao Yeh and Sabre Kais Citation: The Journal of Chemical Physics 141, 234105 (2014); doi: 10.1063/1.4903546 View online: http://dx.doi.org/10.1063/1.4903546 View Table of Contents: http://scitation.aip.org/content/aip/journal/jcp/141/23?ver=pdfcov Published by the AIP Publishing Articles you may be interested in Two-dimensional electronic spectroscopy of bacteriochlorophyll a in solution: Elucidating the coherence dynamics of the Fenna-Matthews-Olson complex using its chromophore as a control J. Chem. Phys. 137, 125101 (2012); 10.1063/1.4752107 Two-dimensional patterning of bacterial light-harvesting 2 complexes on lipid-modified gold surface Appl. Phys. Lett. 100, 233701 (2012); 10.1063/1.4726105 Signatures of correlated excitonic dynamics in two-dimensional spectroscopy of the Fenna-Matthew-Olson photosynthetic complex J. Chem. Phys. 136, 104505 (2012); 10.1063/1.3690498 Simulation of the two-dimensional electronic spectra of the Fenna-Matthews-Olson complex using the hierarchical equations of motion method J. Chem. Phys. 134, 194508 (2011); 10.1063/1.3589982 Efficient and accurate simulations of two-dimensional electronic photon-echo signals: Illustration for a simple model of the Fenna–Matthews–Olson complex J. Chem. Phys. 132, 014501 (2010); 10.1063/1.3268705
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 132.206.27.24 On: Tue, 10 Feb 2015 11:52:33
THE JOURNAL OF CHEMICAL PHYSICS 141, 234105 (2014)
Simulated two-dimensional electronic spectroscopy of the eight-bacteriochlorophyll FMO complex Shu-Hao Yeh1 and Sabre Kais1,2,a) 1
Department of Chemistry and Birck Nanotechnology Center, Purdue University, West Lafayette, Indiana 47907, USA 2 Qatar Environment and Energy Research Institute, Qatar Foundation, Doha, Qatar
(Received 13 August 2014; accepted 25 November 2014; published online 15 December 2014) The Fenna-Matthews-Olson (FMO) protein-pigment complex acts as a molecular wire conducting energy between the outer antenna system and the reaction center; it is an important photosynthetic system to study the transfer of excitonic energy. Recent crystallographic studies report the existence of an additional (eighth) bacteriochlorophyll a (BChl a) in some of the FMO monomers. To understand the functionality of this eighth BChl, we simulated the two-dimensional electronic spectra of both the 7-site (apo form) and the 8-site (holo form) variant of the FMO complex from green sulfur bacteria, Prosthecochloris aestuarii. By comparing the spectrum, it was found that the eighth BChl can affect two different excitonic energy transfer pathways: (1) it is directly involved in the first apo form pathway (6 → 3 → 1) by passing the excitonic energy to exciton 6; and (2) it facilitates an increase in the excitonic wave function overlap between excitons 4 and 5 in the second pathway (7 → 4,5 → 2 → 1) and thus increases the possible downward sampling routes across the BChls. © 2014 AIP Publishing LLC. [http://dx.doi.org/10.1063/1.4903546] I. INTRODUCTION
Solar energy is one of the most abundant energy source that can be utilized by living organisms. In plants and bacteria, the conversion of this solar electromagnetic energy into chemical energy is called photosynthesis.1–3 One of the most studied systems within photosynthesis is the FennaMatthews-Olson (FMO) protein-pigment complex in green sulfur bacteria,4–6 which funnels the excitonic energy absorbed from the outer antenna system (chlorosome/baseplate) to the reaction center. Excitonic energy transfer (EET) in such photosynthetic protein-pigment complexes has been studied by various nonlinear optical spectroscopy methods such as pump-probe, hole burning, and photon-echo.1, 7 Recent experimental studies using 2D electronic photon-echo spectroscopy reported the existence of beating patterns with oscillation periods that approximately correspond to the difference of the eigenenergies. This phenomenon has been found in several photosynthetic protein-pigment complexes at both cryogenic temperatures8–10 and room temperature,11, 12 which reveals the importance of quantum coherence in excitonic energy transfer. For the FMO complex, the oscillation time is on the order of hundreds of femtoseconds at ambient temperature and may reach picoseconds at 77 K,13 and it has been suggested that quantum coherence may be a crucial component to the reported high-efficiency energy transfer.14–16 However, the origin of this long coherence time is still an open question. Recently, several theoretical and experimental work reported that vibronic coupling is a possible explanation for such long-lasting beating.17–22 The FMO complex has a homotrimer structure where each subunit comprises seven BChl a. Recent crystallographic a) Author to whom correspondence should be addressed. Electronic mail:
[email protected] 0021-9606/2014/141(23)/234105/8/$30.00
studies have found the presence of an additional (eighth) bacteriochlorophyll a (BChl a) in some of the monomers; the eighth BChl sits in the cleft at the protein surface close to the baseplate.23, 24 Busch et al. calculated the site energies of BChls in the FMO complex based on a new crystal structure from the green sulfur bacteria, Prosthecochloris aestuarii.25 In their calculation, it was found that BChl 8 has the highest site energy, and after modeling the baseplate and the reaction center (RC) by two low dielectric constant layers ( = 4) they also found that BChl 8 was red-shifted comparing to a homogeneous high dielectric constant environment = 80 and the overlap with the baseplate fluorescence spectrum was increased. Furthermore, due to the greater vicinity of BChl 8 to the baseplate (BChl 8 being closer than any other BChl of the FMO) it has also been proposed that BChl 8 functions as a linker between the baseplate and the FMO complex. Experimental findings of quantum coherence within the photosynthetic protein-pigment complexes have also instigated great interest in the use of theoretical models for describing photosynthetic EET. Theoretical description of some particularly challenging systems – such as the FMO complex and the B850 ring of the light harvesting complex II (LH2) found in purple bacteria – have their difficulties due to the magnitude of the electronic coupling and the exciton-phonon coupling are of same order and therefore neither of them can be treated perturbatively. To tackle this problem, a methodology referred to as hierarchical equation of motion (HEOM) has been developed by Tanimura and Kubo;26, 27 a high temperature approximation (HTA) of HEOM has also been developed by Ishizaki and Fleming to describe the EET dynamics in a model dimer system.28, 29 The HEOM has already been implemented in calculations for the excitonic transfer dynamics of photosynthetic protein-pigment complexes30–34 and more recently it is also used for simulations of both the
141, 234105-1
© 2014 AIP Publishing LLC
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 132.206.27.24 On: Tue, 10 Feb 2015 11:52:33
234105-2
S.-H. Yeh and S. Kais
J. Chem. Phys. 141, 234105 (2014)
linear and the 2D electronic spectra of the above mentioned systems.35–37 The aim of this study is to understand the role of the BChl 8 in the EET pathways of the FMO complex. To achieve this, we first analyse the difference of excitonic wave function overlap among BChls between the apo (7-sites) and holo (8sites) forms and propose corresponding EET pathways. Then the HEOM was applied for simulations of the 2D photon-echo spectra for both forms by using a model Hamiltonian formed in a previous study (PDBID: 3EOJ).24 Meanwhile, we have also explored the effects of both temperature and the presence of static disorder to the EET. A. Model Hamiltonian
fined by the coupling between the non-sharing pairs. The final system Hamiltonian used for linear absorption and 2D spectra calculation are both in a block diagonal form and can be written as 0 0 , (4) HS,1D = 0 H1ex ⎞ ⎛ 0 0 0 ⎟ ⎜ (5) HS,2D = ⎝ 0 H1ex 0 ⎠, 0
(1)
where HS , HB , and HSB represent the Hamiltonian of the system, environment and system-environment coupling, respectively. We employed the Frenkel exciton model to describe the system with each site treated as a two-level system. For the system Hamiltonian, HS , both the 7-sites and 8-sites singleexciton Hamiltonian of the FMO BChls are obtained from a structural-based theoretical study by Busch et al.25 This single exciton Hamiltonian is given as ⎞ ⎛ N ⎝(εj + λj )|j j | + Jj k (|j k| + |kj |)⎠ , H1ex = k=j
j =1
(2) where the state |j represents the QY excitation of the jth BChl with site energy εj , and Jjk corresponds to the excitonic dipole-dipole coupling between BChl sites j and k. λj is the reorganization energy which is introduced to counteract the energy shift due to the system-bath coupling. It is assumed that this reorganization energy is identical for each site, therefore λj = λ. In order to calculate the two-dimensional electronic spectra, it is necessary to include both the ground state and the two-exciton manifold into consideration. An approximation made herein is that each BChl is modeled as a two-level system. In the two-exciton manifold, the system is capable of exciting two different BChls simultaneously and its Hamiltonian can be constructed based on information of the single-exciton manifold
j k|H2ex |j k = δjj δkk (εj + εk ) + δjj (1 − δkk )Jkk + δj k (1 − δkj )Jkj + δkj (1 − δj k )Jj k + δkk (1 − δjj )Jjj .
(3)
|jk represents a two-exciton state describing the simultaneous excitation of sites j and k, where for counting 1 ≤ j < k ≤ N. The first term defines the diagonal elements by summing the site energy of the two involved sites. The remaining terms describe the off-diagonal elements under the stipulation that the excitonic transition is permitted if the initial and the final states share only one common site, where the strength is de-
H2ex
where the energy of the ground state is taken to be zero. The environment Hamiltonian, HB , describes the nuclear degrees of freedom by a phonon bath,
The total Hamiltonian both for linear absorption and for 2D spectra calculations can be described through H = HS + HB + HSB ,
0
HB =
N pˆ j2ξ j =1
ξ
2mj ξ
1 + mj ξ ωj2ξ xˆj2ξ , 2
(6)
where mjξ , ωjξ , xˆj ξ , pˆ j ξ are the mass, frequency, position, and momentum operators of the ξ th bath oscillator associate with the jth BChl, respectively. For the system-environment coupling Hamiltonian, HSB , it is assumed that each BChl is linearly coupled to the vibrational modes of its own bath and therefore induces site energy fluctuations independently. This interaction Hamiltonian is given by HSB =
N j =1
cj ξ |j j |xˆj ξ =
N
Pˆj Bˆ j ,
(7)
j =1
ξ
where Pˆj = |j j | is the excitonic projection operator and Bˆ j = ξ cj ξ xˆj ξ is defined as the collective bath operator. cjξ represents the electron-phonon coupling between the jth BChl and the ξ th phonon mode. B. Optical response functions
The interaction between an electric dipole and a weak electromagnetic field can be simulated by using optical response functions. These response functions describe how a weakly perturbed nonequilibrium system deviates from its equilibrium behaviors. The linear absorption spectra can be simulated first by calculating the first-order optical response function38, 39 R (1) (t) =
i [μ(t), ˆ μ], ˆ ¯
(8)
ˆ −iH t/¯ is the dipole operator within the where μ(t) ˆ ≡ eiH t/¯ μe Heisenberg picture and · · · ≡ Tr{· · · ρˆeq }. The 2D spectra can be simulated similarly by first calculating the third-order optical response function, which takes the form R (3) (t3 , t2 , t1 ) =−
i [[[μ(t ˆ 1 + t2 + t3 ), μ(t ˆ 1 + t2 )], μ(t ˆ 1 )], μ] ˆ . ¯3
(9)
Because the excitation energy is much larger than the scale of the thermal energy (kB T) for our system of interest, it can be assumed that the system is initially equilibrated in its ground state. The initial equilibrium state can thus be written as
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 132.206.27.24 On: Tue, 10 Feb 2015 11:52:33
234105-3
S.-H. Yeh and S. Kais
J. Chem. Phys. 141, 234105 (2014)
ρˆeq = |00| ⊗ e−βHB /TrB {e−βHB }, where β refrains its standard thermodynamic definition as 1/kB T and |0 represents the ground state of the system. μˆ is the total transition dipole op† erator defined as μˆ = μˆ − + μˆ + , where μˆ + = μˆ − and μˆ −,1D =
N
μj |0j |,
(10)
j =1
μˆ −,2D =
N
μj |0j | +
j =1
N−1
N
μj |kj k| + μk |j j k|.
j =1 k=j +1
(11) Once we obtain the first-order response function, the linear absorption spectra I(ω) can be calculated as ∞ I (ω) ∝ Im dt eiωt R (1) (t) . (12) 0
√ ±1/φ, φ), where φ = (1 + 5)/2, these being the same vectors used by Hein et al.37 Armed with these laser polarization vectors one can then specify the transition dipole moment, μj , used in Eqs. (10) and (11) with μj = d j · l, where d j is the transition dipole pointing from the NB to ND atom of the jth BChl. C. The scaled HEOM and high temperature approximation
The correlation function of the collective bath operator can be written as 1 Tr {e−βHB eiHB t/¯ Bˆ j e−iHB t/¯ Bˆ j } Cj (t) = TrB {e−βHB } B 1 ∞ e−iωt = dωDj (ω) , (17) π −∞ 1 − e−β¯ω
To include the random orientation of the samples when its absorption spectra is measured, the rotational average of the signal must be considered. To calculate the rotational average, instead of sampling over the whole unit sphere it can be done by averaging over the signal from three orthogonal laser polarizations.37 In 2D photon echo experiments, three laser pulses with wave vectors k1 , k2 , and k3 are sequentially applied to the analyte so to generate a third-order polarization, and this nonlinear response can be detected along a phase-matched direction.40–43 The rephasing signal is detected along the direction kRP = −k1 + k2 + k3 and its contribution to the third order nonlinear response function can be represented as
where the spectral distribution function is given by Dj (ω) = ¯−1 ξ [cj2ξ /2mj ξ ωj ξ ]δ(ω − ωj ξ ). The spectral density of the bath is assumed to be identical for each site, and we employed the Drude model for spectral density in this study,
i ˆ )μˆ × ˆ ˆ Tr μˆ − G(t ˆ× ˆ× + G(t2 )μ + G(t1 )μ − ρˆeq , 3 ¯3 (13) by using the superoperator notation Oˆ × fˆ ≡ Oˆ fˆ − fˆOˆ for any operator Oˆ and operand operator fˆ. The nonrephasing signal is detected along the direction kNR = k1 − k2 + k3 and its contribution in response function can be written as
where v0 = γ and vk = 2kπ/β¯ (k ≥ 1) are known as the Matsubara frequencies. The constants ck are given by β¯γ λγ c0 = cot −i (20) ¯ 2
(3) (t3 , t2 , t1 ) = − RRP
i ˆ )μˆ × ˆ ˆ = − 3 Tr μˆ − G(t ˆ× ˆ× + G(t2 )μ − G(t1 )μ + ρˆeq . 3 ¯ (14) After employing a double Fourier transform of Eqs. (13) and (14) for temporal variables t1 and t3 , the absorptive part of the 2D electronic spectra then can be calculated as S(ω3 , t2 , ω1 ) ≡ SRP (ω3 , t2 , ω1 ) + SNR (ω3 , t2 , ω1 ), where the rephasing contribution is ∞ ∞ (3) (3) (ω3 , t2 , ω1 ) ≡ Im dt1 dt3 ei(−ω1 t1 +ω3 t3 ) RRP (t3 , t2 , t1 ) , SRP (3) (t3 , t2 , t1 ) RNR
0
0
(15) and the nonrephasing contribution is ∞ ∞ (3) (3) SNR (ω3 , t2 , ω1 ) ≡ Im dt1 dt3 ei(ω1 t1 +ω3 t3 ) RNR (t3 , t2 , t1 ) . 0
0
(16) Because the third-order nonlinear response function possesses inversion symmetry along the laser polarization vector, the rotational average of the 2D spectra is calculated by sampling ten laser polarization vectors on vertices of a dodecahedron in half-space. The polarization vectors, l, used in this study are (±1, ±1, 1), (±φ, 0, 1/φ), (±1/φ, φ, 0), and (0,
Dj (ω) =
2λ ωγ , ¯ ω2 + γ 2
(18)
where λ = λj = ξ cj2ξ /2mj ξ ωj2ξ is the reorganization energy and γ is the Drude decay constant. The correlation function now can be written as ∞ ck e−vk t , (19) Cj (t > 0) = k=0
and ck =
vk 4λγ , for k 1 . 2 2 β¯ vk − γ 2
(21)
Using the scaled HEOM approach developed by Shi and co-workers,44 we apply the Ishizaki-Tanimura scheme27, 45 which truncates the hierarchy level at the Kth Matsubara frequency by treating higher phonon frequencies (vk>K ) with the Markovian approximation: vk e−vk t δ(t). The timeevolution of the density operator can be described by the following set of hierarchically coupled equations of motion: N K i d ρˆn = − [HS , ρˆn ] − nj k vk ρˆn dt ¯ j =1 k=0
−i
N K j =1 k=0
(nj k + 1)|ck | [Pˆj , ρˆn+j k ]
∞ cm [Pˆj , [Pˆj , ρˆn ]] − v j =1 m=K+1 m N K nj k (ck Pˆj ρˆn−j k − ck∗ ρˆn−j k Pˆj ) . −i |c | k j =1 k=0 N
(22)
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 132.206.27.24 On: Tue, 10 Feb 2015 11:52:33
234105-4
S.-H. Yeh and S. Kais
J. Chem. Phys. 141, 234105 (2014)
Within the above n is defined as a set comprised of non-negative integers, n ≡ {n1 , n2 , · · · , nN } n±j k = {{n10 , n11 · · · , n1K }, · · · , {nN0, nN1 · · · , nNK }}. refers to change in the value of njk to njk ± 1 in the global index n. The density operator with all indices equal to 0 is the system’s reduced density operator (RDO) while all other operators are the auxiliary density operators (ADOs). Here, another truncation level, NC , is introduced and is defined as NC = j,k nj k , which limits the total number of density
1 2
operators to Unfortunately, simulating 2D electronic spectra by the HEOM method requires an enormous amount of computational resource. Chen et al. developed a modified version of the HTA to HEOM by applying the Ishizaki-Tanimura truncating scheme to all the Matsubara frequencies.35 The modified HTA is capable of providing similar peak shapes and amplitude evolution to the original HEOM method, and the equations of motion can be described as N N d i ρˆ = − [HS , ρˆn ] − γ nj ρˆn − i [Pˆj , ρˆn+j ] dt n ¯ j =1 j =1
N 2λ λ β¯γ − cot [Pˆj , [Pˆj , ρˆn ]] − 2γ β¯ ¯ 2 j =1 −i
N j =1
nj (c0 Pˆj ρˆn−j − c0∗ ρˆn−j Pˆj ) .
(23)
II. RESULTS
In all simulations, the reorganization energy, λ, is set to be 35 cm−1 ; this value has been employed in many previous works.46–48 The bath response time, γ −1 , is assumed to be 100 fs−1 . Both parameters are chosen to fit the one-phonon profile from the experimental spectral density in Ref. 49. Because it is insufficient to use a single Lorentz-Drude peak to reproduce the experimental spectral density, the parameters are tuned to fit the one-phonon profile around 100– 150 cm−1 , which is about the range of energy difference between excitonic states. Arcsinh scaling is applied to all 2D spectra by applying linear scaling the maximum signal Imax = 10; the final √ signal is then obtained by calculating arcsinh(I) = ln(I + 1 + I 2 )8 . To simulate the effect of the static disorder, a Gaussian distributed noise with a standard deviation of 25 cm−1 is added to every site energy in the single-exciton Hamiltonian Eq. (2). A. Exciton delocalization comparison between the apo- and holo-FMO
Because the transfer pathway of the exciton energy is highly dependent on the spatial overlap between excitonic wave functions on each site, a comparison between the apo- and holo-FMO delocalized excitons has been completed by calculating the site (BChl) occupation probabilities associated with each exciton. These probabilities are obtained by solving for the eigenvectors of the single-exciton Hamilto-
5
7 3
2
3
5 4
7 3
E5
3
0.8
4 E4
0.6
0.4 5
7 3
E6
5
6
2
4
6 7
E3
1 5
2
4
3
6
2
5
7
4
1
1
6
2
E2
1
7 3
5
7
4
6
1
6
2
E1
1
N +(K+1)N . CN C C
1
6
0.2
4 E7
0
FIG. 1. The probability density distribution of excitonic states on each BChl in apo-FMO. The thick lines represent the excitonic couplings that are greater than 50% of the largest excitonic coupling. The thin lines and the dotted lines represent the couplings that are 30%-50% and 10%-30% of the largest coupling, respectively.
nian, Eq. (2), for each form. The probability density distribution for each excitonic state of apo- and holo-FMO are shown in Figs. 1 and 2, respectively. The calculated values for the probability density of each exciton are provided in Tables I and II for apo- and holo-FMO, respectively. Site population coefficients for each excitonic state for the 7-site system. En /Sm : The state of the nth exciton/mth site in site basis; should be read as, the degree to which the nth exciton is present on the mth site. The sites with larger populations in each exciton are shown in bold. Site population coefficients for each excitonic state for the 8-site system. En /Sm : The state of the nth exciton/mth site in site basis; should be read as, the degree to which the nth exciton is present on the mth site. The sites with larger populations in each exciton are shown in bold. From the excitonic wave function overlap of apo-FMO, two major EET pathways can be identified. The first pathway traverses through excitons 6 and 3 which are both delocalized over BChl 1 and 2. The energy will then be funnelled downstream to exciton 1, which is mostly localized on the energy sink of FMO complex–BChl 3. The second pathway features the excitonic wave packet starting from exciton 7, flowing through excitons 5 and/or 4, then passing through exciton 2 and finally to exciton 1. This result is similar to previous 1 8
8 1 2
8 1
6 5
7 3
2
4
3
8
1
6 5
7
E1
2
5
3
4 E5
1
6 5
7
4
3
2
5
3
4 E6
6 5
7
4
3
4
2
0.4 1
6 5
7 3
4 E7
0.8
0.6
E4
8 1
6 7
2
E3
8 1
6 7
2
E2
8 1
8
2
6 5
7 3
0.2
4 E8
0
FIG. 2. The probability density distribution of excitonic states on each BChl in holo-FMO. The thick lines represent the excitonic couplings that are greater than 50% of the largest excitonic coupling. The thin lines and the dotted lines represent the couplings that are 30%-50% and 10%-30% of the largest coupling, respectively.
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 132.206.27.24 On: Tue, 10 Feb 2015 11:52:33
234105-5
S.-H. Yeh and S. Kais
J. Chem. Phys. 141, 234105 (2014)
TABLE I. Squares of the eigenvector elements of apo-FMO.
Experiment apo−FMO holo−FMO
1.0
E2
E3
E4
E5
E6
E7
0.00098 0.01660 0.87520 0.10106 0.00243 0.00081 0.00292
0.00936 0.00484 0.09882 0.61940 0.08566 0.02144 0.16049
0.37792 0.56177 0.00786 0.00977 0.01492 0.02775 0.00002
0.00182 0.02734 0.00936 0.15324 0.19006 0.57384 0.04435
0.00371 0.02476 0.00457 0.05440 0.16516 0.01061 0.73681
0.60617 0.36325 0.00372 0.00003 0.00516 0.00581 0.01587
0.00004 0.00145 0.00047 0.06212 0.53661 0.35975 0.03956
0.8 I(ω) (a.u.)
S1 S2 S3 S4 S5 S6 S7
E1
0.6 0.4 0.2 0 1.2
studies using another green sulfur bacteria species, Chlorobium tepidum.48 The only major difference in our result is the initial excitonic state in each pathway. In the work of Cho et al.,48 the two pathways they obtained are: 7 → 3 → 1 and 6 → 5,4 → 2 → 1, where we find: 6 → 3 → 1 and 7 → 5,4 → 2 → 1. For the holo-FMO case, the BChl 8 modifies the EET in two ways. Due to the vicinity of BChl 8 to BChl 1, an analogous path to the first pathway discussed in the previous paragraph for the apo-FMO system (8 → 6 → 3 → 1) is now permitted. Furthermore, the second pathway is similarly affected by the enhanced delocalization of excitons 4 and 5 over BChl sites 4, 5, 6, and 7 and thus the excitonic wave function overlap between excitons 4 and 5 is increased. This larger spacial overlap may facilitate the energy transfer through the second pathway: 7 → 5,4 → 2 → 1. B. Linear absorption spectra
The linear absorption spectra of apo- and holo-FMO are simulated at 77 K. The data points are collected by evolving the system for 2048 fs with a time step of 0.25 fs. The simulated spectra are compared with the experimental data in Fig. 3. We observed that for the simulated spectra, the lowest exciton peak is underestimated in reference to the experimental peak. Also, instead of a single-peak ∼12 450 cm−1 as seen in experimental spectra, our calculation shows a double-peak structure around 12 500 cm−1 with underestimated intensity. Finally, in the holo-FMO simulated spectrum there is a significant absorption at 12 700 cm−1 , which is mostly contributed to by BChl 8. However, it is worth noting that in the crystal structure of holo-FMO BChl 8 does not occur with full occupancy.24 This means that a spectrum which includes the occupancy of BChl 8 will be bounded by the apo-FMO and
1.22
1.24 1.26 ω (cm−1)
1.28
1.3 4
x 10
FIG. 3. Simulated linear absorption spectra (with no static disorder) of apoFMO (red line) and holo-FMO (green line). The experimental spectra of P. aestuarii from Ref. 50 are shown in black circles. All the spectra are at 77 K.
the holo-FMO spectra and may better resemble the experimental data around 12 700 cm−1 . C. 2D electronic spectra without static disorder
Under the condition of disregarded static disorder, we have simulated multiple 2D electronic spectra with different waiting times, t2 , from 0 fs to 1024 fs, using an interval of 32 fs. The 2D photon-echo spectra of apo- and holo-FMO at 77 K are shown in Figs. 4 and 5, respectively. The white lines represent the eigenenergies of the single-exciton Hamiltonian. It is worth noting here that for the Hamiltonian used in this study, Eq. (2), the eigenenergies of excitons 2 and 3 are nearly degenerate to each other in both forms of FMO; in apo, holo form the energy differences are 19.33 cm−1 , 19.36 cm−1 , respectively. Fig. 4 shows a simulated 2D electronic spectrum of the apo-FMO, in which four diagonal peaks can be consistently identified throughout the simulation, these being excitons 1, 2-3, 4, and 6. The exciton 7 diagonal peak is only distinguishable within the first 300 fs of simulation. As time evolves, the EET of excitons 4,6 → 2 becomes increasingly significant; the same trend is also found in the exciton 7 → 3 transfer. Compared with the apo form, the holo form shown in Fig. 5 has a long-lived diagonal peak at exciton 7-8 region. The 2D spectra data, when compared with the previous excitonic state delocalization analysis, informs that the apo-FMO spectrum shows a broader peak around excitonic state 4; this implies an enhanced transition from excitons 5 and 4. Additionally,
TABLE II. Squares of the eigenvector elements of holo-FMO.
S1 S2 S3 S4 S5 S6 S7 S8
E1
E2
E3
E4
E5
E6
E7
E8
0.00094 0.01794 0.87897 0.09624 0.00187 0.00093 0.00312 0
0.00990 0.00787 0.09301 0.59312 0.06057 0.02093 0.21459 0.00001
0.30699 0.62873 0.00836 0.01395 0.01128 0.02835 0.00064 0.00170
0.00104 0.01932 0.01196 0.17855 0.09129 0.57981 0.11751 0.00051
0.00056 0.01984 0.00402 0.06248 0.17284 0.09980 0.63738 0.00307
0.62309 0.29432 0.00314 0 0.00053 0.00501 0.00418 0.06973
0.00038 0.00153 0.00043 0.05558 0.66160 0.26023 0.01845 0.00179
0.05709 0.01044 0.00012 0.00006 0.00001 0.00494 0.00414 0.92318
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 132.206.27.24 On: Tue, 10 Feb 2015 11:52:33
234105-6
S.-H. Yeh and S. Kais
(a)
J. Chem. Phys. 141, 234105 (2014) 2.5
(b)
(a)
2.5
(b)
2 1.25 1.5 1.23
−1 ω3 (10 4 cm )
1.27
4
−1
ω3 (10 cm )
1.27
2 1.25 1.5 1.23
1
1
1.21
1.21 0.5 (c)
0.5
(d)
(c)
1.27
(d)
1.27
1.25 −0.5 1.23 −1
0 4 −1 ω3 (10 cm )
−1
ω3 (10 4 cm )
0
1.21
1.25 −0.5 1.23 −1 1.21
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
−1.5
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
−1.5
FIG. 4. Simulated full 2D electronic spectrum of apo-FMO complex from P. aestuarii at 77 K without applied static disorder; a waiting time, t2 , of (a) 0 fs, (b) 192 fs, (c) 384 fs, and (d) 1024 fs was used.
FIG. 6. Simulated full 2D electronic spectrum of apo-FMO complex from P. aestuarii at 150 K without applied static disorder; a waiting time, t2 , of (a) 0 fs, (b) 192 fs, (c) 384 fs, and (d) 1024 fs was used.
at the lower diagonal there exists a significant cross-peak between exciton 7-8 and 2-3, which is very likely to be contributed to by the energy transfer from exciton 8 to 3. There also exist upper diagonal cross-peaks between excitons 8, 6, and 3 which is consistent with the first pathway of EET (8 → 6 → 3 → 1) in our previous analysis. In both forms, the frequencies of diagonal peak in the 2D spectra at 77 K have good correspondence to the peak frequencies calculated in the linear absorption spectra shown in Fig. 3. To investigate temperature effects on the EET, the 2D spectra of both forms of FMO have also been simulated at 150 K; the results are shown in Figs. 6 and 7. Compared to the 77 K spectrum, peaks in the 150 K spectrum are broader
due to faster dephasing between the ground state and the oneexciton states. The enhanced EET from exciton 5 to 4 in the holo-FMO is still very apparent at 150 K, and a similar trend exists for the exciton 8 → 6 → 3 transfer pathway can also be seen from the lower diagonal of the same spectrum, Fig. 7. However, it is worth noting that the upper diagonal portion of the 8 → 6 → 3 EET no longer exists at higher temperatures, as compared to 77 K spectrum of the holo-FMO. The dephasing rate associated with cross-peak (2-3, 1) were analysed by fitting the oscillation of the corresponding cross-peak (lower diagonal) from the 2D spectra. The analysis was performed by first integrating the cross-peak amplitude over squares of length 20 cm−1 , then applying a
(a)
2.5
(b)
(a)
2.5
(b)
1.27 2
1.25 1.5 1.23
−1 ω3 (10 4 cm )
−1 ω3 (10 4 cm )
1.27
2 1.25 1.5 1.23
1
1
1.21
1.21 0.5 (c)
0.5
(d)
(c)
1.27
(d)
1.27
−0.5 1.23 −1 1.21
0 4 −1 ω3 (10 cm )
4 −1 ω3 (10 cm )
0 1.25
1.25 −0.5 1.23 −1 1.21
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
−1.5
FIG. 5. Simulated full 2D electronic spectrum of holo-FMO complex from P. aestuarii at 77 K without applied static disorder; a waiting time, t2 , of (a) 0 fs, (b) 192 fs, (c) 384 fs, and (d) 1024 fs was used.
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
−1.5
FIG. 7. Simulated full 2D electronic spectrum of holo-FMO complex from P. aestuarii at 150 K without applied static disorder; a waiting time, t2 , of (a) 0 fs, (b) 192 fs, (c) 384 fs, and (d) 1024 fs was used.
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 132.206.27.24 On: Tue, 10 Feb 2015 11:52:33
(a)
S.-H. Yeh and S. Kais
Waiting time (fs) 0 150
450
J. Chem. Phys. 141, 234105 (2014)
(b)
750
1050
0.08 0.04 0 −0.04 −0.08
77 K
0.08 0.04 0 −0.04 −0.08
125 K
0.08 0.04 0 −0.04 −0.08
150 K
Waiting time (fs) 0 150
450
(a)
750
0.08 0.04 0 −0.04 −0.08
77 K
0.08 0.04 0 −0.04 −0.08
125 K
0.08 0.04 0 −0.04 −0.08
150 K
2.5
(b)
1.27
1050 −1 ω3 (10 4 cm )
234105-7
2 1.25 1.5 1.23
50
100
Temperature (K)
150
4 −1 ω3 (10 cm )
0
50 γ31(T)/T = 0.24 ± 0.03 cm /K 40 30 20 10 50
100
150
Temperature (K)
FIG. 8. Temperature dependence of the dephasing rates from (a) apo-FMO, and (b) holo-FMO. The beating signal is obtained from integrated cross-peak amplitudes after subtraction of population transfer modeled by decaying exponential function from 77 K, 125 K, and 150 K (colored lines with circles). The sum of two exponentially decaying functions were used for the fitting of the beating signal, shown in black solid lines. The dephasing rate for the 3 → 1 transition, γ 31 , and its linearly fitted temperature dependence are shown at the bottom (dashed black line).
multi-exponential fit to the signal in an effort to model the population transfer. After subtracting the population transfer portion of the signal, the remainder is fitted to the sum of two exponentially decaying sinusoidal functions (Fig. 8). Simulating results with a waiting time, T, less than 96 fs are not included to avoid the pulse overlap effects found in experimental studies.11 We have found in all cases that the dephasing rates increase linearly with rising temperature (77 K, 125 K, and 150 K are the temperatures used in this study). For the apo-FMO, the temperature dependence of the dephasing rates are γ 21 (T)/T = 0.33 ± 0.04 and γ 31 (T)/T = 0.31 ± 0.08, whereas in holo-FMO the values are γ 21 (T)/T = 0.23 ± 0.02 and γ 31 (T)/T = 0.24 ± 0.03. This linear dependence of the dephasing rate on temperature is similar to the experimental results in Ref. 11; the difference being that our values are smaller. This suggests that under the condition of experiments11 the bath response time, γ −1 , should be less than 100 fs.
D. 2D electronic spectra with static disorder
To investigate how static disorder affects the EET, a Gaussian function noise distribution with standard deviation of 25 cm−1 is applied to every site energy, j , of Eq. (2). The simulated spectrum are averaged over 1000 samples of different laser polarizations along with static disorder. Figs. 9 and 10 show the result of apo- and holo-FMO, respectively, at 77 K with various waiting times as described in the captions.
1.25 −0.5 1.23 −1 1.21
−1
0
(d)
1.27
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
−1.5
FIG. 9. Simulated full 2D electronic spectrum of apo-FMO complex from P. aestuarii at 77 K with an applied static disorder of 25 cm−1 . The waiting time, t2 , are taken at (a) 0 fs, (b) 192 fs, (c) 384 fs, and (d) 1024 fs.
The spectroscopic features in Figs. 9 and 10 are very similar to those features within the 2D spectra neglecting static disorder (Figs. 4 and 5), apart from the increased inhomogeneous broadening. Furthermore, in the upper diagonal of Fig. 10 the cross-peaks of excitons (8, 6) and (8, 3) still exists, implying that the disappearance of these cross-peaks in Fig. 7 may not merely be due to an increase in static disorder. III. DISCUSSIONS AND CONCLUSIONS
In this study, we generated the 2D photon-echo electronic spectra of apo- and holo-FMO at different temperatures by application of a high temperature approximation of
(a)
2.5
(b)
1.27 −1 ω3 (10 4 cm )
0
0.5 (c)
2 1.25 1.5 1.23 1 1.21 0.5 (c)
(d)
1.27 0 4 −1 ω3 (10 cm )
50 γ31(T)/T = 0.31 ± 0.08 cm /K 40 30 20 10
−1
−1
Dephasing rate (cm )
−1
Dephasing rate (cm )
1 1.21
1.25 −0.5 1.23 −1 1.21 1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
−1.5
FIG. 10. Simulated full 2D electronic spectrum of holo-FMO complex from P. aestuarii at 77 K with an applied static disorder of 25 cm−1 . The waiting time, t2 , are taken at (a) 0 fs, (b) 192 fs, (c) 384 fs, and (d) 1024 fs.
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 132.206.27.24 On: Tue, 10 Feb 2015 11:52:33
234105-8
S.-H. Yeh and S. Kais
hierarchical equation of motion. We have found that the 8th BChl can affect the excitonic energy transfer through both available pathways. From our simulated 2D spectrum, the cross-peak of exciton (8,3) shows that the 8th BChl is able to enhance the excitonic energy transfer through the 8 → 6 → 3 → 1 pathway due to an enhancement of the coupling, although the excitonic wave function overlap between exciton 3 and 1 is quite low and the transfer rate is comparatively slower than the other pathway. Along the 7 → 5,4 → 2 → 1 pathway, the 8th BChl increases the excitonic wave function overlap between exciton 4 and 5 and hence facilitates the downward energy transfer towards the energy sink–BChl 3. Introducing static disorder into the system increases the levels of inhomogeneous broadening; yet, this static interference did not affect the existence of the upper diagonal cross-peaks of exciton (8,6) and (8,3), suggesting these peaks are more sensitive to temperature than to static disorder. In both forms of the FMO complex, cross-peak oscillations–which last at minimum ∼600 fs–were observed. To compare apo- and holo-FMO, we analysed the dephasing rate of the important 3 → 1 excitonic transition and found that the 8th BChl decreases the dephasing rate slightly. This reduction in dephasing alleviates excitonic energy loss to the environment. Furthermore, we have also shown that this rate is linearly dependent to the temperature, which is in agreement with findings in other experimental and theoretical studies.11, 37, 51 ACKNOWLEDGMENTS
The authors would like to thank the ITaP Research Computing of Purdue University and support from National Science Foundation (NSF) Centers for Chemical Innovation (Grant No. CHE-1037992). 1 Photosynthetic Excitons, edited by H. van Amerongen, L. Valkunas, and R.
van Grondelle (World Scientific, Singapore, 2000). Mechanisms of Photosynthesis, edited by R. E. Blankenship (Blackwell Science, Oxford, 2002). 3 T. Renger, V. May, and O. Kuhn, Phys. Rep. 343, 137 (2001). 4 R. E. Fenna and B. W. Matthews, Nature (London) 258, 573 (1975). 5 C. Sybesma and J. M. Olson, Proc. Natl. Acad. Sci. U.S.A. 49, 248 (1963). 6 A. Ishizaki, T. R. Calhoun, G. S. Schlau-Cohen, and G. R. Fleming, Phys. Chem. Chem. Phys. 12, 7319 (2010). 7 R. van Grondelle and V. I. Novoderezhkin, Phys. Chem. Chem. Phys. 8, 793 (2006). 8 G. S. Engel, T. R. Calhoun, E. L. Read, T. K. Ahn, T. Mancal, Y. C. Cheng, R. E. Blankenship, and G. R. Fleming, Nature (London) 446, 782 (2007). 9 H. Lee, Y. C. Cheng, and G. R. Fleming, Science 316, 1462 (2007). 10 T. R. Calhoun, N. S. Ginsberg, G. S. Schlau-Cohen, Y.-C. Cheng, M. Ballottari, R. Bassi, and G. R. Fleming, J. Phys. Chem. B 113, 16291 (2009). 11 G. Panitchayangkoon, D. Hayes, K. A. Fransted, J. R. Caram, E. Harel, J. Wen, R. E. Blankenship, and G. S. Engel, Proc. Natl. Acad. Sci. U.S.A. 107, 12766 (2010). 12 E. Collini, C. Y. Wong, K. E. Wilk, P. M. G. Curmi, P. Brumer, and G. D. Scholes, Nature (London) 463, 644 (2010). 2 Molecular
J. Chem. Phys. 141, 234105 (2014) 13 K.
Huang and A. Rhys, Proc. R. Soc. London, Ser. A 204, 406 (1950).
14 M. Mohseni, P. Rebentrost, S. Lloyd, and A. Aspuru-Guzik, J. Chem. Phys.
129, 174106 (2008). Cao and R. J. Silbey, J. Phys. Chem. A 113, 13825 (2009). 16 F. Caruso, A. W. Chin, A. Datta, S. F. Huelga, and M. B. Plenio, Phys. Rev. A 81, 062346 (2010). 17 N. Christensson, H. F. Kauffmann, T. Pullerits, and T. Mancal, J. Phys. Chem. B 116, 7449 (2012). 18 J. Schulze, M. Torbjornsson, O. Kuehn, and T. Pullerits, New J. Phys. 16, 045010 (2014). 19 A. W. Chin, J. Prior, R. Rosenbach, F. Caycedo-Soler, S. F. Huelga, and M. B. Plenio, Nat. Phys. 9, 113 (2013). 20 F. D. Fuller, J. Pan, A. Gelzinis, V. Butkus, S. S. Senlik, D. E. Wilcox, C. F. Yocum, L. Valkunas, D. Abramavicius, and J. P. Ogilvie, Nat. Chem. 6, 706 (published online). 21 A. Halpin, P. J. M. Johnson, R. Tempelaar, R. S. Murphy, J. Knoester, T. L. C. Jansen, and R. J. D. Miller, Nat. Chem. 6, 196 (2014). 22 V. Tiwari, W. K. Peters, and D. M. Jonas, Proc. Natl. Acad. Sci. U.S.A. 110, 1203 (2013). 23 A. Ben-Shem, F. Frolow, and N. Nelson, FEBS Lett. 564, 274 (2004). 24 D. E. Tronrud, J. Z. Wen, L. Gay, and R. E. Blankenship, Photosynth. Res. 100, 79 (2009). 25 M. S. A. Busch, F. Mueh, M. E.-A. Madjet, and T. Renger, J. Phys. Chem. Lett. 2, 93 (2011). 26 Y. Tanimura and R. Kubo, J. Phys. Soc. Jpn. 58, 101 (1989). 27 Y. Tanimura, J. Phys. Soc. Jpn. 75, 082001 (2006). 28 A. Ishizaki and G. R. Fleming, J. Chem. Phys. 130, 234110 (2009). 29 A. Ishizaki and G. R. Fleming, J. Chem. Phys. 130, 234111 (2009). 30 J. Strumpfer and K. Schulten, J. Chem. Phys. 131, 225101 (2009). 31 J. Strumpfer and K. Schulten, J. Chem. Phys. 134, 095102 (2011). 32 J. Zhu, S. Kais, P. Rebentrost, and A. Aspuru-Guzik, J. Phys. Chem. B 115, 1531 (2011). 33 J. Zhu, S. Kais, A. Aspuru-Guzik, S. Rodriques, B. Brock, and P. J. Love, J. Chem. Phys. 137, 074112 (2012). 34 S.-H. Yeh, J. Zhu, and S. Kais, J. Chem. Phys. 137, 084110 (2012). 35 L. Chen, R. Zheng, Q. Shi, and Y. Yan, J. Chem. Phys. 132, 024505 (2010). 36 Y. Jing, L. Chen, S. Bai, and Q. Shi, J. Chem. Phys. 138, 045101 (2013). 37 B. Hein, C. Kreisbeck, T. Kramer, and M. Rodriguez, New J. Phys. 14, 023018 (2012). 38 S. Mukamel, Principles of Nonlinear Optical Spectroscopy (Oxford University Press, New York, 1995). 39 P. Hamm and M. T. Zanni, Concepts and Methods of 2D Infrared Spectroscopy (Cambridge University Press, 2011). 40 M. Khalil, N. Demirdoven, and A. Tokmakoff, J. Phys. Chem. A 107, 5258 (2003). 41 M. Cho, J. Chem. Phys. 115, 4424 (2001). 42 M. Cho, Chem. Rev. 108, 1331 (2008). 43 D. Abramavicius, B. Palmieri, D. V. Voronine, F. Sanda, and S. Mukamel, Chem. Rev. 109, 2350 (2009). 44 Q. Shi, L. P. Chen, G. J. Nan, R. X. Xu, and Y. J. Yan, J. Chem. Phys. 130, 084105 (2009). 45 A. Ishizaki and Y. Tanimura, J. Phys. Soc. Jpn. 74, 3131 (2005). 46 L. Chen, R. Zheng, Y. Jing, and Q. Shi, J. Chem. Phys. 134, 194508 (2011). 47 A. Ishizaki and G. R. Fleming, Proc. Natl. Acad. Sci. U.S.A. 106, 17255 (2009). 48 M. H. Cho, H. M. Vaswani, T. Brixner, J. Stenger, and G. R. Fleming, J. Phys. Chem. B 109, 10542 (2005). 49 M. Wendling, T. Pullerits, M. Przyjalgowski, S. Vulto, T. Aartsma, R. van Grondelle, and H. van Amerongen, J. Phys. Chem. B 104, 5825 (2000). 50 M. Wendling, M. Przyjalgowski, D. Gulen, S. Vulto, T. Aartsma, R. van Grondelle, and H. van Amerongen, Photosynth. Res. 71, 99 (2002). 51 P. Rebentrost, M. Mohseni, I. Kassal, S. Lloyd, and A. Aspuru-Guzik, New J. Phys. 11, 033003 (2009). 15 J.
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 132.206.27.24 On: Tue, 10 Feb 2015 11:52:33
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 132.206.27.24 On: Tue, 10 Feb 2015 11:52:33
THE JOURNAL OF CHEMICAL PHYSICS 141, 234105 (2014)
Simulated two-dimensional electronic spectroscopy of the eight-bacteriochlorophyll FMO complex Shu-Hao Yeh1 and Sabre Kais1,2,a) 1
Department of Chemistry and Birck Nanotechnology Center, Purdue University, West Lafayette, Indiana 47907, USA 2 Qatar Environment and Energy Research Institute, Qatar Foundation, Doha, Qatar
(Received 13 August 2014; accepted 25 November 2014; published online 15 December 2014) The Fenna-Matthews-Olson (FMO) protein-pigment complex acts as a molecular wire conducting energy between the outer antenna system and the reaction center; it is an important photosynthetic system to study the transfer of excitonic energy. Recent crystallographic studies report the existence of an additional (eighth) bacteriochlorophyll a (BChl a) in some of the FMO monomers. To understand the functionality of this eighth BChl, we simulated the two-dimensional electronic spectra of both the 7-site (apo form) and the 8-site (holo form) variant of the FMO complex from green sulfur bacteria, Prosthecochloris aestuarii. By comparing the spectrum, it was found that the eighth BChl can affect two different excitonic energy transfer pathways: (1) it is directly involved in the first apo form pathway (6 → 3 → 1) by passing the excitonic energy to exciton 6; and (2) it facilitates an increase in the excitonic wave function overlap between excitons 4 and 5 in the second pathway (7 → 4,5 → 2 → 1) and thus increases the possible downward sampling routes across the BChls. © 2014 AIP Publishing LLC. [http://dx.doi.org/10.1063/1.4903546] I. INTRODUCTION
Solar energy is one of the most abundant energy source that can be utilized by living organisms. In plants and bacteria, the conversion of this solar electromagnetic energy into chemical energy is called photosynthesis.1–3 One of the most studied systems within photosynthesis is the FennaMatthews-Olson (FMO) protein-pigment complex in green sulfur bacteria,4–6 which funnels the excitonic energy absorbed from the outer antenna system (chlorosome/baseplate) to the reaction center. Excitonic energy transfer (EET) in such photosynthetic protein-pigment complexes has been studied by various nonlinear optical spectroscopy methods such as pump-probe, hole burning, and photon-echo.1, 7 Recent experimental studies using 2D electronic photon-echo spectroscopy reported the existence of beating patterns with oscillation periods that approximately correspond to the difference of the eigenenergies. This phenomenon has been found in several photosynthetic protein-pigment complexes at both cryogenic temperatures8–10 and room temperature,11, 12 which reveals the importance of quantum coherence in excitonic energy transfer. For the FMO complex, the oscillation time is on the order of hundreds of femtoseconds at ambient temperature and may reach picoseconds at 77 K,13 and it has been suggested that quantum coherence may be a crucial component to the reported high-efficiency energy transfer.14–16 However, the origin of this long coherence time is still an open question. Recently, several theoretical and experimental work reported that vibronic coupling is a possible explanation for such long-lasting beating.17–22 The FMO complex has a homotrimer structure where each subunit comprises seven BChl a. Recent crystallographic a) Author to whom correspondence should be addressed. Electronic mail:
[email protected] 0021-9606/2014/141(23)/234105/8/$30.00
studies have found the presence of an additional (eighth) bacteriochlorophyll a (BChl a) in some of the monomers; the eighth BChl sits in the cleft at the protein surface close to the baseplate.23, 24 Busch et al. calculated the site energies of BChls in the FMO complex based on a new crystal structure from the green sulfur bacteria, Prosthecochloris aestuarii.25 In their calculation, it was found that BChl 8 has the highest site energy, and after modeling the baseplate and the reaction center (RC) by two low dielectric constant layers ( = 4) they also found that BChl 8 was red-shifted comparing to a homogeneous high dielectric constant environment = 80 and the overlap with the baseplate fluorescence spectrum was increased. Furthermore, due to the greater vicinity of BChl 8 to the baseplate (BChl 8 being closer than any other BChl of the FMO) it has also been proposed that BChl 8 functions as a linker between the baseplate and the FMO complex. Experimental findings of quantum coherence within the photosynthetic protein-pigment complexes have also instigated great interest in the use of theoretical models for describing photosynthetic EET. Theoretical description of some particularly challenging systems – such as the FMO complex and the B850 ring of the light harvesting complex II (LH2) found in purple bacteria – have their difficulties due to the magnitude of the electronic coupling and the exciton-phonon coupling are of same order and therefore neither of them can be treated perturbatively. To tackle this problem, a methodology referred to as hierarchical equation of motion (HEOM) has been developed by Tanimura and Kubo;26, 27 a high temperature approximation (HTA) of HEOM has also been developed by Ishizaki and Fleming to describe the EET dynamics in a model dimer system.28, 29 The HEOM has already been implemented in calculations for the excitonic transfer dynamics of photosynthetic protein-pigment complexes30–34 and more recently it is also used for simulations of both the
141, 234105-1
© 2014 AIP Publishing LLC
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 132.206.27.24 On: Tue, 10 Feb 2015 11:52:33
234105-2
S.-H. Yeh and S. Kais
J. Chem. Phys. 141, 234105 (2014)
linear and the 2D electronic spectra of the above mentioned systems.35–37 The aim of this study is to understand the role of the BChl 8 in the EET pathways of the FMO complex. To achieve this, we first analyse the difference of excitonic wave function overlap among BChls between the apo (7-sites) and holo (8sites) forms and propose corresponding EET pathways. Then the HEOM was applied for simulations of the 2D photon-echo spectra for both forms by using a model Hamiltonian formed in a previous study (PDBID: 3EOJ).24 Meanwhile, we have also explored the effects of both temperature and the presence of static disorder to the EET. A. Model Hamiltonian
fined by the coupling between the non-sharing pairs. The final system Hamiltonian used for linear absorption and 2D spectra calculation are both in a block diagonal form and can be written as 0 0 , (4) HS,1D = 0 H1ex ⎞ ⎛ 0 0 0 ⎟ ⎜ (5) HS,2D = ⎝ 0 H1ex 0 ⎠, 0
(1)
where HS , HB , and HSB represent the Hamiltonian of the system, environment and system-environment coupling, respectively. We employed the Frenkel exciton model to describe the system with each site treated as a two-level system. For the system Hamiltonian, HS , both the 7-sites and 8-sites singleexciton Hamiltonian of the FMO BChls are obtained from a structural-based theoretical study by Busch et al.25 This single exciton Hamiltonian is given as ⎞ ⎛ N ⎝(εj + λj )|j j | + Jj k (|j k| + |kj |)⎠ , H1ex = k=j
j =1
(2) where the state |j represents the QY excitation of the jth BChl with site energy εj , and Jjk corresponds to the excitonic dipole-dipole coupling between BChl sites j and k. λj is the reorganization energy which is introduced to counteract the energy shift due to the system-bath coupling. It is assumed that this reorganization energy is identical for each site, therefore λj = λ. In order to calculate the two-dimensional electronic spectra, it is necessary to include both the ground state and the two-exciton manifold into consideration. An approximation made herein is that each BChl is modeled as a two-level system. In the two-exciton manifold, the system is capable of exciting two different BChls simultaneously and its Hamiltonian can be constructed based on information of the single-exciton manifold
j k|H2ex |j k = δjj δkk (εj + εk ) + δjj (1 − δkk )Jkk + δj k (1 − δkj )Jkj + δkj (1 − δj k )Jj k + δkk (1 − δjj )Jjj .
(3)
|jk represents a two-exciton state describing the simultaneous excitation of sites j and k, where for counting 1 ≤ j < k ≤ N. The first term defines the diagonal elements by summing the site energy of the two involved sites. The remaining terms describe the off-diagonal elements under the stipulation that the excitonic transition is permitted if the initial and the final states share only one common site, where the strength is de-
H2ex
where the energy of the ground state is taken to be zero. The environment Hamiltonian, HB , describes the nuclear degrees of freedom by a phonon bath,
The total Hamiltonian both for linear absorption and for 2D spectra calculations can be described through H = HS + HB + HSB ,
0
HB =
N pˆ j2ξ j =1
ξ
2mj ξ
1 + mj ξ ωj2ξ xˆj2ξ , 2
(6)
where mjξ , ωjξ , xˆj ξ , pˆ j ξ are the mass, frequency, position, and momentum operators of the ξ th bath oscillator associate with the jth BChl, respectively. For the system-environment coupling Hamiltonian, HSB , it is assumed that each BChl is linearly coupled to the vibrational modes of its own bath and therefore induces site energy fluctuations independently. This interaction Hamiltonian is given by HSB =
N j =1
cj ξ |j j |xˆj ξ =
N
Pˆj Bˆ j ,
(7)
j =1
ξ
where Pˆj = |j j | is the excitonic projection operator and Bˆ j = ξ cj ξ xˆj ξ is defined as the collective bath operator. cjξ represents the electron-phonon coupling between the jth BChl and the ξ th phonon mode. B. Optical response functions
The interaction between an electric dipole and a weak electromagnetic field can be simulated by using optical response functions. These response functions describe how a weakly perturbed nonequilibrium system deviates from its equilibrium behaviors. The linear absorption spectra can be simulated first by calculating the first-order optical response function38, 39 R (1) (t) =
i [μ(t), ˆ μ], ˆ ¯
(8)
ˆ −iH t/¯ is the dipole operator within the where μ(t) ˆ ≡ eiH t/¯ μe Heisenberg picture and · · · ≡ Tr{· · · ρˆeq }. The 2D spectra can be simulated similarly by first calculating the third-order optical response function, which takes the form R (3) (t3 , t2 , t1 ) =−
i [[[μ(t ˆ 1 + t2 + t3 ), μ(t ˆ 1 + t2 )], μ(t ˆ 1 )], μ] ˆ . ¯3
(9)
Because the excitation energy is much larger than the scale of the thermal energy (kB T) for our system of interest, it can be assumed that the system is initially equilibrated in its ground state. The initial equilibrium state can thus be written as
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 132.206.27.24 On: Tue, 10 Feb 2015 11:52:33
234105-3
S.-H. Yeh and S. Kais
J. Chem. Phys. 141, 234105 (2014)
ρˆeq = |00| ⊗ e−βHB /TrB {e−βHB }, where β refrains its standard thermodynamic definition as 1/kB T and |0 represents the ground state of the system. μˆ is the total transition dipole op† erator defined as μˆ = μˆ − + μˆ + , where μˆ + = μˆ − and μˆ −,1D =
N
μj |0j |,
(10)
j =1
μˆ −,2D =
N
μj |0j | +
j =1
N−1
N
μj |kj k| + μk |j j k|.
j =1 k=j +1
(11) Once we obtain the first-order response function, the linear absorption spectra I(ω) can be calculated as ∞ I (ω) ∝ Im dt eiωt R (1) (t) . (12) 0
√ ±1/φ, φ), where φ = (1 + 5)/2, these being the same vectors used by Hein et al.37 Armed with these laser polarization vectors one can then specify the transition dipole moment, μj , used in Eqs. (10) and (11) with μj = d j · l, where d j is the transition dipole pointing from the NB to ND atom of the jth BChl. C. The scaled HEOM and high temperature approximation
The correlation function of the collective bath operator can be written as 1 Tr {e−βHB eiHB t/¯ Bˆ j e−iHB t/¯ Bˆ j } Cj (t) = TrB {e−βHB } B 1 ∞ e−iωt = dωDj (ω) , (17) π −∞ 1 − e−β¯ω
To include the random orientation of the samples when its absorption spectra is measured, the rotational average of the signal must be considered. To calculate the rotational average, instead of sampling over the whole unit sphere it can be done by averaging over the signal from three orthogonal laser polarizations.37 In 2D photon echo experiments, three laser pulses with wave vectors k1 , k2 , and k3 are sequentially applied to the analyte so to generate a third-order polarization, and this nonlinear response can be detected along a phase-matched direction.40–43 The rephasing signal is detected along the direction kRP = −k1 + k2 + k3 and its contribution to the third order nonlinear response function can be represented as
where the spectral distribution function is given by Dj (ω) = ¯−1 ξ [cj2ξ /2mj ξ ωj ξ ]δ(ω − ωj ξ ). The spectral density of the bath is assumed to be identical for each site, and we employed the Drude model for spectral density in this study,
i ˆ )μˆ × ˆ ˆ Tr μˆ − G(t ˆ× ˆ× + G(t2 )μ + G(t1 )μ − ρˆeq , 3 ¯3 (13) by using the superoperator notation Oˆ × fˆ ≡ Oˆ fˆ − fˆOˆ for any operator Oˆ and operand operator fˆ. The nonrephasing signal is detected along the direction kNR = k1 − k2 + k3 and its contribution in response function can be written as
where v0 = γ and vk = 2kπ/β¯ (k ≥ 1) are known as the Matsubara frequencies. The constants ck are given by β¯γ λγ c0 = cot −i (20) ¯ 2
(3) (t3 , t2 , t1 ) = − RRP
i ˆ )μˆ × ˆ ˆ = − 3 Tr μˆ − G(t ˆ× ˆ× + G(t2 )μ − G(t1 )μ + ρˆeq . 3 ¯ (14) After employing a double Fourier transform of Eqs. (13) and (14) for temporal variables t1 and t3 , the absorptive part of the 2D electronic spectra then can be calculated as S(ω3 , t2 , ω1 ) ≡ SRP (ω3 , t2 , ω1 ) + SNR (ω3 , t2 , ω1 ), where the rephasing contribution is ∞ ∞ (3) (3) (ω3 , t2 , ω1 ) ≡ Im dt1 dt3 ei(−ω1 t1 +ω3 t3 ) RRP (t3 , t2 , t1 ) , SRP (3) (t3 , t2 , t1 ) RNR
0
0
(15) and the nonrephasing contribution is ∞ ∞ (3) (3) SNR (ω3 , t2 , ω1 ) ≡ Im dt1 dt3 ei(ω1 t1 +ω3 t3 ) RNR (t3 , t2 , t1 ) . 0
0
(16) Because the third-order nonlinear response function possesses inversion symmetry along the laser polarization vector, the rotational average of the 2D spectra is calculated by sampling ten laser polarization vectors on vertices of a dodecahedron in half-space. The polarization vectors, l, used in this study are (±1, ±1, 1), (±φ, 0, 1/φ), (±1/φ, φ, 0), and (0,
Dj (ω) =
2λ ωγ , ¯ ω2 + γ 2
(18)
where λ = λj = ξ cj2ξ /2mj ξ ωj2ξ is the reorganization energy and γ is the Drude decay constant. The correlation function now can be written as ∞ ck e−vk t , (19) Cj (t > 0) = k=0
and ck =
vk 4λγ , for k 1 . 2 2 β¯ vk − γ 2
(21)
Using the scaled HEOM approach developed by Shi and co-workers,44 we apply the Ishizaki-Tanimura scheme27, 45 which truncates the hierarchy level at the Kth Matsubara frequency by treating higher phonon frequencies (vk>K ) with the Markovian approximation: vk e−vk t δ(t). The timeevolution of the density operator can be described by the following set of hierarchically coupled equations of motion: N K i d ρˆn = − [HS , ρˆn ] − nj k vk ρˆn dt ¯ j =1 k=0
−i
N K j =1 k=0
(nj k + 1)|ck | [Pˆj , ρˆn+j k ]
∞ cm [Pˆj , [Pˆj , ρˆn ]] − v j =1 m=K+1 m N K nj k (ck Pˆj ρˆn−j k − ck∗ ρˆn−j k Pˆj ) . −i |c | k j =1 k=0 N
(22)
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 132.206.27.24 On: Tue, 10 Feb 2015 11:52:33
234105-4
S.-H. Yeh and S. Kais
J. Chem. Phys. 141, 234105 (2014)
Within the above n is defined as a set comprised of non-negative integers, n ≡ {n1 , n2 , · · · , nN } n±j k = {{n10 , n11 · · · , n1K }, · · · , {nN0, nN1 · · · , nNK }}. refers to change in the value of njk to njk ± 1 in the global index n. The density operator with all indices equal to 0 is the system’s reduced density operator (RDO) while all other operators are the auxiliary density operators (ADOs). Here, another truncation level, NC , is introduced and is defined as NC = j,k nj k , which limits the total number of density
1 2
operators to Unfortunately, simulating 2D electronic spectra by the HEOM method requires an enormous amount of computational resource. Chen et al. developed a modified version of the HTA to HEOM by applying the Ishizaki-Tanimura truncating scheme to all the Matsubara frequencies.35 The modified HTA is capable of providing similar peak shapes and amplitude evolution to the original HEOM method, and the equations of motion can be described as N N d i ρˆ = − [HS , ρˆn ] − γ nj ρˆn − i [Pˆj , ρˆn+j ] dt n ¯ j =1 j =1
N 2λ λ β¯γ − cot [Pˆj , [Pˆj , ρˆn ]] − 2γ β¯ ¯ 2 j =1 −i
N j =1
nj (c0 Pˆj ρˆn−j − c0∗ ρˆn−j Pˆj ) .
(23)
II. RESULTS
In all simulations, the reorganization energy, λ, is set to be 35 cm−1 ; this value has been employed in many previous works.46–48 The bath response time, γ −1 , is assumed to be 100 fs−1 . Both parameters are chosen to fit the one-phonon profile from the experimental spectral density in Ref. 49. Because it is insufficient to use a single Lorentz-Drude peak to reproduce the experimental spectral density, the parameters are tuned to fit the one-phonon profile around 100– 150 cm−1 , which is about the range of energy difference between excitonic states. Arcsinh scaling is applied to all 2D spectra by applying linear scaling the maximum signal Imax = 10; the final √ signal is then obtained by calculating arcsinh(I) = ln(I + 1 + I 2 )8 . To simulate the effect of the static disorder, a Gaussian distributed noise with a standard deviation of 25 cm−1 is added to every site energy in the single-exciton Hamiltonian Eq. (2). A. Exciton delocalization comparison between the apo- and holo-FMO
Because the transfer pathway of the exciton energy is highly dependent on the spatial overlap between excitonic wave functions on each site, a comparison between the apo- and holo-FMO delocalized excitons has been completed by calculating the site (BChl) occupation probabilities associated with each exciton. These probabilities are obtained by solving for the eigenvectors of the single-exciton Hamilto-
5
7 3
2
3
5 4
7 3
E5
3
0.8
4 E4
0.6
0.4 5
7 3
E6
5
6
2
4
6 7
E3
1 5
2
4
3
6
2
5
7
4
1
1
6
2
E2
1
7 3
5
7
4
6
1
6
2
E1
1
N +(K+1)N . CN C C
1
6
0.2
4 E7
0
FIG. 1. The probability density distribution of excitonic states on each BChl in apo-FMO. The thick lines represent the excitonic couplings that are greater than 50% of the largest excitonic coupling. The thin lines and the dotted lines represent the couplings that are 30%-50% and 10%-30% of the largest coupling, respectively.
nian, Eq. (2), for each form. The probability density distribution for each excitonic state of apo- and holo-FMO are shown in Figs. 1 and 2, respectively. The calculated values for the probability density of each exciton are provided in Tables I and II for apo- and holo-FMO, respectively. Site population coefficients for each excitonic state for the 7-site system. En /Sm : The state of the nth exciton/mth site in site basis; should be read as, the degree to which the nth exciton is present on the mth site. The sites with larger populations in each exciton are shown in bold. Site population coefficients for each excitonic state for the 8-site system. En /Sm : The state of the nth exciton/mth site in site basis; should be read as, the degree to which the nth exciton is present on the mth site. The sites with larger populations in each exciton are shown in bold. From the excitonic wave function overlap of apo-FMO, two major EET pathways can be identified. The first pathway traverses through excitons 6 and 3 which are both delocalized over BChl 1 and 2. The energy will then be funnelled downstream to exciton 1, which is mostly localized on the energy sink of FMO complex–BChl 3. The second pathway features the excitonic wave packet starting from exciton 7, flowing through excitons 5 and/or 4, then passing through exciton 2 and finally to exciton 1. This result is similar to previous 1 8
8 1 2
8 1
6 5
7 3
2
4
3
8
1
6 5
7
E1
2
5
3
4 E5
1
6 5
7
4
3
2
5
3
4 E6
6 5
7
4
3
4
2
0.4 1
6 5
7 3
4 E7
0.8
0.6
E4
8 1
6 7
2
E3
8 1
6 7
2
E2
8 1
8
2
6 5
7 3
0.2
4 E8
0
FIG. 2. The probability density distribution of excitonic states on each BChl in holo-FMO. The thick lines represent the excitonic couplings that are greater than 50% of the largest excitonic coupling. The thin lines and the dotted lines represent the couplings that are 30%-50% and 10%-30% of the largest coupling, respectively.
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 132.206.27.24 On: Tue, 10 Feb 2015 11:52:33
234105-5
S.-H. Yeh and S. Kais
J. Chem. Phys. 141, 234105 (2014)
TABLE I. Squares of the eigenvector elements of apo-FMO.
Experiment apo−FMO holo−FMO
1.0
E2
E3
E4
E5
E6
E7
0.00098 0.01660 0.87520 0.10106 0.00243 0.00081 0.00292
0.00936 0.00484 0.09882 0.61940 0.08566 0.02144 0.16049
0.37792 0.56177 0.00786 0.00977 0.01492 0.02775 0.00002
0.00182 0.02734 0.00936 0.15324 0.19006 0.57384 0.04435
0.00371 0.02476 0.00457 0.05440 0.16516 0.01061 0.73681
0.60617 0.36325 0.00372 0.00003 0.00516 0.00581 0.01587
0.00004 0.00145 0.00047 0.06212 0.53661 0.35975 0.03956
0.8 I(ω) (a.u.)
S1 S2 S3 S4 S5 S6 S7
E1
0.6 0.4 0.2 0 1.2
studies using another green sulfur bacteria species, Chlorobium tepidum.48 The only major difference in our result is the initial excitonic state in each pathway. In the work of Cho et al.,48 the two pathways they obtained are: 7 → 3 → 1 and 6 → 5,4 → 2 → 1, where we find: 6 → 3 → 1 and 7 → 5,4 → 2 → 1. For the holo-FMO case, the BChl 8 modifies the EET in two ways. Due to the vicinity of BChl 8 to BChl 1, an analogous path to the first pathway discussed in the previous paragraph for the apo-FMO system (8 → 6 → 3 → 1) is now permitted. Furthermore, the second pathway is similarly affected by the enhanced delocalization of excitons 4 and 5 over BChl sites 4, 5, 6, and 7 and thus the excitonic wave function overlap between excitons 4 and 5 is increased. This larger spacial overlap may facilitate the energy transfer through the second pathway: 7 → 5,4 → 2 → 1. B. Linear absorption spectra
The linear absorption spectra of apo- and holo-FMO are simulated at 77 K. The data points are collected by evolving the system for 2048 fs with a time step of 0.25 fs. The simulated spectra are compared with the experimental data in Fig. 3. We observed that for the simulated spectra, the lowest exciton peak is underestimated in reference to the experimental peak. Also, instead of a single-peak ∼12 450 cm−1 as seen in experimental spectra, our calculation shows a double-peak structure around 12 500 cm−1 with underestimated intensity. Finally, in the holo-FMO simulated spectrum there is a significant absorption at 12 700 cm−1 , which is mostly contributed to by BChl 8. However, it is worth noting that in the crystal structure of holo-FMO BChl 8 does not occur with full occupancy.24 This means that a spectrum which includes the occupancy of BChl 8 will be bounded by the apo-FMO and
1.22
1.24 1.26 ω (cm−1)
1.28
1.3 4
x 10
FIG. 3. Simulated linear absorption spectra (with no static disorder) of apoFMO (red line) and holo-FMO (green line). The experimental spectra of P. aestuarii from Ref. 50 are shown in black circles. All the spectra are at 77 K.
the holo-FMO spectra and may better resemble the experimental data around 12 700 cm−1 . C. 2D electronic spectra without static disorder
Under the condition of disregarded static disorder, we have simulated multiple 2D electronic spectra with different waiting times, t2 , from 0 fs to 1024 fs, using an interval of 32 fs. The 2D photon-echo spectra of apo- and holo-FMO at 77 K are shown in Figs. 4 and 5, respectively. The white lines represent the eigenenergies of the single-exciton Hamiltonian. It is worth noting here that for the Hamiltonian used in this study, Eq. (2), the eigenenergies of excitons 2 and 3 are nearly degenerate to each other in both forms of FMO; in apo, holo form the energy differences are 19.33 cm−1 , 19.36 cm−1 , respectively. Fig. 4 shows a simulated 2D electronic spectrum of the apo-FMO, in which four diagonal peaks can be consistently identified throughout the simulation, these being excitons 1, 2-3, 4, and 6. The exciton 7 diagonal peak is only distinguishable within the first 300 fs of simulation. As time evolves, the EET of excitons 4,6 → 2 becomes increasingly significant; the same trend is also found in the exciton 7 → 3 transfer. Compared with the apo form, the holo form shown in Fig. 5 has a long-lived diagonal peak at exciton 7-8 region. The 2D spectra data, when compared with the previous excitonic state delocalization analysis, informs that the apo-FMO spectrum shows a broader peak around excitonic state 4; this implies an enhanced transition from excitons 5 and 4. Additionally,
TABLE II. Squares of the eigenvector elements of holo-FMO.
S1 S2 S3 S4 S5 S6 S7 S8
E1
E2
E3
E4
E5
E6
E7
E8
0.00094 0.01794 0.87897 0.09624 0.00187 0.00093 0.00312 0
0.00990 0.00787 0.09301 0.59312 0.06057 0.02093 0.21459 0.00001
0.30699 0.62873 0.00836 0.01395 0.01128 0.02835 0.00064 0.00170
0.00104 0.01932 0.01196 0.17855 0.09129 0.57981 0.11751 0.00051
0.00056 0.01984 0.00402 0.06248 0.17284 0.09980 0.63738 0.00307
0.62309 0.29432 0.00314 0 0.00053 0.00501 0.00418 0.06973
0.00038 0.00153 0.00043 0.05558 0.66160 0.26023 0.01845 0.00179
0.05709 0.01044 0.00012 0.00006 0.00001 0.00494 0.00414 0.92318
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 132.206.27.24 On: Tue, 10 Feb 2015 11:52:33
234105-6
S.-H. Yeh and S. Kais
(a)
J. Chem. Phys. 141, 234105 (2014) 2.5
(b)
(a)
2.5
(b)
2 1.25 1.5 1.23
−1 ω3 (10 4 cm )
1.27
4
−1
ω3 (10 cm )
1.27
2 1.25 1.5 1.23
1
1
1.21
1.21 0.5 (c)
0.5
(d)
(c)
1.27
(d)
1.27
1.25 −0.5 1.23 −1
0 4 −1 ω3 (10 cm )
−1
ω3 (10 4 cm )
0
1.21
1.25 −0.5 1.23 −1 1.21
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
−1.5
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
−1.5
FIG. 4. Simulated full 2D electronic spectrum of apo-FMO complex from P. aestuarii at 77 K without applied static disorder; a waiting time, t2 , of (a) 0 fs, (b) 192 fs, (c) 384 fs, and (d) 1024 fs was used.
FIG. 6. Simulated full 2D electronic spectrum of apo-FMO complex from P. aestuarii at 150 K without applied static disorder; a waiting time, t2 , of (a) 0 fs, (b) 192 fs, (c) 384 fs, and (d) 1024 fs was used.
at the lower diagonal there exists a significant cross-peak between exciton 7-8 and 2-3, which is very likely to be contributed to by the energy transfer from exciton 8 to 3. There also exist upper diagonal cross-peaks between excitons 8, 6, and 3 which is consistent with the first pathway of EET (8 → 6 → 3 → 1) in our previous analysis. In both forms, the frequencies of diagonal peak in the 2D spectra at 77 K have good correspondence to the peak frequencies calculated in the linear absorption spectra shown in Fig. 3. To investigate temperature effects on the EET, the 2D spectra of both forms of FMO have also been simulated at 150 K; the results are shown in Figs. 6 and 7. Compared to the 77 K spectrum, peaks in the 150 K spectrum are broader
due to faster dephasing between the ground state and the oneexciton states. The enhanced EET from exciton 5 to 4 in the holo-FMO is still very apparent at 150 K, and a similar trend exists for the exciton 8 → 6 → 3 transfer pathway can also be seen from the lower diagonal of the same spectrum, Fig. 7. However, it is worth noting that the upper diagonal portion of the 8 → 6 → 3 EET no longer exists at higher temperatures, as compared to 77 K spectrum of the holo-FMO. The dephasing rate associated with cross-peak (2-3, 1) were analysed by fitting the oscillation of the corresponding cross-peak (lower diagonal) from the 2D spectra. The analysis was performed by first integrating the cross-peak amplitude over squares of length 20 cm−1 , then applying a
(a)
2.5
(b)
(a)
2.5
(b)
1.27 2
1.25 1.5 1.23
−1 ω3 (10 4 cm )
−1 ω3 (10 4 cm )
1.27
2 1.25 1.5 1.23
1
1
1.21
1.21 0.5 (c)
0.5
(d)
(c)
1.27
(d)
1.27
−0.5 1.23 −1 1.21
0 4 −1 ω3 (10 cm )
4 −1 ω3 (10 cm )
0 1.25
1.25 −0.5 1.23 −1 1.21
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
−1.5
FIG. 5. Simulated full 2D electronic spectrum of holo-FMO complex from P. aestuarii at 77 K without applied static disorder; a waiting time, t2 , of (a) 0 fs, (b) 192 fs, (c) 384 fs, and (d) 1024 fs was used.
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
−1.5
FIG. 7. Simulated full 2D electronic spectrum of holo-FMO complex from P. aestuarii at 150 K without applied static disorder; a waiting time, t2 , of (a) 0 fs, (b) 192 fs, (c) 384 fs, and (d) 1024 fs was used.
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 132.206.27.24 On: Tue, 10 Feb 2015 11:52:33
(a)
S.-H. Yeh and S. Kais
Waiting time (fs) 0 150
450
J. Chem. Phys. 141, 234105 (2014)
(b)
750
1050
0.08 0.04 0 −0.04 −0.08
77 K
0.08 0.04 0 −0.04 −0.08
125 K
0.08 0.04 0 −0.04 −0.08
150 K
Waiting time (fs) 0 150
450
(a)
750
0.08 0.04 0 −0.04 −0.08
77 K
0.08 0.04 0 −0.04 −0.08
125 K
0.08 0.04 0 −0.04 −0.08
150 K
2.5
(b)
1.27
1050 −1 ω3 (10 4 cm )
234105-7
2 1.25 1.5 1.23
50
100
Temperature (K)
150
4 −1 ω3 (10 cm )
0
50 γ31(T)/T = 0.24 ± 0.03 cm /K 40 30 20 10 50
100
150
Temperature (K)
FIG. 8. Temperature dependence of the dephasing rates from (a) apo-FMO, and (b) holo-FMO. The beating signal is obtained from integrated cross-peak amplitudes after subtraction of population transfer modeled by decaying exponential function from 77 K, 125 K, and 150 K (colored lines with circles). The sum of two exponentially decaying functions were used for the fitting of the beating signal, shown in black solid lines. The dephasing rate for the 3 → 1 transition, γ 31 , and its linearly fitted temperature dependence are shown at the bottom (dashed black line).
multi-exponential fit to the signal in an effort to model the population transfer. After subtracting the population transfer portion of the signal, the remainder is fitted to the sum of two exponentially decaying sinusoidal functions (Fig. 8). Simulating results with a waiting time, T, less than 96 fs are not included to avoid the pulse overlap effects found in experimental studies.11 We have found in all cases that the dephasing rates increase linearly with rising temperature (77 K, 125 K, and 150 K are the temperatures used in this study). For the apo-FMO, the temperature dependence of the dephasing rates are γ 21 (T)/T = 0.33 ± 0.04 and γ 31 (T)/T = 0.31 ± 0.08, whereas in holo-FMO the values are γ 21 (T)/T = 0.23 ± 0.02 and γ 31 (T)/T = 0.24 ± 0.03. This linear dependence of the dephasing rate on temperature is similar to the experimental results in Ref. 11; the difference being that our values are smaller. This suggests that under the condition of experiments11 the bath response time, γ −1 , should be less than 100 fs.
D. 2D electronic spectra with static disorder
To investigate how static disorder affects the EET, a Gaussian function noise distribution with standard deviation of 25 cm−1 is applied to every site energy, j , of Eq. (2). The simulated spectrum are averaged over 1000 samples of different laser polarizations along with static disorder. Figs. 9 and 10 show the result of apo- and holo-FMO, respectively, at 77 K with various waiting times as described in the captions.
1.25 −0.5 1.23 −1 1.21
−1
0
(d)
1.27
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
−1.5
FIG. 9. Simulated full 2D electronic spectrum of apo-FMO complex from P. aestuarii at 77 K with an applied static disorder of 25 cm−1 . The waiting time, t2 , are taken at (a) 0 fs, (b) 192 fs, (c) 384 fs, and (d) 1024 fs.
The spectroscopic features in Figs. 9 and 10 are very similar to those features within the 2D spectra neglecting static disorder (Figs. 4 and 5), apart from the increased inhomogeneous broadening. Furthermore, in the upper diagonal of Fig. 10 the cross-peaks of excitons (8, 6) and (8, 3) still exists, implying that the disappearance of these cross-peaks in Fig. 7 may not merely be due to an increase in static disorder. III. DISCUSSIONS AND CONCLUSIONS
In this study, we generated the 2D photon-echo electronic spectra of apo- and holo-FMO at different temperatures by application of a high temperature approximation of
(a)
2.5
(b)
1.27 −1 ω3 (10 4 cm )
0
0.5 (c)
2 1.25 1.5 1.23 1 1.21 0.5 (c)
(d)
1.27 0 4 −1 ω3 (10 cm )
50 γ31(T)/T = 0.31 ± 0.08 cm /K 40 30 20 10
−1
−1
Dephasing rate (cm )
−1
Dephasing rate (cm )
1 1.21
1.25 −0.5 1.23 −1 1.21 1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
1.21
1.23 1.25 4 −1 ω1 (10 cm )
1.27
−1.5
FIG. 10. Simulated full 2D electronic spectrum of holo-FMO complex from P. aestuarii at 77 K with an applied static disorder of 25 cm−1 . The waiting time, t2 , are taken at (a) 0 fs, (b) 192 fs, (c) 384 fs, and (d) 1024 fs.
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 132.206.27.24 On: Tue, 10 Feb 2015 11:52:33
234105-8
S.-H. Yeh and S. Kais
hierarchical equation of motion. We have found that the 8th BChl can affect the excitonic energy transfer through both available pathways. From our simulated 2D spectrum, the cross-peak of exciton (8,3) shows that the 8th BChl is able to enhance the excitonic energy transfer through the 8 → 6 → 3 → 1 pathway due to an enhancement of the coupling, although the excitonic wave function overlap between exciton 3 and 1 is quite low and the transfer rate is comparatively slower than the other pathway. Along the 7 → 5,4 → 2 → 1 pathway, the 8th BChl increases the excitonic wave function overlap between exciton 4 and 5 and hence facilitates the downward energy transfer towards the energy sink–BChl 3. Introducing static disorder into the system increases the levels of inhomogeneous broadening; yet, this static interference did not affect the existence of the upper diagonal cross-peaks of exciton (8,6) and (8,3), suggesting these peaks are more sensitive to temperature than to static disorder. In both forms of the FMO complex, cross-peak oscillations–which last at minimum ∼600 fs–were observed. To compare apo- and holo-FMO, we analysed the dephasing rate of the important 3 → 1 excitonic transition and found that the 8th BChl decreases the dephasing rate slightly. This reduction in dephasing alleviates excitonic energy loss to the environment. Furthermore, we have also shown that this rate is linearly dependent to the temperature, which is in agreement with findings in other experimental and theoretical studies.11, 37, 51 ACKNOWLEDGMENTS
The authors would like to thank the ITaP Research Computing of Purdue University and support from National Science Foundation (NSF) Centers for Chemical Innovation (Grant No. CHE-1037992). 1 Photosynthetic Excitons, edited by H. van Amerongen, L. Valkunas, and R.
van Grondelle (World Scientific, Singapore, 2000). Mechanisms of Photosynthesis, edited by R. E. Blankenship (Blackwell Science, Oxford, 2002). 3 T. Renger, V. May, and O. Kuhn, Phys. Rep. 343, 137 (2001). 4 R. E. Fenna and B. W. Matthews, Nature (London) 258, 573 (1975). 5 C. Sybesma and J. M. Olson, Proc. Natl. Acad. Sci. U.S.A. 49, 248 (1963). 6 A. Ishizaki, T. R. Calhoun, G. S. Schlau-Cohen, and G. R. Fleming, Phys. Chem. Chem. Phys. 12, 7319 (2010). 7 R. van Grondelle and V. I. Novoderezhkin, Phys. Chem. Chem. Phys. 8, 793 (2006). 8 G. S. Engel, T. R. Calhoun, E. L. Read, T. K. Ahn, T. Mancal, Y. C. Cheng, R. E. Blankenship, and G. R. Fleming, Nature (London) 446, 782 (2007). 9 H. Lee, Y. C. Cheng, and G. R. Fleming, Science 316, 1462 (2007). 10 T. R. Calhoun, N. S. Ginsberg, G. S. Schlau-Cohen, Y.-C. Cheng, M. Ballottari, R. Bassi, and G. R. Fleming, J. Phys. Chem. B 113, 16291 (2009). 11 G. Panitchayangkoon, D. Hayes, K. A. Fransted, J. R. Caram, E. Harel, J. Wen, R. E. Blankenship, and G. S. Engel, Proc. Natl. Acad. Sci. U.S.A. 107, 12766 (2010). 12 E. Collini, C. Y. Wong, K. E. Wilk, P. M. G. Curmi, P. Brumer, and G. D. Scholes, Nature (London) 463, 644 (2010). 2 Molecular
J. Chem. Phys. 141, 234105 (2014) 13 K.
Huang and A. Rhys, Proc. R. Soc. London, Ser. A 204, 406 (1950).
14 M. Mohseni, P. Rebentrost, S. Lloyd, and A. Aspuru-Guzik, J. Chem. Phys.
129, 174106 (2008). Cao and R. J. Silbey, J. Phys. Chem. A 113, 13825 (2009). 16 F. Caruso, A. W. Chin, A. Datta, S. F. Huelga, and M. B. Plenio, Phys. Rev. A 81, 062346 (2010). 17 N. Christensson, H. F. Kauffmann, T. Pullerits, and T. Mancal, J. Phys. Chem. B 116, 7449 (2012). 18 J. Schulze, M. Torbjornsson, O. Kuehn, and T. Pullerits, New J. Phys. 16, 045010 (2014). 19 A. W. Chin, J. Prior, R. Rosenbach, F. Caycedo-Soler, S. F. Huelga, and M. B. Plenio, Nat. Phys. 9, 113 (2013). 20 F. D. Fuller, J. Pan, A. Gelzinis, V. Butkus, S. S. Senlik, D. E. Wilcox, C. F. Yocum, L. Valkunas, D. Abramavicius, and J. P. Ogilvie, Nat. Chem. 6, 706 (published online). 21 A. Halpin, P. J. M. Johnson, R. Tempelaar, R. S. Murphy, J. Knoester, T. L. C. Jansen, and R. J. D. Miller, Nat. Chem. 6, 196 (2014). 22 V. Tiwari, W. K. Peters, and D. M. Jonas, Proc. Natl. Acad. Sci. U.S.A. 110, 1203 (2013). 23 A. Ben-Shem, F. Frolow, and N. Nelson, FEBS Lett. 564, 274 (2004). 24 D. E. Tronrud, J. Z. Wen, L. Gay, and R. E. Blankenship, Photosynth. Res. 100, 79 (2009). 25 M. S. A. Busch, F. Mueh, M. E.-A. Madjet, and T. Renger, J. Phys. Chem. Lett. 2, 93 (2011). 26 Y. Tanimura and R. Kubo, J. Phys. Soc. Jpn. 58, 101 (1989). 27 Y. Tanimura, J. Phys. Soc. Jpn. 75, 082001 (2006). 28 A. Ishizaki and G. R. Fleming, J. Chem. Phys. 130, 234110 (2009). 29 A. Ishizaki and G. R. Fleming, J. Chem. Phys. 130, 234111 (2009). 30 J. Strumpfer and K. Schulten, J. Chem. Phys. 131, 225101 (2009). 31 J. Strumpfer and K. Schulten, J. Chem. Phys. 134, 095102 (2011). 32 J. Zhu, S. Kais, P. Rebentrost, and A. Aspuru-Guzik, J. Phys. Chem. B 115, 1531 (2011). 33 J. Zhu, S. Kais, A. Aspuru-Guzik, S. Rodriques, B. Brock, and P. J. Love, J. Chem. Phys. 137, 074112 (2012). 34 S.-H. Yeh, J. Zhu, and S. Kais, J. Chem. Phys. 137, 084110 (2012). 35 L. Chen, R. Zheng, Q. Shi, and Y. Yan, J. Chem. Phys. 132, 024505 (2010). 36 Y. Jing, L. Chen, S. Bai, and Q. Shi, J. Chem. Phys. 138, 045101 (2013). 37 B. Hein, C. Kreisbeck, T. Kramer, and M. Rodriguez, New J. Phys. 14, 023018 (2012). 38 S. Mukamel, Principles of Nonlinear Optical Spectroscopy (Oxford University Press, New York, 1995). 39 P. Hamm and M. T. Zanni, Concepts and Methods of 2D Infrared Spectroscopy (Cambridge University Press, 2011). 40 M. Khalil, N. Demirdoven, and A. Tokmakoff, J. Phys. Chem. A 107, 5258 (2003). 41 M. Cho, J. Chem. Phys. 115, 4424 (2001). 42 M. Cho, Chem. Rev. 108, 1331 (2008). 43 D. Abramavicius, B. Palmieri, D. V. Voronine, F. Sanda, and S. Mukamel, Chem. Rev. 109, 2350 (2009). 44 Q. Shi, L. P. Chen, G. J. Nan, R. X. Xu, and Y. J. Yan, J. Chem. Phys. 130, 084105 (2009). 45 A. Ishizaki and Y. Tanimura, J. Phys. Soc. Jpn. 74, 3131 (2005). 46 L. Chen, R. Zheng, Y. Jing, and Q. Shi, J. Chem. Phys. 134, 194508 (2011). 47 A. Ishizaki and G. R. Fleming, Proc. Natl. Acad. Sci. U.S.A. 106, 17255 (2009). 48 M. H. Cho, H. M. Vaswani, T. Brixner, J. Stenger, and G. R. Fleming, J. Phys. Chem. B 109, 10542 (2005). 49 M. Wendling, T. Pullerits, M. Przyjalgowski, S. Vulto, T. Aartsma, R. van Grondelle, and H. van Amerongen, J. Phys. Chem. B 104, 5825 (2000). 50 M. Wendling, M. Przyjalgowski, D. Gulen, S. Vulto, T. Aartsma, R. van Grondelle, and H. van Amerongen, Photosynth. Res. 71, 99 (2002). 51 P. Rebentrost, M. Mohseni, I. Kassal, S. Lloyd, and A. Aspuru-Guzik, New J. Phys. 11, 033003 (2009). 15 J.
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 132.206.27.24 On: Tue, 10 Feb 2015 11:52:33
