Recently, various spectroscopic techniques have been developed, which can measure the 2D response of the inter-molecular degrees of freedom of liquids in the THz regime. By employing hybrid Raman-THz pulse sequences, the inherent experimental problems of 2D-Raman spectroscopy are circumvented completely, culminating in the recent measurement of the 2D-Raman-THz responses of water and aqueous salt solutions. This review article focuses on the possibility to observe echoes in such experiments, which would directly reveal the inhomogeneity of the typically extremely blurred THz bands of liquids, and hence the heterogeneity of local structures that are transiently formed, in particular, in a hydrogen-bonding liquid such as water. The generation mechanisms of echoes in 2D-Raman-THz spectroscopy are explained, which differ from those in “conventional” 2D-IR spectroscopy in a subtle but important manner. Subsequently, the circumstances are discussed, under which echoes are expected, revealing a physical picture of the information content of an echo. That is, the echo decay reflects the lifetime of local structures in the liquid on a length scale that equals the delocalization length of the intermolecular modes. Finally, recent experimental results are reviewed from an echo perspective.

The low-frequency spectrum of liquids, in particular water and aqueous solutions, reflects the intermolecular degrees of freedom of the liquid, and as such carries a wealth of information on the intermolecular forces and dynamics at least in principle. Most importantly, the thermally excited modes are observed directly in the frequency range k B T 200 cm−1, which actually render the liquid a liquid. Unfortunately however, these low-frequency spectra are extremely blurred, both in THz absorption spectra (Fig. 1)1–4 and in Raman or optical Kerr effect (OKE) spectra,2,5–13 owing to the large anharmonicity and ultrafast dynamics of the intermolecular modes and the structural complexity of the liquid. Often, THz spectroscopy is therefore referred to as “blob spectroscopy.”14 

FIG. 1.

Absorption spectrum of liquid water in the spectral range of the intermolecular modes with the assignment of the various spectral features indicated. Adapted from J. Bertie and Z. Lan, Appl. Spectrosc. 50, 1047 (1996). ©1996 SAGE Publications, Ltd.

FIG. 1.

Absorption spectrum of liquid water in the spectral range of the intermolecular modes with the assignment of the various spectral features indicated. Adapted from J. Bertie and Z. Lan, Appl. Spectrosc. 50, 1047 (1996). ©1996 SAGE Publications, Ltd.

Close modal

The intermolecular THz absorption spectrum of pure water consists of three to some extent resolved features: the librational (hindered rotational) mode at ≈600 cm−1, the hydrogen bond stretch vibration at ≈200 cm−1, and the hydrogen bond bend vibration at ≈50 cm−1, the latter of which being evident only as a shoulder (Fig. 1). Adding salts,11,15–17 small organic molecules,12,18 or proteins3,19 to water modulates the low-frequency spectrum somewhat, again in a very blurred fashion, by either changing the hydrogen bond network of water or by introducing a new set of water-solute modes and intra-solute modes. Due to the low spectral resolution, however, the information content, which could be extracted from these spectra, remains rather limited and—even more so—controversial.

The very idea of multidimensional (2D) spectroscopy is to enhance the spectral resolution by thinning out spectra and spreading them into more than one dimension with the help of intelligently designed pulse sequences, and the NMR community has developed that concept to perfection.20,21 With regards to liquids, it will be important that 2D spectroscopy, which according to the definition given below also includes photon echo spectroscopy, can distinguish homogeneous from inhomogeneous line broadening, in contrast to 1D methods, such as spontaneous Raman, coherent anti-Stokes Raman scattering (CARS), and pump-probe spectroscopy.22 Liquids are structurally as heterogeneous as glasses, albeit on a very short time scale, where the heterogeneity facilitates inhomogeneous broadening and the fast dynamics homogeneous broadening. It is a priori not clear which of the two effects wins and whether THz spectra of liquids are so blurred because of large inhomogeneous or large homogeneous broadening. Furthermore, 2D spectroscopy might reveal couplings between various degrees of freedom, e.g., water-water modes with water-solute modes. After establishing the very concept of 2D spectroscopy for NMR,23 it became also feasible for vibrational transitions in the mid-IR24 as well as for electronic transitions in the UV/VIS25 owing to the development of laser-based pulse sources (for general accounts on optical 2D spectroscopy, see Refs. 26 and 27). It is the purpose of the present perspective article to discuss what one might expect to learn about liquids from extending vibrational 2D spectroscopy into a far-IR, or THz, regime ≲200cm−1.

For the purpose of this article, “2D spectroscopy” is defined as a multi-pulse experiment with two experimentally controllable time delays between incident laser pulses that are related to periods, in which the molecular system under study is in a “coherence state,” i.e., during which the bra and the ket of the density matrix are in different molecular levels. The meaning of that statement is expressed by the “energy ladder diagram” shown in Fig. 2 (the original term introduced by Albrecht et al.28), sometimes also called “wave mixing energy level (WMEL) diagrams,” which visualize the “coherence pathways” in response to a pulse sequence applied to the sample. In these WMEL diagrams, a solid arrow depicts a field interaction acting on the bra of the molecule’s density matrix, and a dotted arrow on the ket. The first field interaction will inevitably generate a coherence state with the ket and bra being in | 1 0 | , respectively (see Fig. 2(a)). For subsequent pulses, it depends on the particular coherence pathways whether a coherence or a population state is reached; in the example of Fig. 2(a), the system reaches a population state after the second pulse with the density matrix being | 1 1 | , and again a coherence state | 0 1 | after the third pulse.

FIG. 2.

Examples of rephasing coherence pathways in (a) 2D-THz spectroscopy, in (b) Raman-echo spectroscopy, in (c) 2D-Raman spectroscopy (i.e., the Raman-Raman-Raman pulse sequence, RRR), and in (d) the various hybrid methods with the Raman-THz-THz (RTT), the THz-Raman-THz (TRT), and the THz-THz-Raman (TTR) pulse sequences. During time periods t1 and t2, the system is in a coherence state with the corresponding oscillation frequencies indicated. In the cases of 2D-THz and Raman-echo spectroscopy, an additional population time T exists. Raman interactions are depicted in blue, and THz interactions in red.

FIG. 2.

Examples of rephasing coherence pathways in (a) 2D-THz spectroscopy, in (b) Raman-echo spectroscopy, in (c) 2D-Raman spectroscopy (i.e., the Raman-Raman-Raman pulse sequence, RRR), and in (d) the various hybrid methods with the Raman-THz-THz (RTT), the THz-Raman-THz (TRT), and the THz-THz-Raman (TTR) pulse sequences. During time periods t1 and t2, the system is in a coherence state with the corresponding oscillation frequencies indicated. In the cases of 2D-THz and Raman-echo spectroscopy, an additional population time T exists. Raman interactions are depicted in blue, and THz interactions in red.

Close modal

The prototype 2D experiment in the IR or visible spectral range uses three input pulses (“3-in-1-out”) with the most fundamental coherence pathways shown in Fig. 2(a).26,27 With typical optical parametrical amplifiers (OPAs) for the frequency conversion,29 this type of experiment can nowadays be done routinely down to frequencies of ≈1000 cm−1.30,31 Below that frequency, however, applications of 2D spectroscopy become extremely rare. There are a couple of realisations of experiments,32–37 which work with pulses in the ≈100 cm−1 range, but otherwise are conceptually the same as 2D-IR spectroscopy with three input THz pulses and one emitted THz field. In some cases, the second and third field interactions are coming from one and the same laser pulse, in which case the population time T shown in Fig. 2(a) is forced to zero. Due to the still relatively low intensity of available THz pulses, however, these experiments are so far restricted to low-lying electronic transitions in solids,32–34 or to rotational states of gas-phase molecules.37 To date, no 2D-THz experiment of a liquid sample has been realized.

To circumvent the limited intensities of available THz pulses, one might replace each resonant THz field interaction with two non-resonant field interactions from a laser pulse in the VIS/NIR spectral range, which induces a Raman process. The Raman effect is weak, but that might be overcompensated by the much higher available intensities of VIS/NIR laser pulses. When replacing all THz interactions of a 2D-THz experiment (Fig. 2(a)) with Raman processes, one ends up with what became known as Raman-echo spectroscopy (Fig. 2(b)). It is based on the 7th-order nonlinear response,38 and as such, the signal is extraordinarily weak. Nevertheless, attempts to perform that type of experiment have appeared in literature.39–41 

