Two-photon absorption spectroscopy of stilbene and phenanthrene: Excited-state analysis and comparison with ethylene and toluene Free

Non-linear spectroscopies have revolutionized many areas of materials and life sciences. Techniques based on two-photon absorption (2PA), for example, provide temporal and spatial control of the optical excitation process, allowing localized fluorescence excitation and chemical activation. These 2PA methods enable biological applications ranging from 3D and super-resolution imaging to photodynamic therapy, optogenetics, and targeted cell deactivation,^1–7 often at wavelengths that are long enough for deep tissue penetration. Two-photon activation has similarly been used for non-linear optical (NLO) material applications, including microfabrication, nanolithography, 3D optical data storage, ultrafast electro-optic switching, and a variety of novel nanobiophotonics applications.^8–14 In addition, 2PA spectroscopy is a powerful research tool that provides complementary (and sometimes more detailed) information about the electronic structure of a molecule relative to one-photon absorption (1PA) spectroscopy.^15–17 Thus, the theoretical prediction of 2PA spectra is essential for both the in silico design of better two-photon chromophores and the interpretation of experimental spectra.

Designing new two-photon active materials and interpreting experimental 2PA spectra require a fundamental understanding of the quantum-mechanical principles behind the NLO process. As described by Maria Göppert-Mayer in 1931,¹⁸ 2PA becomes allowed in the second order of perturbation theory and has a quadratic dependence of the absorption strength on the light intensity. Consequently, two-photon transitions are governed by different selection rules than one-photon transitions, and the magnitude of 2PA cross sections is generally much lower than for 1PA. Specific details of the 2PA spectrum are more difficult to obtain because of both the difficulty in measuring experimental 2PA spectra and the inherent challenges in calculating higher-order molecular properties.¹⁹ If asked to predict the main features of the 2PA spectrum of a molecule, little can be said beyond symmetry-imposed selection rules^20,21 (e.g., only gerade-gerade and ungerade-ungerade transitions are allowed for centro-symmetric molecules, in contrast to the 1PA selection rules), possible resonance enhancement of some transitions, and a handful of structural predictions derived from simple empirical (“few states”) models.^11,22 Moreover, the computational protocols for computing 2PA spectra and the reliability of existing methods have not yet been fully established. Although several benchmark studies carefully investigated the effects of various approximations on the computed cross sections by comparing different theoretical methods,^23–26 no comprehensive comparisons between calculated and experimental spectra exist, the set of molecules investigated in previous validation studies is rather small, and the emphasis was on low-lying electronic states.

The goal of this work is to investigate trends in the 2PA spectra of several prototypical molecules (ethylene, toluene, stilbene, and phenanthrene). These compounds represent a series of conjugated hydrocarbons of increasing size, where ethylene and toluene can be viewed as building blocks of the NLO-active stilbene and phenanthrene. The choice of stilbene is motivated by its interesting photochemistry, which includes cis/trans isomerization and electrocyclization pathways that represent the key reaction channels in many photochromic molecular switches.^13,27–34 The 2PA spectroscopy of stilbene has been a subject of theoretical studies for several decades,^15,35,36 however, only a handful of studies have employed correlated many-body electronic structure methods.^37–39 To treat stilbene and larger conjugated systems, many studies employed semi-empirical approaches and time-dependent density functional theory (TDDFT).^40–49 Several important general trends have emerged from these studies, including an increasing 2PA intensity with an increasing oligomer length,^40,44 larger response of trans- relative to cis-polyenes,⁴⁰ and proposed design guidelines for enhanced 2PA cross sections in dendritic molecules.^41–43 Other studies examined steric effects on the 2PA cross sections of substituted stilbenes,^46,47 solvent effects on the 2PA of trans-stilbene and charge-transfer homologues (e.g., D-π-D, D-π-A, and A-π-A),⁵⁰ the photoswitching of 2PA properties,⁴⁵ and enhanced 2PA upon dimerization,⁵¹ among other NLO properties of stilbene and stilbene analogues.

In this paper, we examine the calculated 2PA spectra of ethylene, toluene, stilbene, and phenanthrene using a recently implemented²⁵ method for calculating 2PA cross sections with equation-of-motion coupled-cluster wave functions with single and double substitutions (EOM-CCSD). Analyzing the electronic structure of the states that give rise to the dominant features in the 2PA spectra using wave function analysis tools^52–54 provides insight into the 2PA process and allows further calibration of computational protocols for the 2PA calculations. The structure of the paper is as follows. In Sec. II, we outline the details of theoretical methods and computational protocols, including a detailed description of the 2PA cross section calculation and related symmetry-based properties. We discuss general trends and symmetry-imposed selection rules from group theory and for atomic transitions. We then present the results for the computed 1PA and 2PA spectra and analyze the character of the states that give rise to the most prominent features in the 2PA spectra. The comparison of the computed 2PA spectra with the experimental spectra can be found in Ref. 55, where we also analyzed the performance of the electronic structure methodology.

EOM-CC methods provide an efficient and robust computational approach for describing multiple electronic states, including electronically excited states that have multi-configurational wave functions.^56–62 The EOM-CC methods are similar, or equivalent (in some aspects), to linear response CC.^63–67 The EOM-CC wave function has the following form:

