Probing photoinduced proton coupled electron transfer process by means of two-dimensional resonant electronic-vibrational spectroscopy
PProbing photoinduced proton coupled electron transfer process by means oftwo-dimensional electronic-vibrational spectroscopy
Jiaji Zhang, a) Raffaele Borrelli, b) and Yoshitaka Tanimura c)1) Department of Chemistry, Graduate School of Science, Kyoto University,Kyoto 606-8502, Japan DISAFA, University of Torino, Largo Paolo Braccini 2, I-10095 Grugliasco,Italy (Dated: 16 February 2021)
We develop a detailed theoretical model of photo-induced proton-coupled electron transfer(PPCET) processes, which are at the basis of solar energy harvesting in biological systemsand photovoltaic materials. Our model enables to analyze the dynamics and the efficiencyof a PPCET reaction under the influence of a thermal environment by disentangling thecontribution of the fundamental electron transfer (ET) and proton transfer (PT) steps. Inorder to study quantum dynamics of the PPCET process under an interaction with non-Markovian environment we employ the hierarchical equations of motion (HEOM). We cal-culate transient absorption spectroscopy (TAS) and two-dimensional electronic-vibrationalspectroscopy (2DEVS) signals in order to study the nonequilibrium reaction dynamics.Our results show that different transition pathways can be separated by TAS and 2DEVS.Keywords: Electron coupled proton transfer, Two-dimensional electronic-vibrational spec-troscopies, Hierarchical Equations of motion a) Electronic mail: [email protected] b) Electronic mail: [email protected] c) Electronic mail: [email protected] a r X i v : . [ phy s i c s . c h e m - ph ] F e b . INTRODUCTION The simultaneous transfer of protons and electrons play an important role in many natural andartificial energy conversion processes. A typical example is the oxygen evolving complex (OEC)of natural photosynthetic system, where the oxygen generation consists of four stepwise proton-coupled electron transfer (PCET) catalyzed reactions.
The specific pathways taken by electronsand protons, can lead to step-wise (consecutive) or concerted type reactions (CEPT). Unravellingthe detailed mechanistic aspects of the PCET process is fundamental for the design of artificialsolar energy utilization systems, for example, dye sensitized photo-electrochemical cell (DS-PEC)and many other bio-mimetic systems which have been developed for solar energy utilization andhydrogen reduction.
A variety of approximate quantum dynamical studies and analytical expressions of reaction rateconstant have been derived for isolated systems on the basis of the golden rule expression, linear re-sponse theory, and on Marcus’s theory of electron-transfer (ET) processes.
Their applicationshave been extended to condensed phase systems by further assuming a perturbative system-bathinteraction and a classical treatment of an environment representing, for example, solvent.
Rate constant for several PCET systems in thermal equilibrium conditions have also been com-puted with the aids of molecular dynamics simulations and quantum chemistry calculations.
Yet, the sole computation of reaction rates does not provide enough information to fully dis-entangle different ET and PT pathways and can hide important information about the role ofthe environment. Ultrafast nonlinear spectroscopy can be a powerful tool for unravelling themechanistic aspects of PPCET reactions and of photosynthesis in general. For example, in-frared (IR) transient absorption spectroscopy (TAS) has been applied to excited-state protontransfer and chemical bond cleavage, and can help to determine the relaxation mechanism afterinitial photoexcitation.
Luminescence TAS has indicated that the quantum effect of donor-acceptor (D-A) vibrations on PPCET is important for a full quantum treatment of the total reactionsystem.
These spectroscopic techniques have also been extended to multi-dimensional cases. Two-dimensional (2D) vibrational spectroscopy (2DVS) and 2D electronic spectroscopy (2DES)have been applied to condense phase transition and succeeded in investigating the electronic ex-citation dynamics and a structural change of molecules.
Their combination, 2D electronic-vibrational spectroscopy (2DEVS), has also been successfully applied to the photo-isomerization2eactions, mental-to-ligands transitions, conical intersection wavepacket dynamics, and ultrafastexcitonic photosynthetic energy transfer reactions.
By utilizing the UV-vis and IR pulses, weare now able to measure the coupling strength and coherence between the electronic and vibra-tional transitions as the off-diagonal peaks of 2D spectroscopy.
These features are useful forthe investigation of PPCET reaction dynamics.In this paper, we present a model of a PPCET reaction and provide a detailed analysis of itsdynamics by computing TAS and 2DEVS signals. We describe the coupled proton-electron dy-namics using two-dimensional potential energy surfaces, and complex system-bath interactions tosimulate a system in realistic conditions. We employ the numerically “exact” hierarchical equa-tions of motion (HEOM) approach to study the reduced system dynamics under non-perturbativeand non-Markovian system-bath interactions at finite temperature.
The paper is organized asfollows. In Sec. II, we derive a system-bath model for a prototypical PPCET process and introducethe HEOM approach for numerical simulation. The theory of non-linear response function is alsobriefly sketched in this section. In Sec. III, we present the calculated TAS and 2DEVS results andanalyze their profiles. Section IV is devoted to concluding remarks.
II. THEORYA. Model Hamiltonian
The system considered in the present work is depicted in figure 1. In the ground electronic statethe proton is localized at bond distance from donor D , and the D - H moiety is hydrogen bondedto the acceptor A . The x coordinate describes the position of the proton between D and A , while D H AQx
FIG. 1. Model PPCET system with a hydrogen bridge. Here, D ( A ) is the donor (acceptor), H is thetransferring proton, x is the proton coordinate describing its distance from the center of the D and A units,and Q is the distance between the heavy atoms connected by the hydrogen bond. Q is the distance between the heavy atoms which is also referred to as the reaction promoting3ode. We wish to describe the dynamics of the system resulting from the photo-excitation of D ,which is followed by a coupled transfer of an electron and a proton to the A moiety. As a result ofthe process an hydrogen atom is transferred from D - H to A , i.e. A is reduced to A - H and D - H isoxidized D .In order to model the coupled PT and ET processes we consider an electronic active spacecomprising the ground state | φ g (cid:105) of the system, the localized excited state | φ LE (cid:105) , in which onlythe moiety D is in the first excited electronic state, while A is in the ground electronic state,and the charge-transfer state | φ CT (cid:105) , in which D has transferred an electron to A . The diabaticrepresentation of system is shown in Fig. 2. The motion along the x and Q coordinates is describedemploying realistic two-dimensional potential energy surfaces. Furthermore, we assume that thesystem interacts with a condensed phase environment which can be either a solvent or a proteinscaffold. The overall Hamiltonian can therefore be expressed asˆ H = ∑ i ˆ H i ( ˆ x , ˆ Q , { q j } ) (cid:12)(cid:12) i (cid:105)(cid:104) i (cid:12)(cid:12) + ∑ i (cid:54) = j ∆ i j (cid:12)(cid:12) i (cid:105)(cid:104) j (cid:12)(cid:12) + ˆ H B , (1)where ˆ H i ( ˆ x , ˆ Q , { q j } ) is the Hamiltonian for the electronic states i = g , CT , and LE , ∆ i j are theelectronic couplings among different electronic states, and ˆ H B is the Hamiltonian of the thermalbath which is modeled as a collection of harmonic oscillatorsˆ H B = ∑ j (cid:32) ˆ p j m j + m j ω j ˆ q j (cid:33) , (2)where ˆ p j , ˆ q j , m j and ω j are the momentum, position, mass and frequency of j th bath oscillator,respectively.The operators ˆ H i can be explicitly written in the formˆ H i ( ˆ x , ˆ Q , { q j } ) = ˆ p x m x + ˆ p Q m Q + U i ( ˆ x , ˆ Q , { q j } ) + ε i . (3)Here, ˆ x , ˆ p x , m x , and ˆ Q , ˆ p Q , m Q are the coordinate, momentum and mass of the proton and of theD-A vibration, respectively, and ε i is the the energy of the electronic state for i . The potential U i ( ˆ x , ˆ Q , { q j } ) describes the variation of the electronic energy as a function of the coordinates x , Q and, furthermore, explicitly includes the interaction between these coordinates and the bath modes { q j } . Following previous work, we use an asymmetric double well Morse potential for the4 FIG. 2. Diabatic representation of the PES for the reduced system. The black curve represent the groundstate | φ g (cid:105) . The red, and blue curves represent the local excited state | φ LE (cid:105) , and charge transfer state | φ CT (cid:105) ,respectively. proton mode, and a harmonic potential for the D-A mode U i ( ˆ x , ˆ Q , { q j } ) = D li (cid:104) − e − α ( ˆ Q / + ˆ x − x e ( { q j } )) (cid:105) + D ri (cid:104) − e − α ( ˆ Q / − ˆ x − x e ( { q j } )) (cid:105) + D k ( ˆ Q − Q e ) , (4)where D li and D ri are the dissociation energy of donor (left well) and acceptor (right well), x e and Q e are the equilibrium distance of the proton and D-A vibrations, α represents the curvature ofthe Morse potentials, and D k is the force constant of the D-A vibration. In our model the roleof bath modes { q j } is to dynamically perturb the equilibrium position of the proton via a linearinteraction, that is x e ( { q j } ) = x − ∑ j g j ˆ q j , (5)where x is the equilibrium distance without the heat bath, and { g j } are coupling strength param-eters. Finally, we simplify this potential by expanding Eq. (4) in terms of the collective coordinateˆ X = ∑ j g j q j up to the first-order, which is similar to the reaction surface approach. The systemPES and the resulting exponential-linear (EL) system-bath interaction are then expressed as U i ( ˆ x , ˆ Q , X ) = U i ( ˆ x , ˆ Q ) + ˆ V i ( ˆ x , ˆ Q ) ˆ X , (6)5here U i ( ˆ x , ˆ Q ) = D li (cid:104) − e − α ( ˆ Q / + ˆ x − x ) (cid:105) + D ri (cid:104) − e − α ( ˆ Q / − ˆ x − x ) (cid:105) + D K ( ˆ Q − Q ) . (7)The operator V i depends solely on system variables and can be explicitly written as ˆ V i ( ˆ x , ˆ Q ) = α D li (cid:104) − e − α ( ˆ Q / + ˆ x − x ) (cid:105) e − α ( ˆ Q / + ˆ x − x ) (8) + α D ri (cid:104) − e − α ( ˆ Q / − ˆ x − x ) (cid:105) e − α ( ˆ Q / − ˆ x − x ) . (9)The structure of this rather complex form of system operator can be easily understood once weexpand it in terms of ˆ x and ˆ Q asˆ V i ( ˆ x , ˆ Q ) = V ( ) i + V ( ) i , x ˆ x + V ( ) i , Q ˆ Q + + V ( ) i , x ˆ x + V ( ) i , xQ ˆ x ˆ Q + ..., (10)where the V ( ) i , V ( ) i , x , etc. are constants whose analytical expressions are given in Appendix A.Hence it is clear that the coupling of Eq. (5) introduces linear interactions with the electronicsubsystem via the constant V ( ) i , and with the nuclear coordinates ˆ x and ˆ Q via V ( ) i , x and V ( ) i , Q . Asdiscussed in Ref. 47, the linear-linear (LL) interaction, such as V ( ) i , x ˆ xq j contributes mainly toenergy relaxation, while the square-linear (SL) system-bath interaction, such as V ( ) i , x ˆ x q j leads tovibrational dephasing in the slow modulation case, due to the frequency fluctuation of the systemvibrations. Finally we note that the LL contribution in the proton mode vanishes for symmetricdouble well potential, i.e. D li = D ri . For simplicity, we further assume that all of the electronic states are coupled to the same heatbath. Then, Eq. (1) can be rewritten asˆ H tot = ∑ i ˆ H i (cid:12)(cid:12) i (cid:11)(cid:10) i (cid:12)(cid:12) + ∑ i (cid:54) = j ∆ i j (cid:12)(cid:12) i (cid:11)(cid:10) j (cid:12)(cid:12) + ˆ H B + I , (11)where ˆ H i = ˆ p x m x + ˆ p Q m Q + U i ( ˆ x , ˆ Q ) (12)and ˆ H B + I = ∑ j ˆ p x , j m j + m j ω j (cid:34) ˆ q j − g j (cid:0) ∑ i (cid:12)(cid:12) i (cid:11)(cid:10) i (cid:12)(cid:12) ˆ V i (cid:1) m j ω j (cid:35) . (13)We also include the counter-term in the definition of ˆ H B + I in order to maintain the translationalsymmetry of the system. . Hierarchical Equations of Motion Approach Next, we briefly introduce the hierarchical equations of motion (HEOM) approach, which isemployed to investigate quantum dynamics of the PCET system in a numerically rigorous way.
We can also employ the multistate quantum hierarchical Fokker-Planck equations (MQHFPE),which has been applied to both optical and nonadiabatic transition problems described by complexPES.
However, here we choose the regular HEOM in the energy eigenstate representation forboth electronic and vibrational modes. This is because the proton motion is well confined in thePES, and the computational cost for using MQHFPE is much higher that regular HEOM.The heat bath is described by the spectral density function (SDF) J ( ω ) = π ∑ j g j m j ω j δ ( ω − ω j ) (14)and the inverse temperature, β = / k B T , where k B is the Boltzmann constant. The overall noiseeffect on the system is characterized by the correlation function C ( t ) = ¯ h (cid:90) ∞ d ω J ( ω ) (cid:20) coth (cid:18) β ¯ h ω (cid:19) cos ( ω t ) − i sin ( ω t ) (cid:21) . (15)In this paper, we use a Drude formed SDF, J ( ω ) = ζ π ωγ γ + ω , (16)where ζ represents the coupling strength, and γ is the reciprocal of the noise correlation time,representing the width of the spectral distribution. Then, Eq. (15) can be expressed in terms of acombination of linear exponential functions and of the δ ( t ) function, as C ( t ) = K ∑ k = ( c (cid:48) k + ic (cid:48)(cid:48) k ) γ k e − γ k t + c δ · δ ( t ) , (17)where c (cid:48) k , c (cid:48)(cid:48) k , γ k and c δ are constants determined by the chosen decomposition method. Here weemploy the Pad´e decomposition method which is known to enhance the efficiency of numer-ical calculations. By introducing the auxiliary density operators (ADO) ˆ ρ (cid:126) n , the HEOM can bederived as ∂∂ t ˆ ρ (cid:126) n ( t ) = − (cid:34) i ¯ h ˆ L S + ∑ k n k γ k + c δ ˆ Φ (cid:35) ˆ ρ (cid:126) n ( t ) − ∑ k ˆ Φ ˆ ρ (cid:126) n + (cid:126) e k ( t ) − ∑ k n k ˆ Θ k ˆ ρ (cid:126) n − (cid:126) e k ( t ) , (18)7here the superoperators are defined asˆ L S ˆ A ≡ [ ˆ H S , ˆ A ] , ˆ Φ u ˆ A ≡ i ¯ h (cid:2) ˆ V u , ˆ A (cid:3) , ˆ Ψ u ˆ A ≡ h { ˆ V u , ˆ A } , (19)for any physical operator ˆ A , and ˆ Θ k ≡ c (cid:48) k ˆ Φ − c (cid:48)(cid:48) k ˆ Ψ . The components of multi-index vector (cid:126) n =( ..., n k , ... ) are all non-negative integers, and (cid:126) e k is the k th unit vector. In the HEOM formalism,only the first element, (cid:126) n = ( , ..., ) , has a physical meaning, corresponding to the reduced densityoperator of system. The others are served as the treatment of non-perturbative and non-Markovianheat bath effect. Although Eq. (18) consists of infinite equations, we can truncate it at aproperly chosen large N value, for N = ∑ k n k . In order to reduce the computational cost for thetime integration, we rescale the ADOs as ˆ ρ (cid:126) n = ˆ ρ (cid:126) n / ∏ u , k (cid:112) n u , k !. Then, Eqs. (18) are rewrittenas ∂∂ t ˆ ρ (cid:126) n ( t ) = − (cid:34) i ¯ h ˆ L S + ∑ k n k γ k + c δ ˆ Φ (cid:35) ˆ ρ (cid:126) n ( t ) − ∑ k (cid:112) n k + Φ ˆ ρ (cid:126) n + (cid:126) e k ( t ) − ∑ k √ n k ˆ Θ k ˆ ρ (cid:126) n − (cid:126) e k ( t ) . (20) C. Projection operators for PT and ET states
In order to analyze the PCET process, next we introduce a set of projection operators definedas ˆ θ li = (cid:12)(cid:12) i (cid:105)(cid:104) i (cid:12)(cid:12) ˆ h ( − x ) , ˆ θ ri = (cid:12)(cid:12) i (cid:105)(cid:104) i (cid:12)(cid:12) ˆ h ( x ) , (21)where i = CT and LE , ˆ h ( x ) is the Heaviside step function for proton coordinate, and the symbols l and r represent the proton localized in the left (donor) and right (acceptor) well, respectively.The corresponding population of the superposition is P α i ( t ) = Tr (cid:8) ˆ θ α i ˆ ρ ( t ) (cid:9) for α = l or r . Thepopulations of | φ LE (cid:105) and | φ CT (cid:105) is then separated as P i ( t ) = P li ( t ) + P ri ( t ) , whereas that in the leftand right well is expressed as P α ( t ) = P α LE ( t ) + P α CT ( t ) .As shown in Fig. 2, the superposition (cid:12)(cid:12) φ lLE (cid:11) represents the configuration D ∗ − H · · · A , and | φ rLE (cid:105) represents D ∗− · · · H − A + . Similarly, (cid:12)(cid:12) φ lCT (cid:11) represents D + − H · · · A − , and (cid:12)(cid:12) φ rCT (cid:11) representsD · · · H − A. Thus, the pure PT process corresponds to the transitions (cid:12)(cid:12) φ lLE (cid:11) ↔ | φ rLE (cid:105) and (cid:12)(cid:12) φ lCT (cid:11) ↔ (cid:12)(cid:12) φ rCT (cid:11) . The pure ET process corresponds to the transitions (cid:12)(cid:12) φ lLE (cid:11) ↔ (cid:12)(cid:12) φ lCT (cid:11) and | φ rLE (cid:105) ↔ (cid:12)(cid:12) φ lCT (cid:11) .The CEPT process corresponds to the transitions (cid:12)(cid:12) φ lLE (cid:11) ↔ (cid:12)(cid:12) φ rCT (cid:11) and | φ rLE (cid:105) ↔ (cid:12)(cid:12) φ lCT (cid:11) .8 . Nonlinear Response Function The nonlinear response functions can be calculated within the framework of the HEOMformalism.
The third-order response function can be expressed as R ( ) ( t , t , t ) = (cid:18) i ¯ h (cid:19) Tr (cid:8) ˆ µ G ( t ) ˆ µ × G ( t ) ˆ µ × G ( t ) ˆ µ × ˆ ρ eq (cid:9) , (22)where ˆ µ i is the dipole operator of the i -th laser interaction, G ( t ) is the Green’s function of the totalHamiltonian without laser interactions, and ˆ ρ eq is the initial state density operator. In the HEOMapproach, the density matrix is replaced by a reduced one, and G ( t ) is evaluated from Eq. (18)(or Eq.(20)). The operator ˆ µ × i is the commutator of the dipole operator ˆ µ i . The right-hand sideof Eq. (22) can be evaluated as follows: The system is first in the initial equilibrium state ˆ ρ eq ,and is excited by the first interaction ˆ µ × at t =
0. The time evolution is computed by numericallyintegrating Eq. (18) up to a chosen time t . Then, the system is excited by the second and thirdinteractions ˆ µ × and ˆ µ × in a similar way. The final signal is computed by the expectation value ofˆ µ . We compute R ( ) ( t , t , t ) for a set of values of t , t , and t .Here we assume that the PES of | φ g (cid:105) and | φ LE (cid:105) have the same equilibrium positions, theirenergy difference is large and the population relaxation in the excited states is small. The directexcitation from (cid:12)(cid:12) φ g (cid:11) to | φ CT (cid:105) is also prohibited. Thus, the initial state is described by the thermalequilibrium distribution of the | φ g (cid:105) as ˆ ρ eq = ˆ ρ eqg . Assuming the PPCET reaction is initializedby a pair of impulsive pump pluses that excite the system from (cid:12)(cid:12) φ g (cid:11) to | φ LE (cid:105) , we set the initialconditions as ˆ ρ ( ) ( ) = − ˆ µ × ˆ µ × ˆ ρ eqg / ¯ h in the following response function analysis. With previousassumption, we can further set ˆ ρ ( ) ( ) = ˆ ρ eqLE , where ˆ ρ eqLE is evaluated as the steady state solution ofthe HEOM for the | φ LE (cid:105) state without non-adiabatic coupling with the | φ CT (cid:105) . Thus, our discussionin the following only considers the dynamics between | φ LE (cid:105) and | φ CT (cid:105) .The transient absorption response function can be evaluated from Eq. (22) by keeping t = R TA ( t , t (cid:48) ) = i ¯ h Tr (cid:110) ˆ µ G ( t ) ˆ µ × G ( t (cid:48) ) ˆ ρ ( ) ( ) (cid:111) . (23)Transient absorption spectrum (TAS) at different t (cid:48) is evaluated as I TA ( ω , t (cid:48) ) ≡ ω Im (cid:90) ∞ dte i ω t R TA ( t , t (cid:48) ) , (24)which also corresponds to linear absorption spectrum for non-equilibrium initial conditions. Thefifth-order transient 2DEVS is defined in a similar way as R ( ) ( t , t , t ) = (cid:18) i ¯ h (cid:19) Tr (cid:110) ˆ µ G ( t ) ˆ µ × G ( t ) ˆ µ × G ( t ) ˆ µ × ˆ ρ ( ) ( ) (cid:111) , (25)9here ˆ ρ ( ) ( ) is the same as TAS. The transient 2D correlation spectroscopy signals are thenevaluated as I Corr ( ω , t , ω ) = I ( NR ) ( ω , t (cid:48) , ω ) + I ( R ) ( ω , t , ω ) , (26)where the non-rephasing and rephrasing parts of the signal are expressed as I NR ( ω , t , ω ) = Im (cid:90) ∞ dt (cid:90) ∞ dt e i ω t e i ω t R ( ) ( t , t , t ) (27)and I R ( ω , t , ω ) = Im (cid:90) ∞ dt (cid:90) ∞ dt e i ω t e − i ω t R ( ) ( t , t , t ) . (28)The dipole operators ˆ µ j for j ≥ µ e = | φ LE (cid:105)(cid:104) φ CT | , or protonpart ˆ µ p = ˆ x · ( | φ LE (cid:105)(cid:104) φ LE | + | φ CT (cid:105)(cid:104) φ CT | ) , or their summation ˆ µ e + ˆ µ p . TABLE I. System parameters α − x Q D k − ˚A − ∆
50 cm − D lLE − D rLE − D lCT − D rCT − III. NUMERICAL RESULTS
The system parameters chosen to simulate our PPCET model are listed in Table. I, based ona typical PT system. The determination of the electronic couplings ∆ is a critical point of anyPCCET reaction, in that it provides the major contribution to the discrimination between adia-batic and non-adiabatic mechanisms. Here, we choose to study the system under moderate non-adiabatic conditions, and set ∆ =
50 cm − , which is close to previously reported studies. The10
ABLE II. The lowest 10 energy eigenvalues of each electronic state as a unit of ω .Eigen numbers ( m , n ) (cid:12)(cid:12)(cid:12) φ ( m , n ) LE (cid:69) (cid:12)(cid:12)(cid:12) φ ( m , n ) CT (cid:69) (0, 0) 0.00 -0.02(0, 1) 0.82 0.81(0, 2) 1.63 1.64(1, 0) 1.94 1.95(0, 3) 2.44 2.47(1, 1) 2.78 2.74(2, 0) 3.06 3.06(0, 4) 3.26 3.31(1, 2) 3.64 3.55(2, 1) 3.94 3.97 energy eigenstates of the system (cid:12)(cid:12) φ ( m , n ) i (cid:11) are obtained by diagonalizing the matrix representationof the system Hamiltonian. The lowest several energy eigenvalues are presented in Table. II,and a schematic view of them are illustrated in Fig. 3. Here, the m and n represent the quantumnumbers of the proton and D-A modes that are determined from the number of nodes along the x and Q directions. While, (cid:12)(cid:12) φ ( m , n ) i (cid:11) with m = (cid:12)(cid:12) φ lLE (cid:11) and (cid:12)(cid:12) φ rCT (cid:11) , those with m = | φ rLE (cid:105) and (cid:12)(cid:12) φ lCT (cid:11) , respectively. The states for m (cid:61) ζ . The time integration is carried out using the low-storage fourth-order Runge-Kutta (LSRK4)method. The time step is chosen as δ t = . ω − , where ω is a characteristic frequency takenas the unit for all the other physical variables. Here, we choose ω = − . We also fix theinverse correlation time as γ = . ω and the bath temperature as β ¯ h ω = . N =
10 and K =
5. In the following,we investigate the effects of the environment on the PCET mechanism as a function of ζ on thebasis of both the population dynamics and the nonlinear optical signals.11 (0, 0)(0, 1)(0, 2) (1, 0)(2, 0)(3, 0) (0, 0)(0, 1)(0, 2)(1, 0)(1, 1) ECG F AB DHI FIG. 3. A schematic view for | φ LE (cid:105) (red curve) and | φ CT (cid:105) (blue curve) in the diabatic representation alongˆ x at the minimum of ˆ Q . The lowest several eigenstates for each PES are also plotted. The labeled or-ange and green arrows represent the corresponding proton and electron transitions appearing in nonlinearspectroscopy. See main text for the meaning of the labels. A. Population dynamics
First, we illustrate the time evolution of the population states for various values of ζ . Theelectron and proton transfer rates can be estimated from P CT ( t ) and P l ( t ) . The calculated resultsare depicted in Fig. 4, for a (a) weak ( ζ = . ω ) , (b) moderate ( ζ = . ω ) , and (c) strong ( ζ = . ω ) coupling cases. Note that, as illustrated in our PT investigation, the effectivecoupling strength on the present EL system-bath coupling model is different from the conventionalLL coupling model. The strength of the coupling parameter is determined on the basis of therelaxation dynamics of the populations and spectral line shape of TAS, as we will show below.For the weak coupling case in Fig. 4(a), coherent recursive oscillations of state populationsare observed. These oscillations do not affect the equilibrium distribution and do not contributeto the population transfer rates. Although the contribution is minor, the population exchangebetween P CT ( t ) and P l ( t ) suggests the presence of a charge transfer process. For the moderate and12trong coupling cases in Fig. 4(b) and (c), the linear term of ˆ V i causes the population relaxationsuppressing the coherent oscillations. In the ˆ x direction (proton mode), the nonlinear terms ofˆ V i also lead to a decrease of the energy barrier so that proton transfer is promoted. In the ˆ Q direction, the linear term of ˆ V i leads to a decrease of proton distance from the heavy atoms andincrease the PT efficiency. A constant term ˆ V i ( , ) (see Eq. (30)) corresponding to the interactionbetween the electronic states and the heat bath also play a role; for larger ζ , both the electron andthe proton are equally distributed in the two wells because of the symmetric PES. Note that wecannot disentangle the contribution of the CEPT, ET, and PT processes only from the analysis ofpopulation dynamics, because their contributions are mixed in the population states. B. Transient absorption spectroscopy (TAS)
Next we present the results of TAS analysis. Although TAS has a capability to analyze thepopulations in the ET and PT states separately on the basis of the peak positions of absorptionspectrum, this is not easy in the present case, because the excitation energies of the ET and PTprocesses are similar, and the abortion peaks are often overlapped. Hence, here we calculate TASfor the electronic and vibrational modes separately to help the analysis of 2DEVS. In TAS, the P o p u l a t i o n (a), = 0.0005 Time / (b), = 0.005 P lLE P rLE P lCT P rCT P l P CT (c), = 0.03 FIG. 4. The population dynamics that represents proton and electron localization. The red, yellow, greenand blue curves represent P lLE ( t ) , P rLE ( t ) , P lCT ( t ) , and P rCT ( t ) , respectively. The configurations of modelsystem is illustrated in Fig. 2. The population P l ( t ) and P CT ( t ) are also presented as the dashed purple andblack curves. / ABCDEFG (a), = 0.0005 t / (b), = 0.005 (c), = 0.03 FIG. 5. The contour map of TAS calculated for the vibrational excitation ( ˆ µ = ˆ µ p ) in the cases of the (a)weak, (b) moderate, and (c) strong system-bath interactions. All the peaks are normalized with respect tomaximum intensity of t (cid:48) = . ω − . The contour lines are drawn from -0.5 to 0.5. The red and blue areasrepresent the positive absorption and negative emission, respectively. charge transition rate can be evaluated from the intensity of corresponding transition peaks, whilecoherent oscillation appears as a δ -function like peak. The characteristic time scale of varioustransitions can also be evaluated as a function of t (cid:48) .In Fig. 5 we present TAS for the vibrational excitation of the proton mode for a waiting timeup to t (cid:48) = . ω − with setting ˆ µ = ˆ µ = ˆ µ p . In each figure, the negative and positive peaksrepresent the emission and absorption, respectively. Note that, although the energy eigenvaluesof the (cid:12)(cid:12) φ ( m , n ) LE (cid:11) and (cid:12)(cid:12) φ ( m , n ) CT (cid:11) in the diabatic representation are degenerate, those in the adiabaticrepresentation are separated by the frequency ∆ because of the diabatic coupling.In the weak coupling case, Fig. 5(a), the peak “A” ( . ω ) predominantly arises from theCPET, (cid:12)(cid:12) φ ( m , n ) i (cid:11) → (cid:12)(cid:12) φ ( m , n ) j (cid:11) . This transition always occurs due to the large overlap between twoelectronic potential surfaces. The peak “B” ( . ω ) and “E” ( . ω ) arise from the pure PT withand without the participation of the D-A mode, where “E” represents (cid:12)(cid:12) φ ( , n ) i (cid:11) → (cid:12)(cid:12) φ ( , n ) i (cid:11) , and “B”represents (cid:12)(cid:12) φ ( , n ) i (cid:11) → (cid:12)(cid:12) φ ( , n − ) i (cid:11) . The peak “C” ( . ω ) represents the excitation of the D-A mode, (cid:12)(cid:12) φ ( m , n ) i (cid:11) → (cid:12)(cid:12) φ ( m , n + ) i (cid:11) , which arises because the proton and the D-A mode are strongly coupled.The proton distribution varies as a function of the quantum number n in the D-A mode, even whenthe quantum number of the proton mode m is unchanged. The other three peaks represent thedelocalization of the proton in the higher energy states ( m (cid:53) ( . ω ) , “F” ( . ω ) , and “G” ( . ω ) represent (cid:12)(cid:12) φ ( , n ) i (cid:11) → (cid:12)(cid:12) φ ( , n ) i (cid:11) , (cid:12)(cid:12) φ ( , n ) i (cid:11) → (cid:12)(cid:12) φ ( , n ) i (cid:11) , and (cid:12)(cid:12) φ ( , n ) i (cid:11) → (cid:12)(cid:12) φ ( , n ) i (cid:11) , respectively. Most of these peaks consist of severalsmall peaks because of the participation of the D-A mode excited states ( n > ) . A schematic viewof all the transitions is illustrated in Fig. 3.We then analyze the effects of ζ through the peak intensities as a function of t (cid:48) . In the weakcoupling case, Fig. 5(a), most of the peaks are unchanged regardless of t (cid:48) except for the peak “A”,whose sign of the intensity changes near t (cid:48) = .
0. In the moderate and strong coupling cases, Fig.5(b) and 5(c), the peak intensity of “A” changes from almost 0 to a positive value in the initialtime period. This indicates that the CEPT process is promoted by the system-bath interaction, andoccurs in a relatively short time period. This promotion effect can be explained by the linear term V ( ) i , Q in the ˆ Q direction, which reduces the transfer distance that enhances the vibronic coupling.The intensity of the peak “C” changes from positive to negative values near t (cid:48) = . t (cid:48) = . / HI (a), = 0.0005 t / (b), = 0.005 (c), = 0.03 FIG. 6. The contour maps of TAS for the electronic excitation ( ˆ µ = ˆ µ e ) in the cases of the (a) weak, (b)moderate, and (c) strong system-bath interactions. The figure is depicted in the same way with Fig. 5. Finally we present TAS for electronic excitation case in Fig. 6, which is computed by settingˆ µ = ˆ µ = ˆ µ e . In the weak coupling case in Fig. 6(a), the peak labeled by “H” corresponds tothe transition (cid:12)(cid:12) φ ( , n ) i (cid:11) → (cid:12)(cid:12) φ ( , n ) j (cid:11) , and “I” corresponds to (cid:12)(cid:12) φ ( , n ) i (cid:11) → (cid:12)(cid:12) φ ( , n + ) j (cid:11) . These two peaksrepresent the pure ET; a possible n increases for larger ζ x . In the moderate and strong coupling15 (a)-(i), t = 0.0 (a)-(ii), t = 1.0 (a)-(iii), t = 10.0 / a1a2 / a3 b1 b2 b3 b4 b5 FIG. 7. The contour maps of 2DEVS for a weak coupling case. The intensities are normalized with respectto the maximum value of each case in order to see the peak profile with the contour lines are drawn from-0.5 to 0.5. The red and blue curves represent the absorption and emission, respectively. The boxes areplotted in different cases in order to have a clear understanding. We also plot the peaks along the diagonalline in the outside above, i.e. I ( ω , t , ω ) . cases, 6(b) and 6(c), peaks “H” and “I” are significantly broadened and enhanced because of theconstant term V ( ) i , which introduces linear interactions between the electron subsystem and theheat bath. Furthermore, several additional peaks appear in the range of 0 . ≤ ω ≤ . ω . Thesepeaks arise from the electronic transition, and the corresponding transitions are the same withproton mode, as shown in Fig. 5. This result can be ascribed to the strong correlation between theelectronic subsystem and the vibrational coordinates. The increase of ζ has a promotion effect onthe proton transfer, which in turn opens additional transition pathways of the electron transfer.We also find that most of the peaks are unchanged regardless of t (cid:48) even in the strong couplingcase. Thus, the pure ET process is not favored in all the ζ cases because of the pretty smallelectronic coupling strength ∆ . With regard to the CEPT peak near ω = . ω , the sign ofthe peak intensity changes in both weak and strong coupling cases, as can be clearly seen fromFigs. 6(a) and 6(c). Such variation becomes prominent in the moderate coupling case around t (cid:48) = . ω − , if Fig. 6(b), which indicates that a turn-over feature under a strong enough interactionoccurs. According to the above results, we find that the CEPT is the predominant process mostly,because of the exact resonance conditions between initial and final states.16 . Two-dimensional electronic-vibrational spectroscopy (2DEVS) Next, we describe the 2DEVS signals as computed from Eq. (25). We consider an experimentalcondition in which the dipole operators is ˆ µ = ˆ µ p + ˆ µ e . The calculated signals involves the con-tribution from 2DES, 2DVS and 2DEVS. The contour maps of the 2D correlation spectroscopy inthe weak, moderate, and strong coupling cases are shown in Figs. 7, 8, and 9, in which we keepthe system parameters the same used to obtain the TAS signals. Note that most of the peaks alongthe diagonal line are relatively weak and not clearly visible in contour maps. Therefore, we plotthese peaks outside above as I Diag ( ω , t ) = I Corr ( ω , t , ω ) .Using the information obtained from TAS, we classify all the observed peaks into three parts.The cross peaks of 2DEVS in the red boxes arise from the coherent transition between the ETand PT processes. The diagonal and cross peaks in the blue box arise from the CEPT. The otherpeaks are the cross peaks of 2DVS, representing the coherent transition of the proton and D-Amodes. The cross peaks of 2DES for the ET process are not visible in the weak coupling case,because of the small electronic coupling ∆ . For each peak, the positive intensity arises from thestimulated emission (SE) or the ground state bleaching (GSB), and the negative intensity arisesfrom the excited state absorption (ESA) for n > ω position asFig. 6, and at the same ω position as Fig. 5, representing the corresponding ET-PT transitions.Here, we only concentrate on the peaks “a1”, “a2” and “a3” that do not appear in TAS. In the ω direction, the peak “a1” represents the transition (cid:12)(cid:12) φ ( , n ) i (cid:11) → (cid:12)(cid:12) φ ( , n + ) i (cid:11) that arises from to the SLkind of interaction in ˆ V i in the ˆ Q direction. The peaks “a2” and “a3” represent (cid:12)(cid:12) φ ( , ) i (cid:11) → (cid:12)(cid:12) φ ( , ) i (cid:11) that arise from the back PT process with a participation of the D-A mode. These results indicatethat we can analyze the combination of the ET and PT transitions from the cross peaks in 2DEVS,while these contributions are mixed and appear as a single peak in TAS.Next, we discuss the peaks in the blue box. The diagonal peak “b1” arises from the CEPTtransition denoted as “A” in Fig. 3. The other cross peaks represent the combination of CEPT-PT and CEPT-ET, where “b2” and “b5” correspond to “A”-“D” and “A”-“G”, and “b3” and “b4”correspond to “A”-“H” and “A”-“I”, respectively. The peaks associated with the CEPT processin blue box appear at symmetric positions with the diagonal line. The other 2DIR cross peaksmainly occur between the proton delocalization and D-A vibration, representing a reduction of17nergy barrier by the D-A motion. Because of the weak heat-bath effect, the position and profilesof the ET-PT and CEPT peaks do not change very much with t . By contrast, the intensities ofthe 2DIR cross peaks decrease and finally vanish when the excited proton reaches the equilibriumdistribution due to the linear interaction V ( ) i , x . (b)-(i), t = 0.0 (b)-(ii), t = 1.0 (b)-(iii), t = 10.0 / / FIG. 8. The contour maps of 2DEVS for a moderate coupling case. The contour lines are drawn from -1.0to 1.0, while the other parameters are the same with Fig. 7.
In the moderate and strong coupling cases shown in Figs. 8 and 9, most of the peaks that arerelated to the proton and D-A transitions are broadened either in the ω or ω direction. For theET-PT cross peaks in the red box, the peak positions in the ω direction are almost unchanged,which indicates that the system-bath interaction has a minor effect on the pure ET, as observedin TAS. In the moderate coupling case, the intensities of the ET-PT peaks increase in Figs. 8(b)-(ii) and 8(b)-(iii) in comparison with Fig. 8(b)-(i), which indicates that the time scale of ET-PTtransition is relatively short. Most of these peaks are also enhanced and broadened in either the ω or ω direction if they are related to proton or D-A mode, compared with the case in Figs. 7.However, most of them almost disappear in the strong coupling case, due to the suppression on theET-PT coherence. This turn-over feature with ζ is also observed in TAS for electronic excitation,and is more clear in the 2DEVS.For the CEPT peaks in blue box, the intensities increase with t , as evident from Fig. 8, andbecome more apparent in Figs. 9 for a larger ζ . This result indicates the existence of a bathinduced CEPT process, which is also observed in TAS. For the 2DIR cross peaks, most of them18o not change till t = .
0, as illustrated in Figs. 8(b)-(ii) and 9(c)-(ii), and almost vanish after along t time as illustrated in Figs. 8(b)-(iii) and 9(c)-(iii). Moreover, the intensity of the twistedpositive and negative cross peak around ( ω , ω ) = ( . ω , . ω ) is reversed at t = .
0. Thispeak mainly arises from the combination of the “B” and “C” transitions (see Fig. 3), and thereverse transition indicates the relaxation of excited proton. Thus, both proton and D-A motionhave relatively longer time scales compared to CPET and ETPT, and they are always mixed. (c)-(i), t = 0.0 (c)-(ii), t = 1.0 (c)-(iii), t = 10.0 / / FIG. 9. The contour map of 2DEVS for a strong coupling case. The contour lines are drawn from -1.0 to1.0, while the other parameters are the same with Fig. 7.
Finally, we concentrate on the peaks along the diagonal line, I Diag ( ω , t ) , which represent theadiabatic transitions. For all the coupling cases, these peaks occur in the same position found inthe TAS signals, representing the corresponding transitions. Among different t cases, the peaksrepresenting the CEPT transition “A” play a major role. The proton and D-A mode vibrationsare only visible after t = . V. CONCLUSION
In this paper we introduce a system-bath model in a multi-state two-dimensional configurationspace to describe the dynamics of PPCET process. Using the HEOM in the eigenstate representa-tion of the system, it is possible to investigate the environment effects under a realistic system-bathinteraction that causes not only fluctuation and relaxation, but also vibrational dephasing. Our re-sults of population dynamics and TAS indicate that CEPT is the predominant process and has ashorter time scale when resonance conditions between initial and final states occurr. Pure ET andPT processes also take place at much longer time. The overall reaction would be a summation ofboth concerted and sequential reaction mechanism. It is shown that 2DEVS provides a wealth ofinformation due to the coherence among the excitation and detection periods. With the aids of theoff-diagonal peaks, we could detect the pathway of sequential ET-PT and PT-ET transition, andconcerted CEPT transition separately, whereas the diagonal peaks could reproduce the results ofTAS.Although calculating nonlinear spectra is numerically intensive, 2DEVS with TAS provides avaluable framework for studying PPCET processes. Since we use the eigenstate representationof the system, it is also possible to improve the description of the reacting system by increas-ing the dimension of its configuration space, and by introducing a more complex and structuredsystem-bath interaction, for example, with the help of machine learning approaches.
This pro-vides a powerful tool to analyze the non-equilibrium reaction dynamics for rather complex PPCETreactions.
ACKNOWLEDGMENTS
The financial support from The Kyoto University Foundation is acknowledged. RB acknowl-edges the support of the University of Torino for the local research funding Grant No. BORR-RILO-18-01.
DATA AVAILABILITY
The data that support the findings of this study are available from the corresponding authorupon reasonable request. 20 ppendix A EXPANSION OF ˆ V i In this Appendix, we expand the interaction function, ˆ V i ( ˆ x , ˆ Q ) , with respect to ˆ x and ˆ Q up tosecond order as,ˆ V i ( ˆ x , ˆ Q ) = V ( ) i + V ( ) i , x ˆ x + V ( ) i , Q ˆ Q + V ( ) i , x ˆ x + V ( ) i , Q ˆ Q + V ( ) i , xQ ˆ x ˆ Q + ... (29)where V ( ) i = ˆ V i ( , ) = α ( D li + D ri ) (cid:0) α e α x − α e α x (cid:1) , (30) V ( ) i , x = ∂ ˆ V i ( ˆ x , ˆ Q ) ∂ ˆ x (cid:12)(cid:12)(cid:12)(cid:12) ( , ) = α ( D li − D ri ) (cid:0) α e α x − α e α x (cid:1) , (31) V ( ) i , Q = ∂ ˆ V i ( ˆ x , ˆ Q ) ∂ Q (cid:12)(cid:12)(cid:12)(cid:12) ( , ) = α ( D li + D ri ) (cid:0) α e α x − α e α x (cid:1) , (32) V ( ) i , x = ∂ ˆ V i ( x , Q ) ∂ ˆ x (cid:12)(cid:12)(cid:12)(cid:12) ( , ) = α ( D li + D ri ) (cid:0) α e α x − α e α x (cid:1) , (33) V ( ) i , Q = ∂ ˆ V i ( ˆ x , ˆ Q ) ∂ Q (cid:12)(cid:12)(cid:12)(cid:12) ( , ) = α ( D li + D ri ) (cid:0) α e α x − α e α x (cid:1) , (34) V ( ) i , xQ = ∂ ˆ V i ( ˆ x , ˆ Q ) ∂ ˆ x ˆ Q (cid:12)(cid:12)(cid:12)(cid:12) ( , ) = α ( D li − D ri ) (cid:0) − α e α x + α e α x (cid:1) . (35) REFERENCES A. Migliore, N. F. Polizzi, M. J. Therien, and D. N. Beratan, “Biochemistry and theory ofproton-coupled electron transfer,” Chemical Reviews , 3381–3465 (2014), pMID: 24684625,https://doi.org/10.1021/cr4006654. J. D. Megiatto Jr, D. D. M´endez-Hern´andez, M. E. Tejeda-Ferrari, A.-L. Teillout, M. J. Llansola-Portol´es, G. Kodis, O. G. Poluektov, T. Rajh, V. Mujica, T. L. Groy, D. Gust, T. A. Moore, andA. L. Moore, “A bioinspired redox relay that mimics radical interactions of the tyr–his pairs ofphotosystem ii,” Nature Chemistry , 423–428 (2014). M. T. Huynh, S. J. Mora, M. Villalba, M. E. Tejeda-Ferrari, P. A. Liddell, B. R. Cherry,A.-L. Teillout, C. W. Machan, C. P. Kubiak, D. Gust, T. A. Moore, S. Hammes-Schiffer,and A. L. Moore, “Concerted one-electron two-proton transfer processes in models in-spired by the tyr-his couple of photosystem ii,” ACS Central Science , 372–380 (2017),https://doi.org/10.1021/acscentsci.7b00125. 21 A. R. Offenbacher and B. A. Barry, “A proton wire mediates proton coupled electron trans-fer from hydroxyurea and other hydroxamic acids to tyrosyl radical in class ia ribonucleotidereductase,” The Journal of Physical Chemistry B , 345–354 (2020), pMID: 31904962,https://doi.org/10.1021/acs.jpcb.9b08587. C. J. Gagliardi, B. C. Westlake, C. A. Kent, J. J. Paul, J. M. Papanikolas, and T. J. Meyer,“Integrating proton coupled electron transfer (pcet) and excited states,” Coordination ChemistryReviews , 2459 – 2471 (2010), 18th International Symposium on the Photochemistry andPhotophysics of Coordination Compounds Sapporo, 2009. R. Becker, S. Amirjalayer, P. Li, S. Woutersen, and J. N. H. Reek, “Aniron-iron hydrogenase mimic with appended electron reservoir for efficient protonreduction in aqueous media,” Science Advances (2016), 10.1126/sciadv.1501014,https://advances.sciencemag.org/content/2/1/e1501014.full.pdf. A. Yamaguchi, R. Inuzuka, T. Takashima, T. Hayashi, K. Hashimoto, and R. Nakamura, “Reg-ulating proton-coupled electron transfer for efficient water splitting by manganese oxides atneutral ph,” Nature Communications , 4256 (2014). J. F. Allen, “Photosynthesis of atp& , 273–276 (2002). S. Papa, G. Capitanio, and P. Luca Martino, “Concerted involvement of cooperative pro-ton–electron linkage and water production in the proton pump of cytochrome c oxidase,”Biochimica et Biophysica Acta (BBA) - Bioenergetics , 1133 – 1143 (2006), mitochon-dria: from Molecular Insight to Physiology and Pathology. R. I. Cukier, “Proton-coupled electron transfer reactions: evaluation of rate constants,” TheJournal of Physical Chemistry , 15428–15443 (1996), https://doi.org/10.1021/jp961025g. Z. K. Goldsmith, A. V. Soudackov, and S. Hammes-Schiffer, “Theoretical analysis of the in-verted region in photoinduced proton-coupled electron transfer,” Faraday Discuss. , 363–378(2019). J. M. Mayer, “Simple marcus-theory-type model for hydrogen-atom transfer/proton-coupledelectron transfer,” The Journal of Physical Chemistry Letters , 1481–1489 (2011), pMID:21686056, https://doi.org/10.1021/jz200021y. A. V. Soudackov and S. Hammes-Schiffer, “Proton-coupled electron transfer reactions: analyt-ical rate constants and case study of kinetic isotope effects in lipoxygenase,” Faraday Discuss. , 171–189 (2016). 22 A. Soudackov and S. Hammes-Schiffer, “Derivation of rate expressions for nonadiabatic proton-coupled electron transfer reactions in solution,” The Journal of Chemical Physics , 2385–2396 (2000), https://doi.org/10.1063/1.482053. A. V. Soudackov and S. Hammes-Schiffer, “Nonadiabatic rate constants for proton trans-fer and proton-coupled electron transfer reactions in solution: Effects of quadratic term inthe vibronic coupling expansion,” The Journal of Chemical Physics , 194101 (2015),https://doi.org/10.1063/1.4935045. R. I. Cukier, “Mechanism for proton-coupled electron-transfer reactions,” The Journal of Physi-cal Chemistry , 2377–2381 (1994), https://doi.org/10.1021/j100060a027. A. Hazra, A. V. Soudackov, and S. Hammes-Schiffer, “Role of solvent dynamics in ultrafastphotoinduced proton-coupled electron transfer reactions in solution,” The Journal of PhysicalChemistry B , 12319–12332 (2010), pMID: 20809583, https://doi.org/10.1021/jp1051547. J. S. Kretchmer and T. F. Miller, “Tipping the balance between concerted versus sequen-tial proton-coupled electron transfer,” Inorganic Chemistry , 1022–1031 (2016), pMID:26440812, https://doi.org/10.1021/acs.inorgchem.5b01821. M. N. Kobrak and S. Hammes-Schiffer, “Molecular dynamics simulation of proton-coupledelectron transfer in solution,” The Journal of Physical Chemistry B , 10435–10445 (2001),https://doi.org/10.1021/jp012102x. J. Grimminger and W. Schmickler, “Dynamics of combined electron- and proton transfer atmetal electrodes,” Chemical Physics , 8 – 17 (2007). B. Auer, L. E. Fernandez, and S. Hammes-Schiffer, “Theoretical analysis of proton relays inelectrochemical proton-coupled electron transfer,” Journal of the American Chemical Society , 8282–8292 (2011), pMID: 21524104, https://doi.org/10.1021/ja201560v. A. Peluso, M. Di Donato, and G. A. A. Saracino, “An alternative way of thinking about electrontransfer in proteins: Proton assisted electron transfer between the primary and the secondaryquinones in photosynthetic reaction centers,” The Journal of Chemical Physics , 3212–3218(2000), https://doi.org/10.1063/1.1286918. M. Cho,
Coherent Multidimensional Spectroscopy (Springer, Singapore, 2019). M. Pfeiffer, C. Chudoba, A. Lau, K. Lenz, and T. Elsaesser, “Coherent motion of low frequencyvibrations in ultrafast excited state proton transfer,” Laser Chemistry , 101–103 (1999). D. G. Hogle, A. R. Cunningham, and M. J. Tucker, “Equilibrium versus nonequilibrium peptidedynamics: Insights into transient 2d ir spectroscopy,” The Journal of Physical Chemistry B ,23783–8795 (2018), pMID: 30040900, https://doi.org/10.1021/acs.jpcb.8b05063. M. Kuss-Petermann, H. Wolf, D. Stalke, and O. S. Wenger, “Influence of donor–acceptordistance variation on photoinduced electron and proton transfer in rhenium(i)–phenol dyads,”Journal of the American Chemical Society , 12844–12854 (2012), pMID: 22809316,https://doi.org/10.1021/ja3053046. Y. Giret, P. Guo, L.-F. Wang, and J. Cheng, “Theoretical study of kinetics of proton cou-pled electron transfer in photocatalysis,” The Journal of Chemical Physics , 124705 (2020),https://doi.org/10.1063/5.0001825. R. Zheng, Y. Jing, L. Chen, and Q. Shi, “Theory of proton coupled electron transfer reactions:Assessing the born–oppenheimer approximation for the proton motion using an analyticallysolvable model,” Chemical Physics , 39 – 45 (2011). K. Song and Q. Shi, “Theoretical study of photoinduced proton coupled electron transfer reac-tion using the non-perturbative hierarchical equations of motion method,” The Journal of Chem-ical Physics , 184108 (2017), https://doi.org/10.1063/1.4982928. P. F. Barbara, G. C. Walker, and T. P. Smith, “Vibrational modes and the dy-namic solvent effect in electron and proton transfer,” Science , 975–981 (1992),https://science.sciencemag.org/content/256/5059/975.full.pdf. Y. Tanimura and S. Mukamel, “Two-dimensional femtosecond vibrational spec-troscopy of liquids,” The Journal of Chemical Physics , 9496–9511 (1993),https://doi.org/10.1063/1.465484. S. Mukamel,
Principles of Nonlinear Optical Spectroscopy (Oxford University Press, New York,1995). P. Hamm and M. Zanni,
Concepts and Methods of 2D Infrared Spectroscopy (Cambridge Uni-versity Press, 2011). J. D. Gaynor, T. L. Courtney, M. Balasubramanian, and M. Khalil, “Fourier transform two-dimensional electronic-vibrational spectroscopy using an octave-spanning mid-ir probe,” Opt.Lett. , 2895–2898 (2016). T. A. A. Oliver, “Recent advances in multidimensional ultra-fast spectroscopy,” Royal Society Open Science , 171425 (2018),https://royalsocietypublishing.org/doi/pdf/10.1098/rsos.171425. Z. W. Fox, T. J. Blair, and M. Khalil, “Determining the orientation and vibronic couplingsbetween electronic and vibrational coordinates with polarization-selective two-dimensional24ibrational-electronic spectroscopy,” The Journal of Physical Chemistry Letters , 1558–1563(2020), pMID: 32004009, https://doi.org/10.1021/acs.jpclett.9b03752. S. Mukamel, “Multidimensional femtosecond correlation spectroscopies of electronic and vi-brational excitations,” Annual Review of Physical Chemistry , 691–729 (2000), pMID:11031297, https://doi.org/10.1146/annurev.physchem.51.1.691. H. Dong, N. H. C. Lewis, T. A. A. Oliver, and G. R. Fleming, “Determining the static electronicand vibrational energy correlations via two-dimensional electronic-vibrational spectroscopy,”The Journal of Chemical Physics , 174201 (2015), https://doi.org/10.1063/1.4919684. E. C. Wu, Q. Ge, E. A. Arsenault, N. H. C. Lewis, N. L. Gruenke, M. J. Head-Gordon, and G. R.Fleming, “Two-dimensional electronic-vibrational spectroscopic study of conical intersectiondynamics: an experimental and electronic structure study,” Phys. Chem. Chem. Phys. , 14153–14163 (2019). T. Ikeda and Y. Tanimura, “Phase-space wavepacket dynamics of internal conversion via coni-cal intersection: Multi-state quantum fokker-planck equation approach,” Chemical Physics ,203 – 213 (2018), ultrafast Photoinduced Processes in Polyatomic Molecules:Electronic Struc-ture, Dynamics and Spectroscopy (Dedicated to Wolfgang Domcke on the occasion of his 70thbirthday). N. H. C. Lewis and G. R. Fleming, “Two-dimensional electronic-vibrational spectroscopy ofchlorophyll a and b,” The Journal of Physical Chemistry Letters , 831–837 (2016), pMID:26894783, https://doi.org/10.1021/acs.jpclett.6b00037. T. L. Courtney, Z. W. Fox, K. M. Slenkamp, and M. Khalil, “Two-dimensionalvibrational-electronic spectroscopy,” The Journal of Chemical Physics , 154201 (2015),https://doi.org/10.1063/1.4932983. E. C. Wu, E. A. Arsenault, P. Bhattacharyya, N. H. C. Lewis, and G. R. Fleming, “Two-dimensional electronic vibrational spectroscopy and ultrafast excitonic and vibronic photosyn-thetic energy transfer,” Faraday Discuss. , 116–132 (2019). Y. Tanimura and R. Kubo, “Time evolution of a quantum system in contact with a nearlygaussian-markoffian noise bath,” Journal of the Physical Society of Japan , 101–114 (1989),https://doi.org/10.1143/JPSJ.58.101. Y. Tanimura, “Nonperturbative expansion method for a quantum system coupled to a harmonic-oscillator bath,” Phys. Rev. A , 6676–6687 (1990).25 A. Ishizaki and Y. Tanimura, “Quantum dynamics of system strongly coupled to low-temperaturecolored noise bath: Reduced hierarchy equations approach,” Journal of the Physical Society ofJapan , 3131–3134 (2005), https://doi.org/10.1143/JPSJ.74.3131. Y. Tanimura, “Stochastic liouville, langevin, fokker–planck, and master equation approachesto quantum dissipative systems,” Journal of the Physical Society of Japan , 082001 (2006),https://doi.org/10.1143/JPSJ.75.082001. Y. Tanimura, “Reduced hierarchical equations of motion in real and imaginary time: Corre-lated initial states and thermodynamic quantities,” The Journal of Chemical Physics , 044114(2014), https://doi.org/10.1063/1.4890441. Y. Tanimura, “Real-time and imaginary-time quantum hierarchal fokker-planck equations,” TheJournal of Chemical Physics , 144110 (2015), https://doi.org/10.1063/1.4916647. Y. Tanimura, “Numerically “exact” approach to open quantum dynamics: The hierarchi-cal equations of motion (heom),” The Journal of Chemical Physics , 020901 (2020),https://doi.org/10.1063/5.0011599. J. Zhang, R. Borrelli, and Y. Tanimura, “Proton tunneling in a two-dimensional potential energysurface with a non-linear system–bath interaction: Thermal suppression of reaction rate,” TheJournal of Chemical Physics , 214114 (2020), https://doi.org/10.1063/5.0010580. D. P. Tew, N. C. Handy, and S. Carter, “A reaction surface hamiltonian study of malonaldehyde,”The Journal of Chemical Physics , 084313 (2006), https://doi.org/10.1063/1.2338891. Y. Tanimura and Y. Maruyama, “Gaussian–markovian quantum fokker–planck approach tononlinear spectroscopy of a displaced morse potentials system: Dissociation, predissocia-tion, and optical stark effects,” The Journal of Chemical Physics , 1779–1793 (1997),https://doi.org/10.1063/1.474531. T. Ikeda and Y. Tanimura, “Probing photoisomerization processes by means of multi-dimensional electronic spectroscopy: The multi-state quantum hierarchical fokker-planck equation approach,” The Journal of Chemical Physics , 014102 (2017),https://doi.org/10.1063/1.4989537. T. Ikeda, A. G. Dijkstra, and Y. Tanimura, “Modeling and analyzing a photo-driven molecu-lar motor system: Ratchet dynamics and non-linear optical spectra,” The Journal of ChemicalPhysics , 114103 (2019), https://doi.org/10.1063/1.5086948. J. Hu, M. Luo, F. Jiang, R.-X. Xu, and Y. Yan, “Pad´e spectrum decompositions ofquantum distribution functions and optimal hierarchical equations of motion construc-26ion for quantum open systems,” The Journal of Chemical Physics , 244106 (2011),https://doi.org/10.1063/1.3602466. J.-J. Ding, R.-X. Xu, and Y. Yan, “Optimizing hierarchical equations of motion for quantumdissipation and quantifying quantum bath effects on quantum transfer mechanisms,” The Journalof Chemical Physics , 224103 (2012), https://doi.org/10.1063/1.4724193. Q. Shi, L. Chen, G. Nan, R.-X. Xu, and Y. Yan, “Efficient hierarchical liouville space propa-gator to quantum dissipative dynamics,” The Journal of Chemical Physics , 084105 (2009),https://doi.org/10.1063/1.3077918. T. Ikeda and Y. Tanimura, “Low-temperature quantum fokker–planck and smoluchowski equa-tions and their extension to multistate systems,” Journal of Chemical Theory and Computation , 2517–2534 (2019), pMID: 30776312, https://doi.org/10.1021/acs.jctc.8b01195. N. Sato and S. Iwata, “Promotion of the proton transfer reaction by the intermolecu-lar stretching mode: Application of the two-dimensional finite element method to thenuclear schr ¨odinger equation,” The Journal of Chemical Physics , 2932–2937 (1988),https://doi.org/10.1063/1.454998. S. Hammes-Schiffer and A. A. Stuchebrukhov, “Theory of coupled electron and pro-ton transfer reactions,” Chemical Reviews , 6939–6960 (2010), pMID: 21049940,https://doi.org/10.1021/cr1001436. S. Ueno and Y. Tanimura, “Modeling intermolecular and intramolecular modes of liquid waterusing multiple heat baths: Machine learning approach,” J. Chem. Theory Comput. , 2099–2108 (2020). S. Ueno and Y. Tanimura, “Modeling and simulating excited-state dynamics of a system incondensed phases: Machine learning approach,” e-print arxiv.org/abs/2102.02427 xxxx