In a seminal paper, Tanimura and Mukamel proposed a simplified 2D-Raman experiment with one Raman process less (“2-in-1-out”), which thus relies on only the 5th-order nonlinear response (Fig. 2(c)).42 We will refer to it as the Raman-Raman-Raman (RRR) pulse sequence. Much attention was given to that paper from both the theoretical43–56 and the experimental point of view.57–63 However, it was realized very quickly that the desired 5th-order Raman signal is contaminated by cascaded 3rd-order signals48,58,61,64 (the same is probably also true for the 7th-order Raman-echo experiment39–41). After learning how to suppress these cascades, 2D-Raman spectroscopy became feasible for certain liquids such as CS2 or formamide,59–63 but has to date not been possible for water because of its small Raman cross section. Newer developments in this regard include new detection schemes to suppress cascading65 as well as 5th-order 2D-resonance-Raman spectroscopy.66,67

A pulse sequence with coherence pathways that are the same as for 5th-order 2D-Raman spectroscopy, but replacing all Raman interactions by direct THz interactions (i.e., “2-THz-in-1-THz-out”), would be related to a 2nd-order (i.e., even order) nonlinear response, which is forbidden in isotropic samples. But if one replaces only two of the three Raman processes in 5th-order 2D-Raman spectroscopy by direct THz interactions, it is an allowed 3rd-order (odd-order) technique, which leads to the 2D-Raman-THz hybrid methods initially proposed by Cho68 with three possible time-orderings: the Raman-THz-THz (RTT), the THz-Raman-THz (TRT) and the THz-THz-Raman (TTR) pulse sequences (Fig. 2(d)). The theory of 2D-Raman-THz spectroscopy has been worked out,69–78 and the experiment became feasible experimentally even for water in the RTT and TRT pulse sequences,79,80 while another experimental realisation in the TTR pulse sequence studied bromoform (CHBr3) and related liquids.78 

Owing to the only two input pulses of 5th-order 2D-Raman spectroscopy and 2D-Raman-THz spectroscopy, inevitably one interaction needs to induce a two-quantum transition; for the coherence pathways shown in Figs. 2(c) and 2(d), it is the 0 2 -transition induced by the second input-pulse. Such two-quantum transitions would be forbidden in the harmonic case; however, one might assume that this is hardly an issue for the intermolecular degrees of freedom of a liquid due to their large anharmonicity. Furthermore, the lack of a population time T in a 2-in-1-out-pulse sequence does not allow one to measure spectral diffusion or exchange processes, which became one of the most common targets of 2D-IR spectroscopy.27 

A 2D response can be presented either in the time-domain or in the frequency-domain, the latter after a 2D Fourier-transformation with respect to the two coherence times t1 and t2. The frequency-domain representation is the standard nowadays in the mid-IR or the UV/VIS spectral range, while the time-domain representation is more common in the THz regime. The principle information content of both representations is of course the same, and the choice is made mostly based on which of the two representations is more intuitive. In that regard, the most rewarding, and also most intuitive, possible outcome from a 2D experiment in the THz regime is the observation of an “echo,” which is the time-domain variant of a response that can distinguish homogeneous from inhomogeneous line broadening.

The very concept of an echo in spectroscopy has been introduced a long time ago first for NMR,81 and later for both electronic82 and vibrational transitions.83–85 In all these cases, the explanation of the echo is essentially the same: A first pulse interacting with the sample “generates a coherence,” i.e., makes the spectroscopic transitions of all molecules in the measured ensemble initially oscillate “in phase.” The oscillators then run out of phase, i.e., the “coherence dephases” as a function of time, since molecules sitting in different environments oscillate with different frequencies. A second pulse interacting with the sample at a time t1 after the first pulse (or in some cases a sequence of two pulses separated by a population time T, see Figs. 2(a) and 2(b)) is then used to “invert the coherence.” That is, molecules whose oscillation was leading just before that inversion of coherence will lag behind just afterwards, and vice versa. If the ensemble of molecules is “inhomogeneous,” for example, because they are embedded in a glass, the oscillation frequency of a given molecule will stay essentially the same during the whole pulse sequence. Hence, those molecules, which happen to be faster and hence lag behind just after the inversion of coherence, will catch up with the slower ones that lead the ensemble, and the coherence will “rephase.” Rephasing will peak after a time t2 = t1, which is why it is called an “echo” (in higher-order spectroscopic techniques that involve multi-quantum excitations, echoes might also appear at different times47,86).

Electronic or vibrational spectroscopy typically works in weak-field regime, where the interaction of the light field with the sample is treated perturbatively.87 In that case, the inversion of coherence is achieved with a term changing from e i ω t 1 to e + i ω t 2 , i.e., a flip in sign of the oscillation frequency. At first sight, the picture used in NMR appears to be different. That is, the rotation sense of a precessing spin is always the same and therefore, the oscillation frequency cannot flip sign. Rather, the action of the second π -pulse in an NMR echo experiment is to mirror the magnetisation vector, which is equivalent to taking the conjugate complex of the corresponding density matrix. This, in turn, is mathematically identical to a frequency of opposite sign during period t1.

A 3-in-1-out-pulse sequence, such as 2D-THz32–37 or 7th-order Raman-echo spectroscopy,39–41 contains certain symmetries with respect to the coherence pathways. That is, there are rephasing pathways with identical but oppositely signed oscillation frequencies ω 01 and + ω 01 during the two coherence times t1 and t2, both originating from one and the same level pair | 0 and | 1 (Figs. 2(a) and 2(b)). These rephasing pathways reveal an echo if there is an inhomogeneous distribution of frequencies ω 01 . Furthermore, a symmetric set of non-rephasing pathways always exists with oscillation frequencies + ω 01 and + ω 01 during the two coherence times, which can be obtained from Figs. 2(a) and 2(b) by interchanging the time-ordering of the first two field interactions that act either on the bra or the ket of the density matrix.

There are subtle but important differences when discussing echoes in a 2-in-1-out-pulse sequence, such as 5th-order 2D-Raman or 2D-Raman-THz spectroscopy. First, due to the even number of input pulses, there is no intrinsic rephasing pathway that oscillates between the same level pair during both coherence times. Coherence pathways of the sort shown in Figs. 2(c) and 2(d) have been classified as rephasing pathways,38,78 in the sense that the frequency ω 01 during period t1 is signed oppositely to ω 12 during period t2. However, in most accounts of 2D-Raman or 2D Raman-THz spectroscopy, e.g., in Refs. 38 and 78, state | 0 , | 1 , and | 2 are considered to be just any vibrational state of a molecule. In that case, there is no mechanism that would assure ω 01 ω 12 and more importantly, there is no mechanism that would render the environment-induced fluctuations Δ ω 01 and Δ ω 12 of these frequencies correlated. Both conditions are a prerequisite for an echo to appear at t2 = t1. The only situation when both conditions are actually fulfilled is when states | 0 , | 1 , and | 2 are the ground, first excited, and second excited states of one and the same vibrational mode, which must not deviate too much from harmonic. Furthermore, as long as one is staying within one particular normal mode, there is no way to construct complementary non-rephasing pathways that would oscillate with the same but and equally signed frequencies + ω 01 and + ω 12 as the related rephasing pathways.

The concept of echoes is much more general and does not only exist in spectroscopy. A particularly easy to understand example is temperature quench echoes.88–90 In that numerical experiment, the instantaneous temperature of a MD simulation is quenched to zero at two instances of time (at 0 ps and 0.35 ps in Fig. 3), i.e., the momenta (velocities) of all particles in the simulation box are set to zero. Right after such a quenching event, the various normal modes of the system will still be deflected, and they will continue to oscillate as the simulation continues to run, thereby converting potential energy back into kinetic energy. Hence, right after the first quenching event, all normal modes of the system will oscillate in phase with random amplitudes that reflect the random phases just before the first quenching event. With the second quenching event, the modes are synchronized again, but now with oscillation amplitudes that reflect the phases ω i t 1 just before the second quenching event (where ω i is the frequency of mode i). This is how the information on t1 is imprinted into the ensemble of modes. Monitoring the instantaneous temperature as a function of time, an echo is observed in the form of a dip at time t2 after the second quenching event, which equals the time interval t1 between both quenching events (Fig. 3).