where the linear EOM operator $R^$ acts on the reference CC wave function, $eT^|Φ0⟩$⁠. The operator $T^$ is an excitation operator satisfying the reference-state CC equations

where $Φμ$ are $μ$-tuply excited (with respect to $Φ0$⁠) determinants, $H¯$ is the similarity-transformed Hamiltonian (⁠$H¯=e−T^H^eT^$⁠), and the reference CC energy is

The EOM amplitudes R and the corresponding energies are found by diagonalizing $H¯$ in the space of target configurations defined by the choice of operator $R^$ and reference $Φ0$ (see Ref. 59 for examples of different EOM models),

In EOM-EE-CCSD (EOM-CCSD for excitation energies), the CC and EOM operators are truncated as follows:

where only the single and double excitation operators (1h1p and 2h2p) are retained in $T^$ and $R^$⁠. Since $H¯$ is a non-Hermitian operator, its left and right eigenstates, $Lk†$ and R^k, are not Hermitian conjugates of each other. Rather, they form a biorthogonal set

where $L^=L^1+L^2$⁠.

For exact wave functions, the two-photon transition moments at frequencies $ω1$ and $ω2$ between states $|0⟩$ and $|k⟩$ are given by the following expressions:

where $μ^$ denotes the dipole moment operator, $Ωn0$ is the transition energy between states $|0⟩$ and $|n⟩$⁠, and the sums run over all electronic states of the system.

In the expectation-value approach to molecular properties,¹⁹ this expression serves as a starting point for deriving two-photon moments for approximate wave functions. This strategy was used in the implementation of 2PA cross sections within the EOM-EE-CCSD and ADC (algebraic diagrammatic construction) frameworks.^25,68 In this approach, the states and the transition energies in Eqs. (9) and (10) are replaced by the EOM-CCSD (or ADC) wave functions and transition energies. Then the expressions are converted into a computationally efficient form by using the resolvent expression and by introducing auxiliary vectors that are solutions of response-like equations.²⁵ Although the programmable expressions have different form from the above sum-over-state equations, they are identical, both formally and numerically. Thus, for the purposes of this paper, we will use Eqs. (9) and (10) for the analysis of the computed spectra.

From the two-photon transition moments, the components of the transition strength matrix, S_ab,cd, are computed as the products of the left and right two-photon transition amplitudes

Using atomic units, the rotationally averaged 2PA strength, $⟨δ2PA⟩$⁠, reads^21,69

where F, G, and H are integer constants, which depend on the polarization of the incident light. F = G = H = 2 for parallel linearly polarized light, whereas F = −1, G = 4, H = −1 for perpendicular linearly polarized light. We caution that these coefficients, which were originally derived by Monson and McClain,⁶⁹ are sometimes reported incorrectly in the literature (i.e., a sign error for perpendicular polarization⁷⁰). The macroscopic 2PA cross section, $σ2PA$⁠, is measured in GM (Göppert-Mayer) units. It is related to the molecular transition strength in the degenerate case²⁶ (⁠$ω1=ω2=Eexc2$⁠) through the following expression:

where $α$ is the fine structure constant, a₀ is the Bohr radius, $ω$ is the photon energy for the excitation of $2ω$⁠, c is the speed of light, and $S(2ω,ω0,Γ)$ is the line-shape function describing spectral broadening effects. In the non-degenerate case, i.e., when two different frequencies are used, the macroscopic cross section reads

where $ωσ=ω1+ω2$⁠. Note that this expression includes an additional factor of two compared with Eq. (13) to account for the lower photon dissipation rate when the two photons are absorbed. The photon dissipation rate of each laser equals the 2PA rate for excitation with two beams, whereas the dissipation rate is twice the 2PA rate when both photons come from the same laser beam.²⁶ 

The line-shape function introduces a phenomenological broadening for the computed stick spectrum. It includes an empirical damping parameter Γ, which describes the spectral broadening due to rovibrational excitations and collisional dynamics. The Lorentzian line-shape function is usually used for gas-phase calculations

where $ω0$ is the excitation energy and Γ is the half-width at half-maximum (HWHM). This function has a maximum at $ωσ=ω0$⁠,

In the condensed phase, a Gaussian line shape is commonly used to describe inhomogeneous broadening

The maximum of the Gaussian is

In this paper, we use a relatively small damping factor (Γ = 0.1 eV) for the 2PA spectra so that the dominant transitions do not spill over the entire spectrum. Using a different broadening factor affects not only the width of each peak but also the maximum amplitude. The uncertainty in the choice of damping factors can limit the predictive power of the calculations.

For calibration purposes, we also calculated 1PA spectra. For ethylene and toluene, 1PA spectra were computed with a Gaussian broadening of Γ = 0.4 eV applied to each transition (see Section 1 of the supplementary material). The calculations presented here also neglect Franck-Condon factors (arising due to the structural relaxation in the excited states), non-Condon electronic effects,⁷¹ and solvent contributions, all of which can broaden or shift the 2PA bands. Additionally, we report only the 2PA spectra for degenerate excitation: $ω1$ = $ω2$ = E_ex/2.

From group theory, one-photon transitions between states I and F are allowed when

