Symmetry breaking fluctuations split the porphyrin $Q$ bands

Porphyrins are ubiquitous in biology and chemistry and play crucial roles in, for example, oxygen transport and storage in blood, charge transfer in cellular respiration, and photosynthesis in plantsMERCER-SMITH and Mauzerall (1984); Amanullah et al. (2022). Their synthetic tunabilityTanaka and Osuka (2015) make porphyrins ideal in advancing bioinspired charge and energy transport and storage technologiesMathew et al. (2014); Swierk et al. (2015); Zhang et al. (2023); Kumagai et al. (2014), especially as they offer a convenient chemical scaffold that can be functionalized with substituentsNing et al. (2020); Roy et al. (2022a), paired with transition metalsMoisanu et al. (2023), or wreathed together in polymer chainsGotfredsen et al. (2022); Kosugi et al. (2023); Magnera et al. (2023); Chen et al. (2024) to tailor and exploit their optoelectronic properties. Numerous studies have turned to sophisticated multidimensional spectroscopies to disentangle the palette of spectral features generated by porphyrin and porphyrin-like systemsNiedringhaus et al. (2018); Policht, Niedringhaus, and Ogilvie (2018); Song et al. (2019); Arsenault et al. (2018); Jana et al. (2023); Zhao et al. (2023); Petropoulos et al. (2024). Yet, predicting and interpreting even the linear electronic spectra of small porphyrin monomers remains challenging as their photophysics continues to elude a conclusive microscopic interpretation.

In photosynthesis, energy transfer between carotenoids and porphyrin-based chlorophyll $Q$ bands facilitate photoprotectionMa et al. (2003); Holt et al. (2005); Ruban et al. (2007); Polívka and Frank (2010); Bode et al. (2009); Son et al. (2020). These $Q$ bands in free-base porphyrins present four peaks corresponding to two distinct electronic transitions, in contrast to their metalloporphyrin counterparts which present only two peaks. The widely accepted group-theoretic description of porphyrin spectraGouterman, Wagnière, and Snyder (1963) interprets these transitions as having mutually orthogonal transition dipoles in the plane of the tetrapyrrole ring: the lower-energy half of the $Q$ bands ($\rm S_{0}\to S_{1}$) denoted as the $Q_{x}$ band and the higher-energy half ($\rm S_{0}\to S_{2}$) denoted as the $Q_{y}$ band, with transition dipoles aligned parallel and orthogonal to the central proton axis, respectively. Yet, even the pioneering work of GoutermanGouterman (1961); Gouterman, Wagnière, and Snyder (1963), which explained the competition of intensities and center frequencies of the Soret and $Q$ bands in metalloporphyrins, does not explain why the $Q$ band peaks split from two to four upon removal of the metal center. Crucially, the microscopic mechanism underlying this spectral splitting continues to spark disagreementGouterman (1961); Roy et al. (2022b); Yu, Baskin, and Zewail (2002); Schalk et al. (2008); He et al. (2012); Kullmann et al. (2012). Indeed, this doubling (four peaks instead of two) of the $Q$ bands has been variously ascribed to a Condon-allowed vibronic progressionGouterman (1961); Yu, Baskin, and Zewail (2002), where the transition dipole is assumed to be insensitive to nuclear motions, or non-Condon effects via symmetry-breaking nuclear motions that tune the brightness, or selection rules, of transitionsRoy et al. (2022b); Neville and Schuurman (2024). Given this controversy, it is not surprising that this critical problem in porphyrin photophysics has been the subject of intense theoretical focus.

Yet, state-of-the-art theoretical spectroscopy tools can only tackle isolated aspects of this question. For example, one can capture non-Condon effects under the microscopic harmonic approximation using the Franck-Condon Herzberg-Teller (FCHT) methodSantoro et al. (2008) but only in the gas phase and for fairly small and rigid molecules Zuehlsdorff et al. (2019); address non-Condon fluctuations with atomistic anharmonicities in the condensed phase but only in the static limit (Ensemble)Bergsma et al. (1984); or treat the dynamics and atomistic anharmonicities in the condensed phase but only in the Condon limit (2nd order cumulant)Mukamel (1985); Zuehlsdorff et al. (2019). Thus, a more general approach is required.

