New Research In
Physical Sciences
Social Sciences
Featured Portals
Articles by Topic
Biological Sciences
Featured Portals
Articles by Topic
 Agricultural Sciences
 Anthropology
 Applied Biological Sciences
 Biochemistry
 Biophysics and Computational Biology
 Cell Biology
 Developmental Biology
 Ecology
 Environmental Sciences
 Evolution
 Genetics
 Immunology and Inflammation
 Medical Sciences
 Microbiology
 Neuroscience
 Pharmacology
 Physiology
 Plant Biology
 Population Biology
 Psychological and Cognitive Sciences
 Sustainability Science
 Systems Biology
Doublequantum resonances and excitonscattering in coherent 2D spectroscopy of photosynthetic complexes

Communicated by Peter M. Rentzepis, University of California, Irvine, CA, April 9, 2008 (received for review December 10, 2007)
Abstract
A simulation study demonstrates how the nonlinear optical response of the Fenna–Matthews–Olson photosynthetic lightharvesting complex may be explored by a sequence of laser pulses specifically designed to probe the correlated dynamics of double excitations. Cross peaks in the 2D correlation plots of the spectra reveal projections of the doubleexciton wavefunctions onto a basis of direct products of single excitons. An alternative physical interpretation of these signals in terms of quasiparticle scattering is developed.
The photosynthetic apparatus depends on lightharvesting complexes, which absorb photons and funnel their energy to reaction centers where it is converted and stored as chemical energy (1, 2). The Fenna–Matthews–Olson (FMO) complex is the prototype photosynthetic antenna (3, 4). The complex (Fig. 1) is a trimer of identical units, each containing seven chlorophyll a chromophores embedded in the protein matrix. The structure and properties of the complex have been studied extensively over the past decade (3–7).
Electronically excited FMO complexes prepared by the absorption of a single photon are well understood, and their properties are described by the Frenkel exciton model (2, 8, 9). The elaborate excitonrelaxation pattern can be monitored by using multidimensional coherent optical spectroscopy (10–12): peak redistribution on the picosecond time scale reflects excitedstate population relaxation (7), whereas femtosecond oscillations indicate electronic coherences (13). Excitedstate lifetimes, intraband excitonrelaxation pathways (7, 14), and longlived electronic quantum coherences (13) have been reported, and the Hamiltonian parameters were refined to simulate these measurements.
In the native environment, under the intense flux of sunlight, photosynthetic complexes have multiple electronic excitations, the interactions of which cause dissipation of the excess energy (1, 2, 9, 15, 16). Biological complexes have developed various protective mechanisms for excess energy discharge to avoid overheating and damage (17, 18). Understanding the coherent manyexciton dynamics, which precedes the incoherent relaxation, is necessary for revealing the initial steps in excitation dynamics. Information about the twoexciton manifold is also important for the applications of the coherent control of excitedstate dynamics of photosynthetic complexes (19). Doubleexciton resonances are much more complicated and less studied than the single excitations. Exciton annihilation, which depends on incoherent multiexciton dynamical properties, often complicates the analysis of nonlinear optical measurements.
In this article we present simulations of an impulsive thirdorder 2D coherent spectroscopic (2DCS) technique aimed at directly probing doubleexciton features in an FMO complex. In most commonly used fourwavemixing techniques, such as pumpprobe, threepulse peak shift and photon echo, doubleexciton information is convoluted with singleexciton resonances, which complicates the analysis (15, 16, 20–23). The present technique (20, 24), analogous to doublequantum coherence techniques in multidimensional NMR (25), has shown high sensitivity to coupling patterns and high spectral resolution in vibrational excitons (11, 26–28).
When the excitons form a set of independent quasiparticles, the doubleexciton wavefunctions are given by simple products of pairs of single excitons. Because of exciton interactions, the actual wavefunctions should be represented as superpositions of such products. We show that the trails of peaks in the 2DCS spectra reflect the form of the doubleexciton wavefunction and are sensitive to its delocalized projections into the space of singleexciton products.
In the following sections we initially describe the properties of single and doubleexciton states in the FMO complex. The 2DCS pulse sequence designed for probing doubleexciton resonances and dynamics in chromophore aggregates is introduced next. We then present the simulated signals for the FMO complex. Finally, these signals are analyzed by using an alternative (quasiparticlescattering) description of doubleexciton dynamics.
Single and DoubleExciton Manifolds of the FMO Complex
We shall describe the electronic excitations of the FMO complex by using the Frenkel exciton Hamiltonian representing N twolevel chromophores: Here, h_{m} is the excitation energy of the mth chromophore, J_{mn} is the excitonic coupling between chromophores m and n, and B̂_{m}^{†} (B̂_{m}) is the creation (annihilation) operator for an exciton on the mth chromophore. These operators satisfy the Pauli commutation relations [B̂_{m}, B̂_{n}^{†}] = δ_{mn}(1 − 2B̂_{n}^{†}B̂_{n}). The dipole interaction between the complex and the optical electric field E is given by where P̂ = Σ_{m}μ_{m}(B̂_{m}^{†} + B̂_{m}), μ_{m} is the transition dipole of chromophore m.
This Hamiltonian is commonly used for describing the optical responses of coupled chromophores in aggregates (1, 2, 9, 22, 29). The eigenstates of this Hamiltonian form independent manifolds that can be classified by the number of excitations. Each manifold is obtained by diagonalizing a given block of the Hamiltonian. Thirdorder 2DCS signals only depend on the single and doubleexciton states. These states will be the focus of our study.
We denote the state where chromophore m is excited by ∣m〉 ≡ B̂_{m}^{†}∣0〉. The N singleexciton eigenstates ∣e〉 are related to ∣m〉 by the transformation matrix φ_{me}: The eigenstate creation operator is similarly given by ê^{†} = Σ_{m}φ_{me}B̂_{m}^{†} and its energy is ε_{e} = Σ_{mn} J_{nm}φ_{ne}φ_{me} (we use J_{mm} ≡ h_{m}).
The doubleexciton states ∣f〉 will be described by using a basis set of direct product of realspace excitations (PRSE) ∣mn〉 with m ≥ n. This set has M = N(N + 1)/2 elements. We then have where ν≡ mn with m ≥ n. Thus, Φ is a (M × M) transformation matrix. The eigen energies are ε_{f} = Σ_{νν′}J_{νν′}^{(2)}Φ_{ν,f}Φ_{ν′,f}, where J_{νν′}^{(2)} ≡ J_{mn,m′n′}^{(2)} = δ_{mm′}J_{nn′}+δ_{nn′}J_{mm′} + δ_{mn′}J_{nm′} + δ_{m′n}J_{mn′}. Note that for our model of hard core bosons two excitations cannot reside on the same chromophore, so the states ∣mm〉 = 2^{−1/2}B̂_{m}^{†2}∣0〉 should be excluded. To simplify the notation we include it in the basis set but require that Φ_{mm,f} ≡ 0. The Φ matrix is obtained by diagonalizing the doubleexciton block of the Hamiltonian.
Doubleexciton states may be alternatively expressed in the basis of products of singleexciton eigenstate space excitations (PESE). To that end, we introduce the boson operators [ê,ê′^{†}] = δ_{ee′}. The doubleexciton basis is ∣ee′〉 = ζ_{ee′}ê^{†}ê′^{†}∣0〉, where ζ_{ee′} = 1 + δ_{ee′}(2^{−1/2} − 1) (we note that for bosons ê^{†}∣e〉 =
The FMO complex has n = 7 singleexciton states and M = 21 doubleexciton states. The Hamiltonian parameters were acquired from previous simulations (7, 8, 30). In a recent study we examined singleexciton properties and their coherent versus incoherent dynamics (30). Here we focus on doubleexciton coherent dynamics and show how nonlinear signals can be designed to resolve doubleexciton wavefunction, localization, and scattering. The singleexciton state energies are shown in Fig. 2. The participation ratios vary between 1 and 3, indicating that the single excitons are essentially localized. The entire set of doubleexciton states is given in Fig. 2 as well. Their participation ratios indicate that fewer doubleexciton states are delocalized in the PESE than in the PRSE basis, which implies that the PESE basis is better suited for describing double excitons.
Coherent DoubleQuantum Spectroscopy of Excitons
The proposed 2DCS technique is performed with four temporally wellseparated laser pulses (Fig. 1). The optical electric field is given by The first three pulses generate a nonlinear polarization in the complex, which is heterodynedetected with the fourth pulse. We shall focus on the signal generated along the phasematching direction k_{III} ≡ +k_{1} + k_{2} − k_{3} (see ref. 31 for this notation). The signal recorded versus the three delay times between pulses t_{1}, t_{2}, and t_{3} will be denoted S(t_{3}, t_{2}, t_{1}). As is common in resonant spectroscopies, we shall invoke the rotating wave approximation (RWA) and only retain the dominant contributions to S, where all interactions are resonant. For our exciton model, there are only two contributions to the k_{III} signal. These are represented by the Feynman diagrams shown in Fig. 1. The two diagrams represent the same evolution during the first two time intervals: during t_{1} the density matrix oscillates with frequency ω_{eg} = ε_{e} − ε_{g}, and during t_{2} the density matrix oscillates with frequency ω_{fg} = ε_{f} − ε_{g}. During t_{3} the diagrams are different: the oscillation frequency is either ω_{fe′} (Fig. 1 Right, diagram B) or ω_{e′g} (Fig. 1 Right, diagram A). S(t_{3}, t_{2}, t_{1}) constitutes a 3D signal. It can be represented conveniently in the frequency domain by a triple Fourier transform with respect to the delay times: This signal is given by (32): The two terms correspond to Fig. 1 Right, diagrams B and A, respectively. Here G_{ab}(Ω) = (Ω − ω_{ab} + iγ_{ab})^{−1} is the frequency domain Green's function for density matrix coherence ∣a〉〈b∣, and γ_{ab} is the dephasing rate. In the time domain we have G_{ab}(t) = θ(t)exp(−iω_{ab}t − γ_{ab}t); angular brackets denote averaging over fluctuations caused by other degrees of freedom (e.g., phonons, solvent). ω_{j} are the carrier frequencies and ℰ_{j}(ω) are the pulse envelopes centered at ω = 0. Note that ω_{4} = ω_{1} + ω_{2} − ω_{3} is required by phase matching. We shall display some 2D sections of the complete 3D signal. This can be done either in the time or frequency domain or in a mixed representation where we replace any of the Green's function G(Ω_{j}) by its Fourier transform G(t_{j}). We shall focus on two signals. The first is S_{21} ≡ S(t_{3}; Ω_{2},Ω_{1}) displayed in (Ω_{2},Ω_{1}) space for various values of t_{3}. This signal vanishes for t_{3} = 0. The second choice will be S_{32} ≡ S(Ω_{3},Ω_{2}; t_{1}). This signal will be displayed in (Ω_{3},Ω_{2}) space. The delay time t_{1} in S_{32} induces phase rotation and does not change the peak amplitudes. Thus, we set t_{1} = 0.
DoubleExciton Resonances of the FMO Complex: 2D Signals
The signals (Eq. 11) were calculated by using the cumulant expansion for Gaussian fluctuations as implemented in the SPECTRON package, which incorporates correlated bath fluctuations (33–36). Each chromophore is coupled to its own, statistically independent bath; the fluctuation statistics of all chromophore frequencies is identical and described by the overdamped Brownian oscillator spectral density. By transforming the bath fluctuation parameters to the eigenstate basis we obtain the correlated statistical properties of fluctuations of eigenstates. All parameters are the same as those described in ref. 30 except one: numerical averaging over static disorder with 20 cm^{−1} variance (inhomogeneous linewidth) made no noticeable difference on the signal and was eliminated. The calculated homogeneous linewidth is approximately 70 cm^{−1}.
S_{21} is displayed in Fig. 3 for several values of t_{3}. Only singleexciton resonances appear along Ω_{1}, and doubleexciton states show up along Ω_{2}. The signal vanishes for t_{3} = 0 and quickly grows until t_{3} ≈ 50 fs. The subsequent variation of the signal with t_{3} reflects the evolution in the doubleexciton manifold. We present the absolute value of the signal (A) as well as its real (R) and imaginary (I) parts. The I (absorptive) part clearly shows oscillations of peaks with t_{3}, whereas A helps to identify the peaks and relate them with the exciton states. R is shown for completeness. To reveal how each doubleexciton resonance is connected to a specific set of singleexciton states, we present on the bottom row in Fig. 3 the same signal at which all linewidths were reduced by a factor of 20. Each Ω_{2} selects a given doubleexciton state, and the corresponding peaks along Ω_{1} show its projection onto the various singleexciton states, the resonant frequencies of which are marked by red arrows. The peak positions along Ω_{1} agree with the singleexciton resonances. By using the states listed in Fig. 2, we find that the spectrum is dominated by three doubleexciton states: 1, 7, and 18. Three other states (9, 16, and 17) make a weaker contribution. The participation ratios of these states (see Fig. 2) indicate that these are the most delocalized states in the PESE basis.
The contributions of various singleexciton states to each doubleexciton state can be rationalized by examining the transformation matrix φ. Ψ_{ee′,f}^{2} is the probability that the system in the doubleexciton state f be found in the pair of states ee′. These probabilities for states 1, 7, and 18 are displayed in Fig. 4[the probability distributions for all 21 states are shown in Supporting Information (SI) Figs. S1 and S2]. As suggested by the participation ratios of states 1, 7, and 18 given in Fig. 2, these states are delocalized in the PESE representation, indicating that 2D signals are sensitive to exciton delocalization. We further note that the singleexciton states contributing to the specific doubleexciton state directly correspond to the series of peaks along Ω_{1} for a fixed Ω_{2}. Thus, the S_{21} signal directly reflects the doubleexciton wavefunction in the PESE space, and the peaks for each Ω_{2} resonance reflect its projections onto the singleexciton basis.
The S_{32} signals are depicted in Fig. 5A. The amplitude spectrum clearly shows a contribution from three doubleexciton states. On the Ω_{3} axis there are now two types of resonances: ω_{fe′} and ω_{e′g}. This spectrum, therefore, carries more details than S_{21} (Fig. 3), which only has ω_{eg} resonances along Ω_{1}; however, interference between the ω_{fe′} and ω_{e′g} peaks complicates the peak assignment. We also present the same signal at which the linewidth was reduced by a factor of 20. Each peak can now be assigned to a specific set of eigenstates. We find that the doubleexciton states 1, 9, and 18 with the largest PESE participation ratios (Fig. 2) are dominant.
QuasiparticleScattering Picture of Double Excitations and the MeanField Approximation
So far, our analysis was based on the properties of single and doubleexciton wavefunctions. Single excitons carry information on the couplings between chromophores and their interaction with the environment. The doubleexciton states reflect the manybody properties: exciton–exciton interactions.
An alternative physical picture for excitons in aggregates is provided by the quasiparticle approach. Rather than computing doubleexciton wavefunctions, we view the excitons as interacting quasiparticles. All relevant properties can then be viewed in terms of their scattering. The excitonscattering matrix, rather than twoexciton wavefunction, then plays a central role in the analysis. The quasiparticle approach provides a highly intuitive physical picture of exciton dynamics. Powerful approximation schemes, stemming from the shortrange nature of exciton–exciton interactions, make the calculations much easier compared with the eigenstate calculations and provide deep insights into the nature of multiexciton dynamics.
In the quasiparticle approach, the optical signals are calculated by solving equations of motion for relevant exciton variables 〈B̂〉 and 〈B̂B̂〉. These are known as the nonlinear exciton equations (10, 12, 33, 37–39). In the singleexciton eigenstate representation of quasiparticle scattering, the response function for k_{III} technique is (12, 39, 40): where I_{e}(ω) = i(ω − ε_{e} + iγ)^{−1} and 𝒢_{ee′}(ω) = i(ω − ε_{e} − ε_{e′} + 2iγ)^{−1} are the singleexciton and noninteracting doubleexciton Green's functions, and γ is the dephasing rate.
Γ(ω) is the excitonscattering matrix: where V and ℙ are the tetradic matrices. In real space for twolevel chromophores V_{mn,kl} = 2J_{ml}(δ_{nk}δ_{kl} − δ_{mn}δ_{nk}) and ℙ_{mn,kl} = δ_{mn}δ_{nk}δ_{kl}; 1_{mnkl} ≡ δ_{mk}δ_{nl} (39). This is equivalent to sumovereigenstates (SOS) expression (11). The two approaches only differ by the dephasing model.
A meanfield approximation (MFA) is often used to simplify the excitonscattering picture (41, 42). It assumes 〈B̂_{m}B̂_{n}〉 = 〈B̂_{m}〉〈B̂_{n}〉. The MFA response is also given by Eq. 12, but the scattering matrix is much simpler: This approximation ignores the effects of exciton statistics in doubly excited resonances. Therefore, doubly excitedstate energies are simply given by sums of singleexciton excitation energies. Fewer peaks will survive the interference of the two Feynman diagrams in this case (38). Exciton scattering takes place during the intervals t_{2} and t_{3} where the particles interact. The MFA neglects pairwise interactions of excitons but includes the perturbation of an exciton pair by a third exciton.
The MFA S_{21} and S_{32} signals are displayed in Fig. 5B. As expected, they contain fewer peaks compared with complete SOS calculations (for S_{21}, see Fig. 3). The peak positions are different in the MFA and SOS; however, the doubleexciton resonances still appear between 24,500 and 25,000 cm^{−1}. Thus, the fine structure of doubleexciton states of the FMO complex in MFA is distorted, whereas the average resonance energy is the same. The peak positions along Ω_{1} can be associated with the singleexciton states, because the MFA does not affect the singleexciton dynamics. The peaks along Ω_{2} are related to combinations ε_{e} + ε_{e′}. Compared with the full calculation (Eq. 13), the MFA (Eq. 14), thus, does not miss frequency shifts caused by exciton–exciton interactions.
By expanding Eq. 13 to first order in V we get Γ(ω) ≈ V[1 − 𝒢(ω)ℙ𝒢(ω)^{−1}] − ℙ𝒢(ω)^{−1}. This result is similar but not identical to the MFA. Only when the particle statistics is a weak perturbation, such that 𝒢ℙ𝒢^{−1} ℙ, do we recover the MFA. The two are identical for bosons, where ℙ = 0 and the scattering matrix is Γ_{(b)}^{(MFA)} = V. The MFA, therefore, is an approximation in terms of coupling strength and nonbosonic nature of excitons, whereas the weak coupling limit still retains nonbosonic nature.
Finally, we note that a different strategy for treating exciton scattering based on bosonization may be used, as well. With bosonization, particle statistics is greatly simplified, and its effects are incorporated in a modified exciton Hamiltonian. Paulion operators can be exactly mapped into a boson representation (43). An approximate empirical bosonization has also been used (20, 37, 44) for calculation of multidimensional signals in the UV region of polypeptides (44, 45).
Conclusions
The various levels of description of 2DCS signals employed here provide insights into the excited state dynamics of photosynthetic complexes. Weak coupling/weak scattering limits can then be described by comparing the signals obtained by using eigenstate expressions (11), the excitonscattering picture (Eq. 12), the MFA, and the weak coupling limits. Comparison of Figs. 3 and 5A shows that doubleexciton dynamics in the FMO complex can be approximately described by the MFA: the couplings and nonbosonic nature of excitons are strong, and a full calculation of the scattering matrix is preferred.
The localized doubleexciton states in PESE affect the 2D signal only weakly, which can be rationalized by the quasiparticlescattering picture. The localized doubleexciton states ∣f〉 = ∣ee′〉 can be isolated; by energy conservation their energy must be ε_{f} = ε_{e} + ε_{e′}, which implies that scattering of excitons does not occur and such states do not contribute to the nonlinear signal. This argument does not hold in the realspace representation, because realspace single excitons are coupled by excitonic interactions at the singleexcitation level. Accordingly, the PESE basis better connects with the experiment rather than the realspace basis. The doubleexciton state delocalization, χ_{E}, is a good indicator for the strength of the doubleexciton features in the 2D spectrum: doubleexciton states with χ_{E} = 1 do not contribute to the signal. These contributions increase with χ_{E}. The requirement Φ_{mm,f} ≡ 0 reflects Pauli statistics and, therefore, is responsible for exciton scattering; two excitations cannot reside on the same molecule: the molecules are hardcore exciton scatterers. Because singleexciton eigenstates are delocalized over a range of chromophores, Ψ_{ee,f} ≠ 0; thus, double excitations of singleexciton eigenstates are possible. Their energy, however, is different than 2ε_{e}. Therefore, they are softcore exciton scatterers. We note that the system–bath interaction induces damping of density matrix coherences and have no significance for our conclusions.
The (Ω_{2},Ω_{1}) signal presented in Fig. 3 shows interesting spectral dynamics with the delay t_{3}. Technically, this time dependence results from the interference of two Feynman diagrams (see Fig. 1 and Eq. 11). After the summation over exciton states e′ and f and in the absence of dephasing, the dynamics along t_{3} will show modulations of the signal peaks. When ω_{e′g} ≈ ω_{fe′} the modulation will be very slow; thus, the characteristic modulation time scale τ_{m} = ∣ω_{e′g} − ω_{fe′}∣^{−1} determines the modulation period. Note that ω_{e′g} = ω_{fe′}, with τ = ∞ corresponding to a harmonic system (when the entire signal vanishes). In the paragraph above we related doubleexciton delocalization in the PESE basis with the signal amplitude. It follows that shorter τ_{m} implies more delocalized doubleexciton states. Both τ_{m} and the delocalization length χ_{E} depend on the nonlinear (anharmonic) part of the Hamiltonian and are signatures of exciton scattering.
The power of the signals proposed here and the quasiparticle analysis will become even more pronounced for large complexes such as PS1 (46), because the density of eigenstates in the energy interval grows rapidly with aggregate size. The localized nature of exciton scattering then becomes crucial for simulating their nonlinear optical response.
Acknowledgments
This research was supported by National Institutes of Health Grant GM59230 and National Science Foundation Grant CHE0745892.
Footnotes
 ^{‡}To whom correspondence should be addressed. Email: smukamel{at}uci.edu

Author contributions: D.A. and S.M. designed research; D.A. and D.V.V. performed research; D.A., D.V.V., and S.M. analyzed data; and D.A. wrote the paper.

↵*Present address: Physikalisches Institut, Universität Würzburg, Am Hubland, 97074 Würzburg, Germany.

The authors declare no conflict of interest.

This article contains supporting information online at www.pnas.org/cgi/content/full/0802926105/DCSupplemental.
 © 2008 by The National Academy of Sciences of the USA
References
 ↵
 Andrews DL,
 Demidov AA
 ↵
 van Amerogen H,
 Valkunas L,
 van Grondelle R
 ↵
 Olson J
 ↵
 ↵
 ↵
 Savikhin S,
 Buck DR,
 Struve WS
 ↵
 ↵
 Vulto SIE,
 et al.
 ↵
 ↵
 ↵
 Abramavicius D,
 Mukamel S
 ↵
 Abramavicius D,
 Mukamel S
 ↵
 ↵
 ↵
 ↵
 Brüggemann B
 ↵
 Holt NE,
 et al.
 ↵
 Berera R,
 et al.
 ↵
 Brüggemann B,
 May V
 ↵
 ↵
 ↵
 Knoester J,
 Agranovich VM
 ↵
 ↵
 ↵
 Ernst RR,
 Bodenhausen G,
 Wokaun A
 ↵
 ↵
 ↵
 Zhuang W,
 Abramavicius D,
 Mukamel S
 ↵
 Davydov A
 ↵
 ↵
 ↵
 Schweigert I,
 Mukamel S
 ↵
 Mukamel S,
 Abramavicius D
 ↵
 ↵
 ↵
 Abramavicius D,
 Valkunas L,
 Mukamel S
 ↵
 Leegwater JA,
 Mukamel S
 ↵
 ↵
 Mukamel S,
 Oszwaldowski R,
 Abramavicius D
 ↵
 ↵
 Mukamel S
 ↵
 Mukamel S,
 Oszwaldowski R,
 Yang L
 ↵
 ↵
 Abramavicius D,
 Zhuang W,
 Mukamel S
 ↵
 ↵
Citation Manager Formats
Sign up for Article Alerts
Jump to section
You May Also be Interested in
More Articles of This Classification
Physical Sciences
Chemistry
Related Content
 No related articles found.