## Abstract

Two types of diagrammatic approaches for the design and simulation of nonlinear optical experiments (closed-time path loops based on the wave function and double-sided Feynman diagrams for the density matrix) are presented and compared. We give guidelines for the assignment of relevant pathways and provide rules for the interpretation of existing nonlinear experiments in carotenoids.

## 1. Introduction

Multi-dimensional (MD) spectroscopy provides unique insights into the photophysical or photochemical dynamics occurring in biological complexes, by monitoring the evolution of coherences and populations during controlled time-delay periods. Proceeding from atoms through molecules and to complex biological systems with increasing size, the interpretation of these signals becomes more challenging because the molecular states become congested and many processes take place over a broad range of time scales. Overlapping signals of different origin result in complex line shapes that hamper the assignment of individual spectroscopic features to appropriate degrees of freedom and often prevent the differentiation of underlying processes. Intra- and intermolecular couplings lead to shortening of lifetimes and fast dephasing, which restrict the window of observation to the femtosecond regime. The temporal resolution required to resolve the elementary processes and intermediate structures is achieved by ultrashort laser pulses [1].

Optical signals are given by superpositions of elementary molecular pathways [2]. By the selection of phase-matching geometries, number of pulses, times and excitation frequencies or different polarization schemes of pulse sequences, interference effects between these pathways can be used to isolate groups of pathways. Additionally, by the introduction of several controlled time delays, higher dimensional spectroscopy offers the possibility to resolve correlations of molecular degrees of freedom accessible by a broad bandwidth [2,3]. Recently, it has been shown that under certain conditions, the results of nonlinear optical experiments can be used to determine the time-dependant wave function [4,5], the nuclear potential energy surfaces [6] or the time-evolution operator [7].

For the investigation of biological complexes, the large variety of nonlinear femtosecond spectroscopy and its high temporal resolution is complementary to X-ray crystallography or nuclear magnetic resonance spectroscopy. This has been demonstrated in the recent spectroscopic studies of photosynthetic-light-harvesting complexes and protein folding [8–11].

The interpretation of nonlinear experiments therefore requires a concerted theoretical and experimental effort. The enormous size of biological systems inhibits a full first-principles treatment. The design of nonlinear spectroscopy experiments necessitates a careful analysis of the induced dynamics, including population and phase relations between molecular states along the relevant pathways.

Time-dependent perturbative expansion in the external fields offers an adequate tool for the classification and simulation of nonlinear spectra. Two types of diagrammatic descriptions have been employed, based either on the wave function (closed-time path-loop diagrams) [12,13] or on the density matrix (double-sided Feynman diagrams) [2]. These result in two types of molecular pathways in Hilbert space and Liouville space, respectively. In this article, we compare the two types of descriptions and apply them to several common techniques.

## 2. Diagrammatic representations of the nonlinear response

Optical signals are related to the nonlinear polarization
2.1
where *V* is the dipole operator and |*Ψ*(*t*)〉 is the wave function of the optically driven molecule.

In the wave-function formulation, the relative time ordering between field–matter interactions with the bra and ket in equation (2.1) is not maintained. In this case, the signal can be interpreted as a loop, first proposed by Schwinger [14], and loop time flows in a counter-clockwise direction from the left branch to the right branch. The left branch of the diagram corresponds to forward evolution of the ket, while the right branch corresponds to backward evolution coming from the bra.

We illustrate this representation by considering one contribution to the signal generated in the *k*_{s}=*k*_{1}−*k*_{2}+*k*_{3} phase-matching direction depicted in figure 1. In this experiment, three laser beams, with wavevectors *k*_{1},*k*_{2} and *k*_{3}, impinge upon the molecule. The electric fields of these pulses are decomposed into positive and negative frequency components as
2.2
The molecule–field interaction Hamiltonian in the rotating wave approximation is then given by
2.3
where *V*^{†} and *V* are components of the dipole that induce excitation and deexcitation, respectively, in the molecule. The three interactions create a third-order polarization, which radiates a background-free signal in the *k*_{s} direction. In an ideal time-domain experiment, this polarization depends parametrically on the delays between laser pulses, while in a pure frequency-domain experiment, the signal is parametrized by the laser beam frequencies.