where $Γn$ denote irreducible representations (irreps; not to be confused with linewidths from Sec. II B) and $μi$ represents the x, y, or z polarized component of the incident light. In centrosymmetric molecules, for example, only transitions between states of opposite parity (g to u or u to g) are allowed.²⁰ For two-photon transitions, the dipole moment operator acts twice in the sum-over-states expressions, Eqs. (9) and (10). Consequently, two-photon transitions between states I and F are allowed when

where $Γint$ is the irreducible representation (irrep) of an intermediate state in the sum-over-states expression, and $μi$ and $μj$ describe the polarizations of the two photons. As the 2PA transition moments are obtained by summing over all states of the system, there should exist at least one intermediate state that is dipole-allowed for both the initial and final states in order for the transition to be 2PA-allowed. For example, in the centrosymmetric molecules, the dipole-allowed transitions are g to u; thus, in 2PA, only g to g and u to u are allowed. Furthermore, one may expect larger cross sections for 2PA transitions between the states that are connected by a larger number of intermediate states. Using these arguments, we illustrate below that the largest 2PA cross sections are expected for transitions where I and F belong to the same irrep. For example, for closed-shell molecules, the excited states from the fully symmetric irrep will feature the brightest 2PA transitions.

To illustrate this point, let us consider two specific examples: the C_2υ point group and the case of atomic transitions. Figure 1 shows the allowed two-photon transition pathways as well as the corresponding two-photon transition matrices for molecules in the C_2υ point group, which has four irreps (A₁, A₂, B₁, and B₂). For 1A₁ → nA₁ transitions, the intermediate states can include states from three irreps (B₁, B₂, and A₁), each providing one contribution to the two-photon transition matrix (M_xx, M_yy, and M_zz, respectively). For example, the B₁ intermediate states only contribute to $MxxnA1←1A1$⁠,

Similarly, the $MyynA1←1A1$ and $MzznA1←1A1$ elements include contributions from all B₂ and A₁ intermediate states, respectively. Note that only the diagonal elements of the matrix are nonzero in this case due to symmetry requirements for the totally symmetric transition. For 1A₁ → nB₁ transitions, on the other hand, only the xz and zx off-diagonal components of the two-photon transition matrix are nonzero (see Figure 1).

The diagonal matrix elements contribute to all three terms in the orientationally averaged parallel 2PA cross section, Eq. (12), giving a total of 15 terms for 1A₁ → nA₁ transitions

In contrast, the two off-diagonal components for the 1A₁ → nB₁ transitions give only four terms in the orientationally averaged 2PA cross section

Similar results are obtained for transitions to the nA₂ and nB₂ states, for which there are only two off-diagonal elements in the two-photon transition matrix, leading to four terms in the 2PA cross section. Thus, the 2PA spectrum of a C_2υ molecule will be dominated by totally symmetric two-photon transitions based on the larger number of pathways to the A₁ states as well as the larger contribution of the diagonal elements in calculating the 2PA cross sections (i.e., 15 terms compared with only four). In general, the diagonal tensor elements of the two-photon transition matrix belong to the totally symmetric irreps for all point groups; therefore, the totally symmetric states will dominate parallel 2PA spectra under non-resonant conditions (i.e., when there is no resonant one-photon transition to an intermediate state). We note that a similar analysis of the structure of the 2PA tensor has been presented in Ref. 20.

A similar analysis of the atomic selection rules helps to rationalize the cross sections for transitions involving Rydberg states of polyatomic molecules. Since only $Δl=±1$ transitions are dipole-allowed, only $Δl=0,±2$ transitions are allowed in 2PA. Figure 2 illustrates this point in more detail; it also shows intermediate states for selected transitions. Based on the number of intermediate states available for a particular transition, one may expect $Δl=0$ transitions (e.g., $s→s$ and $p→p$⁠) featuring, in general, larger 2PA cross sections than the $Δl=±2$ transitions (e.g., $s→d$ and $p→f$⁠). For example, the top row of the figure shows that there are three possible pathways for $s→s$ transitions, each contributing a diagonal element to the two-photon transition matrix, $Maa(n+1)s←ns$⁠. In contrast, the $s→d$ transitions have only two off-diagonal pathways to each of the d_xy, d_xz, and d_yz orbitals and one diagonal pathway to $dz2$⁠. Given larger contributions from the diagonal elements to Eq. (12), the $s→s$ transition will dominate the 2PA spectrum relative to the $s→$ d_xy, d_xz, d_yz, or $dz2$ transitions. Two-photon transitions to $dx2−y2$ are not symmetry-allowed.

The last two rows of Figure 2 illustrate a similar result for transitions originating from the np_z orbital, for which there are two off-diagonal pathways to each of the p_x and p_y final orbitals (i.e., via the s and d_xz or d_yz intermediate states), compared with four diagonal pathways to p_z (one passing through s, one through $dz2$⁠, one through d_xz, and one through d_yz). The key result is that the totally symmetric $s→s$ and $pz→pz$ atomic transitions have larger cross sections than the non-totally symmetric transitions based on simple symmetry arguments.

