Torsionally broken symmetry assists infrared excitation of biomimetic charge-coupled nuclear motions in the electronic ground state

The concerted interplay between reactive nuclear and electronic motions in molecules actuates chemistry. Here, we demonstrate that out-of-plane torsional deformation and vibrational excitation of stretching motions in the electronic ground state modulate the charge-density distribution in a donor-bridge-acceptor molecule in solution. The vibrationally-induced change, visualised by transient absorption spectroscopy with a mid-infrared pump and a visible probe, is mechanistically resolved by ab initio molecular dynamics simulations. Mapping the potential energy landscape attributes the observed charge-coupled coherent nuclear motions to the population of the initial segment of a double-bond isomerization channel, also seen in biological molecules. Our results illustrate the pivotal role of pre-twisted molecular geometries in enhancing the transfer of vibrational energy to specific molecular modes, prior to thermal redistribution. This motivates the search for synthetic strategies towards achieving potentially new infrared-mediated chemistry.


A2. Generating the vibrationally-excited ground-state population
The thermalised population in S0 is subsequently projected onto S1, while conserving positions and momenta. Following S1 relaxation for 15 fs, the population is projected back onto S0. This is because for linear π-delocalised systems like conjugated polyenes, it is known [9] that the bond-length alternation (BLA) stretching mode is activated (with a period of ~15 fs) during the S1 relaxation from the Franck-Condon region. Consequently, a snapshot of the S1 population after 15 fs, when projected onto S0, is expected to result in a vibrationally-excited S0 population with energy funneled into the BLA stretching mode. This population represents the initial conditions for studying the S0 vibrational dynamics for the following 2 ps (Fig. 4), which constitutes the very premise of our key observations.

B. Survey of the initial conditions
Here, we enumerate the various calculations performed to examine the role of BLA at every intermediate simulation step leading up to the initial conditions -namely, the creation of a vibrationallyexcited BLA-activated S0 population. Our first aim is to identify the S0 normal mode associated with the pronounced BLA motion localised at the central ethylenic bridge. This requires S0 optimisation of a selected geometry. The evolution of this optimised geometry may then be traced and visually inspected at every simulation step. Notably, the geometry optimisation erases all thermal information. This treatment is markedly different from the trajectory calculations, which include temperature-related effects originating from the creation of the initial thermal S0 population, as described above. However, a comparison between the evolution of the optimised geometry and the trajectories representing the entire population illustrates quantitative differences in a qualitatively unanimous and consistent picture. This justifies the methodology adopted here, which facilitates visual inspection of the evolved geometries and provides valuable insight into the induction of BLA character.
In order to select a geometry for the S0 optimisation, the geometry corresponding to the peak of the simulated absorption spectrum for the thermal S0 population (see Fig. S7A for further details) is considered, since it is representative of the average of the total population. Following the geometry optimisation at the 2-root-state-average CASSCF/6-31G* (including the effect of the solvent-shell configuration), a vibrational frequency analysis reveals an intense spectral peak localized at 1796 cm -1 (Fig. S8A). The corresponding normal mode is characterised by a distinct BLA motion of the skeletal structure, particularly in the region surrounding the central ethylenic bridge (Fig. S8C). Incorporating the correctional factor of 0.9 associated with the intrinsic overestimation of the computed frequencies due to the simulation framework employed, the frequency of this normal mode translates to 1616 cm -1 . This is close to the experimental value of 1577 cm -1 in the steady-state infrared absorption spectrum of DASPMI in dichloromethane (Fig. 1B), representing the S0 vibrational mode excited by the mid-infrared pump laser in the experiment.