The time intervals between successive interactions ordered along the loop (not to be confused with the externally controlled interpulse delays in real time) are denoted as *s*_{i}, *i*=1,2,3…. The *n*th-order polarization induced in the signal direction at time *t* is written as a convolution of the fields with a (*n*+1)-point correlation function of the dipole operator, and can be calculated from the diagram by the application of the following rules.

— The diagram is read clockwise along the loop from the bottom left to the bottom right. The left branch corresponds to the ket, and the right branch to the bra in equation (2.1).

— The loop contains a series of interactions separated by free evolution periods

*s*_{i}.*s*_{i}are positive and vary from 0 to . Interactions are ordered along the loop, not in real (physical) time.— Each interaction with the field is represented by an arrow. Arrows pointing to the right represent coupling with with wavevector

*k*_{i}, and arrows pointing to the leftrepresent coupling with with wavevector −*k*_{i}.— As is expected intuitively, arrows pointing into the loop (ii and iii) correspond to photon absorption, and arrows pointing out of the loop (i and iv) correspond to photon emission.

— The chronologically last interaction is on the left branch of the diagram.

— The corresponding correlation function is constructed from right to left, following the loop in a clockwise fashion. Each interaction is accompanied by a dipole operator

*V*.— Between the dipole operators, we write a Green function that propagates the wave function. Forward propagation in the left branch is indicated by the retarded Green function , where

*H*_{0}is the field-free molecular Hamiltonian. Backward propagation along the right branch is given by the advanced Green function .— Each interaction with

*k*_{i}is accompanied by the negative frequency component of the*i*th pulse, with a time argument corresponding to the interaction time (defined as the observation time*t*minus the relevant time intervals). For interactions with −*k*_{i}, the positive frequency component is needed.— Diagrams representing (

*n*+1)-wave mixing have a prefactor of*i*^{n}(−1)^{NR}, where*N*_{R}is the number of interactions occurring on the right-hand side of the diagram.

Using these rules, we can immediately write the expression for the polarization corresponding to the diagram of figure 1,
2.4
where *V*_{ab} is the dipole matrix element between states *a* and *b*, *G*_{a} is Green's function representing unitary evolution of the wave function in the *a* state and 〈*x*〉 stands for 〈*ψ*_{0}|*x*|*ψ*_{0}〉 for an initial state . The molecule thus undergoes a path in Hilbert space where it is found in states *e*, *g*′ and *e*′ during the three intervals, respectively.

In the alternative (density matrix) description, the diagram represents the joint forward evolution of the ket and bra. Time increases from bottom to top. Timing between interactions with the bra and ket is fully accounted for, and therefore more diagrams are needed to describe the same process: the loop diagram from figure 1 is now split into the three double-sided diagrams shown in figure 2.

The polarization is calculated in Liouville space, where the left and right superoperators are defined by their action on ordinary Hilbert space operators by 2.5 We further define linear combinations of L/R operators by 2.6 The time evolution of the density matrix is indicated by the retarded Liouville space Green function, 2.7 where . Only forward evolution is required in this representation.

The Liouville space expression for the *n*th-order polarization induced in the signal direction at time *t* can be constructed directly from the diagram by the application of the following rules.

— Time flows from top to bottom. The left branch represents the ket, and the right branch corresponds to the bra side.

— The diagram consists of a series of interactions with periods of free evolution in between. The time intervals between interactions

*t*_{i}are positive and vary from 0 to .*t*_{i}represent intervals in real (physical) time.— Each interaction with the field is represented by an arrow. Arrows pointing to the right are labelled with a

*k*_{i}, and arrows pointing to the left are labelled with a −*k*_{i}. Arrows pointing into the diagram (vi and vii) correspond to photon absorption, and arrows pointing out of the diagram (v and viii) correspond to photon emission.— The chronologically last interaction is on the ket side of the diagram.