The symmetry arguments in Sec. II C are most appropriate for excitation conditions in which none of the intermediate one-photon accessible states dominate the sum-over-states expressions in Eqs. (9) and (10). This condition (i.e., no resonance enhancement) is achieved when all one-photon transition energies $Ωn0$ are very different from the photon energies, $ω1$ and $ω2$⁠. However, it is also instructive to consider a special case in which a single intermediate state dominates the 2PA transition. For example, intermediate state $|i⟩$ will dominate the sum-over-states expressions if it is low-lying and is strongly dipole-coupled to both the initial and the final states. In this case, the 2PA cross section for degenerate excitation (⁠$ω1$ = $ω2$ = $Ωk02$⁠) and parallel polarization takes the form

or, in the case of the non-hermitian EOM theory

Such a three-state model has been used to rationalize trends in the non-linear properties of a variety of conjugated chromophores.^47,72–74 Within the three-state model, we can limit our attention to the ground, intermediate, and target excited states. As follows from Eq. (24), large 2PA cross sections can be attained by maximizing the brightness of the 1PA transitions from the ground to the intermediate state and from the intermediate state to the target state.

The one-photon transition dipole moment between two states is

where $γIF$ is the one-particle transition density matrix connecting two states, $ΨI$ and $ΨF$⁠. The transition density matrix is defined as

where $p^†$ and $q^$ are the creation and annihilation operators corresponding to orbitals $ϕp$ and $ϕq$⁠. As one can see from Eqs. (26) and (27) (and as discussed in Ref. 47), the transition dipole moment is large when the corresponding one-particle transition density is delocalized. This can be achieved when the electronic transition involves significant charge-resonance character over a large area. One can further simplify the analysis by considering a pair of states that differ by the state of one electron, such as HOMO → LUMO excitation from the closed-shell ground state. Then the expression for the transition dipole moment can be written in terms of a single matrix element between the orbitals involved in the transition

As described below, the wave function analysis based on the one-particle transition density matrix allows one to describe electronic transitions between correlated many-electron states in terms of natural transition orbitals, so Eq. (28) can be used to rationalize the brightness of transitions between states described by general multi-configurational wave functions.

Eq. (28) shows that large one-photon transition dipoles can be obtained for transitions between bonding and antibonding orbitals, such as $π→π∗$ transitions in ethylene. The delocalization of the two orbitals due to increased conjugation in larger molecules will result in even larger dipole moments. A large matrix element can also be obtained when one orbital is a delocalized valence orbital (such as an extended π orbital in a conjugated molecule) and the second is a delocalized Rydberg orbital. In order to maximize both M_0i and M_ik, the target state should be a one-electron transition from the intermediate state and the respective molecular orbitals should be delocalized. Assuming that the intermediate bright state is accessed through a bonding-antibonding transition, a large value of M_ik can be achieved for two different types of target states: (i) a doubly excited valence state, such as a hypothetical $π∗π∗$ state of ethylene (note that such doubly excited configuration belongs to the fully symmetric irrep) or (ii) a singly excited state of Rydberg character, such as the $π→Ry(3p)$ state of ethylene. We note that diffuse character of the target state can result in larger 2PA cross sections even beyond the three-state model because a delocalized orbital can lead to larger transition dipole moments for many intermediate states in the sum-over-states expression.

Let us now discuss the molecular orbital character of electronic states giving rise to large 2PA cross sections. The character of an excited state can be assigned by analyzing the leading electronic configuration of the underlying wave function, that is, by identifying the leading amplitude and the respective MOs. A more rigorous (and orbital-invariant) way^{52–54,75–77} to analyze excited states is based on the one-particle transition density matrix, Eq. (27). $γIF$ contains the information needed to compute one-electron transition properties such as transition dipole moments (i.e., Eq. (26)). It also contains a compact description of the changes in electronic structure between $ΨI$ and $ΨF$⁠. The norm of $γIF$ measures the degree of one-electron character in the $ΨI→ΨF$ transition. For example, $||γ||$ = 1 for those transitions which are of purely one-electron character,⁷⁸ e.g., for the transitions between the Hartree-Fock ground state and target configuration interaction single (CIS) states. The norms of the EOM amplitudes R₁ and R₂ from Eq. (5) are $||R1||≈$ 1 and $||R2||≈$ 0 when the transition has one-electron character.

The one-particle transition density matrix can be decomposed using the singular value decomposition procedure giving rise to the most compact description of the electronic transition

where U and V are the unitary matrices and $𝚲$ is a diagonal matrix containing non-negative numbers

Usually, only a small number of $λi$ are significant. For example, for the CIS transition of pure HOMO-LUMO character, there is only one non-zero $λi$ and it equals one (for spinless $γ$ represented in the basis of spatial MOs). For each significant $λi$⁠, one can define a pair of natural transition orbitals (NTOs); they can be interpreted as hole and particle orbitals for the transition

To quantify the number of significant contributions to the specific transition, one can define the NTO participation ratio

For the CIS transition of pure HOMO → LUMO character, PR_NTO = 1, whereas for the state which has equal weights of the HOMO → LUMO and HOMO−1 → LUMO+1 configurations, PR_NTO = 2 (for example, PR_NTO equals 1.11 for the $π→π∗$ transition of ethylene and 2.06 for toluene).

Using the NTO basis, one can define the exciton wave function as follows:

The square of the exciton wave function gives the probability that the hole is located at the position r_h and the electron at the position r_e. One can define the average positions of the hole and the electron as the centroids of the respective NTOs