FIG. 3.

Temperature quench echo in the MD simulation of a small protein, BPTI. Two temperature quenching events occur at times 0 ps and 0.35 ps, and an echo is observed at time t2 = t1 in the form of a dip of the instantaneous temperature. Adapted from O. M. Becker and M. Karplus, Phys. Rev. Lett. 70, 3514 (1993). Copyright 1993 The American Physical Society.

FIG. 3.

Temperature quench echo in the MD simulation of a small protein, BPTI. Two temperature quenching events occur at times 0 ps and 0.35 ps, and an echo is observed at time t2 = t1 in the form of a dip of the instantaneous temperature. Adapted from O. M. Becker and M. Karplus, Phys. Rev. Lett. 70, 3514 (1993). Copyright 1993 The American Physical Society.

Close modal

Two lessons can be learnt from that numerical experiment. First, the term “inversion of coherence,” which is the common language in the description of echoes in spectroscopy, does not appear in the theory of temperature quench echoes.88–90 A classical harmonic oscillator rotates in phase space in the same way as a spin precesses in real space. If all momenta pi in the MD simulation would be replaced by −pi upon the second perturbation, they would be mirrored in phase space precisely like a π -pulse acts on the spins in NMR. But rather, the momenta pi are set to zero by quenching the temperature to zero. An “inversion of coherence” is not required to reveal an echo; a “synchronisation” of modes is sufficient. Any external perturbation of a molecular system will synchronize modes to a certain extent; hence, a sequence of two such perturbations can potentially generate an echo. With its two perturbations, the temperature-quench echo is thus conceptually very similar to a laser induced echo using any one of the “2-in-1-out” pulse sequences (apart from “selection rules,” the temperature quench echo acts on all modes equally and thus does not have any selection rule).

For the second lesson to be learnt, we note that the MD simulation of Fig. 3 has been performed on a single protein,88,90 and the inhomogeneity leading to the echo is not related to an inhomogeneous ensemble of proteins, which would be the typical language used for NMR or IR echoes. Rather, the protein provides a quasi-continuum of vibrational states in particular in the low-frequency range due to its large size and structural complexity,91,92 which serves as an “inhomogeneous ensemble of modes” within one and the same protein. Furthermore, since the protein stays folded with about the same structure on the time scale of the experiment, the set of modes during period t1 will in essence be the same as that during period t2. Hence, an echo is in fact expected in any case for a protein. The same can be said for the Lennard-Jones glass studied in Ref. 89; with respect to the low-frequency vibrations, glasses and proteins share many common properties.91,92

We intend to study intermolecular modes of liquids with 2D-Raman-THz spectroscopy. A snapshot of the structure taken from a MD simulation of a liquid will be very similar to that of a glass; hence, their instantaneous normal modes will resemble each other as well. In both cases, they form a continuum of states, which should reveal an echo according to the discussion above. However, the structures of a liquid, of course, interconvert very quickly, in contrast to that of a glass. In essence, an echo will occur when the set of modes during period t2 is still the same as that during period t1. The decay of an echo therefore reflects the lifetime of intermolecular modes, or, in other words, the lifetime of local structures on a length scale that equals the delocalization length of the intermolecular modes.

Another important issue concerns the fact that the inter-molecular degrees of freedom of a liquid presumably are quite anharmonic. We have seen in the discussion of Figs. 2(c) and 2(d) that ground | 0 , first excited | 1 , and second excited states | 2 of a mode are needed with ω 01 ω 12 in order to form an echo; in other words, a little bit of harmonicity needs to survive. The classical analogue of that argument is that coherence can exist for a system with many coupled degrees of freedom only if they are harmonic to some degrees of approximation. To what extent the low frequency intermolecular degrees of freedom of a liquid are sufficiently harmonic for an echo to be seen is an open and exciting question.

The first demonstration of an echo in simulation work of the 2D Raman response (i.e., the RRR sequence) of water has been reported by Saito et al. in Ref. 49, see Fig. 4(a). It reveals a spike close to t1 = t2 = 0 together with a weak but distinct ridge along the diagonal t1 = t2, the latter of which has been attributed to an echo. From the frequency of the oscillations in the antidiagonal direction, it can be concluded that the echo originates from the librational mode of water at around 600 cm−1, which indeed has been shown to be coherent to a certain extent.93–95 It must however be understood that this signal has been calculated from quenched normal modes together with diagrammatic techniques. That is, snapshots of the structures of a MD simulation of liquid water have been quenched to 0 K, and harmonic normal modes have been calculated, which then entered a sum over all possible coherence pathways, one of which is exemplified in Fig. 2(c). This procedure in essence generates a glassy state of water, which inevitably reveals an echo, as we have seen above. The echo nevertheless remains relatively weak since many of the possible coherence pathways mix two normal modes, i.e., coherence pathways, in which states | 0 , | 1 , and | 2 in Fig. 2(c) are not the ground, first, and second excited states of one and the same mode, but rather | 2 is the first excited state of another higher frequency mode, for example. As discussed above, such coherence pathways will not generate any echo.

FIG. 4.

Simulation of (a) the RRR signal (plotting the intensity) based on an harmonic approximation,49 (b) the RRR signal (plotting the field) obtained directly from a MD simulation, which includes dynamical anharmonicity.52 The various 2D-Raman-THz sequences are shown in (c), which again include dynamical anharmonicity.73 Echoes from the librational mode at 600 cm−1 are indicated. Adapted from S. Saito and I. Ohmine, J. Chem. Phys. 108, 240 (1998); S. Saito and I. Ohmine, J. Chem. Phys. 125, 084506 (2006); and H. Ito, J. Y. Jo, and Y. Tanimura, Struct. Dyn. 2, 054102 (2015). Copyright 2015 AIP Publishing.

FIG. 4.

Simulation of (a) the RRR signal (plotting the intensity) based on an harmonic approximation,49 (b) the RRR signal (plotting the field) obtained directly from a MD simulation, which includes dynamical anharmonicity.52 The various 2D-Raman-THz sequences are shown in (c), which again include dynamical anharmonicity.73 Echoes from the librational mode at 600 cm−1 are indicated. Adapted from S. Saito and I. Ohmine, J. Chem. Phys. 108, 240 (1998); S. Saito and I. Ohmine, J. Chem. Phys. 125, 084506 (2006); and H. Ito, J. Y. Jo, and Y. Tanimura, Struct. Dyn. 2, 054102 (2015). Copyright 2015 AIP Publishing.

Close modal

Follow-up work by Saito and Ohmine,50,52 in which the signal has been calculated directly from a MD simulation and thus implicitly includes all anharmonic and dynamical aspects of the liquid, revealed only the initial spike but the echo essentially disappeared (Fig. 4(b)). The same has been observed for liquid Xe with the echo vanishing once dynamical anharmonicity is considered.46 

More recently, our group as well as that of Tanimura69–73,76,77 extended the theory to the hybrid 2D-Raman-THz pulse sequences, which is conceptually very similar to that of 2D-Raman spectroscopy. One example from Ref. 73 is shown in Fig. 4(c) with an echo of the 600 cm−1 mode reappearing in the TRT sequence. Note that in contrast to Fig. 4(a), that simulation does in fact include all dynamical and anharmonic effects of liquid water; hence, the distinct echo is a non-trivial effect in light of the discussion above. Interestingly, the echo appears only in the TRT sequence,70,73,74,77 which can be understood considering that the 600 cm−1 mode is a librational mode, i.e., a hindered rotation. In the limiting case of a free rotor (linear or spherical) in the gas phase, a Δ J = 1 selection rule would apply for a direct THz interactions and a Δ J = 2 selection rule for a Raman interaction; hence, the rephasing diagram shown in Fig. 2(d) for the TRT sequence would in fact be fully allowed, while all others (RRR, RTT, TTR sequences) are forbidden. These selection rules, which reflect the difference between molecular orientation by a THz pulse versus molecular alignment by a Raman pulse, will be softened for liquid water, but to a certain extent still apply.

