# On the Origin of the Shift Between Vertical Excitation and Band Maximum in Molecular Photoabsorption

Shuming Bai, Ritam Mansour, Ljiljana Stojanović, Josene Toldo, Mario Barbatti

► To cite this version:

Shuming Bai, Ritam Mansour, Ljiljana Stojanović, Josene Toldo, Mario Barbatti. On the Origin of the Shift Between Vertical Excitation and Band Maximum in Molecular Photoabsorption. Journal of Molecular Modeling, 2020, 26 (5), 10.1007/s00894-020-04355-y. http://dx.doi.org/10.1007/s00894-020-04355-y

HAL Id: hal-02612358

https://amu.hal.science/hal-02612358

Submitted on 19 May 2020

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# On the Origin of the Shift Between Vertical Excitation and Band Maximum in Molecular Photoabsorption

Shuming Bai,<sup>1,2,d,*</sup> Ritam Mansour,<sup>1,a</sup> Ljiljana Stojanović,<sup>1,3,c</sup> Josene M. Toldo,<sup>1,b,*</sup> Mario Barbatti<sup>1,e,*</sup>

<sup>1</sup> Aix Marseille University, CNRS, ICR, Marseille, France.

<sup>2</sup> Department of Chemistry, Duke University, Durham, NC 27708.

<sup>3</sup> School of Biological and Chemical Sciences, Queen Mary University of London, Mile End Road, London, United Kingdom, E1 4NS.

<sup>a</sup>ORCID: 0000-0001-8591-2866

<sup>b</sup>ORCID: 0000-0002-8969-6635

<sup>c</sup>ORCID: 0000-0003-2821-5110

<sup>d</sup>ORCID: 0000-0002-4993-8070

<sup>e</sup>ORCID: 0000-0001-9336-6607

* Corresponding authors:

SB: [shuming.bai@duke.edu](mailto:shuming.bai@duke.edu)

JMT: [josene-maria.toldo@univ-amu.fr](mailto:josene-maria.toldo@univ-amu.fr)