The exciton can be characterized by the absolute value of the mean separation vector $|d→h→e|=|⟨r→e−r→h⟩|$ that quantifies the charge separation between the hole and the electron; the hole size $σh=⟨r→h2⟩−⟨r→h⟩2$⁠; the electron size $σe=⟨r→e2⟩−⟨r→e⟩2$⁠; and the electron-hole pair correlation coefficient $Reh=⟨(r→h−⟨r→h⟩)⋅(r→e−⟨r→e⟩)⟩σhσe$⁠. In the analysis below, we provide these quantities as well as the $||γIF||$ and ||R₁||² values, which quantify the degree of single excitation involved in a transition.

Additional insight can be derived from the properties of the electronic states, such as permanent dipole moments and the spatial extent of the wave functions.⁷⁹ The size of the wave function can be characterized by the expectation value of R²,

Since the size of the electronic wave function increases with the system size, we consider the difference between $⟨R2⟩$ of the target EOM-EE-CCSD state and of the CCSD reference state⁷⁹ 

For a valence state, the $⟨R2⟩$ value is similar to the ground-state value, but for Rydberg states, the value is larger. Although the magnitude of the individual Cartesian components of R² depends on the choice of the coordinate system (or molecular orientation), the absolute value (⁠$⟨R2⟩$⁠) is invariant with respect to the molecular orientation, owing to the properties of the trace. In the following, we use this quantity to distinguish between Rydberg and valence states.

The structures of ethylene (D_2h), toluene (C_s), trans-stilbene (C_2h), cis-stilbene (C₂), and phenanthrene (C_2v) were optimized at the B3LYP/cc-pVTZ level of theory. All relevant geometries are given in the supplementary material. The effect of the basis set and of using Cholesky decomposition on the computed properties was investigated using ethylene and toluene. We employed singly and doubly augmented Dunning’s double- and triple-ζ basis sets. We also considered “-df” versions of the bases derived from the parent naug-cc-pVXZ bases by removing d and f functions from the diffuse augmenting sets. The errors due to freezing core electrons and due to using Cholesky decomposition with a threshold of 10⁻² have been also quantified. The results of this analysis are summarized in the supplementary material. We found that (i) at least a doubly augmented Dunning basis set is needed, (ii) when considering low energy states (below ∼9 eV, for ethylene), we can use a double-ζ basis set and remove the d diffuse functions from both augmenting sets, and (iii) using the Cholesky decomposition and freezing the core electrons have a negligible effect on 2PA cross sections (less than 3%). Thus, the d-aug(-d)-cc-pVDZ is the minimal basis set for computing the 2PA cross sections of low energy states. In the production calculations, we employed d-aug-cc-pVTZ basis set for ethylene, d-aug(-df)-cc-pVTZ basis set for toluene, and d-aug(-d)-cc-pVDZ basis set for larger molecules. Except for ethylene, we freeze core electrons and use the Cholesky decomposition with a threshold of 10⁻².

All calculations were performed using the Q-Chem electronic structure program.^80,81 The reported symmetry labels of the electronic states and MOs correspond to the standard molecular orientation used in Q-Chem, which differs from the Mulliken convention.^82,83

The ethylene molecule belongs to the D_2h point group which has 8 irreps, four with ungerade symmetry and four with gerade symmetry. One-photon transitions from the X¹A_g ground state can access only ungerade states with B_1u, B_2u, and B_3u symmetry, but two-photon transitions to all of the gerade states are symmetry-allowed. The totally symmetric transitions to A_g states have three diagonal components of the two-photon transition matrix, whereas each of the other three irreps (B_1g, B_2g, and B_3g) has only two off-diagonal components that are nonzero.⁸⁴ Therefore, based on symmetry considerations alone, transitions to nA_g states are expected to dominate the 2PA spectrum.

Figure 3 shows the computed 1PA and 2PA spectra of ethylene. As expected from the symmetry-imposed selection rules, the two spectra are markedly different. Here, we computed ten EOM-EE-CCSD singlet states per irrep. The lowest ionization energy of ethylene is 10.75 eV (EOM-IP-CCSD/d-aug-cc-pVTZ) and corresponds to ionization from the b_1u π HOMO. Thus, the bound part of each spectrum ends at 10.75 eV and the spectral features beyond this value correspond to autoionizing resonances (the description of these states in our calculations is very approximate).

Applying a phenomenological broadening of Γ = 0.4 eV to each of the transitions in the 1PA spectrum results in several distinct absorption bands, including a single strong band below 9 eV. This low-energy band is dominated by two bright transitions, 1¹B_1u and 1¹B_3u. Similarly, the broadened 2PA spectrum (Γ = 0.1 eV) has a single, strong band below 9 eV and several additional bands at higher energy. The lowest-energy 2PA band contains a dominant transition to the 2¹A_g state at 8.51 eV, with weaker transitions to 1¹B_2g and 1¹B_3g giving rise to a shoulder on the low-energy side of the band.