— The correlation function is read from right to left, starting with the first dipole interaction. Interactions on the ket (bra) side are represented by the dipole superoperator

*V*_{L}(*V*_{R}).— Between the dipoles operators, during the

*t*_{i}, we write a Green function that propagates the density matrix forward in time.— For every interaction with

*k*_{i}, the negative frequency component of the*i*th pulse is needed, with a time argument corresponding to the interaction time (defined as the observation time*t*minus the relevant time intervals). For interactions with −*k*_{i}, the positive frequency component is needed.— Diagrams representing (

*n*+1)-wave mixing have a prefactor of*i*^{n}(−1)^{NR}, where*N*_{R}is the number of interactions occurring on the right-hand side of the diagram.

Using these rules, we can write the contribution to the third-order polarization depicted in figure 2, 2.8 Each term now represents a different Liouville space pathway of the density matrix. In an ideal time-domain experiment, where the pulses are much shorter than their delays, only one of the diagrams in figure 2 will contribute, the one in which the pulses act in the order they arrive at the sample.

When all of the degrees of freedom are treated explicitly, the wave function can adequately describe the dynamics and it is preferable to work in Hilbert space with the loop diagram. The wave-function approach is most useful in the case where the system plus bath eigenstates are known explicitly. When it is necessary to eliminate some bath degrees of freedom and use a reduced density matrix, the Liouville space formulation is needed. The use of a reduced density matrix allows us to incorporate the effect of external degrees of freedom implicitly through a variety of approximate techniques. Pure dephasing can be incorporated in Liouville space phenomenologically. In the wave-function approach, it is recovered upon ensemble averaging.

## 3. Application to typical nonlinear optical signals

### (a) Four-wave mixing in the *k*_{1}−*k*_{2}+*k*_{3} direction

*k*

*k*

*k*

We first examine the various contributions to the signal detected in the *k*_{s}=*k*_{1}−*k*_{2}+*k*_{3} direction, which are depicted in figure 3. Using the rules presented in §2, we can write expressions for the signal contributions presented in figure 3 as
3.1
3.2
3.3
and
3.4

Using homodyne detection, the signal is given as |*P*^{(3)}(*t*)|^{2}. By combining the emitted polarization with a known reference beam, the ‘local oscillator’, we obtain the heterodyned signal ℑ*P*^{(3)}(*t*).

The ground and excited electronic states have been investigated in β-carotene by Hornung *et al.* [16]. Later on, these experiments were expanded to the use of shaped excitation pulses by Hauer *et al.* [15] using this technique. Pulses *k*_{1} and *k*_{2} are shaped, while pulse *k*_{3} is transform limited. The delay time between the two shaped pulses *k*_{1} and *k*_{2} is set to zero, and the signal is detected along *k*_{1}−*k*_{2}+*k*_{3}, and is a function of the delay time *τ* between *k*_{2} and *k*_{3}. Pathways *X*_{i} and *X*_{ii} interrogate ground-state vibrational coherences during the delay time *τ*, while *X*_{iii} and *X*_{iv} probe excited-state vibrational coherences. By tuning the carrier frequency of the laser pulses across to the transition frequency, the relative weight of these contributions can be varied. For the non-resonant case, contributions from *X*_{iii} and *X*_{iv} will be negligible, and the primary effect of the shaped pulse is to create a ground-state vibrational wave packet, via a process known as impulsive stimulate Raman scattering (ISRS) [17], which is then interrogated by pulse *k*_{3}. When the pulses are resonant, all contributions must be accounted for.

By modifying the spectral profile of the first two pulses, it is possible to drive vibrational coherences in certain modes, and to observe the mode-specific dephasing time scales. This was reported for β-carotene by Hornung *et al.* [16]. In the resonant case, they report a strong detection-wavelength dependence of fast oscillations that can be assigned to *S*_{2} modes, while in the non-resonant case, this behaviour was not observed. According to this, it has been argued [15] that by the substitution of the two first pulses of the pulse sequence, a signal enhancement is obtained only in the resonant case, implying that the pulse-shaping effect is restricted to the excited state wavepackets appearing in the two pathways *X*_{iii} and *X*_{iv} in figure 3.