A completely different aspect of 2D-Raman-THz spectroscopy concerns the polarizability of water. The dipole moment of a water molecule in the gas-phase is 1.85 D, while that in bulk water is significantly larger with ≈2.5 D due to charge fluxes within and between water molecules in the polar environment of hydrogen bond networks. Point charge models of water such as TIP4P/200596 take into account the polarization of water in a mean-field sense, which is probably quite acceptable for bulk water,97 but which is prone to fail when a water molecule encounters the apolar interface such as at the water-air surface98 or of a protein. Nevertheless, for reasons of computational efficiency, point-charge water models are still most commonly used in connection with protein force fields such as CHARMM,99 AMBER,100 or OPLS.101 Changing the water model for simulations of bio-macromolecules is not trivial since the water model and protein force field are fitted self-consistently, so a change of the former requires a reparametrization of the latter, which is a formidable task. Nevertheless, in light of the ever-growing computer power and the importance of water polarizability to describe its interaction with a protein, the development of polarizable water models is currently a very active field of research102–117 (see Ref. 118 for a recent review) because they will become a cornerstone for the next generation of biomolecular force fields. However, as to date, no polarizable model of water has emerged as the one standard that would be generally accepted as the best model in terms of a compromise between accuracy and computer time efficiency (in the same way as TIP4P/2005 is now considered to be the best possible non-polarizable point-charge model97).

A realistic description of the 2D-Raman-THz response requires a polarizable water model for two reasons. First, trivially, the spectroscopy includes a Raman interaction, which acts on the polarizability of the system. A second, more subtle reason concerns the hydrogen bond stretch band at ≈200 cm−1 in the THz absorption spectrum (see Fig. 1), which is a major target of the experiment. It is well established that simple point charge models of water, such as TIP4P/2005,96 cannot account for that band, despite the fact that the corresponding vibration does, of course, exist in these models. That is, the hydrogen bond stretch band does not have any transition dipole as long as two neutral water molecules vibrate against each other. The intensity of that band originates from charge fluxes within and between water molecules upon hydrogen bonding.4,119–122 Adding polarizability to a water model, either in an ad hoc manner to a trajectory that has been precalculated with the help of a point-charge model122–125 or explicitly as part of the force field111,113,115 reveals that band in the THz absorption spectrum, albeit, in most cases, with severely underestimated intensity. We believe that this failure reflects the fact that most of the water models considered so far do not include the charge flux between water molecules upon hydrogen bonding, despite the fact that its importance is widely recognized.4,77,119–122,126

In a recent paper,71 we have shown that 2D-Raman-THz spectroscopy is particularly sensitive to the level of description of polarizability of a water model. Tanimura and coworker came to the same conclusion for a different set of water models.77 To that end, Fig. 5 (left column) shows the simulated 2D-Raman-THz response functions for the RTT and the TRT pulse sequence for three comparably simple models of water. First, TIP4P/200596 has been considered, which is the best among all point-charge models of water in terms of its thermodynamic properties.97 It is a non-polarizable water model; hence, it is used to calculate the MD trajectory per se, which is then amended with polarizability for the calculation of the spectroscopy only. In contrast, SWM4-NDP109 is an intrinsically polarizable water model by employing a massless, charged “Drude”-particle that is attached to the oxygen atom via a harmonic spring. Finally, TL4P116 features a Gaussian shaped inducible dipole and consequently avoids the common problem of the polarization catastrophe. Hence, it is a transferable model that describes both the gas and the solution electrostatics of water correctly.

FIG. 5.

Simulated 2D-Raman-THz response for various water models: TIP4P/200596 amended with an isotropic polarizability, SWM4-NDP,109 and TL4P.116 The left column shows the molecular response function with the RTT pulse sequence towards positive t1-values and the TRT-pulse sequence towards negative t1-values. The right column shows the expected signal after convolution of the molecular response function with the laser pulses and including field generation and propagation effects, taking into account the experimental conditions of Ref. 79Positive response is depicted in red, and negative in blue. Adapted from P. Hamm, J. Chem. Phys. 141, 184201 (2014). Copyright 2014 AIP Publishing.

FIG. 5.

Simulated 2D-Raman-THz response for various water models: TIP4P/200596 amended with an isotropic polarizability, SWM4-NDP,109 and TL4P.116 The left column shows the molecular response function with the RTT pulse sequence towards positive t1-values and the TRT-pulse sequence towards negative t1-values. The right column shows the expected signal after convolution of the molecular response function with the laser pulses and including field generation and propagation effects, taking into account the experimental conditions of Ref. 79Positive response is depicted in red, and negative in blue. Adapted from P. Hamm, J. Chem. Phys. 141, 184201 (2014). Copyright 2014 AIP Publishing.

Close modal

In order to compare the molecular response functions (Fig. 5, left column) with experimental results,79 one needs to include the convolution of the response functions with the THz and Raman pulse shapes, where in particular the former has a rather peculiar shape, as well as the process of the generation of an emitted field together with propagation effects from the sample to the detection crystal (see Refs. 71 and 79 for details). The result of these procedures is shown in Fig. 5 (right column), taking into account the experimental parameters from Ref. 79. It can be seen that the responses dramatically differ for the various water models, despite the fact that all models are of about equal complexity and agree with each other regarding their thermodynamic properties and 1D Raman and THz spectra. Given that extraordinary sensitivity of 2D-Raman-THz spectroscopy, it appears that any polarizable water model102–117 should be tested against 2D-Raman-THz spectroscopy for a stringent verification. From all the water models tested in Ref. 71, TL4P116 agrees by far the best with experiment (see Fig. 6(b)).

FIG. 6.

(a) The TRT and RTT pulse sequences with the definitions of the delay times t1 and t2.(b) Corresponding RTT and TRT response of pure water in comparison to (c) the instrument response function. Positive response is depicted in red, negative in blue, and the THz as well as Raman pulses is depicted as green lines on the top and right side, respectively, of panel (b). The dashed lines represent the echo-directions for the RTT pulse sequence (upper-right quadrant) and the TRT pulse sequence (upper triangle of the upper-left quadrant). Adapted with permission from Ref. 79.

FIG. 6.

(a) The TRT and RTT pulse sequences with the definitions of the delay times t1 and t2.(b) Corresponding RTT and TRT response of pure water in comparison to (c) the instrument response function. Positive response is depicted in red, negative in blue, and the THz as well as Raman pulses is depicted as green lines on the top and right side, respectively, of panel (b). The dashed lines represent the echo-directions for the RTT pulse sequence (upper-right quadrant) and the TRT pulse sequence (upper triangle of the upper-left quadrant). Adapted with permission from Ref. 79.

Close modal

Note that none of the simulations shown in Fig. 5 (left column) reveal an echo of the librational mode that would be comparable to the distinct echo shown in Fig. 4(b).77 The echo is missing since all water models considered here used an isotropic polarizability for reasons of computational efficiency, in which case the Raman activity of the librational mode essentially vanishes. Experimentally, the polarizability of water is indeed very close to isotropic, but not completely. For a comparison to experimental results,79 this simplification is not too much of a concern since the librational mode anyway is too high in frequency to be reached with the spectral width of our current THz and Raman pulses.

The RTT and the TRT pulse sequence can be measured with one and the same experimental setup79,80 by interchanging the timing of the Raman and THz pump pulses, since the read-out is a THz field in both cases (conversely, the TTR pulse sequences requires a different experimental setup78). For reasons that will become apparent later on, we focus on the RTT pulse sequence, which appears in the upper-right quadrant of our representation of the data with t 1 > 0 and t 2 > 0 (see Fig. 6). Owing to the resulting arrangement of the optical delay stages, time t2 is the delay between the THz input pulse and the emitted field, and not necessarily between the second interaction and the emitted field, so that t 2 t 1 + t 2 for t 1 < 0 , and TRT pulse sequence appears in the upper triangle of the upper-left quadrant with t 1 < 0 and t 1 + t 2 > 0 (see pulse sequences and the definitions of times t1 and t2 in Fig. 6(a)).

Fig. 6(b) shows the 2D-Raman-THz response of pure water, and Fig. 6(c) the instrument response function (IRF), i.e., the signal one would expect if the molecular response function would be δ -shaped.79 The IRF is quite extended along the anti-diagonal and for the most part reflects the peculiar shape of the THz pulse shown in green on the right side of Fig. 6(c). Significant parts of the water response resemble the IRF and hence originate from contributions of the molecular response function that are significantly faster than our current time resolution. But some features in the water response go beyond the IRF, the most relevant being that along the diagonal with t1 = t2, which is where one would expect an echo to appear. Despite the fact that the signal is very short-lived with an averaged decay time of 65 fs, it does extend a bit farer into the echo-direction than the IRF.

