Experimental Access to Mode-Specific Coupling between Quantum Molecular Vibrations and Classical Bath Modes

The interaction of quantum-mechanical systems with a fluctuating thermal environment (bath) is fundamental to molecular mechanics and energy transport/dissipation. Its complete picture requires mode-specific measurements of this interaction and an understanding of its nature. Here, we present a combined experimental and theoretical study providing detailed insights into the coupling between a high-frequency vibrational two-level system and thermally excited terahertz modes. Experimentally, two-dimensional terahertz-infrared-visible spectroscopy reports directly on the coupling between quantum oscillators represented by CH3 streching vibrations in liquid dimethyl sulfoxide and distinct low-frequency modes. Theoretically, we present a mixed quantum-classical formalism of the sample response to enable the simultaneous quantum description of high-frequency oscillators and a classical description of the bath. We derive the strength and nature of interaction and find different coupling between CH3 stretch and low-frequency modes. This general approach enables quantitative and mode-specific analysis of coupled quantum and classical dynamics in complex chemical systems.


Mixed quantum-classical formalism of response function A. Mixed quantum-classical equation of motion
We consider a quantum system having multiple degrees of freedom.The state of the system is described by the density matrix operator |⟩⟩, and its evolution is described by the Liouville equation where  is Liouville superoperator and  ̂ is Hamiltonian operator.The accent  ̂ denotes a quantum operator in Hilbert space.
In the full system, we distinguish two groups of coordinates represented by the corresponding (multidimensional) operators  ̂ and  ̂.Total Hamiltonian of the system is then composed of three terms We assume that total density matrix of the system is given by the product of the density matrixes of the two subsystems: Thus, the first term on the right-hand side of Eq.S3a The third term in Eq.S3a Terms I and II in this equation describe the evolution of the quantum and classical subsystems without interaction.Term III describes the influence on the quantum subsystem by the classical degrees of freedom.The last term IV reflects the influence on the evolution of classical subsystem by quantum degrees of freedom.Equation S10 is identical to the previously derived mixed quantum-classical equation of motion 1 .
In Eq.S10, we neglect the last Term IV, i.e., we assume that quantum subsystems do not influence the evolution of the classical bath.The same approximation is used in 1D and 2D spectroscopy calculations using the frequency map approach 2 .One can expect a violation of this approximation for systems with fast energy relaxation from quantum degrees of freedom.Nevertheless, simulations using the frequency map approach show very good agreement with experiments for absorption, Raman and 2D IR spectroscopy, even for O-H stretch vibration of liquid water 3 , which has very fast energy dissipation (≈100s of fs 4 ).Hence, this approximation can be reasonable for a wide range of vibrational oscillators.
We can group terms describing the evolution of quantum and classical subsystems and transform Eq.S10 into a system of two coupled equations: In this model, the evolution of the quantum subsystem is affected by interaction with the bath, which evolves independently.
The equation of motion of the bath can be written using the classical Liouville operator 5 : The solution to this equation reads: We can write Eq.S11 using Liouville space notations: The Liouville superoperator describing the interaction of the subsystems parametrically depends on classical coordinates.
We consider a trajectory of classical bath described by functions () and ().For this trajectory, the solution of Eq.S15 is Where |  ( 0 )⟩⟩ is the initial state of the quantum subsystem.Thus, the state of the quantum subsystem at any time depends on the initial state |  ( 0 )⟩⟩ and trajectory () and () of the classical bath.
Because classical motion is deterministic and coordinates  and  at any time define the full trajectory of the system, for a given trajectory, we can consider |  ()⟩⟩ to depend parametrically on the coordinates of the classical bath at the same moment of time: Therefore, for the entire phase space of the classical bath, we can write: where evolution superoperator of quantum subsystem The density matrix of the entire system at time  in partial Wigner representation is given by: ⟨⟨|

B. Mixed quantum-classical response function in 2D TIRV spectroscopy
Using the quantum formalism, the response function of a sample measured in 2D TIRV spectroscopy is given by [6][7][8] : where  = [̂  , … ] is dipole moment superoperator (̂  is dipole moment operator of the system),  is the evolution superoperator and Π ̂ is the polarizability operator.We use identity superoperator  = ∬ |⟩⟩⟨⟨| to write response function in partial Wigner representation 9 : The polarizability is an operator in  ̂-space, which parametrically depends on coordinates  and  of the bath.
By taking the Taylor series for exponents and neglecting terms with ℏ  ,  ≥ 1, we obtain: By substituting Eq.S34 into Eq.S32, we obtain: Because (, , −∞) is an equilibrium distribution function, we can write 10 : Then: Equation S35 takes the form: We assume that the quantum subsystem is initially in the ground state and  ̂(, , −∞) = |⟩⟨|.We

Parameters of THz and mid-infrared vibrational resonances of DMSO in absorption spectrum.
Figure S3 shows the fitting of LFMs and HFMs in DMSO absorption spectrum to a superposition of Lorentzian (LFM2, LFM3, symmetric and asymmetric CH3 stretch) and Gaussian (LFM1) functions.The fitting parameters are summarized in Table S2.

First-order response of low-frequency modes of DMSO
Figure S4 shows the simulated spectrum of the first-order response of DMSO low-frequency modes.
absorbance of the  = 4 mode.These results demonstrate that the parameters of the 2D TIRV model are consistent with the DMSO absorption spectrum in Fig.S1, i.e., absorption of resonance  = 4 is weak and not pronounced.

FIG. S4 .
FIG. S4.The red line shows the imaginary part of the first-order response function of DMSO in the terahertz range.It is calculated using the same transition dipole moments and frequency fluctuations as in the 2D TIRV model (with mechanical coupling).The blue line shows a similar simulation but assumes zero dipole moment for the n=4 mode.

Comparison of lineshapes produced by electrical and mechanical anharmonicity FIG. S6 .
FIG. S6.Comparison of theoretical and experimental spectra of DMSO.Red contour lines show real (a,c) and imaginary (b,d) parts of the 2D TIRV response function calculated for the case of electrical (a,b) and mechanical (c,d) anharmonicity.Black contour lines show the measured spectrum.

TABLE S1 .
The fitting parameters of DMSO Raman spectrum.

TABLE S2 .
The fitting parameters of the DMSO absorption spectrum.