MB: [mario.barbatti@univ-amu.fr](mailto:mario.barbatti@univ-amu.fr); [www.barbatti.org](http://www.barbatti.org)

## Abstract

The analysis of the photoabsorption spectra of molecules shows that the band maximum is usually redshifted in comparison to the vertical excitation. We conducted a throughout analysis of this shift based on low-dimensional analytical and numerical model systems, showing that its origin is rooted in the frequency changing between the ground and the excited states in multidimensional systems. Moreover, we deliver a benchmark of ab initio results for the shift based on a comparison of vertical excitations and band maxima calculated with the nuclear ensemble approach for the 28 organic molecules in the Mülheim molecular dataset. The mean value of the shift calculated over 60 transitions is 0.11 ± 0.08 eV. The mean value of the band width is 0.32 ± 0.14 eV.

## Keywords

Spectrum simulation, Vertical excitations, Excited states, Absorption cross-section
**Introduction**

In the computational analysis of UV/Vis photoabsorption spectra of molecules, several authors have noticed that the band maximum $E_i^{\max}$ and the vertical excitation $E_i^v$ of band $i$, computed at the same electronic structure level, are shifted by up to 0.2 eV, with the band maximum systematically redshifted [1-5]. This shift $\delta_i$ (Figure 1) is a parameter that plays a particularly important role in the comparison between theory and experiments and in the evaluation of the precision of computational methods [4]. Therefore, a direct experimental/computational comparison without taking the shifts into account may lead to systematic errors.

<figure><img src="image_2.png" /><figcaption>Figure 1. Schematic illustration of the absorption band parameters.</figcaption></figure>

In the last years, diverse methods for first-principle spectrum simulations have been proposed [6-9,5,10], leading to high-accurate predictions of band shapes, including their vibrational resolution [11]. Through simulations of vibrationally-resolved absorption spectra of several molecules, Avila Ferrer et al. [3] confirmed the existence of the redshift and suggested to compare vertical excitations not to the band maximum, but to the first moment $M^1$ (center of gravity) of the absorption band. Nevertheless, to the best of our knowledge, the origin of the redshift has not been clarified so far. We know from our experience with spectrum simulations based on the nuclear ensemble approach [12] that the redshift also occurs when the band maximum is simulated with incoherent transition processes. Therefore, its origin should have an underlying semiclassical mechanism, not directly related to quantum effects.

Using a simple analytical model, we show here that when a multidimensional ground-state vibrational wavefunction is projected onto the excited state, the frequency change between
the excited and the ground states causes the maximum spectral intensity to be displaced to lower energies compared to the vertical excitation. We confirmed such dependence of the shift on the frequency changing using low-dimensionality numerical models. Finally, we deliver a benchmark of shift values for the first bands of all 28 organic molecules in the Mülheim molecular dataset [13] (Figure 2). This benchmark was produced comparing ab initio vertical excitations to band maxima simulated with the nuclear ensemble approach at the same electronic structure level.

<figure><img src="image_3.png" /><figcaption>Figure 2. Structures of the 28 molecules in the Mülheim dataset.</figcaption></figure>

## Computational details

The ground state optimization and normal mode analysis of all 28 molecules were done with the resolution-of-the-identity coupled-cluster to approximated second-order (RI-CC2) [14,15] with the cc-pVTZ basis set. The Cartesian coordinates of the ground-state minimum geometries are given in the Supporting Information. Excited states were calculated at RI-CC2 level with the aug-cc-pVDZ basis set [16]. In a few cases with convergence problems, we
alternatively employed the algebraic diagrammatic construction to second-order (ADC(2))
[17].

For the nuclear ensemble simulations [12], the nuclear configurations were sampled
according to a harmonic-oscillator Wigner distribution based on normal modes in the ground
state minimum. An ensemble of 500 nuclear geometries was generated for each molecule of
the dataset. For each nuclear configuration, up to 15 electronic states were computed using RI-
CC2/aug-cc-pVDZ. Some of the molecules presented recurrent convergence problems with
CC2, caused by near-degenerated excited states. In such cases, we employed RI-ADC(2)/aug-
cc-pVDZ in the spectrum simulation. Each spectrum was plotted with a Gaussian linewidth
broadening of 0.05 eV. Details on vertical excitations and simulated spectra are given in the
Supporting Information.

RI-CC2 and RI-ADC(2) calculations were done with TURBOMOLE [18]. Nuclear-
ensemble spectra were done with NEWTON-X [19,20] interfaced with TURBOMOLE. A
program to calculate spectra based on Eq. (13) is freely available for download [21].

**Results and discussion**

Hypotheses to explain the band maximum shift

To understand the shift $\delta_i$ between the vertical excitation and the band maximum, it is useful to analyze it in terms of the nuclear ensemble approach [12]. In this method, the absorption spectrum is built as an incoherent sum of transitions from nuclear geometries representing the vibrational distribution of the electronic ground state. Thus, the nuclear ensemble approach delivers a band envelope without vibrational resolution but with enough information to estimate the shift.

To determine the origin of the shift, we considered three possible causes: (1) the
oscillator strength variation in the sampling (a post-Condon effect); (2) the displacement of the
vibrational mode in the excited state; and (3) the vibrational frequency change between the
ground and excited states. The penetration of the wave function in regions beyond the classical
return point may also cause a shift, but it tends to be negligible for large molecules [22].

The post-Condon effect (1) was discarded after we noticed that nuclear ensemble
simulations using a uniform value for the oscillator strength still yields the shift. We show in
the next sections that the normal mode displacement (2) also does not cause the shift, although it impacts its value. Finally, we found out that the frequency change between the ground and excited states (3) is at the root of the phenomenon. The underlying reason for the shift between the vertical excitation and the band maximum is that for a multidimensional potential energy surface, although most of the nuclear geometries are distributed around the equilibrium geometry, there are even more geometries displaced from equilibrium satisfying the resonance condition $\Delta E_{ge} = E_p$ ($\Delta E_{ge}$ is the energy gap between the ground and the excited states, and $E_p$ is the photon energy).

## Analytical analysis of the shift origin

To understand the origin of the shift, consider the following simple two-modes/two-states analytical model. The two states are separated by a vertical transition energy $E^v$. For each state, the potential energy surface is supposed to be given by two uncoupled harmonic and degenerated potentials in the coordinates $q_1$ and $q_2$, with the ground and excited states sharing the same equilibrium geometry, but with different angular frequencies $\omega_g$ and $\omega_e$. Without losing generality, the energy of the ground-state equilibrium geometry is assumed to be zero. Thus, the ground- and excited-state potential energies are given as [23]

$$ 
\begin{aligned}
E_g &= \frac{1}{2} \mu \omega_g^2 (q_1^2 + q_2^2), \\
E_e &= E^v + \frac{1}{2} \mu \omega_e^2 (q_1^2 + q_2^2).
\end{aligned} 
\quad (1) $$

where $\mu$ is the reduced mass. The resonance condition is

$$ \Delta E_{ge} - E_p = (E^v - E_p) + \frac{1}{2} \mu (\omega_e^2 - \omega_g^2) (q_1^2 + q_2^2) = 0. \quad (2) $$

In a semiclassical approximation, the spectrum is given by

$$ \sigma = k \int \rho(q_1, q_2) g(\Delta E_{eg}(q_1, q_2) - E_p) dq_1 dq_2 \quad (3) $$

where $\rho$ is the distribution of geometries in the ground state, and $g$ is a line-shape function satisfying the resonance condition in Eq. (2). $k$ is a constant, which includes the oscillator strength between the two states (Condon approximation). If $\rho$ is given by a Wigner distribution for the quantum harmonic oscillator,
$$ 
\begin{aligned}
\rho(q_1, q_2) &= \left(\frac{\alpha}{\pi}\right) \exp(-\alpha(q_1^2 + q_2^2)), \\
\alpha &= \frac{\mu\omega_g}{\hbar},
\end{aligned} 
\qquad (4) $$

and the line-shape function is taken as a normalized Gaussian function with arbitrary width w,

$$ g(x) = \frac{1}{w\sqrt{\pi}} \exp\left(-\frac{x^2}{w^2}\right), \qquad (5) $$

we obtain

$$ \sigma = \frac{k\alpha}{w\pi^{3/2}} \int \exp\left(-\alpha(q_1^2 + q_2^2) - \frac{1}{w^2}\left(E^v - E_p\right) + \frac{1}{2}\mu(\omega_e^2 - \omega_g^2)(q_1^2 + q_2^2)\right)^2 dq_1 dq_2. \quad (6) $$

w is assumed to be much narrower than the band width, and in the limit when it tends to zero, g converges to a Dirac delta function.

In polar coordinates $(q_1 = R \cos(\theta), q_2 = R \sin(\theta))$, Eq. (6) becomes

$$ \sigma = \frac{k\alpha}{w\pi^{3/2}} \int \exp\left(-\alpha R^2 - \frac{1}{w^2}\left(E^v - E_p\right) + \frac{1}{2}\mu(\omega_e^2 - \omega_g^2)R^2\right)^2 RdRd\theta, \quad (7) $$

which can be analytically integrated to give

$$ \sigma = -\frac{k}{A} \exp\left(\frac{w^2}{A^2}\left(1 + \frac{2A}{w^2}(E^v - E_p)\right)\right) \left(1 + \operatorname{erf}\left(\frac{w}{A}\left(1 + \frac{A}{w^2}(E^v - E_p)\right)\right)\right), \quad (8) $$

where

$$ A = \frac{\hbar(\omega_e^2 - \omega_g^2)}{\omega_g}. \qquad (9) $$

The maximum of this spectral band happens for the photon energy $E^{\max}$ satisfying $d\sigma / dE_p = 0$. If we expand this derivative to the first order in $E_p$ around $E^v$, we obtain the shift

$$ \delta \equiv E^v - E^{\max} = \frac{\hbar(\omega_g^2 - \omega_e^2)}{2\omega_g}. \qquad (10) $$
Because upon excitation the bond is loosened, it implies that $\omega_e < \omega_g$ and, therefore, $E^v > E^{\max}$. This means that the band maximum is redshifted in relation to the vertical excitation.

<figure><img src="image_4.png" /><figcaption>Figure 3. Illustration of the 2-modes/2-states analytical model. The upper surface shows the potential energy difference $\Delta E_{ge}$ while the bottom surface illustrates the Wigner distribution $\rho$ in the ground state. The dashed curves indicate isoenergetic curves satisfying the resonance condition $\Delta E_{ge} = E_p$.</figcaption></figure>

With this simple model, it is easy to understand the origin of the shift. For a certain photon energy $E_p$, the spectral intensity in Eq. (3) is determined by both the ground state population on the isoenergetic line satisfying the resonance condition and the length of this line. In this model, this line is merely a circle, but it is a hyperline in general multidimensional cases. In one limiting case, the ground state population is maximal ($q_1 = q_2 = 0$) but the line length is too short, yielding a low spectral intensity. In the case of excited and ground states sharing the same equilibrium geometry, this limit happens at $E^v$. In the other limiting case, the curve length tends to infinite ($q_1 = q_2 \to \infty$), but the ground state population is null, also yielding a minimal intensity. Between the two extremes, the band maximum occurs at $E^{\max}$, with optimal length and population values, which is proportional to the difference between the squared frequencies of the excited and ground states (Eq. (10)).

Although there are no analytical results for multidimensional, realistic cases, this dependence of the shift on the difference of the squared frequencies seems to be general. For instance, as mentioned in the Introduction, Avila Ferrer et al. [3] proposed to compare vertical excitations to the first-moment $M^1$ of the absorption band to approximately account for the shift effect. They recall that the $M^1$ can be analytically computed for harmonic ground and
excited states in the Franck-Condon approximation, with the excited state built in the vertical Hessian approach. In this model, the shift between the vertical excitation and M¹ is [24]

$$E^v - M^1 = \sum_r \frac{\hbar (\omega_{g,r}^2 - \omega_{e,r}^2)}{4 \omega_{g,r}} \coth \left( \frac{\hbar \omega_{g,r}^2}{2k_B T} \right), \quad (11)$$

where $k_B$ is the Boltzmann constant, T is the temperature, and r runs over the normal modes. In the limit of low temperature and two degenerated normal modes, this expression yields the shift given by Eq. (10).

Note yet that the shift seems to be a multidimensional phenomenon. We have not been able to find an analytical expression for the shift in a single-dimensional model, but numerical estimates showed that in such a case the shift is a function of the width w, and it tends to zero for small w. In the limit of a Dirac delta function (w = 0), the band of the single-dimensional model is peaked at the vertical excitation [23].

## Analysis of the shift with model systems

To get further numerical insights into the origin of the shift, we simulated the absorption spectrum of a model system with few vibrational modes using the nuclear ensemble approach. The advantage of working with such low-dimensionality systems is that we can systematically control the values of key parameters to verify their impact on the spectrum.

The model systems are composed of two harmonic potentials representing the ground and the first excited states. The simplest model contains a single vibrational mode (1V) with ground-state frequency $\omega_g$ and excited-state frequency given by $\omega_e^2 / \omega_g^2 = \eta$, where $\eta$ is a constant. The potential energies of the ground and excited states are

$$E_g = \frac{1}{2} \mu \omega_g^2 q^2 \\ E_e = E^v + \frac{1}{2} \mu \eta \omega_g^2 (q-d)^2 \quad (12)$$

In the simulations, $E^v = 4.08$ eV (0.15 au), $\mu = 1.0$ amu (1823 au) and $\omega_g = 1500$ cm⁻¹ (0.0068 au). To compute the spectrum, 50,000 geometries were sampled from a Wigner distribution for $E_g$, and the oscillator strength was assumed to be 1.0 for all points.
As expected for the one-dimensional model, the results do not show any shift for small $w$.

We extended the 1V model into three vibrational modes (3V), with ground-state frequencies $\omega_{g1} = 1800$, $\omega_{g2} = 1500$, and $\omega_{g3} = 1000$ cm$^{-1}$; frequency-changing parameters $\eta_1$, $\eta_2$, and $\eta_3$; and displacements $d_1$, $d_2$, and $d_3$. All other parameters were kept at the same values as in the 1V model.

The results for the shift with the 3V model are given in Table 1. Now, we see that displacements without a frequency change lead to a negligible shift. For $\eta_2 = 0.8$, as the displacements of the other two modes increase, the shift increases, too. This result tells that although the displacement does not cause the shift, it modulates the shift value.

<table><caption>**Table 1.** Numerical simulations of the shift (eV) between the vertical excitation and the band maximum for a model with three vibrational modes, considering a shift between the ground and excited-state frequencies (given by $\eta_i$) and a displacement $d_i$ (au) between the modes. All results with $\eta_1 = \eta_3 = 1.0$ and $d_2 = 10$ a.u.</caption><thead><tr><th>η<sub>2</sub></th><th>d<sub>1</sub></th><th>d<sub>3</sub></th><th>Shift /eV</th></tr></thead><tbody><tr><td>1.0</td><td>3.0</td><td>3.0</td><td>-0.01</td></tr><tr><td>0.8</td><td>0.0</td><td>0.0</td><td>0.08</td></tr><tr><td>0.8</td><td>3.0</td><td>3.0</td><td>0.15</td></tr><tr><td>0.8</td><td>8.0</td><td>8.0</td><td>0.27</td></tr><tr><td>0.8</td><td>10.0</td><td>10.0</td><td>0.31</td></tr><tr><td>0.8</td><td>15.0</td><td>15.0</td><td>0.39</td></tr></tbody></table>

## Benchmark of shift values

Based on the vertical excitations and the results from the nuclear ensemble simulations, we investigated the shift between the vertical transition energy $E_{calc}^v$ and the band maximum $E_{calc}^{max}$ in the 28 molecules of the Mülheim dataset. Figure 4 illustrates these simulations for three molecules—cyclopropene, pyridine, and pyrimidine (molecules 5, 13, and 15 in Figure 2).
<figure><img src="image_5.png" /><figcaption>Figure 4. Comparison between the nuclear ensemble spectrum simulation (solid black) and the vertical excitations (sticks) for (top) cyclopropene, (middle) pyridine, and (bottom) pyrimidine. The fitting of the nuclear ensemble result is shown as a dashed red curve. The Gaussian function components of the fitting, corresponding to the sub-bands of each vertical excitation, are shown as well.</figcaption></figure>

For energies in eV and absorption cross section in Å².molecule⁻¹, the spectrum convolution can be written as (see derivation in Appendix A)

$$ \sigma_A(E) = 0.619n \sum_i \frac{f_i}{\Gamma_i} e^{-(E-E_i^v+\delta_i)^2/\Gamma_i^2}, \quad (13) $$

where n is the refractive index, and $E_i^v$, $f_i$, $\delta_i$, and $\Gamma_i$ are, respectively, the vertical excitation, oscillator strength, shift, and width of band i. We used this expression to fit the simulated spectra. In the fitting process, $E_i^v$ and $f_i$ were kept fixed, and we got $\delta_i$, and $\Gamma_i$. In some cases, the oscillator strength was also optimized to improve the fitting (see Supporting Information). The number of Gaussian functions used to fit each spectrum was determined by the number of vertical excitations with considerable oscillator strength in the spectral region of interest. Thus, for some bands, several Gaussian functions (up to five) were needed. This is the case of cyclopropene shown in Figure 4-top (which required three Gaussian functions) and of
pyridine in Figure 4-middle (two Gaussian functions). Other bands, like those in pyrimidine (Figure 4-bottom), required only one Gaussian function.

The shift results for all molecules are collected in Table 2. All fitted parameters are given in the Supporting Information. As shown in this table, for all 28 molecules we observe sub-bands with a redshift between the band maximum and the vertical excitation. The shift varies from zero to 0.4 eV; its mean value is 0.11 eV, and the standard deviation is 0.08 eV. As expected, the shift value does not show any correlation either with $E_i^v$ or with $f_i$.

<table><caption>**Table 2.** Vertical excitation energy ($E_i^v$ in eV) and oscillator strengths ($f_i$) for all molecules 28 molecules in the Mülheim dataset in the fitted spectral region. The peak of each Gaussian sub-band is redshifted by $\delta_i$ (eV) from $E_i^v$. The width of the sub-band is $\Gamma_i$ (eV). Molecules are numbered as in Figure 2.</caption><thead><tr><th>Molecule</th><th>Method</th><th>E<sub>i</sub><sup>v</sup> (eV)</th><th>f<sub>i</sub></th><th>δ<sub>i</sub> (eV)</th><th>Γ<sub>i</sub> (eV)</th></tr></thead><tbody><tr><td rowspan="2">1</td><td rowspan="2">CC2</td><td>7.16</td><td>0.082</td><td>0.07</td><td>0.192</td></tr><tr><td>7.91</td><td>0.381</td><td>0.14</td><td>0.605</td></tr><tr><td>2</td><td>CC2</td><td>6.14</td><td>0.745</td><td>0.13</td><td>0.403</td></tr><tr><td>3</td><td>CC2</td><td>5.18</td><td>1.265</td><td>0.10</td><td>0.400</td></tr><tr><td>4</td><td>CC2</td><td>4.53</td><td>1.778</td><td>0.06</td><td>0.329</td></tr><tr><td rowspan="3">5</td><td rowspan="3">CC2</td><td>6.75</td><td>0.075</td><td>0.02</td><td>0.771</td></tr><tr><td>7.09</td><td>0.031</td><td>0.19</td><td>0.294</td></tr><tr><td>7.40</td><td>0.058</td><td>0.14</td><td>0.287</td></tr><tr><td>6</td><td>CC2</td><td>5.49</td><td>0.105</td><td>0.09</td><td>0.460</td></tr><tr><td rowspan="2">7</td><td rowspan="2">CC2</td><td>5.64</td><td>0.018</td><td>0.07</td><td>0.296</td></tr><tr><td>6.09</td><td>0.018</td><td>0.18</td><td>0.356</td></tr><tr><td rowspan="2">8</td><td rowspan="2">ADC(2)</td><td>6.89</td><td>0.078</td><td>0.06</td><td>0.059</td></tr><tr><td>7.19</td><td>0.770</td><td>0.15</td><td>0.155</td></tr><tr><td>9</td><td>CC2</td><td>4.79</td><td>0.080</td><td>0.13</td><td>0.360</td></tr><tr><td rowspan="2">10</td><td rowspan="2">CC2</td><td>6.39</td><td>0.185</td><td>0.10</td><td>0.328</td></tr><tr><td>6.43</td><td>0.042</td><td>0.00</td><td>0.749</td></tr><tr><td rowspan="3">11</td><td rowspan="3">CC2</td><td>5.77</td><td>0.013</td><td>0.15</td><td>0.400</td></tr><tr><td>5.86</td><td>0.027</td><td>0.21</td><td>0.134</td></tr><tr><td>6.28</td><td>0.189</td><td>0.21</td><td>0.309</td></tr><tr><td rowspan="3">12</td><td rowspan="3">CC2</td><td>6.30</td><td>0.032</td><td>0.03</td><td>0.032</td></tr><tr><td>6.35</td><td>0.027</td><td>0.00</td><td>0.027</td></tr><tr><td>6.44</td><td>0.153</td><td>0.14</td><td>0.153</td></tr><tr><td>13</td><td>CC2</td><td>5.04</td><td>0.005</td><td>0.27</td><td>0.216</td></tr></tbody></table>
<table><tbody><tr><td></td><td></td><td></td><td>5.32</td><td>0.032</td><td>0.05</td><td>0.278</td></tr><tr><td>14</td><td>CC2</td><td></td><td>4.18</td><td>0.006</td><td>0.14</td><td>0.213</td></tr><tr><td></td><td></td><td></td><td>5.14</td><td>0.085</td><td>0.05</td><td>0.229</td></tr><tr><td>15</td><td>CC2</td><td></td><td>4.44</td><td>0.006</td><td>0.07</td><td>0.395</td></tr><tr><td></td><td></td><td></td><td>5.52</td><td>0.030</td><td>0.11</td><td>0.262</td></tr><tr><td></td><td></td><td></td><td>6.11</td><td>0.004</td><td>0.15</td><td>0.323</td></tr><tr><td>16</td><td>CC2</td><td></td><td>3.83</td><td>0.005</td><td>0.08</td><td>0.286</td></tr><tr><td></td><td></td><td></td><td>5.39</td><td>0.018</td><td>0.06</td><td>0.359</td></tr><tr><td>17</td><td>CC2</td><td></td><td>4.80</td><td>0.016</td><td>0.06</td><td>0.422</td></tr><tr><td>18</td><td>CC2</td><td></td><td>5.24</td><td>0.056</td><td>0.09</td><td>0.364</td></tr><tr><td>19</td><td>CC2</td><td></td><td>6.38</td><td>0.022</td><td>0.05</td><td>0.128</td></tr><tr><td>20</td><td>CC2</td><td></td><td>5.74</td><td>0.031</td><td>0.08</td><td>0.226</td></tr><tr><td>21</td><td>ADC(2)</td><td></td><td>5.33</td><td>0.593</td><td>0.14</td><td>0.393</td></tr><tr><td>22</td><td>CC2</td><td></td><td>6.13</td><td>0.028</td><td>0.05</td><td>0.349</td></tr><tr><td></td><td></td><td></td><td>6.55</td><td>0.021</td><td>0.17</td><td>0.215</td></tr><tr><td></td><td></td><td></td><td>6.71</td><td>0.064</td><td>0.01</td><td>0.251</td></tr><tr><td></td><td></td><td></td><td>7.23</td><td>0.029</td><td>0.10</td><td>0.241</td></tr><tr><td></td><td></td><td></td><td>7.53</td><td>0.337</td><td>0.05</td><td>0.500</td></tr><tr><td>23</td><td>CC2</td><td></td><td>5.77</td><td>0.032</td><td>0.00</td><td>0.452</td></tr><tr><td></td><td></td><td></td><td>6.36</td><td>0.019</td><td>0.03</td><td>0.239</td></tr><tr><td></td><td></td><td></td><td>6.88</td><td>0.019</td><td>0.19</td><td>0.241</td></tr><tr><td></td><td></td><td></td><td>7.26</td><td>0.206</td><td>0.09</td><td>0.421</td></tr><tr><td></td><td></td><td></td><td>7.60</td><td>0.056</td><td>0.00</td><td>0.240</td></tr><tr><td>24</td><td>CC2</td><td></td><td>5.77</td><td>0.026</td><td>0.02</td><td>0.399</td></tr><tr><td></td><td></td><td></td><td>6.36</td><td>0.020</td><td>0.08</td><td>0.247</td></tr><tr><td></td><td></td><td></td><td>6.76</td><td>0.014</td><td>0.00</td><td>0.259</td></tr><tr><td></td><td></td><td></td><td>7.19</td><td>0.189</td><td>0.16</td><td>0.394</td></tr><tr><td></td><td></td><td></td><td>7.48</td><td>0.027</td><td>0.00</td><td>0.192</td></tr><tr><td>25</td><td>CC2</td><td></td><td>4.68</td><td>0.050</td><td>0.08</td><td>0.412</td></tr><tr><td></td><td></td><td></td><td>5.57</td><td>0.137</td><td>0.18</td><td>0.391</td></tr><tr><td>26</td><td>ADC(2)</td><td></td><td>5.10</td><td>0.205</td><td>0.18</td><td>0.418</td></tr><tr><td></td><td></td><td></td><td>6.10</td><td>0.028</td><td>0.41</td><td>0.275</td></tr><tr><td></td><td></td><td></td><td>6.19</td><td>0.048</td><td>0.20</td><td>0.254</td></tr><tr><td>27</td><td>CC2</td><td></td><td>5.34</td><td>0.182</td><td>0.17</td><td>0.344</td></tr><tr><td></td><td></td><td></td><td>6.30</td><td>0.068</td><td>0.23</td><td>0.430</td></tr><tr><td>28</td><td>ADC(2)</td><td></td><td>5.13</td><td>0.247</td><td>0.17</td><td>0.357</td></tr><tr><td></td><td></td><td></td><td>5.19</td><td>0.059</td><td>0.00</td><td>0.569</td></tr></tbody></table>

The current results confirm that the redshifts are a common feature shared by many molecules. 80% of the sub-bands presented shifts equal or superior to 0.05 eV. Therefore, it
seems that the standard practice of comparing computed vertical excitations to the experimental band maxima should be reconsidered. Ideally, we should compare the experimental band maximum to the calculated band maximum. If only vertical excitations are known (as it is usually the case), a simple estimate of the band maximum is given by $E_{calc}^{max} = E_i^v - \bar{\delta}$, where $\bar{\delta} = 0.1$ eV is the mean value of our benchmark. This point is further elaborated in the next section, where we discuss a specific example.

The fitting of the simulated spectra yielded the sub-band width $\Gamma_i$ as a side product. Their values are also given in Table 2. The minimum width is 0.03 eV and the maximum 0.77 eV. The mean value is 0.32 eV, and the standard deviation is 0.14 eV.

## Accounting for the shifts in spectrum convolutions

The parameters in the simple spectrum convolution defined in Eq. (13) are illustrated in Figure 1. This equation is often used to make simple spectrum simulations from vertical excitation energies and oscillator strengths calculated at the equilibrium geometry. In such applications, the width is either taken as an arbitrary value common to all bands (based on our benchmark, we recommend to employ $\Gamma = 0.3$ eV) or fitted to experimental values.

Note that the area under each sub-band in Eq. (13) is proportional to the oscillator strength of the corresponding vertical excitation. Nevertheless, it is still common in the literature to find spectrum convolutions taking the oscillator strength as the height of the sub-band. Such a procedure, however, should be avoided, as it does not bear the right functional form $f_i / \Gamma_i$ and delivers spectra in arbitrary units.

Figure 5 shows the absorption cross-section of pyrimidine in the gas phase. Compared to the experimental spectrum from Ref. [25], the simulations deliver satisfactory results. The intensity and position of the two bands are semi-quantitatively described. Neither the spectrum convolution nor the nuclear ensemble can describe the vibrational resolution of the bands. Nevertheless, in the nuclear ensemble, the band envelope is a result of the simulation, while in the simple convolution, it is assumed to be a sum of Gaussian functions with an arbitrary width $\Gamma$. Naturally, there is a large difference in the computational effort invested in the two methods. While the Gaussian convolution only demands the calculation of vertical excitation energies and oscillator strengths for the ground-state equilibrium geometry, the nuclear ensemble needs hundreds of such calculations; in the present case, 500 single points.
<figure><img src="image_6.png" /><figcaption>Figure 5. Absorption cross section for pyrimidine in the gas phase. Spectrum convolutions based on Eq. (13) with the same width value for all sub-bands (0.3 eV). Nuclear ensemble spectra computed with 500 points. Experimental data from Ref. [25].</figcaption></figure>

If now we turn our attention to the spectrum convolutions and the nuclear ensemble
results in Figure 5, we can clearly see the effect of the shift. With $\delta=0$, the convolution is
blueshifted compared to the nuclear ensemble band. A significantly better agreement is
obtained when a shift $\delta=0.1$ eV is uniformly applied to both bands.

Experimentally, the high-energy band is peaked at 5.2 eV [25]. The vertical excitation
energy computed with CC2 is at 5.52 eV, 0.32 eV above the experimental value. However, as
we have been discussing, we should preferentially compare the estimated band maximum to
the experimental result. The band maximum computed with the nuclear ensemble, still with
CC2, is at 5.43 eV, therefore, 0.23 eV above the experimental value. A simpler estimate of the
band maximum is obtained by computing $E_i - \bar{\delta}$, where $\bar{\delta} = 0.1$ eV, the mean shift value
obtained from our benchmark. In this case, the band maximum estimate is 5.42 eV, 0.22 eV
above the experimental value.

**Conclusion**

It is well established in the literature that band maxima of absorption spectra are
redshifted in comparison to the corresponding vertical excitations by values reaching up to $\delta$
0.2 eV (assuming that both are computed with the same electronic-structure method). Using an
analytical model and numerical estimates from low-dimensionality model systems, we
determined that these shifts are primarily caused by the frequency change between the ground
and the excited states in multidimensional systems. The geometry displacement between the ground- and excited-state minima does not cause the shift but may affect its value. Nuclear ensemble simulations for a dataset of 28 small organic molecules completed this analysis by delivering a benchmark of shift values. The mean value of the shift calculated over 60 transitions is 0.11 eV with 0.08 eV standard deviation. The mean value of the band width is 0.32 eV with 0.14 eV standard deviation.

Despite the availability of highly accurate methods for spectrum simulations, it is likely that the workhorse method for routine estimation of the one-photon absorption of molecules will remain the simple spectrum convolution of vertical excitations. Based on the results presented in this paper, we recommend that in such convolutions the $\delta$ shift should be systematically taken into consideration in the comparison between the experimental band maximum and computed vertical excitations. This can be done, for instance, following the suggestion of Ref. [3] of comparing the vertical excitation to the $M^1$ moment of the experimental band. Alternatively, it can be simply done by estimating the band maximum by redshifting the vertical excitation by 0.1 eV, the mean value of our benchmark.

## Appendix A. Derivation of Eq. (13)

Experimentally, the oscillator strength of an absorption band $i$ is given as [26]

$$f_i = \frac{10^3 \ln(10) m_e c}{\pi N_A e^2 n} \int \varepsilon_i(\nu) d\nu, \qquad (14)$$

where $m_e$ is the electron mass, $c$ is the speed of the light, $N_A$ is the Avogadro's number, $e$ is the electron charge, and $n$ is the refractive index (CGS units). $\varepsilon_i$ is the molar absorbance (or extinction coefficient). The integral is done over the absorption frequency $\nu$. Supposing that the band has a Gaussian shape with maximum $\varepsilon_i^{\max}$ at $\nu_i^{\max}$, and width $\gamma_i$, the integral is

$$\int \varepsilon_i(\nu) d\nu = \int \varepsilon_i^{\max} e^{-(\nu-\nu_i^{\max})^2 / l^2} d\nu = \sqrt{\pi} \varepsilon_i^{\max} \gamma_i. \qquad (15)$$

Thus, the oscillator strength is

$$f_i = \frac{10^3 \ln(10) m_e c}{\sqrt{\pi N_A e^2 n}} \varepsilon_i^{\max} \gamma_i. \qquad (16)$$

Following these definitions, the total simulated spectrum can be written as the sum of all bands
$$
\begin{align}
\varepsilon(\nu) &= \sum_i \varepsilon_i(\nu) \nonumber \\
&= \sum_i \varepsilon_i^{\max} e^{-(\nu-\nu_i^{\max})^2 / \gamma_i^2} \tag{17} \\
&= \frac{\sqrt{\pi} N_A e^2 n}{10^3 \ln(10) m_e c} \sum_i \frac{f_i}{\gamma_i} e^{-(\nu-\nu_i^{\max})^2 / \gamma_i^2}. \nonumber
\end{align}
$$

The absorption cross section is [26]

$$
\begin{align}
\sigma(\nu) &= \ln(10) \frac{10^3}{N_A} \varepsilon(\nu) \nonumber \\
&= \frac{\sqrt{\pi} e^2 n}{m_e c} \sum_i \frac{f_i}{\gamma_i} e^{-(\nu-\nu_i^{\max})^2 / \gamma_i^2}. \tag{18}
\end{align}
$$

It is convenient to express the absorption cross section in terms of the absorption energy

$$
\sigma(E) = \frac{2\pi^{3/2}\hbar e^2 n}{m_e c} \sum_i \frac{f_i}{\Gamma_i} e^{-(E-E_i^{\max})^2 / \Gamma_i^2}, \quad (19)
$$

where $\Gamma_i$ is the energy width of the band and $E_i^{\max}$ its energy maximum. As usual, $\hbar$ is the reduced Planck constant.

To account for the shift between the band maximum $E_i^{\max}$ and the vertical excitation $E_i^v$, it is convenient to introduce an energy-shift parameter $\delta_i$ to account for this effect:

$$
E_i^{\max} = E_i^v - \delta_i. \tag{20}
$$

Replacing this equation in Eq. (17) leads to (CGS units)

$$
\sigma(E) = \frac{2\pi^{3/2}\hbar e^2 n}{m_e c} \sum_i \frac{f_i}{\Gamma_i} e^{-(E-E_i^{\max})^2 / \Gamma_i^2}. \quad (21)
$$

For energies in eV and absorption cross section in Å².molecule⁻¹, it results in Eq. (13).

From this result, the molar absorbance is

$$
\varepsilon[\mathrm{M}^{-1}\mathrm{cm}^{-1}] = \frac{\sigma[\mathrm{Å}^2 \cdot \mathrm{molecule}^{-1}]}{3.82353 \times 10^{-5}}. \qquad (22)
$$
Acknowledgments

The authors thank the support of the Excellence Initiative of Aix-Marseille University (A*MIDEX) and the project Equip@Meso (ANR-10-EQPX-29-01), both funded by the French Government "Investissements d'Avenir" program. M.B. and R.M. thank the European Research Council (ERC) Advanced grant SubNano (Grant agreement 832237). M.B. and J.M.T. thank the FetOpen grant BoostCrop (Grant agreement 828753).

References

1. Bomble YJ, Sattelmeyer KW, Stanton JF, Gauss J (2004) On the vertical excitation energy of cyclopentadiene. J Chem Phys 121 (11):5236-5240

2. Barbatti M, Aquino AJA, Lischka H (2010) The UV absorption of nucleobases: semi-classical ab initio spectra simulations. Phys Chem Chem Phys 12 (19):4959-4967

3. Avila Ferrer FJ, Cerezo J, Stendardo E, Improta R, Santoro F (2013) Insights for an Accurate Comparison of Computational Data to Experimental Absorption and Emission Spectra: Beyond the Vertical Transition Approximation. J Chem Theory Comput 9 (4):2072-2082. doi:10.1021/ct301107m

4. Zobel JP, Nogueira JJ, Gonzalez L (2017) The IPEA dilemma in CASPT2. Chem Sci 8 (2):1482-1499. doi:10.1039/C6SC03759C

5. Barton D, König C, Neugebauer J (2014) Vibronic-structure tracking: A shortcut for vibrationally resolved UV/Vis-spectra calculations. J Chem Phys 141 (16):164115. doi:10.1063/1.4898665

6. Reddy CS, Prasad MD (2016) A Gaussian Wave Packet Propagation Approach to Vibrationally Resolved Optical Spectra at Non-Zero Temperatures. J Phys Chem A 120 (16):2583-2590. doi:10.1021/acs.jpca.6b00308

7. Baiardi A, Bloino J, Barone V (2013) General Time Dependent Approach to Vibronic Spectroscopy Including Franck-Condon, Herzberg-Teller, and Duschinsky Effects. J Chem Theory Comput 9 (9):4097-4115. doi:10.1021/ct400450k

8. Segarra-Martí J, Segatta F, Mackenzie TA, Nenov A, Rivalta I, Bearpark MJ, Garavelli M (2020) Modeling multidimensional spectral lineshapes from first principles: application to water-solvated adenine. Faraday Discuss. doi:10.1039/C9FD00072K

9. Reddy CS, Prasad MD (2015) Finite temperature vibronic spectra of harmonic surfaces: a time-dependent coupled cluster approach. Mol Phys 113 (19-20):3023-3030. doi:10.1080/00268976.2015.1070928
10. Santoro F, Improta R, Lami A, Bloino J, Barone V (2007) Effective Method to Compute Franck-Condon Integrals for Optical Spectra of Large Molecules in Solution. J Chem Phys 126 (8):084509. doi:10.1063/1.2437197

11. Charaf-Eddin A, Cauchy T, Felpin F-X, Jacquemin D (2014) Vibronic spectra of organic electronic chromophores. RSC Advances 4 (98):55466-55472. doi:10.1039/C4RA10731D

12. Crespo-Otero R, Barbatti M (2012) Spectrum Simulation and Decomposition with Nuclear Ensemble: Formal Derivation and Application to Benzene, Furan and 2-Phenylfuran. Theor Chem Acc 131 (6):1237. doi:10.1007/s00214-012-1237-4

13. Schreiber M, Silva-Junior MR, Sauer SPA, Thiel W (2008) Benchmarks for electronically excited states: CASPT2, CC2, CCSD, and CC3. J Chem Phys 128 (13):134110

14. Christiansen O, Koch H, Jorgensen P (1995) The 2nd-Order Approximate Coupled-Cluster Singles and Doubles Model CC2. Chem Phys Lett 243 (5-6):409-418

15. Hättig C, Köhn A (2002) Transition Moments and Excited-State First-Order Properties in the Coupled-Cluster Model CC2 Using the Resolution-of-the-Identity Approximation. J Chem Phys 117 (15):6939-6951. doi:10.1063/1.1506918

16. Dunning Jr. TH (1989) Gaussian Basis Sets for Use in Correlated Molecular Calculations. I. The Atoms Boron Through Neon and Hydrogen. J Chem Phys 90 (2):1007-1023. doi:10.1063/1.456153

17. Schirmer J (1982) Beyond the Random-Phase Approximation: A New Approximation Scheme for the Polarization Propagator. Phys Rev A 26 (5):2395-2416

18. Ahlrichs R, Bär M, Häser M, Horn H, Kölmel C (1989) Electronic-Structure Calculations on Workstation Computers - the Program System Turbomole. Chem Phys Lett 162 (3):165-169

19. Barbatti M, Granucci G, Ruckenbauer M, Plasser F, Crespo-Otero R, Pittner J, Persico M, Lischka H (2018) NEWTON-X: A package for Newtonian Dynamics Close to the Crossing Seam (v. 2.2). Available via the Internet at www.newtonx.org.

20. Barbatti M, Ruckenbauer M, Plasser F, Pittner J, Granucci G, Persico M, Lischka H (2014) Newton-X: A Surface-Hopping Program for Nonadiabatic Molecular Dynamics. WIREs: Comp Mol Sci 4 (8):26-33

21. Mansour R, Barbatti M (2019) Spectra Simulation. Available for download at https://github.com/Ritam19/Spectra-Simulation.git.

22. Band YB, Freed KF (1977) Energy distribution in selected fragment vibrations in dissociation processes in polyatomic molecules. J Chem Phys 67 (4):1462-1472. doi:10.1063/1.435020
23. Schinke R (1995) Photodissociation Dynamics: Spectroscopy and Fragmentation of Small Polyatomic Molecules. Cambridge University Press, Cambridge

24. Biczysko M, Bloino J, Santoro F, Barone V (2011) Time-Independent Approaches to Simulate Electronic Spectra Lineshapes: From Small Molecules to Macrosystems. In: Computational Strategies for Spectroscopy. John Wiley & Sons, Inc., pp 361-443. doi:10.1002/9781118008720.ch8

25. da Silva FF, Almeida D, Martins G, Milosavljevic AR, Marinkovic BP, Hoffmann SV, Mason NJ, Nunes Y, Garcia G, Limao-Vieira P (2010) The electronic states of pyrimidine studied by VUV photoabsorption and electron energy-loss spectroscopy. Phys Chem Chem Phys 12 (25):6717-6731. doi:10.1039/b927412j

26. Parson WW (2007) Modern Optical Spectroscopy. Springer-Verlag, Berlin Heidelberg