Fig. 7 shows what happens upon the addition of salts to the solution.80 Depending on the capability of the cation to “structure water,” the decay along the echo-direction either speeds up slightly (e.g., Cs+, see Figs. 7(a) and 7(c)) or slows down (e.g., Mg2+, see Figs. 7(b) and 7(c)). The “structure making” capability of an ion is often characterized in terms of the so-called B-coefficient, i.e., the coefficient of the linear term in the Jones-Dole relationship, which empirically relates the solvent viscosity to the salt concentration.127,128 Fig. 7(c) (red line) shows that the decay in the echo direction correlates extremely well with that B-coefficient.80 

FIG. 7.

RTT and TRT responses of (a) a 2M CsCl and (b) a 2M MgCl2 solution. Positive response is depicted in red, negative in blue. The red line in panel (c) shows a linear fit to the averaged relaxation time of the signal in the echo direction as function of the B-coefficient of the cation (the same anion Cl has been used in all cases). The corresponding averaged relaxation time of neat water is indicated as a blue diamond. The green line shows the same for the signal decay along the t1-axis. Adapted with permission from Ref. 80.

FIG. 7.

RTT and TRT responses of (a) a 2M CsCl and (b) a 2M MgCl2 solution. Positive response is depicted in red, negative in blue. The red line in panel (c) shows a linear fit to the averaged relaxation time of the signal in the echo direction as function of the B-coefficient of the cation (the same anion Cl has been used in all cases). The corresponding averaged relaxation time of neat water is indicated as a blue diamond. The green line shows the same for the signal decay along the t1-axis. Adapted with permission from Ref. 80.

Close modal

A signal in the t1 = t2 direction might be called an echo, when it lives longer than inhomogeneous dephasing. To that end, Fig. 7(c) compares the diagonal decay (red line) with that along the t1-axes (green line). Adopting that rigorous definition of an echo, it is probably wrong to call the diagonal signals of neat water an echo (which is why a question mark is added in Fig. 6(b)), while one is definitely entering that regime in the case of the strongly structure-making cations such as MgCl2. We have provided evidence in Ref. 80 that the echo originates from water-water intermolecular modes and not from the newly introduced water-ion modes.

Arguing in a very hand-waving manner, the addition of salts renders the hydrogen-bond network of water more inhomogeneous and at the same time, the first solvation layer around an ion is stabilized significantly. In particular in the case of Mg2+, the persistent time of the water molecules in the first solvation layer is very long with reported time scales varying from 100’s of ps129,130 to > 1 μ s .131–134 It has been argued in Sec. III A that the decay of an echo reflects the lifetime of local structures on a length scale that equals the delocalization length of the intermolecular modes. In light of that statement, it appears reasonably that both effects of ions, i.e., larger inhomogeneity together with longer lifetimes of hydrogen bond structures, should facilitate the formation of an echo. The only simulation work to date that addressed that question is that by Zhuang and co-workers,75 who saw a slight prolongation of the signal in the echo direction for a 3.5M MgCl2 solution, albeit to a much lesser extent than what has been observed experimentally.80 

In contrast to Refs. 70, 73, 74, and 77 (e.g., Fig. 4(c)), we observe an echo only for the RTT pulse sequence. Two things can be said about that. First, the echo-direction in the RTT pulse sequence is “perpendicular” to the IRF. That is, the IRF extends dominantly along the antidiagonal in our data representation (Fig. 6(c)), along which it reflects the rather long THz pulse, while the time-resolution in the echo direction of the RTT pulse sequence is in essence dictated by the much shorter Raman pump pulse. In the contrary, the echo direction in the TRT pulse sequence lies pretty close to the IRF (see dashed line in the upper triangle of the upper-left quadrant of Figs. 6(b) and 6(c)), and our experiment is much less sensitive to a possible echo in that case. Significantly shorter THz pulses would be needed to possibly resolve an echo for the TRT pulse sequence.

Second, as discussed in Sec. III B, the echo observed in MD simulations for the TRT pulse sequence originates from the librational mode at 600 cm−1, which anyway is outside the spectral window of our current experiment (≲200 cm−1). Different “selection rules” might apply for the other modes of water. In that regard, the results from Blake and coworkers are quite revealing, who measured the TTR response of boroform (CHBr3) and related liquids, see Fig. 8.78 In that case, the 2D response is dominated by two bands originating from intramolecular modes (Fig. 8(b)), which can be unambiguously assigned to the coherence pathways shown in Fig. 8(c), based on the energies of the known vibrational states of that small molecule. Interestingly, the THz pulse can induce two- and even three-quantum transitions ( | 01 | 10 and | 01 | 20 , respectively), while the Raman interaction can only induce a one-quantum transition ( | 20 | 10 ). Apparently, the dipole nonlinearity (i.e., the electrical anharmonicity) is much larger than the polarizability nonlinearity for the intramolecular modes of liquid bromoform. The origin of that effect, and to what extent it also applies for the intermolecular hydrogen-bond modes of water, is currently not clear. Nevertheless, it would potentially explain the appearance of an echo in the RTT sequence, which requires a two-quantum transition induced by the second (THz) interaction for an inversion of coherence to occur (Fig. 2(d)).

FIG. 8.

TTR response of CHBr3 in (a) the time domain and (b) the frequency domain (showing only the non-rephasing quadrant). Panel (c) depicts the WMEL diagrams giving rise to the two peaks observed in the 2D spectra with the corresponding intramolecular states indicated, where | i j refers to i quanta in mode ν 6 and j quanta in mode ν 3 (both are C-Br bending modes). Adapted with permission from Ref. 78.

FIG. 8.

TTR response of CHBr3 in (a) the time domain and (b) the frequency domain (showing only the non-rephasing quadrant). Panel (c) depicts the WMEL diagrams giving rise to the two peaks observed in the 2D spectra with the corresponding intramolecular states indicated, where | i j refers to i quanta in mode ν 6 and j quanta in mode ν 3 (both are C-Br bending modes). Adapted with permission from Ref. 78.

Close modal

In the present perspective article, it has been discussed under which circumstances one might expect to see a photon echo for the low-frequency inter-molecular modes of liquids. While our work, so far, concentrates on water and aqueous solutions,79,80 other liquids should be investigated as well, and it has indeed by suggested by MD simulations that the 2D Raman THz response of formaldehyde, DMSO, or formamide looks significantly different from that of water.72 Due to signal-to-noise limitations, we are currently limited to bulk water or solutions at very high concentrations that do not allow one to distinguish bulk water from a hydration shell. As technology evolves, we hope to be able at some point to also investigate water in more complex situations, such as water in micelles or the hydration layer of a protein.

In any case, the observation of an echo is the probably the most one might expect from a 2D experiment in the THz regime since low-frequency inter-molecular spectra of liquids, in particular water, are so extremely blurred. Being able to look underneath, the potentially inhomogeneously broadened line can enhance the resolution, and consequently the information content that can be extracted from the spectroscopy. Subtle but important differences exist in comparison to the well established concept of echoes in NMR,81 electronic,82 or vibrational spectroscopy.83–85 On the one hand, these differences are related to the fact that any complex molecular system, be it a protein, a glass, or a liquid, contains an “inhomogeneous ensemble” of low-frequency modes; hence, an echo might trivially be expected in any case, even for a single molecule. On the other hand, the “2-in-1-out”-pulse sequences considered here induce different coherence pathways than the “3-in-1-out”-pulse sequences known from NMR, electronic or vibrational spectroscopy. Rephasing and hence an echo are obtained only if the degrees of freedom of the molecular system under study are sufficiently harmonic, and if the persistent time of those modes is sufficiently long. It is a priori not clear to what extent these two conditions are fulfilled for the inter-molecular degrees of freedom of a liquid.

Therefore, the observation of an echo in Ref. 80 for aqueous solutions of strongly “structure-making” salts is in our opinion the first compelling evidence that multi-dimensional THz spectroscopy can indeed enhance the information content extracted from the extremely blurred inter-molecular spectra in the THz range. As such, we believe that Refs. 79 and 80 have introduced the most decisive spectroscopy of water and ion solvation as to date, despite the fact that admittedly the full information cannot be retrieved at this stage. A true understanding of the experimental signals of Figs. 6 and 7 will require extensive further simulation work, most likely including polarizable MD models of ions as well.