This impasse has led to seemingly contradictory explanations of the physical origin of the $Q$ bands. One of the dominant interpretationsYu, Baskin, and Zewail (2002) assigns the first two peaks (1.9-2.2 eV) of the tetraphenylporphyrin (TPP) $Q$ band (see, for example, the experimental spectrum of TPP in $\rm CS_{2}$ in Fig. 1) as the 0-0 and 0-1 features in the vibronic progression of the S₁ ($Q_{x}$) electronic transition, and the latter two (2.2-2.5 eV) as the 0-0 and 0-1 vibronic progression of the S₂ ($Q_{y}$) transition. However, recent studies have used FCHT-type arguments to assign the $Q$ bands of porphine to mixed Condon-allowed (first and third peaks) and non-Condon-enabled electronic transitions (second and fourth)Santoro et al. (2008); Roy et al. (2022b). Yet, these interpretations do not provide a fully predictive theory of the lineshape that unambiguously captures the splitting and broadening of spectral peaks. Hence, it is crucial to employ a theoretical framework that can capture both dynamical vibronic progressions and the non-Condon fluctuations that tune transition brightness as molecules vibrate and jostle in solution.

We bridge this theoretical gap by leveraging our recent theory of non-Condon fluctuations that offers a quantum dynamical means to process the output of classical molecular dynamics (MD) simulations and electronic structure theory to predict linear optical spectra in the condensed phase Wiethorn et al. (2023). Our theory predicts the optical response of solvated chromophores while faithfully incorporating their anharmonic atomistic interactions, arbitrarily strong non-Condon effects (nuclear motion-induced variations of the transition dipoles), and significant environmental coupling. Our expressions encode the chemical details that determine optical spectra through spectral densities that quantify how nuclear motions modulate the energy gap, the transition dipole, and their correlations. Crucially, one can directly calculate these spectral densities from atomistic simulations of a material, offering a direct route to identifying and interpreting the optical response. This ultimately allows us to resolve the long-standing controversy over the physical origin of the $Q$ band splitting of porphyrins by quantifying the influence of the following processes:

1. 1.

   Vibronic features on Condon-allowed transitions.

2. 2.

   Static and dynamic non-Condon fluctuations.

3. 3.

   Energy transfer between $Q_{x}$ and $Q_{y}$ bands.

Our work is the first to accurately reproduce and explain the spectral splitting and intensity distribution of the $Q$ bands of representative porphyrins: porphine, TPP, and tetraphenylporpholactone (TPPL). Our analysis reveals that the $Q$-band splitting arises from extreme timescale separation between the slow torsions of the pendant phenyls (in TPP and TPPL) and collective solvent motions and the fast tetrapyrrole ring vibrations that modulate the transition dipole, making this splitting a primarily non-Condon feature. We further introduce a perturbative treatment of excited-state energy transfer to confirm that this splitting is robust to energy flow across the $Q$ bands.