To illustrate the orbital character of the dominant 1PA and 2PA transitions below 9 eV, Figure 4 shows the natural transition orbitals (NTOs) associated with the 1¹B_1u, 1¹B_3u, 2¹A_g, 1¹B_2g, and 1¹B_3g states of ethylene. The relevant properties derived from the wave function analysis are summarized in Table I. All transitions have participation ratios around 1; thus, each transition can be described by a single pair of hole-particle NTOs. The hole NTO is very similar for all five states and can be described as a valence π orbital nearly identical to the b_1u HOMO. The particle NTO for the 1¹B_3u state has $π∗$ character and is accessible via $π→π∗$ valence excitation, whereas the 1¹B_1u state has more diffuse Ry(s)-like character. The valence and Rydberg-like character of the two 1PA-accessible states is evident from the very different spatial extent of the wave functions, as measured by $Δ⟨R2⟩$⁠, which equals 3.7 and 10.0 Ų for 1¹B_3u and 1¹B_1u, respectively. The other three particles’ NTOs in Fig. 4 represent the excited states responsible for the lowest-energy 2PA band. Their shapes suggest that these states can be described as Rydberg-type states of 2p character. As predicted on the basis of the atomic selection rules (Fig. 2), the two-photon transition to the 2¹A_g state is brighter than the transitions to 1¹B_2g and 1¹B_3g. The transition to 2¹A_g resembles a $pz→pz$ excitation, which can be accessed via the intermediate s, d_z², d_xz, and d_yz states. Each of these four pathways contributes to a diagonal element to the two-photon transition matrix. The 1¹B_2g and 1¹B_3g states resemble p_x and p_y Rydberg orbitals that are accessible through the intermediate s and d_xz or d_yz states, respectively, and therefore contribute only to two off-diagonal elements for each transition (see Fig. 2). As discussed above, the additional pathways and larger contributions from diagonal elements explain the larger relative 2PA cross section for the 2¹A_g state.

Alternatively, one can rationalize the relative 2PA cross sections in the context of the three-state model in Eq. (24). The weight of the leading amplitude in the EOM-CCSD wave function of the 2¹A_g state is 0.95; thus, this is a singly excited state. Singly excited character of this state is further confirmed by the square norm of transition density matrix, $γ$⁠, which is equal to 0.942. In a three-state model, excitation to the 2¹A_g state can be described as $π→π∗→pz$ transition via the brightest 1PA intermediate state, 1¹B_3u. Both the valence $π→π∗$ and Rydberg $π∗→pz$ transitions have large transition dipoles according to Eq. (28), therefore, the overall transition has a large 2PA cross section. In addition, the one-photon accessible 1¹B_1u intermediate state facilitates an alternative $π→s→pz$ pathway to the 2¹A_g state, as well as the $π→s→px,py$ transitions to the 1¹B_3g and 1¹B_2g states, respectively.

The diffuse character of the 2PA active states is confirmed by the quantitative wave function analysis, which gives $Δ⟨R2⟩$ values in the range 13-19 Ų. The size difference between the electron and hole NTOs, $σe−σh$⁠, provides an additional measure of the relative diffuseness of the states, with values ranging from 2.7 to 3.5 Å. We note that for these three states, a more diffuse character correlates with larger 2PA cross sections, as predicted above.