This technique can be extended by adding an electronically resonant shaped pump pulse before the three pulse sequence; the pulse sequence and loop diagrams corresponding to the dominant contributions to the detected polarization are given in figure 4. The signal now depends on two time delays, *τ* and *T*, during which the system is in an electronic population, and transient components to the signals can be attributed purely to vibrational coherences.

Correlations between the vibrational dynamics occurring during the different time delays are thus monitored. In the contributions labelled *S*_{i} and *S*_{ii} in figure 4, the vibrational dynamics probed during both time delays occur on the same excited electronic state. From the rules in §2, we can write expressions for the signal contributions presented in figure 4,
3.5
3.6
and
3.7

Vibrational modes whose period matches, or is a whole multiple of, the interpulse delay will be driven to large amplitude motion by the pulse train [18]. More highly shaped pulse trains, in which the centre frequency and chirp of the sub pulses are further variable, can drive these modes into even larger amplitude motion and limit wave-packet spreading [19]. By selectively driving some vibrational modes in the ground state prior to electronic excitation, an insight is gained on which modes are relevant to the subsequent surface crossings [15].

### (b) Double-quantum-coherence spectroscopy

In nonlinear spectroscopy, where the electric field interacts several times with the molecular system, higher order transitions are possible if electronic states in the appropriate energetic regime exist. Double-quantum spectroscopy is known from nuclear magnetic resonance spectroscopy, where it is used to resolve chemical shifts as well as structural and dynamical properties.

The density matrix pathways contributing to this signal generated along *k*_{1}+*k*_{2}−*k*_{3} are shown in figure 5*a*. The signal contributions depicted in figure 5*a* are ℑ*P*^{(3)}(*t*), with
3.8
When pulses *k*_{1} and *k*_{2} are coincident in time, two more diagrams like figure 5*a* must be accounted for, diagrams in which pulse *k*_{2} acts first.

For comparison, we also write this polarization in Hilbert space, using the loop diagram from figure 5*b*, as
3.9

The first direct evidence for the existence and spectral contribution of the doubly excited state *S*_{2n} above the bright *S*_{2} state in β-carotene was found by Christensson *et al.* [20]. With an experimental setup using phase-stable pulses of 16 fs length with a bandwidth of 1280 cm^{−1} full width at half-maximum, centred at 18 350 cm^{−1} (545 nm) in the folded boxcar geometry, they detected the double-quantum-coherence (DQC) signal of β-carotene. By setting the time delay *t*_{1}=0, it is ensured that only the *S*_{2n} contribution of the double exciton manifold is addressed (see figure 6 for the β-carotene-level scheme). No other contributions by either excitation in the *S*_{1} state or other pathways will contribute to the signal. Further investigations of these DQC signals during the first delay *t*_{1}, where the system evolves in the electronic coherence might give further insights into the process of internal conversion to the dark *S*_{1} state. These spectra allowed detection of the doubly excited state *S*_{2n} and the according transition dipole moment *μ*_{2,2n}, which could be estimated as approximately a third of the ground-state bleach contribution, shows the significant importance of the incorporation of the doubly excited state *S*_{2n} in spectroscopic investigations.

### (c) Coherence-modulated third harmonic generation

In this technique, a non-stationary vibrational wave packet is created on the ground electronic state through interaction with an ultrashort off-resonant laser pulse via impulsive Raman scattering. A second pulse, the probe, then interacts with the system. The signal is the total photon count at 3*ω* where *ω* is the centre frequency of the laser (see figure 7 for the pulse sequence and loop diagrams). The vibrational coherence created by the pump pulse modulates the third harmonic signal, and the technique is deemed coherence modulated third harmonic generation (CM-THG). Using tight focusing conditions, Kupka *et al.* [21] have been able to distinguish bulk and surface contributions to this *χ*^{(5)} signal [21]. By using a shaped pulse train, as in the experiments of Hauer *et al*. [15], one could selectively probe different vibrational modes in the surface and bulk.