Our theory of non-Condon spectroscopy exploits the idea that the macroscopic number of nuclear motions that modulate the energy gap $U(\hat{\mathbf{q}})$ and transition dipole $\bm{\mu}(\hat{\mathbf{q}})$ in the condensed phase are likely to cause Gaussian fluctuations. This enables us to write an effective Gaussian model for which we have derived exact closed-form expressions for the optical response function Wiethorn et al. (2023) that can be parameterized directly from MD simulations (see Methods Section LABEL:Methods_Sec:Theory-GNCT-nature). In contrast to previous works, our Gaussian non-Condon theory (GNCT) does not resort to fitting parameters that restrict its predictive power and can lead to ill-converging spectra when strongly anharmonic nuclear motions render harmonic approximations inapplicable (see FCHT prediction in Fig. 5(c). In addition, it enables us to isolate contributions to the absorption spectrum, $\sigma(\omega)$, arising from transition dipole fluctuations that are static and dynamic relative to electronic dephasing times,

Here, $\alpha(\omega)$ is the linear absorption coefficientMukamel (1995) centered at $\omega=0$ where as $\sigma^{\bm{\mu}{\rm-stat.}}_{\rm GNCT}(\omega)$ and $\sigma^{\mu{\rm-dyn.}}_{\rm GNCT}(\omega)$ are shifted by the average transition energy $\langle\omega_{eg}\rangle$ between the ground ($g$) and excited ($e$) state. Below, we describe the physical meaning of the two components in Eq. (1) within the pure-dephasing limit (i.e., in the absence of energy transfer), and discuss analogous expressions we derive and employ to capture the influence of excited-state energy flow on the spectrum in Sec. IV of the supplementary information (SI).

Our GNCT predicts that the average transition dipole, $\langle\hat{\bm{\mu}}_{ge}\rangle$, determines the magnitude of the static non-Condon contribution to the spectrum,

Here, $\langle\hat{\bm{\mu}}_{ge}\rangle$ arises from the average effect of nuclear motions that break the symmetry of a molecule. When $\sigma^{\mu{\rm-stat.}}_{\rm GNCT}(\omega)$ dominates the spectrum, it captures non-Condon disorder while accounting for the energy gap dynamics, which encode vibronic progressions. When $\langle\hat{\bm{\mu}}_{ge}\rangle$ is insensitive to nuclear configurations, it can be replaced by the transition dipole at any nuclear configuration, recovering the celebrated second-order cumulant within the Condon approximation, $|\hat{\bm{\mu}}_{ge}|^{2}\sigma_{\delta U\delta U}(\omega)$Zuehlsdorff et al. (2019); Mukamel (1985) (see Methods A for the full expression). The energy gap spectral density, $C_{\delta U\delta U}(\omega)=\left[1+\coth(\omega/k_{\rm B}T)\right]C_{\delta U\delta U}^{\prime\prime}(\omega)$Mukamel (1995), encodes the dynamical energy gap fluctuations ($\delta U(\hat{\mathbf{q}})$) that control the lineshape of $\sigma_{\delta U\delta U}(\omega)$ and its ability to exhibit distinct vibronic progressions or its related asymmetry.

Our GNCT quantifies how dynamic non-Condon fluctuations impact the spectrum via $s_{\delta\bm{\mu}\delta\bm{\mu}}(\omega)$ and $\mathbf{s}_{\delta U\delta\bm{\mu}}(\omega)$,

Just as $C_{\delta U\delta U}(\omega)$ encodes how energy gap fluctuations control the prominence of vibronic progressions in $\sigma_{\delta U\delta U}(\omega)$, here, two new spectral densities, $C_{\delta\bm{\mu}\delta\bm{\mu}}(\omega)$ and $\mathbf{C}_{\delta U\delta\bm{\mu}}(\omega)$, encode how nuclear motions cause dynamic transition dipole fluctuations ($\delta\bm{\mu}(\hat{\mathbf{q}};t)$) and determine $s_{\delta\bm{\mu}\delta\bm{\mu}}(\omega)$ and $\mathbf{s}_{\delta U\delta\bm{\mu}}(\omega)$. $C_{\delta\bm{\mu}\delta\bm{\mu}}(\omega)$ quantifies how nuclear motions modulate the brightness of transitions while $\mathbf{C}_{\delta U\delta\bm{\mu}}(\omega)$ measures how nuclear motions jointly tune the energy gap and transition brightness. Together, Eqs. (2) and (3) determine the optical spectrum arising from static and dynamic transition dipole modulations. Our GNCT thus reveals that dynamic non-Condon effects arise from the convolution (denoted by $[f\ast g](\omega)$) of the energy gap lineshape, $\sigma_{\delta U\delta U}(\omega)$, with the contributions arising from the transition dipole fluctuations.

Because our theory provides a transparent protocol to parameterize Eqs. (2) and (3) using atomistic simulations, we can assign spectral signals to specific nuclear motions and assess their effect on both energy gap and transition dipole fluctuations. We now show how our GNCT, unlike previous theories, recapitulates the condensed phase spectra of porphyrins and uncovers the mechanism leading to the splitting of their $Q$ bands as a function of their chemical structure.

We first ask if the observed peak splitting arises purely from Condon-allowed vibronic progressions. The second-order cumulant method in the Condon limit (to which we refer henceforth as “GCT” for Gaussian Condon theory) successfully captures vibronic progressionsMukamel (1985); Loco et al. (2018); Zuehlsdorff et al. (2019), and is analogous to $\sigma_{\rm GNCT}^{\mu{\rm-stat.}}(\omega)$ in our GNCT (Eq. (2)) in the Condon limit (i.e., $\langle\hat{\bm{\mu}}_{ge}\rangle\to\hat{\bm{\mu}}_{ge}$). Figure 1 shows that the predicted GCT spectrum (teal) only recovers the first ($Q_{x}$) and third ($Q_{y}$) peaks, revealing that the $Q$-band splitting does not arise from a Condon-allowed vibronic progression and non-Condon effects must be considered.

Since Condon-allowed transitions do not cause the peak splitting, we turn to static non-Condon effects. The ensemble methodBergsma et al. (1984), widely used to predict condensed phase spectraZuehlsdorff et al. (2021), treats non-Condon fluctuations under the approximation that nuclear motions are slower than electronic dephasing times, thereby capturing the static limit of both non-Condon and energy gap fluctuations. However, similar to the GCT, Fig. 1 shows that the ensemble method (purple) can capture only the first ($Q_{x}$) and third ($Q_{y}$) peaks. Hence, the splitting does not arise solely from static non-Condon effects.

To incorporate both static and dynamic non-Condon fluctuations, we turn to our GNCT. Figure 1 illustrates that our GNCT (dark green) is the only method that correctly reproduces the four peaks of the TPP $Q$ bands, their relative intensities, and even their linewidths. Thus, peak splitting in the TPP $Q$ bands appears to arise from dynamic non-Condon fluctuations.

We are now poised to interrogate the microscopic mechanism leading to $Q$ band splitting. First, Fig. 2(a.1) shows that excluding $\mathbf{s}_{\delta U\delta\bm{\mu}}(\omega)$ (dotted) causes only negligible changes to the intensity distribution of the spectrum (see also SI Fig. S1), demonstrating that the joint modulation of transition dipoles and energy gaps by nuclear motions does not participate in splitting the TPP $Q$ bands. Further, Fig. 2(a.2) shows that zeroing the purely dynamical non-Condon contribution, $\sigma_{\rm GNCT}^{\delta\bm{\mu}\delta\bm{\mu}}(\omega)~{}\equiv~{}[s_{\delta\bm{\mu}\delta\bm{\mu}}\ast\sigma_{\delta U\delta U}](\omega)$, causes the second and fourth peaks in the TPP spectrum to disappear (pink-dashed). Since this is tantamount to neglecting the full dynamical non-Condon contribution, $\sigma^{\mu{\rm-dyn.}}_{\rm GNCT}(\omega)$, one can assign the second ($Q_{x}$) and fourth ($Q_{y}$) peaks in TPP to purely dynamical non-Condon effects. In contrast, zeroing the static non-Condon component, $\sigma^{\mu{\rm-stat.}}_{\rm GNCT}(\omega)$, decimates the first ($Q_{x}$) and third ($Q_{y}$) peaks (blue-dashed), revealing these as largely static non-Condon contributions. This suggests that peak splitting in TPP arises from a static-dynamic divide in the nuclear motion-induced modulation of the selection rules. However, the persistence of small components in the first and third peaks upon neglecting $\sigma^{\mu{\rm-stat.}}_{\rm GNCT}(\omega)$ motivates a more microscopic interrogation of $\sigma^{\delta\bm{\mu}\delta\bm{\mu}}_{\rm GNCT}(\omega)$.

To uncover the origin of the persistent splitting in the surviving dynamical non-Condon contribution to the GNCT lineshape, $\sigma_{\rm GNCT}^{\delta\bm{\mu}\delta\bm{\mu}}(\omega)$, we analyze the transition dipole spectral density that determines its behavior, $C_{\delta\bm{\mu}\delta\bm{\mu}}(\omega)$ (see Fig. 1(b)). $C^{Q_{x}}_{\delta\bm{\mu}\delta\bm{\mu}}(\omega)$ and $C^{Q_{y}}_{\delta\bm{\mu}\delta\bm{\mu}}(\omega)$ display complex structure with a small contribution at low frequency ($<200$ cm^-1) and a region of considerable spectral weight centered around $\sim 1500$ $\rm cm^{-1}$. Restricted normal mode analyses in the gas phase enable us to assign the high-frequency region of the spectral density (region 2 in Fig. 1(b)) to inner-ring C-C and C-N vibrations (see SI Sec. II). The low-frequency contribution below $<200$ cm^-1 arises from collective solvent motions in solution and torsional motions of the phenyl substituents with timescales $\sim 33$ cm^-1 (see SI Sec. III). But to what extent does the structure of these spectral densities matter and how does this apparent timescale separation in nuclear motions contribute to the splitting of the $Q$-bands?

To address the first question, we replace the MD-derived transition dipole spectral densities ($C_{\delta\bm{\mu}\delta\bm{\mu}}(\omega)$) with model counterparts composed of low- and high-frequency contributions, $C^{\rm model}_{\delta\bm{\mu}\delta\bm{\mu}}(\omega)=C^{\rm model}_{\delta\bm{\mu}\delta\bm{\mu}}(\omega_{\rm lf})+C^{\rm model}_{\delta\bm{\mu}\delta\bm{\mu}}(\omega_{\rm hf})$. For $C^{\rm model}_{\delta\bm{\mu}\delta\bm{\mu}}(\omega_{\rm lf})$, we employ an often adopted Debye form,Weiss (2012) and subsume the high-frequency intramolecular vibrations into a Gaussian envelope for $C^{\rm model}_{\delta\bm{\mu}\delta\bm{\mu}}(\omega_{\rm lf})$ (see Methods Sec. LABEL:Methods_Sec:Model-system-spectral-densities-GNCT-nature). Figure 2(b.1)-(b.2) illustrate the MD-based (line) and our simplified (shaded) $C_{\delta\bm{\mu}\delta\bm{\mu}}(\omega)$ for the $Q_{x}$ and $Q_{y}$ bands, with their corresponding absorption spectra in (b.3)-(b.4). Our simplified spectral densities faithfully reproduce the splitting of the $Q$-bands (dashed-peach), agreeing semi-quantitatively with the MD-based spectrum (solid), indicating that the detailed structure of the MD-derived spectral density is not essential to recover the TPP $Q$ band splitting and validates our use of these models in our investigation of timescale separation.

We now assess the influence of timescale-separated fluctuations in $\sigma^{\rm model}_{\delta\bm{\mu}\delta\bm{\mu}}(\omega)$. By turning off the high-frequency Gaussian envelopes in $C^{\rm model}_{\delta\bm{\mu}\delta\bm{\mu}}(\omega)$, one observes that the second ($Q_{x}$) and fourth ($Q_{y}$) peaks vanish in Fig. 2(b.3) and (b.4), leaving only the small contributions to the first ($Q_{x}$) and third ($Q_{y}$) peaks (green-dashed), which we had largely assigned as static non-Condon contributions. What is more, the energy difference between split $Q$ band peaks (0.19 eV) in Fig. 2(b.3)-(b.4) matches the energetic split between the low- versus high-frequency components of the simplified non-Condon spectral densities in Fig. 2(b.1)-(b.2). A consistent picture emerges: the first ($Q_{x}$) and third ($Q_{y}$) peaks arise from static ($\langle\hat{\bm{\mu}}_{ge}\rangle$) and slow (collective solvent motions and phenyl torsions) non-Condon contributions, while second ($Q_{x}$) and fourth ($Q_{y}$) peaks arise from the fast modulation of the transition dipoles arising from intra-tetrapyrrole ring vibrations, revealing that the extreme timescale separation of the nuclear motions that modulate the transition dipole, and thus the optical selection rules, offers the mechanism that splits the TPP $Q$ bands.

While our GNCT addresses the pure dephasing limit where no energy transfer occurs between excited states, TPP is known to undergo internal conversion as energy flows from the $Q_{y}$ to the $Q_{x}$ bandYu, Baskin, and Zewail (2002). Hence, we developed a framework to quantify the impact of $Q_{y}\rightarrow Q_{x}$ energy flow on the spectrum (see SI Sec. IV) by treating the driver of energy flow as a perturbation. We tailor an excited state transformationMedders et al. (2017) to move into the frame where states have spatially orthogonal transition dipoles, consistent with the traditional interpretation of porphyrin photophysics (see SI Sec. IVB). We find that including energy flow across the $Q$ bands has a marginal effect on the computed spectrum (Fig. S10), indicating that our timescale separation mechanism for $Q$ band splitting remains valid even in the presence of internal conversion.

Besides offering mechanistic insights, our GNCT yields design principles to tune optoelectronic properties. Specifically, the spectral splitting mechanism in porphyrins suggests that one can balance the contributions of $\sigma^{\mu\rm-stat.}_{\rm GNCT}(\omega)$ and $\sigma^{\mu\rm-dyn.}_{\rm GNCT}(\omega)$ via simple chemical changes. For instance, we hypothesize that one can suppress the first and third peaks of the porphyrin $Q$ bands by removing the meso-phenyls or replacing them with bulky chemical moieties. In contrast, asymmetrically modifying the tetrapyrrole ring by hydrating one of the pyrroles to form a chlorin or substituting a functional group into the ring can be expected to: (1) change the intensity distribution across the $Q$ band by modulating $|\langle\bm{\hat{\mu}}_{ge}\rangle|^{2}$ via steric or electron conjugation effects that tune the rotation of the pendant phenyls; and (2) modify the individual $Q_{x}$/$Q_{y}$ lineshapes by tuning the energy gap spectral density $C_{\delta U\delta U}(\omega)$ which determines the prominence of vibronic progressions in $\sigma_{\delta U\delta U}(\omega)$. To assess the validity of our design hypotheses, we turn to porphine and TPPL. These molecules offer an ideal testing ground since replacing the pendant phenyls in TPP with hydrogens yields porphine, whereas TPPL retains TPP’s phenyls but replaces an external ($\beta$-periphery) $\rm C_{2}H_{2}$ on one pyrrole orthogonal ($y$) to the proton axis ($x$) with a lactone.

We begin by leveraging the transition dipole spectral densities, $C_{\delta\bm{\mu}\delta\bm{\mu}}(\omega)$, to analyze the effect of the chemical modifications required to transform TPP to porphine and TPPL. We employ the non-Condon factorWiethorn et al. (2023), $\varphi^{\mu}=\frac{C^{\prime\prime}_{\delta\bm{\mu}\delta\bm{\mu}}(\omega)-|\langle\hat{\bm{\mu}}_{ge}\rangle|^{2}}{C^{\prime\prime}_{\delta\bm{\mu}\delta\bm{\mu}}(\omega)+|\langle\hat{\bm{\mu}}_{ge}\rangle|^{2}}$, to quantify how these changes tune the balance between static and dynamic non-Condon effects, with an increase (decrease) in $\varphi^{\mu}$ indicating that dynamical (static) non-Condon effects become more significant. Figure 3 shows $C^{Q_{x}}_{\delta\bm{\mu}\delta\bm{\mu}}(\omega)$ and $C^{Q_{y}}_{\delta\bm{\mu}\delta\bm{\mu}}(\omega)$ for porphine (a) and TPPL (b). In porphine, the absence of phenyl torsions that drive low-frequency transition dipole fluctuations in TPP leads to spectral densities with no low-frequency weight, and an increase in $\varphi_{Q_{x}}^{\mu}$ and $\varphi_{Q_{y}}^{\mu}$ by factors of 1.51 and 1.65, respectively, indicating that the optoelectronic properties of porphine depend more strongly on dynamic non-Condon fluctuations. In TPPL, these spectral densities exhibit significant low-frequency spectral weight from the phenyl torsions, albeit with greater intensity than in TPP. In contrast to TPP, one also observes enhanced intensity centering at 1000 $\rm cm^{-1}$ emerging from dynamic modulations of the transition dipoles by the C-O stretch. Moreover, the enhanced low-frequency motions lead to a reduction in $\varphi_{Q_{x}}^{\mu}$ and $\varphi_{Q_{y}}^{\mu}$ by factors of 2.42 and 2.79, respectively, indicating that static non-Condon effects are increasingly significant to the optoelectronic properties of TPPL as compared to TPP.

We now turn to the optical spectra in Fig. 3(c) and (d). In porphine and TPPL (as in TPP), setting $\mathbf{s}_{\delta U\delta\bm{\mu}}(\omega)=0$ (dotted) reveals that $s_{\delta\bm{\mu}\delta\bm{\mu}}(\omega)$ captures the dynamic contribution to the $Q$ band splitting (see also SI Sec. I). Hence, we analyze the resulting spectra through the lens of static and dynamic non-Condon contributions. Figure 3(c) shows the porphine spectrum and confirms our hypothesis that removing the phenyl groups reduces both the low-frequency dynamical and static contributions arising from non-Condon fluctuations, leading to suppressed first and third peaks in the porphyrin $Q$ bands. TPPL is more complex. Figure 3(d) reveals that the TPPL $Q$ bands contain three features that are absent in porphine and TPP: (1) spectral intensity is greatest at the third peak rather than the fourth; (2) a vibronic shoulder now emerges from the fourth peak and; (3) the energetic separation of $\sim 919.5$ cm^-1 between second (dynamic-$Q_{x}$) and third (static-$Q_{y}$) peaks deviates from the former $\sim 1500$ cm^-1 splitting between subsequent peaks. Two of these changes align with our hypotheses: the shift in spectral intensity across the $Q$ bands from the lactone-induced change to the static-dynamic balance of non-Condon effects, and the larger vibronic progression from the lactone-induced modification of the high-frequency component of the energy gap spectral density that shifts peak prominence. However, the modified splitting between the second and third peaks is surprising, motivating a closer inspection of our hypotheses.

Hypothesis (1) centers on the balance between $\sigma^{\mu-\rm stat.}_{\rm GNCT}(\omega)$ and $\sigma^{\mu-\rm dyn.}_{\rm GNCT}(\omega)$. Figure 4 decomposes our GNCT prediction for the $Q_{x}$ and $Q_{y}$ bands of porphine, TPP, and TPPL, illustrating the relative contributions of $\sigma^{\mu-\rm stat.}_{\rm GNCT}(\omega)$ (red and orange) and $\sigma^{\delta\bm{\mu}\delta\bm{\mu}}_{\rm GNCT}(\omega)$ (pink and peach) to the individual bands (dotted) as a function of the non-Condon factors. This figure shows that the $\sigma^{\mu-\rm stat.}_{\rm GNCT}(\omega)$ contribution is negligible in porphine, becomes significant in TPP, and reaches commensurate values with the dynamic contribution in TPPL, thus aligning with our hypothesis that the phenyl groups continue to emphasize the first and third peaks relative to the second and fourth. Figure 4 further reveals the true complexity of the TPPL spectrum as all four contributions overlap significantly, with each displaying a prominent lactone-induced vibronic progression. This is in stark contrast to the four unimodal peaks that lack energy gap-derived vibronic structure underlying the $Q$ bands of porphine and TPP. This observation also confirms hypothesis (2) as it reveals that the asymmetric lactone functionalization of the inner ring induces vibronic progressions that can redistribute spectral weight.

Yet, it seems that a contradiction lingers in our analysis. Figure 3(d) shows that the third peak in the experimental TPPL spectrum carries greater intensity than the fourth peak, which one may naively ascribe to a greater static non-Condon contribution in the $Q_{y}$ band based on our current analysis. However, Fig. 4(c) shows the opposite is true: the dynamical non-Condon $Q_{y}$ peak is predicted to be greater. This seeming contradiction is a result of considering the static and dynamic non-Condon peaks separately. Fortuitously, the peak separation in the vibronic progression of the dynamic non-Condon $Q_{x}$ feature (light pink solid line) is nearly identical to the peak separation between dynamic $Q_{x}$ ($\sim 2.1$ eV) and static $Q_{y}$ ($\sim 2.2$ eV) in the full spectrum. As a result, the overlapping dynamic $Q_{x}$ satellite peak (vibronic progression) supplies spectral weight to the static non-Condon $Q_{y}$ peak at $\sim 2.2$ eV, causing it to appear to dominate over that of the dynamic $Q_{y}$ at $\sim 2.42$ eV in the total spectrum. This analysis thus demonstrates how the interplay between Condon-allowed vibronic progressions and symmetry-forbidden transitions control the TPPL $Q$-band intensity distribution. Yet, despite our GNCT’s success in predicting the relative peak heights, it appears to overestimate the prominence of a lactone-induced vibronic progression for the static $Q_{x}$ peak that is absent in the experimental spectrum (consider the $\sim 0.17$ eV gap between first and second peaks in Fig. 3(d)). Since energy gap spectral densities encode vibronic progressions, an analysis of $C^{Q_{x}}_{\delta U\delta U}(\omega)$ ought to reveal why our GNCT predicts this spurious effect in the TPPL $Q_{x}$ band.

Figure 5(a) compares $C^{Q_{x}}_{\delta U\delta U}(\omega)$ for porphine, TPP, and TPPL (see SI Figs. S2(a), S3(a), and S4(a) for $C^{Q_{y}}_{\delta U\delta U}(\omega)$ of porphine, TPP, and TPPL). This figure shows that, while $C^{Q_{x}}_{\delta U\delta U}(\omega)$ for porphine and TPP are similar, the spectral density of TPPL shows a greater weight, especially in the region between 600–1400 cm^-1, which roughly corresponds to internal-ring C-O stretching around 1000 cm^-1. Because the lactone group contains a carbonyl, one might expect that the vibration of this C=O should modulate the energy gap, leading to a peak near $1800-1900$ $\rm cm^{-1}$. Yet, this peak is conspicuously absent from the TPPL spectral density, while porphine and TPP display greater spectral weight in this region. This is suspect and can indicate a violation of Gaussian statistics. Indeed, employing the Kullback-Leibler divergence between each MD-derived energy gap distribution and their corresponding normal distributions of identical mean and variance (denoted as $D_{\rm KL}^{{\rm MD}|\mathcal{N}}$) to quantify the extent of non-Gaussian behavior, we find that the TPPL $Q_{x}$ energy gap is more non-Gaussian than both its $Q_{y}$ counterpart and the energy gaps for porphine and TPP (see SI Sec. V and Figure S9). This suggests the need for a closer interrogation of the $Q_{x}$ energy gap.

Since the only change between TPP and TPPL is the lactone substitution, which largely affects high-frequency intramolecular vibrations, one can interrogate the source of non-Gaussian behavior via a microscopic harmonic analysis of the nuclear motions that modulate the energy gap for an isolated propholactone ring (see Sec. SVI of the SI). This analysis reveals that the normal modes on the ground electronic surface mix to yield those on the $Q_{x}$ excited state. Figure 5(b) illustrates this by visualizing the Duschinsky rotation, $\mathbf{J}$,Duschinsky (1937) of the normal modes upon photoexcitation. The extent to which the Duschinsky matrix displays non-diagonal entries quantifies the amount of mode mixing, which causes non-Gaussian behavior and masks high-frequency vibrational coupling as lower-frequency weight in the spectral density.Zuehlsdorff et al. (2019) Although the Duschinsky matrix is diagonally dominant, indicating that mode mixing is an overall weak effect (see SI Fig. S11), Fig. 5(b) shows a region of the Duschinsky matrix for the $Q_{x}$ band of porpholactone that reveals non-negligible mode mixing between ground state (horizontal-axis) and excited state (vertical-axis) normal modes with frequencies $>900$ cm^-1.

Does the spurious vibronic progression predicted by our GNCT for the TPPL $Q_{x}$ band arise from mode mixing-induced non-Gaussianity? To answer this question, we correct the energy gap contribution in our GNCT, $\sigma_{\delta U\delta U}(\omega)$ in Eqs. (2) and (3), with an energy gap contribution that incorporates mode mixing (see SI Sec. VI). Specifically, we calculate the exact response function for the microscopically harmonic model of the porpholactone ring and add the effect of the collective solvent motions and phenyl torsions in solvated TPPL by extracting the broadening one would obtain from fitting a Debye contribution to the low-frequency portion of our atomistically derived spectral density. Figure 5(c) shows the resulting spectrum ($\mathbf{J}^{Q_{x}}_{\rm PL}$-GNCT) compared to the predicted TPPL $Q$ bands from our GNCT result and the FCHT which yields an unphysical negative absorption feature at $\sim 1.8$ eV. Here, we confirm that by accounting for microscopically harmonic mode mixing corrections in the TPPL $Q_{x}$ band, our GNCT (purple line) now accurately captures the splitting of both $Q_{x}$ and $Q_{y}$ bands of TPPL, resolving all four peaks and agreeing well with the experimental lineshape. Thus, our results for porphine and TPPL confirm our hypothesis that our theory’s mechanistic insights enable us to tailor and optimize the optical properties of a molecular system by manipulating nuclear motions that induce statistical transition dipole fluctuations with characteristic timescales.

In summary, we have leveraged our GNCT to successfully simulate, reproduce, and analyze, for the first time, the condensed phase $Q$-band spectra of representative porphyrins: TPP, its parent porphine without phenyl substituents, and the lactone-substituted TPPL. The $Q$ band region plays a key role in the fundamental energy and charge transfer steps in photosynthesis Amanullah et al. (2022); Ma et al. (2003); Holt et al. (2005); Ruban et al. (2007); Polívka and Frank (2010); Bode et al. (2009); Son et al. (2020), highlighting the importance of our analysis, which enabled us to settle the long-standing debate over the physical origin of the splitting of the $Q$ bands observed for metal-free porphyrins. We showed that peak-splitting arises from extreme timescale separation in the nuclear motions that modulate selection rules: low-frequency phenyl torsions and high-frequency inner-ring vibrations. The insights we uncovered can be expected to shed light on the dynamical mechanisms at play in metal-substituted porphyrins that often appear in natural and artificial light-harvesting complexes and photo- and electro-catalysis. Such timescale-separated transition dipole fluctuations are also likely crucial in understanding the linewidths and relative intensities between the $Q_{x}$ and $Q_{y}$ transitions in the optical spectra of metallo-porphyrins even if the respective bands do not appear to split.

We further demonstrated that our theory can reveal design principles to tune porphyrin optoelectronic behavior. To this end, we exploited the atomistic detail of our GNCT to selectively suppress slow transition dipole fluctuations and tune fast fluctuations of the transition dipoles and energy gap via chemical modifications. By simulating the $Q$ bands of porphine, we proved that severing the pendant phenyls on TPP almost entirely eliminates the low-frequency nuclear motions leading to the systematic suppression of the first and third peaks in the $Q$ band region, in agreement with our theory-led hypothesis. In turn, we also confirmed that substituting the most distant carbons in one of the pyrrole groups in the tetrapyrrole ring indeed shifted the peak height ratios of the second and fourth peaks by tuning the fast transition dipole fluctuations and caused a vibronic progression-induced asymmetry of all peaks, especially those corresponding to the $Q_{y}$ transition, through the enhancement of energy gap modulations from fast intra-ring motions.

Our results constitute an exciting breakthrough that opens the door to elucidating the physical mechanisms that produce the versatile optoelectronic properties of the large porphyrin family and beyond. With the promising properties that porphyrins offer as chemically tunable scaffolds to guide energy and charge transfer for energy conversion, catalysis, and even electrochemistry, the ability to employ theory to design chemical structures that exhibit desired functionality is essential. By offering a simple and reliable framework to simulate and interpret the spectral features of chromophores that transcend the standard approximations (e.g., harmonic motions, lack of appreciable non-Condon fluctuations, lack of disorder or thermal effects), our theory sets the stage to provide novel insight, disentangle the interplay between thermally accessible molecular motions and spectral lineshapes, and establish robust design principles for next-generation bioinspired charge and energy transfer technologies in terms of tuning material optoelectronic properties as a function of composite chemical structure.