S2
For the purpose of comparison, the selected geometry is also optimised on S1. The normal mode of interest, with conspicuous BLA character around the central ethylenic bridge (Fig. S8D), now corresponds to an intense spectral peak at 1173 cm -1 (Fig. S8B).
Our next aim is to demonstrate the validity of our hypothesis that the BLA stretching mode is activated during the S1 relaxation from the Franck-Condon region. To this end, ten randomly selected trajectories from the entire population of 400 trajectories are analysed and visually inspected, using the Visual Molecular Dynamics (VMD) software. The trajectories clearly depict temporal oscillations along the BLA coordinate (Fig. S9), thereby confirming that energy is released into the BLA mode. This behaviour is qualitatively reproduced (Fig. S9) by the Franck-Condon relaxation trajectory, arising from the vertical transition of the selected S0-optimised equilibrium geometry to the Franck-Condon point on S1. This facilitates a visual representation of the direction along which the Franck-Condon relaxation occurs (Fig.   S10A), which may immediately be compared with the normal-mode representation in S0 (Fig. S8C) or S1 The final aim is to characterise the initial conditions, particularly in terms of the energy injected into the BLA mode, which actuates the ensuing S0 vibrational dynamics (Fig. 4). Following projection back onto S0, the selected geometry is further relaxed for 20 fs, divided into 5-fs segments in order to compute the displacement vectors. The average displacement vector (Fig. S10B) is then compared with the S0 normal-mode vector (Fig. S8C). The relative orientation of the two vectors (computed by the scalar product) reveals a parallelism of 35%. This is indicative of the percentage of the total vibrational energy released into the BLA mode, where the total vibrational energy is the sum of the vibrational energies in S0 and S1 and is ~36 kcal/mol for this selected geometry. Considering the 35% parallelism, ~13 kcal/mol is therefore estimated to be released into the BLA mode. This is validated by a computation of the BLA energy in the harmonic approximation, which also yields a value of 13 kcal/mol. As discussed previously, this analysis is qualitative, since it does not incorporate any thermal effects. The vibrational energy injected into the experimentally relevant BLA mode is ~7 kcal/mol. This value is computed from the hessian at the turning points of the classical trajectories, and averaged over the entire population of 400 trajectories. The discrepancy with the experimental value (~5 kcal/mol) is within acceptable limits, given the simulation protocol generates classical trajectories, which do not produce quantised vibrational motion. Furthermore, the simulations consider a fully anharmonic model with 105 vibrational degrees of freedom for the DASPMI molecule and hundreds of degrees of freedom for the solvent molecules. On calculating the energy for the first vibrationally-excited level in S0 along this particular BLA mode in the harmonic approximation, it may be surmised that vibrational-ladder climbing up to the first vibrationallyexcited level in S0 may have been achieved in the course of the simulations. The zero-point-energy (ZPE) contributions associated with the trajectory vibrational modes are also evaluated. As documented in Table   S2, the maximum ZPE corresponding to the highest frequency vibrational mode is less than 30% of the vibrational energy initially injected into the experimentally relevant BLA mode. Besides, the lowest α vibrational component (that appears to be the one associated with the experimentally observed oscillations) features a small 0.1 kcal/mol ZPE.

C. Ground-state vibrational dynamics
Following the creation of the vibrationally-excited S0 population with the BLA stretching mode activated, the ensuing vibrational dynamics is simulated at the 2-root-state-average CASSCF level for 2 ps (Fig. 4). The time-evolution of the oscillator strength, f, the energy gap, E, and the charge difference, q, between S0 and S1 show distinct peaks at 1250 cm -1 and in the vicinity of 50 cm -1 in the Fourier domain ( Fig. 4). This may be correlated to the time-evolution of the C2-C3=C4-C5 dihedral, α, the H-C3=C4-H dihedral, β, the π-overlap, τ, and the BLA (Fig. 3).
The vibrationally-excited S0 population reveals that the ethylenic stretching incorporated in the activated BLA stretching mode can effectively couple with the ethylenic torsional mode, as it alters the ethylenic π-bonding character, and therefore the vibronic coupling. Thus, the coupling between BLA and α is expected in a situation where it is not prevented by symmetry. This is independently examined by calculating the sensitivity of the S0-S1 energy gap with respect to all the vibrationally populated modes resulting from the simulation (Table S1). The 1475 cm -1 mode (corresponding to the C-C-H bending of the ethylenic bridge), and the 1300 cm -1 mode (corresponding to the C-N stretching) have little influence on the π-overlap, and therefore the S0-S1 energy gap (6% and 14%, respectively). On the contrary, the experimentally relevant 1796 cm -1 mode (corresponding to the bridge-localised stretching BLA mode) changes the energy gap up to 28%. This assignment appears to be logical while considering the change in the overlap of the orbitals describing the π-bonding character. In contrast with the other modes, both stretching and torsion along the ethylenic bridge effectively change the π-overlap (Table S1).
Furthermore, the sensitivity of the S0-S1 energy gap is also examined with respect to: 1. the geometrical changes along the experimentally relevant 1796 cm -1 mode for a suitably selected geometry representing the thermalised population; 2. the geometrical changes induced by the simulation protocol, relaxing the same representative structure in S1.
The result shows that 1 and 2 generate the same energy gap sensitivity (the change in the S0-S1 energy gap being 46 and 49 kcal/mol, respectively). This supports the fact that our simulation protocol, involving projection onto S1 and subsequent relaxation onto S0, is equivalent to distorting a thermalised geometry in equilibrium along the experimentally relevant vibrational mode.