In that regard, it is important to add that the three-time-point correlation functions needed to calculate the 2D-Raman-THz response converge extremely slowly; typically 150-300 μs of total MD simulation time went into the calculation of any of the response functions shown in Fig. 5.135 Extremely computer-time efficient—and still accurate—MD models are therefore required to simulate the 2D-Raman-THz response. For sure, computer-time efficiency will also be a killing criterion, should a polarizable water model ever become the cornerstone of a newly developed biomolecular force field, since the largest part of the computational effort is needed for the water solvating a bio-macromolecule. Biomolecular simulations will always suffer from a time and length scale problem despite the ever-growing computer power. Together with its extraordinary sensitivity on the level of description of polarizability, 2D-Raman-THz spectroscopy will provide an excellent guideline in the development of polarizable force fields, both in terms of accuracy and in terms of computer-time efficiency.

We are indebted to Janne Savolainen, Saima Ahmed, Gustavo Ciardi, Arian Berger, and David Sidler, who made, and continue to make, the work happen in the lab, on which this perspective article is based. We also thank Yoshitaka Tanimuara, Yuki Nagata, Geoff Blake, Gerhard Stock, Markus Meuwly, and Thomas Feurer for many insightful discussions on the topic. The work is a central project of the NCCR MUST, funded by the Swiss National Science Foundation (SNF), whose generous support is highly acknowledged.