The signal contributions shown in figure 7 are 3.10 and 3.11

### (d) Time-resolved Raman signals following impulsive excitation

The relaxation dynamics of β-carotene were also studied using time-resolved femtosecond stimulated Raman spectroscopy following electronic excitation by Kukura *et al.* [22]. In this technique, an ultrashort resonant pump pulse *k*_{1} excites the system to the first optically bright state, usually denoted *S*_{2}. The vibrational spectrum is then interrogated by a combination of two off-resonant pulses, a narrowband picosecond pulse *k*_{2} and a broadband femtosecond pulse *k*_{3}. The signal is detected in the *k*_{3} direction, and frequency dispersed in a monochromator. They observed that for small time delays, the vibrational spectrum could easily be attributed to the *S*_{2} state, and this spectrum gradually changed over 200–300 fs, consistent with the internal conversion time scale. These measurements provide a unique window into the evolution of the vibrational modes correlating with the electronic dynamics. The polarization for this signal, depicted in figure 8, is given by Mukamel & Biggs [23],
3.12

Two extensions to the experiment are depicted in figure 9. In the first (figure 9*a*), the electronic excitation from figure 9*b* is replaced by an off-resonant ISRS excitation. This signal will reveal the coupling between the low-frequency vibrations excited by the first Raman interaction and the modes detected by the process. Fifth-order experiments of this kind have been carried out by Wilson *et al.* [24]; however, their interpretation is complicated by the presence of cascading third-order processes. Another extension would involve using two different pulses in place of the picosecond Raman pump (pulse *k*_{2} in figure 8*a*). In this scheme, depicted in figure 9*b*, the signal is sought in the *k*_{s}=*k*_{4}+*k*_{3}−*k*_{2} direction. This homodyne-detected signal can then be gated, as described in §3*e*, and the time and frequency resolution can be independently adjusted. The expression for can be obtained from equation (3.12) by changing to , and
3.13

### (e) Fluorescence-detected coherent signals

The phase relations between the incoming laser fields can be used to selectively measure specific pathways for the signal, even when incoherent fluorescence detection is used [25,26]. By repeating the experiment with different phases of the pulses (phase cycling), pathways can be selected, as is done for coherent signals under phase-matching conditions. To keep track of the pulse phase instead of the pulse wavevector, we represent the field in the form
3.14
i.e. to the negative frequency component , we associate the positive phase *ϕ*_{i}. The resulting phase-modulated fluorescence signal in general is a function of the three delays between pulses, and the time between the last pulse and *t* when the signal is observed.

The phase-modulated fluorescence signal-under phase-cycling conditions *Φ*=(*ϕ*_{1}−*ϕ*_{2}+*ϕ*_{3}−*ϕ*_{4}) in response to four impulsive pulses is analogous to four wave mixing in the direction *k*_{s}=*k*_{1}−*k*_{2}+*k*_{3}. The total fluorescence signal bearing the optical phase *Φ*=(+*ϕ*_{1},−*ϕ*_{2},+*ϕ*_{3},−*ϕ*_{4}) depicted in figure 10 is
3.15
where the auxillary functions are defined as
3.16
3.17
3.18
and
3.19
Here, is a projection operator for all fluorescing states. It is possible to apply perfect time-resolved, or perfect frequency-resolved gates or intermediate detection [27]. This technique can be applied to single molecules [28].

### (f) Multi-dimensional photoelectron spectroscopy