Toluene belongs to the C_s symmetry group, which has only two irreps (A′ and A″ both of which are active in both 1PA and 2PA. We calculated 15 states per irrep to obtain the 1PA and 2PA spectra in Fig. 5. The lowest ionization energy of toluene is 9.21 eV (EOM-IP-CCSD/d-aug-cc-pVDZ), well above the energies of the states calculated here. The brightest transitions in the low-energy region of the 1PA spectrum correspond to the 4¹A″ and 7¹A′ states. In the 2PA spectrum, the three states with the largest cross sections are 2¹A″, 5¹A′, and the much stronger 13¹A′. The key parameters characterizing these states are collected in Table II, and the respective NTOs are shown in Fig. 6. The ||R₁||² and $||γ||$ values show that all of these states are singly excited.

The two leading NTOs of the 7¹A′ state, which have the largest 1PA cross section, contribute with weights of 0.49 and 0.33. The first pair of NTOs can be described as a $π→dz2$ Rydberg transition, whereas the second pair of NTOs corresponds to a $π→π∗$ transition. The partial valence character of this transition is evident from the moderately diffuse value of $Δ⟨R2⟩$ = 9.8 Ų. The second brightest 1PA state, 4¹A″, also has mixed Rydberg/valence NTOs (⁠$π→p$ and $π→π∗$⁠), although the Rydberg-like transition carries more weight (0.64 and 0.19, respectively) and therefore gives a final state with slightly more diffuse character of $Δ⟨R2⟩$ = 10.3 Ų.

Each of the three brightest 2PA states has only one dominant pair of NTOs (⁠$λ1$ ≈ 0.85). The particle NTOs for all three states have significant Rydberg character and become increasingly diffuse for higher-lying states. The transition to the 2¹A″ state has $π→s$ Rydberg character (⁠$Δ⟨R2⟩$ = 14.1 Ų), the transition to 5¹A′ has $π→p$ Rydberg character (⁠$Δ⟨R2⟩$ = 20.9 Ų), and the transition to 13¹A′ has the largest 2PA cross section and $π→d$ Rydberg character (⁠$Δ⟨R2⟩$ = 33.4 Ų). As for ethylene, the values of $σe−σh$ correlate perfectly with $Δ⟨R2⟩$ due to the similarity of the hole orbital for all three transitions. Because of the lower symmetry of toluene relative to ethylene, the states in toluene have moderate charge-transfer character, as indicated by the non-zero $|d→h→e|$ values in Table II. We observe the largest $|d→h→e|$ values for the 4¹A″ and 7¹A′ states (the two brightest 1PA states), which have the least diffuse character of the five states discussed here.

Judging from the shapes of the NTOs, the 13¹A′ state resembles a totally symmetric $dyz→dyz$-type transition, whereas 2¹A″ and 5¹A′ resemble $d→s$ and $f→p$ transitions. Thus, on the basis of the atomic selection rules (Fig. 2), one would expect the following trend in the respective 2PA cross sections: 13¹A′ > 5¹A′ > 2¹A″ because of the number of pathways between the respective (atomic orbital-like) NTOs, which is consistent with the computed values. Here we again observe that more diffuse character correlates with larger cross sections.

Within the few-state model, the dominant one-photon transitions to 7¹A′ and 4¹A″ provide a glimpse of the molecular orbital pathways. The leading $Ry(dz2)$-like NTO of 7¹A′ provides a likely intermediate for the $f→dz2→p$ pathway to the 5¹A′ state, and the Ry(p_z)-like NTO for 4¹A″ provides possible pathways for the $d→s$ and $d→d$ transitions to 2¹A″ and 13¹A′, respectively. Although other pathways are also likely to contribute to the sum-over-states 2PA expression, the strong one-photon transitions and partial Rydberg-like character suggest that the 7¹A′ and 4¹A″ states play an important role as intermediates.

The compounds trans-stilbene, cis-stilbene, and phenanthrene have similar chemical structures but belong to different symmetry point groups. trans-stilbene belongs to the C_2h point group, which has two irreps of gerade symmetry that are active in 2PA and two of ungerade symmetry that are active in 1PA. The A_g transitions should dominate the 2PA spectrum because the two-photon transition matrix has five nonzero components including the three diagonal elements, compared with two off-diagonal elements for B_g transitions.⁶⁹ cis-stilbene is in the C₂ point group, with two irreps (A and B) that are active in both 1PA and 2PA. Totally symmetric transitions to the A states should dominate the 2PA spectrum in this point group as well, due to five non-zero matrix elements, compared with only two off-diagonal elements for states with B symmetry. Phenanthrene belongs to the C_2v point group, which has four total irreps that are all active in 2PA but only three of which are active in 1PA. As illustrated in Fig. 1, the states belonging to the totally symmetric irrep of C_2v also should have the largest cross sections because the diagonal elements of the two-photon transition matrix are nonzero for A₁ transitions, whereas the other irreps each have only two off-diagonal elements that are nonzero.

Figure 7 shows the computed 2PA spectra of trans-stilbene, cis-stilbene, and phenanthrene. The spectra were computed using 8 states per irrep. The vertical ionization energies are 7.58 eV, 7.84 eV, and 7.79 eV (EOM-IP-CCSD/d-aug(-d)-cc-pVDZ) for trans-stilbene, cis-stilbene, and phenanthrene, respectively. In the analysis below, we focus on the states with the largest 2PA cross sections, which generally are not the lowest-energy transitions. Table III collects relevant properties of the selected excited states. The ||R₁||² and $||γ||$ of these states reveal singly excited character, and the respective NTOs are shown in Figs. 8–10.

The 2PA spectra of all three compounds are dominated by the fully symmetric (A_g, A, and A₁) transitions, as predicted above. The spectrum of trans-stilbene has a very strong transition to the 5¹A_g state near 6.3 eV. The cross section of that transition is ∼10 times larger than the next most intense transition in the molecule (to the 4¹A_g state) and the strongest of all the 2PA transitions calculated in this paper. The spectrum of cis-stilbene has a strong 2PA band near 6.3 eV that includes transitions to the 6¹A and 9¹A states, as well as contributions from several weaker transitions. Finally, the 2PA spectrum of phenanthrene has two main bands due to the 5¹A₁ and 9¹A₁ states, respectively, with some weaker underlying contributions.

The absolute 2PA cross sections for trans-stilbene, cis-stilbene, and phenanthrene are larger than for ethylene and toluene, which can be explained by the extended delocalization of the respective states. Most notable is the very strong transition to the 5¹A_g state of trans-stilbene, which has two dominant NTO pairs (Fig. 8). The first pair of NTOs resembles a $π→π∗$ transition with some diffuse Ry(p)-like character, but the second pair has more compact $π→π∗$ character. The two pairs of NTOs resemble linear combinations of the 5¹A′ state of toluene, having $π→p$ Rydberg character. The leading pair of NTOs can be described as the symmetric combination of NTOs from toluene, and the second pair as the antisymmetric combination with additional $π∗$ character of the ethylene bridge (resembling the 1¹B_3u particle NTO of ethylene). With $Δ⟨R2⟩$ = 10.2 Ų, the 5¹A_g state of trans-stilbene is more diffuse than a pure $π∗$ valence state but not as diffuse as the Rydberg states of the other compounds due to its mixed Rydberg/valence character.

The strongest 2PA transitions of cis-stilbene (6¹A and 9¹A) also resemble combinations of the toluene NTOs but have more complicated character due to the reduced conjugation of the nonplanar structure. The twisted structure limits the interaction between the two phenyl rings and reduces the electronic coupling across the ethylene bridge. This is evident from the leading NTO of the 6¹A state, which resembles an anti-symmetric combination of the 5¹A′ state (⁠$π→p$⁠) of toluene for each side of the molecule. The deviation from planarity also reduces the intensity of the 1PA transitions and explains the reduced 2PA cross sections in cis-stilbene relative to trans-stilbene.

The lower cross section for cis-stilbene compared with trans-stilbene can be rationalized based on the molecular orbital characteristics and Eqs. (24)–(27). The delocalization effect is quantified in Table IV using a three-state model to compare the calculated 2PA cross sections for a single intermediate state. Despite the very similar 1PA intensities for initial to intermediate states, and for intermediate to target states, the overall 2PA cross section of trans-stilbene via the low-energy 1¹B_u state is an order of magnitude higher than any of the other transitions in the table.

Finally, for phenanthrene, the 5¹A₁ and 9¹A₁ states have the largest 2PA cross sections although they are less intense than the two strongest transitions in cis-stilbene. The 5¹A₁ state is related to the 5¹A_g state of trans-stilbene, with the leading pair of NTOs resembling a $π→p$-like Rydberg transition. The transition to the higher-energy 9¹A₁ state has three main pairs of NTOs: the first and third pairs resembling two different $π→d$-like Rydberg transitions and the second pair corresponding to $π→π∗$-like transition. The mixed valence/Rydberg character of the two phenanthrene states is evident from the moderately diffuse orbitals (⁠$Δ⟨R2⟩$ = 11-14 Ų), compared with predominantly Rydberg transitions in cis-stilbene (⁠$Δ⟨R2⟩>$ 18 Ų). The reduced conjugation of cis-stilbene may be partially compensated by the more delocalized character of the dominant states in that spectrum compared with the other two compounds (see Table III). The NTOs of trans-stilbene and phenanthrene are not as diffuse (⁠$σ≤$ 5 Å), possibly due to the partial valence character of the excited states. We note that the relationship between 2PA cross section and electronic delocalization, as well as the reduced cross sections for twisted structures, has been documented in previous studies.^40,47,51 A more detailed discussion of the valence and Rydberg character of 2PA transitions is provided in Ref. 55. The distinction is especially important for comparison with experimental spectra in solution, where Rydberg states might be pushed up to higher energies.

In this contribution, we investigated 2PA spectra of several prototypical molecules (ethylene, toluene, cis- and trans-stilbene, and phenanthrene) and analyzed the electronic structure of the states that give rise to the dominant spectral features. By employing the best available methodology²⁵ for calculating 2PA cross sections for molecules of this size, which is based on the EOM-CCSD formalism, our calculations provide an important high-quality reference data for future theoretical and experimental studies. We also use an orbital-invariant density-matrix based approach⁵³ for analyzing the characters of the states in question. This analysis, which is more rigorous than qualitative description of the characters of MOs involved in the transition, provides a valuable contribution towards understanding photophysics of stilbene and phenanthrene.

We reviewed the symmetry-imposed selection rules and showed that, since the diagonal tensor elements of the two-photon transition matrix belong to the totally symmetric irreps, those transitions typically dominate the 2PA spectra for parallel polarization. Similar symmetry arguments based on the atomic selection rules predict that $Δl$ = 0 transitions have the largest 2PA cross sections among Rydberg-like transitions, which further helps to rationalize the calculated 2PA spectra. We also note that for the same molecule, a more diffuse character of the target state corresponds to a larger 2PA cross section, which can be explained by considering transition dipole moments between the molecular orbitals involved in the brightest valence transitions and the diffuse Rydberg orbitals.

trans-stilbene features the largest 2PA cross sections among the molecules investigated here owing to the most extensive delocalization of the relevant states. The 5¹A_g state has the largest cross section. The first pair of the corresponding NTOs can be described as $π→p$-like Rydberg transition, where the diffuse NTO of the electron resembles the bonding linear combination of p-like Rydberg MOs on each phenyl ring. In cis-stilbene, due to the nonplanar geometry, the exciton delocalization is reduced, resulting in lower 2PA cross sections (⁠$π→$ 6¹A and 9¹A) relative to the 5¹A_g state of trans-stilbene. In phenanthrene, the largest 2PA cross sections correspond to p and d-like Rydberg target states.

The general trends that affect the 2PA transition strengths derived on the basis of simple symmetry arguments and atomic selection rules provide an important foundation for predicting and rationalizing 2PA spectra of increasingly complex chromophores. Such fundamental insight is essential for the development of new two-photon active compounds and NLO materials, as well as interpreting the spectra of existing compounds. For example, the larger 2PA cross sections for Rydberg states compared with valence states may play an important role in determining the most efficient intermediate states for resonance-enhanced multi-photon ionization (REMPI) schemes involving a two-photon transition.

See supplementary material for additional details and relevant Cartesian geometries.

A.I.K. acknowledges support by the U.S. National Science Foundation (Grant No. CHE-1264018). C.G.E. acknowledges support from the U.S. National Science Foundation through a CAREER Award (No. CHE-1151555). We are grateful to Dr. Kaushik Nanda for valuable discussions and his help with the calculations. We thank Dr. Thomas Jagau for his feedback about the manuscript.