1.
J.
Bertie
and
Z.
Lan
,
Appl. Spectrosc.
50
,
1047
(
1996
).
2.
T.
Fukasawa
,
T.
Sato
,
J.
Watanabe
,
Y.
Hama
,
W.
Kunz
, and
R.
Buchner
,
Phys. Rev. Lett.
95
,
197802
(
2005
).
3.
S.
Ebbinghaus
,
S. J.
Kim
,
M.
Heyden
,
X.
Yu
,
U.
Heugen
,
M.
Gruebele
,
D. M.
Leitner
, and
M.
Havenith
,
Proc. Natl. Acad. Sci U. S. A.
104
,
20749
(
2007
).
4.
M.
Heyden
,
J.
Sun
,
S.
Funkner
,
G.
Mathias
,
H.
Forbert
,
M.
Havenith
, and
D.
Marx
,
Proc. Natl. Acad. Sci. U. S. A.
107
,
12068
(
2010
).
5.
G. E.
Walrafen
,
M. R.
Fischer
,
M. S.
Hokmabadi
, and
W. H.
Yang
,
J. Chem. Phys.
85
,
6970
(
1986
).
6.
E. W.
Castner
, Jr.
,
Y. J.
Chang
,
Y. C.
Chu
, and
G. E.
Walrafen
,
J. Chem. Phys.
102
,
653
(
1995
).
7.
C. J.
Fecko
,
J. D.
Eaves
, and
A.
Tokmakoff
,
J. Chem. Phys.
117
,
1139
(
2002
).
8.
T.
Torre
,
P.
Bartolini
,
R.
Righini
,
R.
Torre
,
P.
Bartolini
, and
R.
Righini
,
Nature
428
,
296
(
2004
).
9.
N. T.
Hunt
,
L.
Kattner
,
R. P.
Shanks
, and
K.
Wynne
,
J. Am. Chem. Soc.
129
,
3168
(
2007
).
10.
N. T.
Hunt
,
A. A.
Jaye
, and
S. R.
Meech
,
Phys. Chem. Chem. Phys.
9
,
2167
(
2007
).
11.
I. A.
Heisler
and
S. R.
Meech
,
Science
305
,
857
(
2010
).
12.
K.
Mazur
,
I. A.
Heisler
, and
S. R.
Meech
,
J. Phys. Chem. B
115
,
2563
(
2011
).
13.
A.
Taschin
,
P.
Bartolini
,
E.
Eramo
,
R.
Righini
, and
R.
Torre
,
Nat. Commun.
4
,
2401
(
2013
).
14.
D. A.
Turton
,
A. R.
Turner
,
N. T.
Hunt
,
G. H.
Welsh
, and
K.
Wynne
, "
An experimental and numerical study of hydrogen-bonding in aqueous salts and methanol
," in
Ultrafast Phenom. XV
, edited by
P.
Corkum
,
D. M.
Jonas
,
R. J. D.
Miller
, and
A. M.
Weiner
(
Springer
,
Berlin, Heidelberg
,
2007
), pp.
427
429
.
15.
D. A.
Schmidt
,
S.
Funkner
,
B. P.
Born
,
R.
Gnanasekaran
,
G. W.
Schwaab
,
D. M.
Leitner
, and
M.
Havenith
,
J. Am. Chem. Soc.
131
,
18512
(
2009
).
16.
S.
Funkner
,
G.
Niehues
,
D. A.
Schmidt
,
M.
Heyden
,
G.
Schwaab
,
K. M.
Callahan
,
D. J.
Tobias
, and
M.
Havenith
,
J. Am. Chem. Soc.
134
,
1030
(
2012
).
17.
I. A.
Heisler
,
K.
Mazur
, and
S. R.
Meech
,
J. Phys. Chem. B
115
,
1863
(
2011
).
18.
M.
Heyden
,
E.
Bründermann
,
U.
Heugen
,
G.
Niehues
,
D. M.
Leitner
, and
M.
Havenith
,
J. Am. Chem. Soc.
130
,
5773
(
2008
).
19.
D. A.
Turton
,
H. M.
Senn
,
T.
Harwood
,
A. J.
Lapthorn
,
E. M.
Ellis
, and
K.
Wynne
,
Nat. Commun.
5
,
3999
(
2014
).
20.
K.
Wüthrich
,
NMR of Proteins and Nucleic Acids
(
Wiley Interscience
,
New York
,
1986
).
21.
R. R.
Ernst
,
G.
Bodenhausen
, and
A.
Wokaun
,
Principles of Nuclear Magnetic Resonance in One and Two Dimensions
(
Clarendon Press
,
Oxford
,
1987
).
22.
R. F.
Loring
and
S.
Mukamel
,
J. Chem. Phys.
83
,
2116
(
1985
).
23.
W. P.
Aue
,
E.
Bartholdi
, and
R. R.
Ernst
,
J. Chem. Phys.
64
,
2229
(
1976
).
24.
P.
Hamm
,
M.
Lim
, and
R. M.
Hochstrasser
,
J. Phys. Chem. B
102
,
6123
(
1998
).
25.
S. M.
Gallagher Faeder
,
A. W.
Albrecht
,
J. D.
Hybl
,
B. L.
Landin
,
B.
Rajaram
, and
D. M.
Jonas
,
J. Opt. Soc. Am. B
15
,
2338
(
1998
).
26.
M.
Cho
,
Two-Dimensional Optical Spectroscopy
(
CRC Press
,
Boca Raton
,
2009
).
27.
P.
Hamm
and
M. T.
Zanni
,
Concepts and Methods of 2D Infrared Spectroscopy
(
Cambridge University Press
,
Cambridge
,
2011
).
28.
D.
Lee
and
A. C.
Albrecht
, in
Advances in Infrared and Raman Spectroscopy
, edited by
R. J. H.
Clark
and
H. R. E.
(
Wiley Heyden
,
1985
), pp.
179
213
.
29.
P.
Hamm
,
R. A.
Kaindl
, and
J.
Stenger
,
Opt. Lett.
25
,
1798
(
2000
).
30.
P. A.
Cazade
,
H.
Tran
,
T.
Bereau
,
A. K.
Das
,
F.
Kläsi
,
P.
Hamm
, and
M.
Meuwly
,
J. Chem. Phys.
142
,
212415
(
2015
).
31.
R.
Costard
,
T.
Tyborski
,
B.P.
Fingerhut
, and
T.
Elsaesser
,
J. Chem. Phys.
142
,
212406
(
2015
).
32.
W.
Kuehn
,
K.
Reimann
,
M.
Woerner
, and
T.
Elsaesser
,
J. Chem. Phys.
130
,
164503
(
2009
).
33.
W.
Kuehn
,
K.
Reimann
,
M.
Woerner
,
T.
Elsaesser
,
R.
Hey
, and
U.
Schade
,
Phys. Rev. Lett.
107
,
67401
(
2011
).
34.
W.
Kuehn
,
K.
Reimann
,
M.
Woerner
,
T.
Elsaesser
, and
R.
Hey
,
J. Phys. Chem. B
115
,
5448
(
2011
).
35.
C.
Somma
,
G.
Folpini
,
K.
Reimann
,
M.
Woerner
, and
T.
Elsaesser
,
J. Chem. Phys.
144
,
184202
(
2015
).
36.
T.
Elsaesser
,
K.
Reimann
, and
M.
Woerner
,
J. Chem. Phys.
142
,
212301
(
2015
).
37.
J.
Lu
,
Y.
Zhang
,
H. Y.
Hwang
,
B. K.
Ofori-Okai
,
S.
Fleischer
, and
K. A.
Nelson
,
Proc. Natl. Acad. Sci. U. S. A.
113
,
11800
(
2016
).
38.
J. T.
Fourkas
,
Adv. Chem. Phys.
117
,
235
(
2001
).
39.
D.
Vanden Bout
,
L. J.
Muller
, and
M.
Berg
,
Phys. Rev. Lett.
67
,
3700
(
1991
).
40.
L. J.
Muller
,
D.
Vanden Bout
, and
M.
Berg
,
J. Chem. Phys.
99
,
810
(
1993
).
41.
R.
Inaba
,
K.
Tominaga
,
M.
Tasumi
,
K. A.
Nelson
, and
K.
Yoshihara
,
Chem. Phys. Lett.
211
,
183
(
1993
).
42.
Y.
Tanimura
and
S.
Mukamel
,
J. Chem. Phys.
99
,
9496
(
1993
).
43.
S.
Palese
,
J. T.
Buontempo
,
L.
Schilling
,
W. T.
Lotshaw
,
Y.
Tanimura
,
S.
Mukamel
, and
R. J. D.
Miller
,
J. Phys. Chem.
98
,
12466
(
1994
).
44.
K.
Okumura
and
Y.
Tanimura
,
J. Chem. Phys.
106
,
1687
(
1997
).
45.
K.
Okumura
and
Y.
Tanimura
,
J. Chem. Phys.
107
,
2267
(
1997
).
46.
A.
Ma
and
R. M.
Stratt
,
Phys. Rev. Lett.
85
,
1004
(
2000
).
47.
Y.
Tanimura
and
T.
Steffen
,
J. Phys. Soc. Jpn.
69
,
4095
(
2000
).
48.
T. L.
Jansen
,
J. G.
Snijders
, and
K.
Duppen
,
J. Chem. Phys.
113
,
307
(
2000
).
49.
S.
Saito
and
I.
Ohmine
,
J. Chem. Phys.
108
,
240
(
1998
).
50.
S.
Saito
and
I.
Ohmine
,
Phys. Rev. Lett.
88
,
207401
(
2002
).
51.
S.
Saito
and
I.
Ohmine
,
J. Chem. Phys.
119
,
9073
(
2003
).
52.
S.
Saito
and
I.
Ohmine
,
J. Chem. Phys.
125
,
084506
(
2006
).
53.
Y.
Nagata
and
Y.
Tanimura
,
J. Chem. Phys.
124
,
24508
(
2006
).
54.
R.
DeVane
,
C.
Kasprzyk
,
B.
Space
, and
T.
Keyes
,
J. Phys. Chem. B
110
,
3773
(
2006
).
55.
T.
Hasegawa
and
Y.
Tanimura
,
J. Chem. Phys.
125
,
074512
(
2006
).
56.
T.
Yagasaki
and
S.
Saito
,
Acc. Chem. Res.
42
,
1250
(
2009
).
57.
A.
Tokmakoff
,
M. J.
Lang
,
D. S.
Larsen
,
G. R.
Fleming
,
V.
Chernyak
, and
S.
Mukamel
,
Phys. Rev. Lett.
79
,
2702
(
1997
).
58.
D. A.
Blank
,
L. J.
Kaufman
, and
G. R.
Fleming
,
J. Chem. Phys.
111
,
3105
(
1999
).
59.
D. A.
Blank
,
L. J.
Kaufman
, and
G. R.
Fleming
,
J. Chem. Phys.
113
,
771
(
2000
).
60.
L. J.
Kaufman
,
J.
Heo
,
L. D.
Ziegler
, and
G. R.
Fleming
,
Phys. Rev. Lett.
88
,
207402
(
2002
).
61.
K. J.
Kubarych
,
C. J.
Milne
, and
R. J. D.
Miller
,
Int. Rev. Phys. Chem.
22
,
497
(
2003
).
62.
O.
Golonzka
,
N.
Demirdöven
,
M.
Khalil
, and
A.
Tokmakoff
,
J. Chem. Phys.
113
,
9893
(
2000
).
63.
Y. L.
Li
,
L.
Huang
,
R. J. D.
Miller
,
T.
Hasegawa
, and
Y.
Tanimura
,
J. Chem. Phys.
128
,
234507
(
2008
).
64.
D. J.
Ulness
,
J. C.
Kirkwood
, and
A. C.
Albrecht
,
J. Chem. Phys.
108
,
3897
(
1998
).
65.
H.
Frostig
,
T.
Bayer
,
N.
Dudovic
,
Y. C.
Eldar
, and
Y.
Silberberg
,
Nat. Photonics
9
,
339
(
2015
).
66.
B. P.
Molesky
,
P. G.
Giokas
,
Z.
Guo
, and
A. M.
Moran
,
J. Chem. Phys.
141
,
114202
(
2014
).
67.
B. P.
Molesky
,
Z.
Guo
,
T. P.
Cheshire
, and
A. M.
Moran
,
J. Chem. Phys.
145
,
180901
(
2016
).
68.
69.
P.
Hamm
and
J.
Savolainen
,
J. Chem. Phys.
136
,
094516
(
2012
).
70.
P.
Hamm
,
J.
Savolainen
,
J.
Ono
, and
Y.
Tanimura
,
J. Chem. Phys.
136
,
236101
(
2012
).
71.
72.
H.
Ito
,
T.
Hasegawa
, and
Y.
Tanimura
,
J. Chem. Phys.
141
,
124503
(
2014
).
73.
H.
Ito
,
J. Y.
Jo
, and
Y.
Tanimura
,
Struct. Dyn.
2
,
054102
(
2015
).
74.
T.
Ikeda
,
H.
Ito
, and
Y.
Tanimura
,
J. Chem. Phys.
142
,
212421
(
2015
).
75.
Z.
Pan
,
T.
Wu
,
T.
Jin
,
Y.
Liu
,
Y.
Nagata
,
R.
Zhang
, and
W.
Zhuang
,
J. Chem. Phys.
142
,
212419
(
2015
).
76.
H.
Ito
and
Y.
Tanimura
,
J. Chem. Phys.
144
,
074201
(
2016
).
77.
H.
Ito
,
T.
Hasegawa
, and
Y.
Tanimura
,
J. Phys. Chem. Lett.
7
,
4147
(
2016
).
78.
I. A.
Finneran
,
R.
Welsch
,
M. A.
Allodi
,
T. F.
Miller
, and
G. A.
Blake
,
Proc. Natl. Acad. Sci. U. S. A.
113
,
6857
(
2016
).
79.
J.
Savolainen
,
S.
Ahmed
, and
P.
Hamm
,
Proc. Natl. Acad. Sci. U. S. A.
110
,
20402
(
2013
).
80.
A.
Shalit
,
S.
Ahmed
,
J.
Savolainen
, and
P.
Hamm
,
Nat. Chem.
9
,
273
(
2017
).
82.
W. P.
de Boeij
,
M. S.
Pshenichnikov
, and
D. A.
Wiersma
,
Annu. Rev. Phys. Chem.
49
,
99
(
1998
).
83.
A.
Tokmakoff
and
M. D.
Fayer
,
J. Chem. Phys.
103
,
2810
(
1995
).
84.
P.
Hamm
,
M. H.
Lim
, and
R. M.
Hochstrasser
,
Phys. Rev. Lett.
81
,
5326
(
1998
).
86.
E. C.
Fulmer
,
F.
Ding
, and
M. T.
Zanni
,
J. Chem. Phys.
122
,
034302
(
2005
).
87.
S.
Mukamel
,
Principles of Nonlinear Optical Spectroscopy
(
Oxford University Press
,
Oxford
,
1995
).
88.
O. M.
Becker
and
M.
Karplus
,
Phys. Rev. Lett.
70
,
3514
(
1993
).
89.
G. S.
Grest
,
S. R.
Nagel
, and
A.
Rahman
,
Solid State Commun.
36
,
875
(
1980
).
90.
D.
Xu
,
K.
Schulten
,
O. M.
Becker
, and
M.
Karplus
,
J. Chem. Phys.
103
,
3112
(
1995
).
91.
X.
Yu
and
D. M.
Leitner
,
J. Chem. Phys.
122
,
054902
(
2005
).
93.
C. P.
Lawrence
and
J. L.
Skinner
,
J. Chem. Phys.
118
,
264
(
2003
).
94.
D.
Laage
and
J. T.
Hynes
,
Chem. Phys. Lett.
433
,
80
(
2006
).
95.
D. E.
Moilanen
,
E. E.
Fenn
,
Y.-S.
Lin
,
J. L.
Skinner
,
B.
Bagchi
, and
M. D.
Fayer
,
Proc. Natl. Acad. Sci. U. S. A.
105
,
5295
(
2008
).
96.
J. L.
Abascal
and
C.
Vega
,
J. Chem. Phys.
123
,
234505
(
2005
).
97.
C.
Vega
and
J. L. F.
Abascal
,
Phys. Chem. Chem. Phys.
13
,
19663
(
2011
).
98.
P. A.
Pieniazek
,
C. J.
Tainter
, and
J. L.
Skinner
,
J. Chem. Phys.
135
,
44701
(
2011
).
99.
A. D.
Mackerell
,
D.
Bashford
,
M.
Bellott
,
R. L.
Dunbrack
,
J. D.
Evanseck
,
M. J.
Field
,
S.
Fischer
,
J.
Gao
,
H.
Guo
,
S.
Ha
,
D.
Joseph-McCarthy
,
L.
Kuchnir
,
K.
Kuczera
,
F. T. K.
Lau
,
C.
Mattos
,
S.
Michnick
,
T.
Ngo
,
D. T.
Nguyen
,
B.
Prodhom
,
W. E.
Reiher
,
B.
Roux
,
M.
Schlenkrich
,
J. C.
Smith
,
R.
Stote
,
J.
Straub
,
M.
Watanabe
,
J.
Wiórkiewicz-Kuczera
,
D.
Yin
, and
M.
Karplus
,
J. Phys. Chem. B
102
,
3586
(
1998
).
100.
K. T.
Debiec
,
D. S.
Cerutti
,
L. R.
Baker
,
A. M.
Gronenborn
,
D. A.
Case
, and
L. T.
Chong
,
J. Chem. Theory Comput.
12
,
3926
(
2016
).
101.
M. J.
Robertson
,
J.
Tirado-Rives
, and
W. L.
Jorgensen
,
J. Chem. Theory Comput.
11
,
3499
(
2015
).
102.
P.
Ahlström
,
A.
Wallqvist
,
S.
Engström
, and
B.
Jönsson
,
Mol. Phys.
68
,
563
(
1989
).
103.
S. W.
Rick
,
S. J.
Stuart
, and
B. J.
Berne
,
J. Chem. Phys.
101
,
6141
(
1994
).
104.
H. A.
Stern
,
F.
Rittner
,
B. J.
Berne
, and
R. A.
Friesner
,
J. Chem. Phys.
115
,
2237
(
2001
).
105.
P.
Ren
and
J. W.
Ponder
,
J. Phys. Chem. B
107
,
5933
(
2003
).
106.
G.
Lamoureux
and
B.
Roux
,
J. Chem. Phys.
119
,
3025
(
2003
).
107.
G.
Lamoureux
,
A. D.
MacKerrell
, Jr.
, and
B.
Roux
,
J. Chem. Phys.
119
,
5185
(
2003
).
108.
P.
Paricaud
,
M.
Predota
,
A. A.
Chialvo
, and
P. T.
Cummings
,
J. Chem. Phys.
122
,
244511
(
2005
).
109.
G.
Lamoureux
,
E.
Harder
,
I. V.
Vorobyov
,
B.
Roux
, and
A. D.
MacKerrell
, Jr.
,
Chem. Phys. Lett.
418
,
245
(
2006
).
110.
G. S.
Fanourgakis
and
S. S.
Xantheas
,
J. Chem. Phys.
128
,
074506
(
2008
).
111.
J.
Liu
,
W. H.
Miller
,
G. S.
Fanourgakis
,
S. S.
Xantheas
,
S.
Imoto
, and
S.
Saito
,
J. Chem. Phys.
135
,
244503
(
2009
).
112.
R.
Kumar
,
F.-F.
Wang
,
G. R.
Jenness
, and
K. D.
Jordan
,
J. Chem. Phys.
132
,
014309
(
2010
).
113.
T.
Hasegawa
and
Y.
Tanimura
,
J. Phys. Chem. B
115
,
5545
(
2011
).
114.
W.
Yu
,
P. E. M.
Lopes
,
B.
Roux
, and
A. D.
MacKerrell
, Jr.
,
J. Chem. Phys.
138
,
034508
(
2013
).
115.
L.-P.
Wang
,
T.
Head-Gordon
,
J. W.
Ponder
,
P.
Ren
,
J. D.
Chodera
,
P. K.
Eastman
,
T. J.
Martinez
, and
V. S.
Pande
,
J. Phys. Chem. B
117
,
9956
(
2013
).
116.
P.
Tröster
,
K.
Lorenzen
,
M.
Schwörer
, and
P.
Tavan
,
J. Phys. Chem. B
117
,
9486
(
2013
).
117.
C. J.
Tainter
,
L.
Shi
, and
J. L.
Skinner
,
J. Chem. Theory Comput.
11
,
2268
(
2015
).
118.
P. E. M.
Lopes
,
B.
Roux
, and
A. D.
MacKerrell
, Jr.
,
Theor. Chem. Acc.
124
,
11
(
2009
).
119.
A.
Pasquarello
and
R.
Resta
,
Phys. Rev. B
68
,
174302
(
2003
).
120.
M.
Sharma
,
R.
Resta
, and
R.
Car
,
Phys. Rev. Lett.
95
,
187401
(
2005
).
121.
122.
H.
Torii
,
J. Chem. Theory Comput.
10
,
1219
(
2014
).
123.
P. A.
Madden
and
R. W.
Impey
,
Chem. Phys. Lett.
123
,
502
(
1986
).
124.
B.
Guillot
,
J. Chem. Phys.
95
,
1543
(
1991
).
125.
M.
Souaille
and
J. C.
Smith
,
Mol. Phys.
87
,
1333
(
1996
).
126.
A. J.
Lee
and
S. W.
Rick
,
J. Chem. Phys.
134
,
184507
(
2011
).
127.
G.
Jones
and
M.
Dole
,
J. Am. Chem. Soc.
51
,
2950
(
1929
).
128.
H. D. B.
Jenkins
and
Y.
Marcus
,
Chem. Rev.
95
,
2695
(
1995
).
129.
S.
Obst
and
H.
Bradaczek
,
J. Phys. Chem.
100
,
15677
(
1996
).
130.
D.
Jiao
,
C.
King
,
A.
Grossfield
,
T. A.
Darden
, and
P.
Ren
,
J. Phys. Chem. B
110
,
18553
(
2006
).
131.
J.
Neely
and
R.
Connick
,
J. Am. Chem. Soc.
92
,
3476
(
1970
).
132.
H.
Ohtaki
and
T.
Radnai
,
Chem. Rev.
93
,
1157
(
1993
).
133.
J. P.
Larentzos
and
L. J.
Criscenti
,
J. Phys. Chem. B
112
,
14243
(
2008
).
134.
K. M.
Callahan
,
N. N.
Casillas-Ituarte
,
M.
Roeselova
,
H. C.
Allen
, and
D. J.
Tobias
,
J. Phys. Chem. A
114
,
5141
(
2010
).
135.

The computational effort grows exponentially with the time ranges of t1 and t2 since classical (chaotic) trajectories diverge exponentially after a perturbation. A time range of ≈1 ps is needed for t1 and t2 (i.e., significantly longer than what is shown in Figs. 4 and 5, left column) in order to be able to perform the convolution with the THz pulse, whose wings are that wide.