Time-resolved photoelectron spectroscopy (TRPES) has proved to be a valuable tool for probing ultrafast electronic dynamics [29–33]. In this technique, the system is subjected to two ultrashort laser pulses with a well-defined delay. The first pulse *k*_{1} brings the system to an excited electronic state and initiates non-stationary dynamics, e.g. internal conversions through one or more conical intersections [31,33] or photo-induced dissociation reactions [30,32,34,35]. The second pulse *k*_{d} then ionizes the molecule, and the ejected electrons are resolved in terms of their kinetic energy formally. The signal is closely related to fluorescence, except that electrons are detected instead of photons. The electron kinetic energy and the interpulse delay represent the two dimensions of this signal, and show the correlation between the initiated dynamics and the detected final state of the system. The pulse configuration and loop diagram for the TRPES experiment are shown in figure 11. The signal is
3.20
where is a projection operator for the final ionic molecular state corresponding to the detected electron with kinetic energy ** ε**.

An extension to TRPES, MD photoelectron spectroscopy (MDPES), was recently proposed by Rahav & Mukamel [36], who also include the effects of arbitrarily shaped pulses and quantum optical fields. Here, we look at 3DPES, where two laser pulses *k*_{1} and *k*_{2} interact with the system before the arrival of the detection pulse *k*_{d}, as in figure 12*a*. One implementation of this experiment, shown in figure 12*b*, involves keeping pulses *k*_{1} and *k*_{2} resonant with an electronic excitation. Using phase cycling or phase modulation, it should be possible to collect only pathways bearing the optical phase *ϕ*_{2}−*ϕ*_{1}, as is commonly done with fluorescence-detected measurements [37–39]. This signal will reveal the correlation between electronic coherences during the first time delay with the reaction probed during the second time delay, as well as correlations between these coherences and the final detected ionic state. Tailoring the spectral profile of the first two pulses will allow us to probe the influences of wave packet motion on the subsequent electronic dynamics. The expression corresponding to figure 12*b* can be obtained from equation (3.17) by changing to .

Another implementation of 3DPES is obtained when pulse *k*_{1} is electronically off-resonant, pulse *k*_{2} is resonant and pulse *k*_{d} is still the ionizing detection pulse, as in figure 12*c*. This gives
3.21
Here, the first interaction is an ISRS process, creating a non-stationary nuclear wave packet on the ground electronic state. During the first time delay, the wave packet moves, and this motion will influence the excited-state dynamics that occur during the second time delay. By shaping the Raman pulse to address certain vibrational modes, it should be possible to determine, for example, which modes are important in non-adiabatic surface-crossing dynamics that occur during the second time delay. Of course, the ISRS process can only drive motion in those low-frequency modes whose period of motion is less than the pulse duration [17].

## 4. Conclusions

Nonlinear spectroscopy is a valuable tool for elucidating the underlying processes in biological complexes. By extending the dimensionality of the experiments, the strength and nature of couplings between various excited states can be revealed. Also, it allows for the probing of multiple excitations that do not contribute to signals from linear spectroscopies. In this paper, we have compared various nonlinear spectroscopy techniques, and shown how a unified diagrammatic approach can aid in the interpretation of existing experiments and suggest new directions for future studies.

When all molecular degrees of freedom are considered explicitly, the wave function may be used to describe the system dynamics. In this case, the resulting signals can be described by the closed-time path-loop diagrams [13], which result in fewer contributing terms for a given technique. This representation treats ket evolution as forward propagation of the wave function, while the evolution of the bra is described by backward propagation of the wave function. When the number of degrees of freedom is large, it is advantageous to use an approximate reduced density matrix formulation. Here, the nonlinear signal must be calculated in Liouville space, and the various pathways contributing to the signal are represented by double-sided Feynman diagrams [2].

## Acknowledgements

We gratefully acknowledge financial support by the National Science Foundation (grant no. CHE-1058791), the National Institutes of Health (grant no. GM-59230), the Defense Advanced Research Projects Agency (grant no. UOFT-49606) and the Chemical Sciences, Geosciences and Biosciences Division, Office of Basic Energy Science, US Department of Energy (grant no. DE-FG02-04ER15571).

## Footnotes

One contribution of 14 to a Theo Murphy Meeting Issue ‘Quantum-coherent energy transfer: implications for biology and new energy technologies’.

- This journal is © 2012 The Royal Society