D. Potential-energy surface mapping
In order to locate the conical intersection, the Franck-Condon geometry, arising from the vertical transition of the selected S0-optimised equilibrium geometry, is gradually twisted in the course of several constrained optimisations. The torsional angle, α, representing the C2-C3=C4-C5 dihedral along the central ethylenic bridge, is decreased from 177 0 for the Franck-Condon geometry to 130 0 . This is followed by S4 structural relaxation until convergence to a geometry with near-zero energy difference between S0 and S1.
From the optimised conical-intersection geometry (Fig. S11A), the local vectors, g and h, spanning the branching plane around the conical intersection point are extracted and identified as BLA and α, respectively. These parameters serve as the coordinates along which the S0 (Fig. S12A) and S1 (Fig.   S12B) potential energy surfaces are mapped in the vicinity of the conical intersection point (Fig. S12D).
An analysis of the charge-transfer character is subsequently performed on S0 in order to distinguish between the diradical (or covalent) region and the charge-transfer region, separated by the conical intersection point (Fig. S12C). The S0-optimised geometry is characterised by a positive charge localised on the acceptor, which is diminished in the S0 charge-transfer region. The charge-transfer character, quantified by the partial charge, ρ, localised on the acceptor is attenuated from a value close to unity in the diradical region to a value close to zero in the charge-transfer region along the BLA coordinate (Fig. S12C). The corresponding transition-state geometries, TSdir (Fig. S11B) and TSCT (Fig. 6C), are determined initially (as a first approximation) by constrained optimisations, fixing the parameter α and varying the BLA, followed by subsequent TS optimisation (without constraints) to find the final TS geometry. Notably, TSCT provides a lower energy barrier (by approximately 10 kcal/mol) compared to TSdir.
To establish the minimum-energy path, two geometries with Δα = ±5 0 , compared to α = 91 0 for TSCT, are selected as starting points to administer paths towards the reactant trans (Fig. 6A) and product cis (Fig. 6D) isomer states. The minimum-energy path along with the corresponding variation in ρ may then be traced in the parameter space defined by BLA and α ( Fig. S13 and Fig. 5B). The substantial change in the Mulliken charge distribution is also shown (Fig. S14).  (Fig. 1B). This may be attributed to the effect of the counter-ion (iodide) present in the experimental sample solution, which is not included in the simulation model and may cause a red-shifted absorption, in addition to the intrinsic error introduced by the CASPT2//CASSCF protocol (usually ~30-nm blue-shift) [10]. The corresponding oscillator strength, averaged over the entire GSth population, is 0.86 (denoted by the dashed line). For GSBLA, the theoretical absorption maximum for the S0-S1 energy gap is red-shifted to 572 nm (50 kcal/mol) and the average oscillator strength increases to 1.11. The increase in the average S0-S1 oscillator strength following the creation of the vibrationally-excited, BLA-activated S0 population may be correlated to the appearance of the transient absorption signal observed experimentally ( Fig. 2) with the electronic probe following the vibrational pump. The normal modes depicted correspond to frequencies of (C) 1796 cm -1 and (D) 1173 cm -1 , respectively. S13 Fig. S9. Excited-state BLA oscillations. S1 relaxation dynamics for ten randomly selected trajectories from the entire population of 400 trajectories. This S1 population is created by the projection of the thermal S0 population onto S1, conserving positions and momenta. The temporal oscillations along the BLA coordinate, defined as the difference between the average single and double bond-lengths, clearly indicates that energy is released into the BLA mode. The average trajectory (red) shows a maximum at ~9 fs, approximately a quarter of the time-period of the oscillation. For comparison, the Franck-Condon (FC) relaxation trajectory, arising from the vertical transition of the selected S0-optimised equilibrium geometry to the Franck-Condon point on S1, is also shown and exhibits a qualitatively similar behaviour.  Table S1: Sensitivity of the S0-S1 energy gap to the modes (Fig. S8A) excited by the simulation protocol