Classification of doubly excited molecular electronic states

Abstract

Electronic states with partial or complete doubly excited character play a crucial role in many areas, such as singlet fission and non-linear optical spectroscopy. Although doubly excited states have been studied in polyenes and related systems for many years, the assignment as singly vs. doubly excited, even in the simplest case of butadiene, has sparked controversies. So far, no well-defined framework for classifying doubly excited states has been developed, and even more, there is not even a well-accepted definition of doubly excited character as such. Here, we present a solution: a physically motivated definition of doubly excited character based on operator expectation values and density matrices, which works independently of the underlying orbital representation, avoiding ambiguities that have plagued earlier studies. Furthermore, we propose a classification scheme to differentiate three cases: (i) two single excitations occurring within two independent pairs of orbitals leaving four open shells (D_OS), (ii) the promotion of both electrons to the same orbital, producing a closed-shell determinant (D_CS), and (iii) a mixture of singly and doubly excited configurations not aligning with either one of the previous cases (D_mix). We highlight their differences in underlying energy terms and explain their signatures in practical computations. The three cases are illustrated through various high-level computational methods using dimers for D_OS, polyenes for D_mix, and cyclobutane and tetrazine for D_CS. The conversion between D_OS and D_CS is investigated using a well-known photochemical reaction, the photodimerization of ethylene. This work provides a deeper understanding of doubly excited states and may guide more rigorous discussions toward improving their computational description while also giving insight into their fundamental photophysics.

---

1. Introduction

Electronic states with doubly excited characters have aroused interest and generated lively debate in recent years. Their involvement in singlet fission^1–3 can provide a promising route towards highly efficient photovoltaic devices, but they are also of particular interest in other technological applications, such as non-linear optical spectroscopy,^4–6 and thermally activated delayed fluorescence.⁷ In the photochemistry of polyenes^8–10 and derived systems, such as carotenoids,^11,12 states with partially doubly excited character play a crucial role and have been investigated for over 50 years, yet, still inciting fiery discussions. Despite this substantial interest, there is no well-defined framework for classifying doubly excited states or even a well-accepted definition of doubly excited character. Different authors use different definitions; consequently, even the simple example of the A_g state of butadiene, and its assignment as singly or doubly excited, has recently sparked considerable controversy.^13,14

Indeed, the description of doubly excited states is still a challenge for computational chemistry.¹⁵ Many commonly used methods, such as linear-response time-dependent density functional theory (TDDFT),¹⁶ second-order approximate singles and doubles coupled cluster (CC2), or the second-order algebraic diagrammatic construction (ADC(2))¹⁷ fail in their description. Describing double excitations within a single-reference framework requires going up hierarchies and using more intricate and expensive methods. For example, the CC3 (ref. 18) and ADC(4) methods¹⁹ are considered appropriate, whereas even ADC(3)²⁰ can be problematic.^19,21 Moreover, in a multireference framework,²² one can produce accurate descriptions of doubly excited states. However, these methods are accompanied by the ever-present problems of choosing an appropriate active space and related parameters. Spin-flip methods^23,24 and state-specific orbital-optimized DFT^25,26 present themselves as interesting alternatives but require special care in their applications as well. The problem is exacerbated by the fact that there is no clear rule to indicate when such more involved methods are required, and related questions are heavily contested in the literature. A typical example is an ongoing discussion of which methods are suitable to describe the lowest A_g state of butadiene.^8,13,14 Therefore, having a well-defined and method-independent quantifier for double excitation character could greatly help for issues of this type.

Aside from methodological questions, it is also desirable to gain a deeper understanding of the underlying physics of doubly excited states with the eventual goal of designing optimised molecules for specific tasks. At this point, it is particularly interesting to contrast doubly excited states on individual molecules with intermolecular doubly excited states in terms of their energies and wave function properties. However, no rigorous and transferable classification scheme exists that would allow comparing these cases meaningfully. Note that discussions of classification schemes for doubly excited states in the literature are restricted to two-electron atoms where notably different physics is at play.^27,28

To obtain a basic definition of singly or doubly excited character, one might sum over the weights of all singly excited configurations (denoted %T₁ henceforth). However, such an assignment is only meaningful within a given wave function model and set of reference orbitals. As a consequence, the assignment may vary if a different computational method is chosen, and it is not even immediately clear whether the concept of a doubly excited state possesses intrinsic physical meaning at all (cf. ref. 29). More specifically, it is unclear whether any given doubly substituted Slater determinant should be interpreted as contributing to correlation or as an actual double excitation.^13,14 Furthermore, one should realise from a fundamental physical viewpoint that a double excitation is represented by four correlated particles—two excitation holes and two excited electrons. Formally, such two-body processes should not be described by orbitals but by geminals.³⁰ Therefore, significantly enhanced complexity can be expected compared to singly excited states. As a consequence, doubly excited states have remained quite elusive and ambiguous in the discussions so far.

We propose solving the first problem, the assignment of doubly excited states, by using density matrices. They are well-defined independently of the wave function model, thus, allowing us to extract molecular orbital pictures and numerical descriptors from correlated wave functions.^31–33 More specifically, we base our analysis on transition and difference density matrices, cancelling out correlation contributions also present in the ground state. The use of density matrices provides meaning to these descriptors via their connection to physical observables. Specifically, we elaborate on the viewpoint that a doubly (or higher) excited state is a state that cannot be coupled to the ground state with any conceivable one-electron operator (cf. ref. 31 and 34). Secondly, to deal with the enhanced complexity of doubly excited states, we combine and contrast several analysis methods to obtain a well-defined and simple yet comprehensive picture. These methods amount to the squared norm of the 1-electron transition density matrix (1TDM),^34,35 the promotion number³⁶ based on the attachment and detachment densities, and the occupation of the natural orbitals (NO), collectively, through the number of unpaired electrons³⁷ or, individually, through the occupation of the lowest unoccupied NOs (LUNO and LUNO+1).³⁸ In addition, we apply an extension of the excitation number as defined by Barca et al.¹⁴

This work aims to comprehensively describe molecular doubly excited states in quantum chemistry computations. We start with a Theory section presenting various definitions of doubly excited character, contrasting different limiting cases for doubly excited states, and discussing the underlying energy contributions. Three illustrative examples follow (Fig. 1). We use the formaldehyde dimer to illustrate the limiting case of a doubly excited state involving two independent pairs of open-shell orbitals (denoted the D_OS case). Subsequently, we study polyenes highlighting the complexity of their A_g excited states involving a mix of partial doubly and singly excited character (denoted D_mix). To examine the interconversion between the open- (D_OS) and closed-shell (D_CS) limiting cases, we investigate the dimerization of ethylene. Finally, the three different archetypes of doubly excited states are reviewed in a more extended set of molecules.

Fig. 1 Structures of the molecules studied within this work: the series of all-trans-polyenes up to 6 alternating double bounds, ethylene dimer and formaldehyde dimer.

2. Theory

2.1 Definition of doubly excited character via transition density matrices

Double excitations are traditionally defined via the %T₁ values, which reflect the total weight of single excitations. However, the challenge in using %T₁ values is that they are only defined within a given computational method, and it is unclear, for example, how to compare results from single- and multireference computations. Therefore, we choose a different route here. We start with a method-independent definition of doubly excited states based only on physical observables without any explicit reference to orbitals or wave functions. Subsequently, we use this starting point to derive the squared 1TDM norm Ω as a rigorous and method-independent substitute for %T₁.

Within the 1TDM picture, we define a state as being doubly (or higher) excited via the condition that it is impossible to couple it to the ground state with any conceivable one-electron operator. We can turn this definition into a practical rule if we first realize that an arbitrary transition property of a one-electron operator between wave functions Ψ_i and Ψ_f is given as

where D^if_pq is the 1TDM, A_pq is the matrix representation of the operator, and both are given for a molecular orbital (MO) basis {ϕ_p}. The 1TDM, in turn, is defined as

D^if_pq = 〈Ψ_i|p^†q|Ψ_f〉 (2)

where p^† and q are the creation and annihilation operators related to the MOs ϕ_p and ϕ_q. Note that Ψ_i is allowed to be a general, correlated wave function here, meaning that separating these MOs into occupied and virtual orbitals is impossible. Applying the Cauchy–Schwarz inequality to eqn (1),³⁴ we obtainwhere the symbol Ω denotes the squared Frobenius norm of the 1TDM, that is,

Importantly, we find that the transition property a_if necessarily vanishes if Ω vanishes, that is, if all elements of the 1TDM are zero. Conversely, if any 1TDM element is non-zero, there is at least a conceivable one-electron operator with a non-vanishing transition property. An Ω value of zero is equivalent to the statement that the state cannot be coupled via a one-electron operator. Thus, a doubly (or higher) excited state exhibits Ω = 0, whereas a purely singly excited state exhibits Ω = 1.

More generally, Ω can be seen as an effective proportionality factor stating how strongly the transition interacts with one-electron operators. Therefore, a value of Ω between 0 and 1 can be used to represent a partial doubly excited character.³¹ In practice, the Ω value is consistent with the fraction of singly excited amplitudes (%T₁) alluded to above³⁵ and presents a natural generalization of this concept. An alternative viewpoint, based on ref. 85, is presented in Section S1.† In Section S2 we discuss the possibility of Ω values larger than 1.†

2.2 Definition of doubly excited character via (difference) density matrices

As an alternative to the 1TDM, it is possible to view doubly excited character via the 1-particle state density matrix (1DM) or difference density matrix (1DDM). The 1DM is defined in analogy to eqn (2) as

D^ff_pq = 〈Ψ_f|p^†q|Ψ_f〉 (5)

Furthermore, the difference density matrix (1DDM) is simply the difference between two state 1DMs

Δ^if = D^ff − Dⁱⁱ (6)

Diagonalization of the 1DDM and separation of the eigenvectors according to their signs provide detachment and attachment densities.^35,36 The sum over all positive or negative eigenvalues of the 1DDM—denoted as the promotion number p—gives the total number of electrons rearranged during the excitation process. In principle, p could lie between zero and the total number of electrons, but in practice, it usually ranges from 1 to 2. The promotion number appears to be a natural measure for defining a multiply excited character. However, since it is not only affected by the electrons taking part in the primary excitation process but also by secondary orbital relaxation,^35,39p has been considered an unsuitable measure.¹⁴ More generally, doubly excited states usually experience increased p values, but increased p values alone are not a sufficient criterion to assign doubly excited character.

As an alternative measure for double excitation character, the excitation number (η) was introduced by Barca et al.¹⁴ in the context of the maximum overlap method (MOM). For two single-determinantal wave functions, Φ_i and Φ_f, η is defined as

where S^if_jk is the overlap between the j-th occupied orbitals of Φ_i and the k-th occupied orbital of Φ_f; n is the total number of electrons in the system. The η value ranges from 0 (when Φ_i is equal to Φ_f) to n (when there is no overlap between any orbitals in Φ_i and Φ_f). A generalization of η to arbitrary wave functions is not trivial, and we discuss this issue in some detail in Section S3.† After several attempts, we suggest using the formula

η = n_eff − tr(DⁱⁱD^ff)/2 (8)

where Dⁱⁱ and D^ff are the spin-traced 1DMs. The value of n_eff in this equation is defined as

n_eff = max(‖Dⁱⁱ‖²,‖D^ff‖²)/2 (9)

This expression reduces to eqn (7) for a single-determinantal wave function with doubly occupied spatial orbitals. Furthermore, just like the original expression, it is invariant to a switch between the initial and final states and vanishes if the initial and final states are the same. Moreover, we have verified that this expression yields the expected result in model systems and a variety of realistic computations.

Finally, we want to point out that p and η are both based on the 1DMs, and thus, they analyze shifts in the electron density rather than probing the actual wave functions. Crucially, if the initial and final state should have the same 1DMs, then p and η would both vanish, formally classifying the state as a “zero-electron transition”. This complicates the assignment in cases of strong ground-state correlation. For example, states where the HOMO and LUMO are singly occupied (HOMO¹LUMO¹) or states constructed as a linear combination of the configurations with doubly occupied HOMO or doubly occupied LUMO (HOMO²−LUMO², HOMO²+LUMO²), all possess the same 1DMs, and, thus, the p and η values between them would vanish.

2.3 Further descriptors

If the state is predominantly singly excited (Ω ≈ 1), then it is meaningful to analyze the 1TDM further and to obtain the natural transition orbitals (NTOs), defined as the singular vectors of the 1TDM.^40,41 The number of independent NTO pairs necessary to describe the transition, called the NTO participation ratio,^42,43 is defined aswhere λ_i are the weights of the NTO pairs. PR_NTO can be used to assess the multiconfigurational character of the transition. A value of 1 means that a single configuration state function can express the excitation; higher values imply that this state has a multiconfigurational nature.

Furthermore, we compute the spin-averaged natural orbitals (NOs), defined as the eigenvectors of the 1DM. Their occupation numbers (n_p) range from 0 (unoccupied) to 2 (doubly occupied). It is common to characterize the 1DM via y₀ and y₁, which correspond to the occupation numbers of the lowest unoccupied natural orbitals (LUNO) and LUNO+1, respectively. For example, (y₀, y₁) equals (0,0), (1,0), and (1,1) correspond to a closed shell, pure diradical, and pure tetra-radical characters, respectively.^38,44 In Section 2.4, we show how these quantities can differentiate between types of doubly excited states. Alternatively, one can also compute the number of effectively unpaired electrons by summing over all NOs of the system. Eqn (12) and (13) show two expressions to obtain the number of unpaired electrons via either

orwhere n_u includes both static and dynamic correlation, while n_u,nl suppresses dynamic correlation, thus, focusing on static contributions.^31,37 In principle, n_u and n_u,nl range from zero to the number of electrons. A value of zero represents a closed shell; a value of two a biradical with two open-shell orbitals, four represents a tetra-radical with four open shells, etc.

2.4 Classification of singly and doubly excited states

The descriptors presented above provide a toolbox for a comprehensive description of electronic excitation processes. This section shows how they can be combined to give a well-defined classification scheme of singly and doubly excited states. For this purpose, we discuss the values of the descriptors for four limiting cases: a single configurational (S_SC) and a multiconfigurational (S_MC) singly excited state, and a closed-shell (D_CS) and an open-shell (D_OS) doubly excited state (Fig. 2). In addition, we will consider the mixed case (D_mix) as a case with notable doubly excited character not conforming with any of the four limiting cases.

Fig. 2 Limiting cases for singly and doubly excited states using a four-orbital four-electron model. Singly excited states distinguish between single configurational (S_SC) and multiconfigurational (S_MC) cases; doubly excited states distinguish between the formal closed-shell (D_CS) and four open-shells (D_OS) cases.

Herein we use Ω as the main characteristic to distinguish between singly (Ω = 1) and doubly (Ω = 0) excited states. Alternatively, p and η can be used where their values correspond to the number of excited electrons. Note, however, that p is also strongly affected by orbital relaxation,^35,39 and it is unclear how η performs in the case of static ground-state correlation.

Within the singly excited states, we distinguish between the single-configurational (S_SC) and multiconfigurational (S_MC) limiting cases. In the single-configurational limiting case (S_SC), the overall excitation can be described as a transition between a single pair of orbitals, e.g., the HOMO → LUMO transition. More generally, we define the S_SC case as a state with only a single contributing NTO pair, leading to a value of PR_NTO = 1. For the S_SC case, (y₀, y₁) is equal to (1, 0) since only one virtual orbital is involved. The multiconfigurational case (S_MC) is obtained if at least one additional pair of NTOs contributes to the state. In the scheme presented in Fig. 2, S_MC is represented by PR_NTO = 2 and (y₀, y₁) = (0.5, 0.5). The significance of these differences is discussed in the literature, for example, in the context of excitons and the L_a/L_b states in aromatic molecules.^43,45–48

Within the doubly excited states, we distinguish between the closed-shell (D_CS) and open-shell (D_OS) limiting cases. In the first case, exemplified by a pure HOMO² → LUMO² transition, two electrons are promoted to the same virtual orbital, and the excited state obtains a closed-shell character (D_CS) similar to the ground state. The second is exemplified by a combined HOMO/HOMO−1 → LUMO/LUMO+1 transition: two electrons are promoted from two different initial orbitals to two different final orbitals, leaving four open-shell orbitals in total. D_CS states can only be realized for singlets due to the Pauli principle, whereas D_OS states can be singlet, triplet, or quintet. Within the presented scheme, the D_CS and D_OS states are resolved via NO occupations. The D_CS limiting case possesses only closed shells and therefore has n_u,nl = 0 and (y₀, y₁) = (0, 0), in analogy to a closed-shell ground state. Here, the LUMO of the ground state becomes a strongly occupied MO of the excited state. The four open-shell orbitals in the D_OS case, on the other hand, are represented by n_u,nl = 4 and (y₀, y₁) = (1, 1).

The idealized D_CS state is single-configurational and behaves like a closed-shell ground state. Such a state would be readily described by a single determinant and would be particularly amenable to the maximum overlap method (MOM).⁴⁹ Furthermore, a CAS(2,2) active space or a single spin flip from a triplet reference would both suffice to describe such a state. On the other hand, a D_OS case would always require a more sophisticated treatment, including at least four active orbitals. Similarly, the D_mix case requires several correlated orbitals to describe its multiconfigurational character (unless the required nondynamic correlation effects can be captured within the exchange–correlation functional employed). We shall explore these issues below in Section 4.2.

2.5 Energies of doubly excited states of monomers and dimers

It is instructive to start by presenting the energies of the various states that can be constructed within a two-orbital two-electron model (TOTEM), as shown in Fig. 3 (see also ref. 8, 16 and 50). For simplicity, we consider transitions from HOMO to LUMO. Four spin-adapted wave functions can be constructed within the TOTEM: the ground state (¹G) with a doubly occupied HOMO, the single configurational singly excited states (¹S_SC/³S_SC) of singlet and triplet multiplicities, and the closed-shell doubly excited state (¹D_CS). The relevant energy terms are the one-electron energies of HOMO and LUMO (h_H, h_L), the three Coulomb integrals (J_HH, J_HL, J_LL), and the exchange integral (K_HL). Fig. 3 indicates the different energy terms contributing to the energy. The energy of the ground state (¹G) is determined by the one-electron energy of the HOMO (h_H) and the Coulomb integral between the two electrons located in the HOMO; the energy of the doubly excited state (¹D_CS) is determined in a completely analogous way only that the HOMO is swapped for the LUMO. Determining the energies of the singlet and triplet S_SC states is slightly more involved since they are described by two interacting configurations, meaning that an exchange term (K_HL) also comes into play. In summary, the energies of the states are given as

E(¹G) = 2h_H + J_HH (13)

E(¹S_SC) = h_H + h_L + J_HL + K_HL (14)

E(³S_SC) = h_H + h_L + J_HL − K_HL (15)

E(¹D_CS) = 2h_L + J_LL (16)

Fig. 3 Excited-state diagrams constructed within a two-orbital two-electron model. Energy terms are shown in red: one-electron energies of HOMO (h_H) and LUMO (h_L), the three Coulomb terms (J_HH, J_HL, J_LL), and the exchange term (K_HL).

Before continuing, we note that the one-electron energies used above include the kinetic energy, the nucleus-electron attraction and possibly the interaction with any other electrons present (treated in the sense of a frozen core) but do not consider any terms involving the HOMO and LUMO. An alternative and equally valid viewpoint is provided in ref. 8 by using the orbital energies (ϵ_H = h_H + J_HH, ϵ_L = h_L + 2J_HL – K_HL) that already include interactions within the HOMO and the LUMO.

Using the TOTEM, we can now examine under what circumstances a closed-shell doubly excited state can be of lower energy than a singly excited state. Solving for E(¹D_CS) < E(¹S_SC) with the definitions given above, we obtain

h_L − h_H < K_HL + J_HL − J_LL ≈ K_HL (17)

where the right-hand side was simplified under the assumption that the J_HL and J_LL Coulomb integrals are of similar magnitude.⁸ In other words, the doubly excited state becomes favorable if the exchange repulsion is large compared to the difference in one-electron energies between HOMO and LUMO.

Note, however, that eqn (17) and the assumption that the Coulomb integrals are of similar magnitude would also imply that E(³S_SC) < E(¹G), i.e. that the triplet lies below the closed shell.⁸K_HL also couples ¹G and ¹D_CS, meaning that the singlet ground state would obtain multiconfigurational character if the exchange interaction were, indeed, of the same order of magnitude as the difference in one-electron energies. This discussion highlights that a large exchange interaction relatively favors doubly excited states by pushing ¹S_SC up in energy. Nevertheless, it also shows that a simple HOMO² → LUMO² transition cannot be the lowest excited state if the ground state is a closed shell. Indeed, the connection between doubly excited character and static electron correlation in the ground state is emphasized in the literature.^13,14

Noting that a double excitation from HOMO to LUMO is not feasible for low-energy excited states, we proceed to an alternative type of doubly excited state. This alternative is present in the case of a dimer where two locally excited states can be combined into one doubly excited state of D_OS type. Such states can be classified according to the spin-multiplicity of the overall state along with the individual transitions. By combining singlet or triplet states on each monomer and considering all possible spin couplings, one obtains the states ¹(TT), ³(TT),⁵(TT), ³(ST),³(TS), and ¹(SS). Here, the ¹(TT) state is crucial for the singlet fission process, where one high-energy singlet excited state can be converted into two low-energy triplets.³ (See also ref. 52 for a discussion of the ¹(TT) state in the context of spin-exchange internal conversion and ref. 1 for a discussion of intramolecular doubly excited states in singlet fission.) We can use the TOTEM model to evaluate the stability of such a state. If we neglect possible biexciton binding effects, the ¹(TT) state will be the lowest state of singlet multiplicity if twice the excitation energy of ³S_SC is lower than the excitation energy of ¹S_SC. More specifically,

h_L − h_H < 3K_HL + J_HH − J_HL ≈ 3K_HL (18)

This condition is certainly easier to satisfy than eqn (17). In a dimer, the exchange repulsion associated with ¹S_SC is avoided, and two exchange integrals are gained due to the two ³S_SC configurations present.

As discussed below, we found that low-energy D_OS type states can be readily constructed in dimers. Conversely, we were unable to find any low-lying D_CS-type states in a variety of investigated molecules. The D_CS states that were indeed found were of σσ* or nπ* character and trivially lay at about twice the energy of the corresponding singly excited state. Importantly, the low-lying A_g ππ* states of polyenes and related systems do not fit either the D_OS or D_CS limitting case. Therefore, we introduce a third class of excited state, D_mix, which is characterized by appreciable double excitation character (as determined by Ω and η) but not conforming to either of the limitting cases. Reviewing the TOTEM, we note that it is a simplified model capturing only the D_OS and D_CS cases but is unable to account for the energies of D_mix type states. Indeed, a more involved model combining a triplet-pair state and charge-transfer exciton has been suggested for the latter.⁵¹

3. Computational details

Ground state geometry optimizations and vibrational frequencies of ethylene, butadiene, hexatriene, octatetraene, decapentaene, dodecahexatriene, s-tetrazine and the tetracene dimer were obtained at density functional theory (DFT) level with CAM-B3LYP functional,⁵³ cc-pVTZ⁵⁴ basis set, and Grimme's D3 dispersion correction.⁵⁵ DFT calculations were done using Gaussian 16 rev a03.⁵⁶ Vertical excitation energies and wave functions obtained at DFT/MRCI level⁵⁷ employed the def2-TZVP⁵⁸ basis set, except for the tetracene dimer where we use def2-SV(P).⁵⁸ In this approach, the CI expansion is built from Kohn–Sham orbitals using the BH-LYP⁵⁹ functional (as implemented in TURBOMOLE 7.5 ⁶⁰) and an effective Hamiltonian. Here, we tested two different parametrizations: the original one proposed by Grimme and Waletzke⁵⁷ and R2018 proposed by Marian et al.⁶¹ As shown in Fig. S1 of the ESI,† R2018 does not reproduce the inversion between states 1¹B_u and 2¹A_g observed in polyenes with the increase of the number of double bonds in the system,⁶² while the original parametrization does.⁶³ However, both parametrizations yield similar trends regarding the wave function analysis (Fig. S2†). In this work, we chose to use the original set of parameters. The initial reference space included configurations obtained from single and double excitations of 10 electrons within 10 orbitals.

TDDFT calculations were performed with the BLYP functional^64,65 and 6-31G** basis set⁶⁶ with Q-Chem 5.3.⁶⁷ Multireference configuration interaction with single and double excitations (MR-CISD) were carried out using COLUMBUS 7.0.^68–70 A complete active reference space (CAS) including all π and π* orbitals was used for polyenes up to four double bonds; polyenes with five and six double bonds were restricted to a CAS(8,8) due to computational cost. A complete active space self-consistent field (CASSCF) with these active spaces was used to construct the orbitals using state-averaging over the first two A_g and the first B_u states. MRCI energies are reported using the Pople extensivity correction (+P); 1s orbitals of all carbon atoms were frozen.⁷¹ Vertical excitation energies were also computed at the third-order algebraic diagrammatic construction method (ADC(3)) for the polarization propagator^20,72 level with the resolution-of-identity approximation and def2-SV(P) basis set, as implemented in Q-Chem.

In ADC(3) calculations, the wave function analysis library (libwfa)^31,32 was used to obtain Ω-values, participation ratio of the natural transition orbitals (PR_NTO), occupation of natural orbitals (y₀ and y₁), number of unpaired electrons (n_u,nl), and promotion number (p). η at ADC(3) and all descriptors at DFT/MRCI and MRCI levels were obtained externally with TheoDORE⁷³ analysis package using a pre-release of version 3.0.

To investigate the cycloaddition of ethylene, we performed a relaxed scan keeping the linear combination of C1–C3 and C2–C4 (Fig. 1) at fixed distances. Ground state geometry optimizations were done at DFT level using B3LYP/cc-pVTZ. Excited-state energies and wave functions were obtained at MRCI+P(8,8)/cc-pVDZ level. All calculations considered D_2h symmetry. CAS(8,8) wave functions with 6 states in the average (four A_g and two B_u) were used as references. MRCI calculations considered 4 frozen core orbitals (belonging to a_g, b_3u, b_2u, and b_1g representations), and 1 orbital in the active space for each of the 8 irreducible representations of the D_2h point group.

4. Results and discussion

4.1 Formaldehyde dimer

First, we discuss the ten lowest excited states of the formaldehyde dimer at 5 Å intermolecular distance calculated at the ADC(3) level (Table 1). The first six excited states are excitonic combinations of the ³nπ*, ¹nπ*, and ³ππ* states on the monomers. They are predominantly singly excited, as seen by their Ω values near or above 0.9 and η values close to 1.0. Their PR_NTO values close to 2 (along with appropriate n_u,nl, y₀, and y₁ values) allow us to classify them as multiconfigurational (S_MC) states, according to Fig. 2. The multiconfigurational nature of excitonic states delocalized between two interacting chromophores has been discussed in detail elsewhere,^43,45 and we shall proceed with the remaining states here. The final four states are unambiguously characterized as doubly excited due to their Ω values of exactly 0.000, along with η values near 2. More specifically, we find that these four states comply with the definition of an open-shell doubly excited state (D_OS) in Fig. 2, and all relevant descriptors are within 0.05 of their idealized values. This discussion highlights that D_OS states are indeed readily constructed in realistic systems and using a high-level method.

Table 1 Excitation energy (ΔE in eV), oscillator strength (f), and wave function descriptors of the lowest excited singlet and triplet states of the formaldehyde dimer at 5 Å intermolecular separation at ADC(3) level

State	Chars	ΔE (eV)	f	Ω	PRNTO	n u,nl	y 0	y 1	p	η
1^3A2	^3nπ*	3.67	—	0.909	2.02	2.46	0.543	0.525	1.15	1.01
1^3B2	^3nπ*	3.67	—	0.909	2.02	2.46	0.543	0.526	1.15	1.01
1^1A2	^1nπ*	4.12	0.000	0.882	2.01	2.51	0.556	0.537	1.15	1.03
1^1B2	^1nπ*	4.12	0.000	0.882	2.01	2.51	0.555	0.537	1.15	1.03
1^3B1	^3ππ*	6.08	—	0.945	2.09	2.44	0.526	0.510	1.16	0.95
1^3A1	^3ππ*	6.08	—	0.945	2.09	2.44	0.523	0.513	1.16	0.95
2^1A1	nπ*:^1(TT)	7.68	0.000	0.000	—	4.05	0.973	0.972	2.00	1.95
2^3B1	nπ*:^3(TT)	7.68	—	0.000	—	4.05	0.973	0.972	2.00	1.95
3^3B1	nπ*:^3(ST)	8.53	—	0.000	—	4.05	0.972	0.971	2.00	1.95
2^3A1	nπ*:^3(ST)	8.53	—	0.000	—	4.05	0.972	0.972	2.00	1.95

---

We find that the excitation energy of the 2¹A₁ state (7.68 eV) is about twice the excitation energy of the singly excited ³nπ* state (3.67 eV), and we, therefore, assign it as the corresponding ¹(TT) state. Note, however, that the agreement is not perfect, and the ¹(TT) state lies about 0.35 eV higher than expected using the monomer energies. In principle, such a difference could derive from (bi)excitonic interaction effects. However, these are probably negligible at 5 Å, as the states come in pairs of the same energy. Conversely, we ascribe the difference to a lack of internal consistency within ADC(3) in terms of describing singly and doubly excited states at exactly the same level (see Section S6.1 in the ESI†). The next doubly excited state is of triplet multiplicity and almost degenerate with ¹(TT); it is assigned as the ³(TT) state. The final two states shown in Table 1 combine singlet and triplet monomer excitations to form a ³(ST) state and are, again, almost degenerate. Their energy (8.53 eV) is significantly higher than the combined ³nπ* and ¹nπ* energies (7.79 eV), which we again attribute to a lack of complete internal consistency within ADC(3).

4.2 Butadiene and larger polyenes

The photophysics of polyenes is usually discussed in terms of two important states close in energy: the 1¹B_u state dominated by the HOMO → LUMO transition and the 2¹A_g state with at least partial admixture of the doubly excited HOMO² → LUMO² configuration.^9,10 Within the following, we study these states in different polyenes using a variety of computational methods. We emphasize that our primary goal is not to obtain an accurate energy gap value between these states but rather to elucidate the nature of the 2¹A_g state.

The characterization and excited state ordering of polyenes is still a subject of discussion in the literature. Experimentally, 2¹A_g becomes the lowest excited state for polyenes with more than four double bonds.^8,74 Computationally, the energy gap and order of states are strongly method-dependent.^20,63,75,76 While MS-CASPT2 calculations correctly predict the state inversion, CC3 predicts that 1¹B_u remains above 2¹A_g, although both methods deliver a small energy gap between those states.⁶³ ADC(2)-x and ADC(3) always predict 2¹A_g state as the lowest excited state, while ADC(2)-s predicts that to be the 1¹B_u state for polyenes up to four double bonds.²⁰ Benchmark studies show that the gap between 2¹A_g and 1¹B_u in polyenes is reproduced correctly by DFT/MRCI, although the excitation energies are lower than the best theoretical estimate.^61,63,75 Due to these discrepancies, here we use four different computational methods to study polyene excited states: ADC(3), MRCI, DFT/MRCI, and TDDFT/BLYP. For DFT/MRCI we use the original parameterisation because it reproduces the inversion between 1¹B_u and 2¹A_g expected in larger polyenes; see ESI Section S4.†

In particular, butadiene has become a paradigmatic case for discussing doubly excited states.^8,13,14 Thus, we first focus on this molecule using the ADC(3) method. To obtain a comprehensive picture, we look at several density descriptors (following Fig. 2) along with the percentage of single excitations (%T₁) within ADC(3). The results are presented in Table 2. We start the discussion with the 1¹B_u state, which at this level of theory is the second excited state lying at 6.72 eV. Ω and η values close to 1 unambiguously assign the state as being singly excited. Furthermore, PR_NTO approximately 1, along with (y₀, y₁) = (0.965, 0.039), allows classifying it as an S_SC single configurational state close to the idealized values presented in Fig. 2.

Table 2 Excitation energy (ΔE in eV), oscillator strength (f), and wave function descriptors of the lowest excited singlet states of butadiene computed at the ADC(3) level

State	ΔE (eV)	f	Ω	PRNTO	n u,nl	y 0	y 1	p	η	%T1
1^1Ag	—	—	—	—	0.186	0.071	0.053	—	—	—
2^1Ag	6.02	0.000	0.305	1.959	2.430	0.890	0.253	1.672	1.537	31.3
1^1Bu	6.72	1.739	0.904	1.095	2.073	0.965	0.039	0.998	0.899	93.5

---

The assignment of 2¹A_g of butadiene, which lies at 6.02 eV, is more involved. An Ω-value of 0.305, a %T₁ value of 31.3%, and an excitation number (η) of 1.537 indicate a partially doubly excited character, which is also supported by a promotion number (p) of 1.672. However, the descriptors are notably different from the idealized case of a doubly excited state (Ω = 0, η = p = 2) shown in Fig. 2. This implies that the admixture of singly excited configurations plays an important role. The dominant contribution to the 2¹A_g state is the HOMO² → LUMO² transition, with a weight of 42%. The HOMO−1 → LUMO (18%) and HOMO → LUMO+1 (13%) transitions come next, followed by many doubly excited configurations, all involving the HOMO−1, HOMO, LUMO, and LUMO+1. Thus, this state is strongly multiconfigurational, with notably different characteristics to the idealized D_CS case. The NO-based characteristics reflect this divergence particularly well. For the idealized D_CS case, n_u,nl, y₀, and y₁ are all zero. However, for the 2¹A_g state of butadiene, these values are n_u,nl = 2.430, and (y₀, y₁) = (0.890, 0.253), which are between the limiting cases shown in Fig. 2. For these reasons, we classify this state as D_mix.

As a next step, we investigate the dependence of the presented results on the electronic structure method and the size of the molecule. For this reason, we computed excitation energies and 1(T)DM descriptors at MRCI, DFT/MRCI, and ADC(3) levels for all-trans-butadiene, all-trans-hexatriene, all-trans-octatetraene, all-trans-decapentaene, and all-trans-dodecahexene (N = 2, 3, 4, 5, 6, respectively, where N is the number of double bonds in the system).

Starting with the excitation energies (Fig. 4A), we find a substantial decrease with increasing N for all methods. However, there is a notable difference in the 2¹A_g excitation energies according to the trend ADC(3) < TDDFT < DFT/MRCI< MRCI+P, spanning a range of up to ≈1 eV. Concerning the Ω values presented in Fig. 4B, all methods aside from TDDFT find substantial double excitation character (Ω < 0.65) for all the molecules. The Ω values generally decrease with increasing N. The only exception is for the largest MRCI computations, which may be affected by size-extensivity problems. Aside from the general trends, there is also a quite notable difference in the obtained Ω values. ADC(3) indicates strong double-excitation character (Ω < 0.30 for N > 2); MRCI delivers intermediate values (0.35 < Ω < 0.5); DFT/MRCI predominantly indicates single excitation character (Ω > 0.5) albeit with notable admixtures of double excitations. TDDFT/BLYP, on the other hand, always predicts a singly excited character (Ω ≈ 1). The low Ω values for 2¹A_g contrast with the high Ω values for 1¹B_u (Fig. S9†), which are consistently above 0.90 for DFT/MRCI and TDDFT and above 0.85 for ADC(3) and MRCI. Thus, the Ω values clearly distinguish between the wave functions of the 2¹A_g and 1¹B_u states, assigning a partial doubly excited character to 2¹A_g.

Fig. 4 Excitation energies ΔE (A) and 1TDM norm Ω (B) of the 2¹Ag state of different polyenes plotted against the number of double bonds (N) calculated at ADC(3), DFT/MRCI (dubbed DFTCI in the insert), MRCI, and TDDFT/BLYP levels.

As outlined in Section 2.1, Ω has a clear physical meaning by acting as an effective proportionality factor specifying how strongly the states are coupled via one-electron operators. Changes in Ω are related to changes in physically observable transition properties. In the present case, the transition dipole moments of the A_g states vanish for symmetry reasons, but the differences should show up via enhanced transition quadrupole moments or angular momenta when computed with methods that produce different Ω values. However, this discussion is left to future work.

The excitation numbers (Fig. 5A) agree with the Ω-values in the sense that they always attribute at least partial doubly excited character to the 2¹A_g state (η > 1.2). Interestingly, the MRCI and DFT/MRCI values are both fairly low, with η ≤ 1.3, whereas ADC(3) provides significantly enhanced doubly excited character (η > 1.5). One can understand this discrepancy by noticing that these descriptors are differently affected by ground-state correlation (see Section 2.2). Generally speaking, it is not clear how the excitation number, initially developed for comparing single determinantal wave functions, should be interpreted in the case of multiconfigurational states. Nonetheless, it is interesting to discuss the excitation numbers of the 1¹B_u excited states (Fig. S10†). For MRCI, these are always below 0.85; for DFT/MRCI and ADC(3), values below 0.95 are obtained. Thus, a clear differentiation between the singly excited 1¹B_u and doubly excited 1¹A_g states is also present when the excitation numbers are considered. Finally, the y₀ values (Fig. 5B) are close to 1 for all methods, whereas the y₁ values are smaller than 0.5. Note that these y₀/y₁ values are inconsistent with any limiting cases presented in Fig. 2. Conversely, they illustrate the multiconfigurational and partial singly and doubly excited nature of the states, in line with a D_mix character. This is, again, markedly different from the singly excited 1¹B_u states (Fig. S13†), which for all methods aside from TDDFT, exhibit the idealized values (y₀ ≈ 1, y₁ ≈ 0) expected for S_SC states.

Fig. 5 Excitation number η (A) and occupations (B) of the lowest unoccupied natural orbital (LUNO, y₀) and LUNO+1 (y₁) of 2¹A_g state of different polyenes plotted against the number of double bonds (N) calculated at ADC(3), DFT/MRCI (dubbed DFTCI in the insert), MRCI, and TDDFT/BLYP levels.

In summary, we find at least a partial double excitation character, classified as D_mix, in the 2¹A_g state of all polyenes investigated, revealed by all methods (except TDDFT) and descriptors in agreement with much of the previous literature.^8–10 Nevertheless, two contradicting viewpoints argue against the double excitation character of butadiene in the literature, which we discuss next.

Shu and Truhlar¹³ have presented butadiene computations at various computational levels to understand the differences between the 2¹A_g and 1¹B_u states and learn which computational methods are suitable for their description. Crucially, they argued that doubly excited states could not be understood separately, but the multireference character of the ground state promotes the contribution of doubly excited configurations in low-lying excited states. This assessment agrees with the data presented in Table 2 and the discussion in Section 2.5. Nonetheless, we emphasise that the unique properties of polyenes cannot be understood by considering the ground state alone. Otherwise, we would observe similar amounts of double excitations in both 2¹A_g and 1¹B_u. In contrast, the 2¹A_g and 1¹B_u states possess distinct characters: only 2¹A_g obtains double excitations, while 1¹B_u retains the singly excited (S_SC) character. Furthermore, Shu and Truhlar have advocated using local functionals to describe the excited states of polyenes.¹³ They argued that local functionals, such as BLYP, revPBE and M06-L, perform well since they minimize the static correlation error included by the Hartree–Fock exchange in nonlocal functionals.^13,77 Considering only the energies, one finds that TDDFT/BLYP does indeed produce results comparable to the wave-function-based methods (Fig. 4A). However, special care has to be taken due to the differences in Ω-values (Fig. 4B), which would predict significantly altered transition properties.

Subsequently, Barca et al.¹⁴ performed DFT/MOM computations on butadiene and analysed them with the help of their original excitation number (η) definition, as shown in eqn (7). A value of η = 1.022 was obtained for 2¹A_g indicating almost perfect singly excited character. The challenge in interpreting these results is that the correlated computations describe the 2¹A_g state as a mixture of the HOMO−1/LUMO, HOMO/LUMO+1, and HOMO²/LUMO² configurations. However, the DFT/MOM method produces only a single open-shell Slater determinant. Thus, it is doubtful whether the MOM method provides a realistic description of this intrinsically multiconfigurational state. Conversely, extending the excitation number to the multiconfigurational case (Fig. 5A) highlights the doubly excited character. Barca et al. argue that it is not clear a priori whether the doubly substituted determinants in CI wave functions account for electron excitation, electron correlation, or orbital relaxation, not allowing an unambiguous assignment.¹⁴ It is precisely for this reason that we have chosen a rigorously defined set of density-matrix-based descriptors to avoid such ambiguities. All the proposed descriptors are invariant to the orbital representation employed. This choice consistently shows the admixture of doubly excited characters in the 2¹A_g state of butadiene. It is important to point out that this contrasts with both Barca et al.'s assignment as a singly excited state, and Shu and Truhlar's argument that the apparent doubly excited character is due to correlation at the ground state.

4.3 Cycloaddition of ethylene

Unable to locate the D_CS case in polyenes, we proceed to a third model, the cycloaddition of ethylene. The dimerization of ethylene to cyclobutane is a typical example illustrating the changes in the wave function character along a chemical reaction. According to the Woodward–Hoffmann rules, the [2 + 2] cycloaddition of ethylene is thermally forbidden due to a change in the ground state electronic configuration, yet it is photochemically allowed.^78,79 The reaction pathway can be explained by the four frontier orbitals shown in Fig. 6. The left represents the case of two isolated ethylene molecules where the frontier orbitals are of π and π* character. As they get closer, the spatial overlap between the orbitals increases until they eventually form the σ and σ* orbitals of cyclobutane shown on the right. Crucially, the HOMO (π₂) of the isolated ethylene molecules corresponds to the LUMO of cyclobutane and vice versa. Thus, the dimerization requires a change in electronic configuration from . This process is thermochemically forbidden but can be facilitated via a doubly excited state. In this section, we assess the involvement of this doubly excited state, focusing on how its character changes along the reaction coordinate.

Fig. 6 Orbital correlation diagram for two ethylene molecules separated by 4 Å (left) and 1.4 Å (right).

The ethylene dimerization is a practical model illustrating the transition from the D_OS to D_CS limits. The closed- or open-shell character is determined by an interplay between LUMO and LUMO+1 energies and the exchange integral, as outlined in Section 2.5. While the exchange term favors D_OS, D_CS is favored if LUMO and LUMO+1 are farther apart in energy. At large separations, HOMO (π₂) and HOMO-1 (π₁), as well as LUMO and LUMO+1 , are pairwise degenerate. Thus, exchange dominates, and the two electrons are promoted to different unoccupied orbitals in a ¹(TT) type state with a configuration. As the distance decreases, the overlap between the ethylene orbitals increases, becoming non-degenerate. Then, the LUMO+1 becomes inaccessible, and a D_CS state arises. Below, we evaluate the validity of this model using ab initio computations.

Generally speaking, three states are relevant in the dimerization process, 1¹A_g, 2¹A_g, and 1¹B_u. Here, the closed shell and doubly excited states—that is for the dimer and for cyclobutane—are always of A_g symmetry. The singly excited states are always B_u. Fig. 7A presents the MRCI+P energies of these states computed along a relaxed scan. The right side shows the case of isolated ethylene molecules, whereas the left side represents the formation of cyclobutane. Starting with the 1¹A_g curve, we find that the dimerization is energetically favorable, but a substantial energy barrier of over 2.5 eV is encountered, making the reaction unfeasible in the ground state, as mentioned. Considering the excited states, the singly excited state (1¹B_u) is generally below any of the A_g states, except in the avoided crossing region, where 2¹A_g becomes lower in energy. The excited states are fairly flat toward the right, whereas a steep increase of the doubly excited state is seen toward the left once cyclobutane is formed. Indeed, the doubly excited state of interest becomes 4¹A_g and reaches an adiabatic energy of 18.5 eV. Fig. 7A suggests a clear mechanism for the photochemical reaction: a photon is initially absorbed by 1¹B_u. Subsequently, the two molecules are attracted, forming an excimer. Near the minimum of the excimer, a crossing with 2¹A_g is encountered. Furthermore, 2¹A_g finally relaxes to the 1¹A_g ground state, forming cyclobutane.

Fig. 7 Relative energies in eV (A), 1TDM norms Ω (B), and number of unpaired electrons (n_u,nl; C) of selected singlet states of the ethylene dimer/cyclobutane system plotted against the intermolecular separation.

To characterize the amount of double excitation character involved, we use the Ω descriptor, as presented in Fig. 7B. The Ω value for 1¹B_u is close to 1 throughout the energy profile, highlighting the singly excited nature of this state. Conversely, the Ω value for 2¹A_g is close to zero on the right side, highlighting that this state is doubly excited in the limit of the separated dimer. On the left side, the doubly excited character is transferred to the 4¹A_g state, as seen by its Ω value close to zero. A large spike in 2¹A_g is observed during the first avoided crossing around 2.5 Å. At this point, the 1¹A_g and 2¹A_g states become multiconfigurational, and a clear definition of doubly excited character becomes more challenging. Note that there is an additional spike in the region around 2.0 Å. It is related to a second avoided crossing involving 2¹A_g, 3¹A_g, and 4¹A_g.

The singly and doubly excited nature of the 1¹B_u and 2¹A_g states for reactant and product can also be represented by the η and p values (Fig. S16 and S17†). Toward the left and right, they represent the singly and doubly excited character, similar to Ω values. Interestingly, both η and p tend towards zero for 1¹B_u and 2¹A_g at the avoided crossing around 2.5 Å, which can be understood following the discussion in the last paragraph of Section 2.2. The and configurations mix and, as a consequence, all states have the same natural orbitals (with singly occupied π₂ and orbitals) and, hence, the same density matrix.

Having verified the overall amount of doubly excited character, we now use n_u,nl to obtain a more detailed classification of the states (see Fig. 7C). As expected, n_u,nl is close to two for the singly excited 1¹B_u state and, aside from the avoided crossing, near zero for the closed-shell 1¹A_g state. For the doubly excited 2¹A_g state, we find that its n_u,nl value is close to 4 for large intermolecular separations, representing the idealized D_OS (¹TT, that is, ) case. The value steadily decreases as the molecules move together (aside from the avoided crossings). Below the second avoided crossing, the doubly excited state becomes 4¹A_g; its n_u,nl value further decreases until reaching a value of 0.26 for the last point probed. This low n_u,nl value, along with the individual (y₀, y₁) = (0.16, 0.08) values shown in Fig. S18–S20,† highlights that the state at this geometry does, indeed, closely conform with the idealized D_CS state hypothesized in Fig. 2. Note, however, that this doubly excited state occurs at an extremely high vertical excitation energy of 18.5 eV. Indeed, its vertical excitation energy is about twice as high as the singly excited 1¹B_u state. Thus, in line with the above discussion (Section 2.5), we can state that a D_CS-type state can only occur at energies significantly higher than the lowest singly excited state.

4.4 Further molecules

Having outlined the different archetypes of doubly excited character in some detail above, we applied the proposed scheme to larger and more complex molecular systems. The selected molecules and a summary of the results are presented in Fig. 8, and more data is presented in the ESI.† We start with a diketopyrrolopyrrole derivative⁷⁶ used in optoelectronics,^80,81 OLEDs⁸² and as singlet-fission chromophores^83,84 and a derivative of its bis-thiophene building block (Fig. 8). A more extended set of derivatives is given in Table S4.† As pointed out previously,⁷⁶ these molecules have an extended polyene backbone with additional functional groups. The lowest doubly excited singlet states (of A_g symmetry) in these molecules lie at similar energies as the lowest singly excited B_u singlet states. In both cases shown and all examples in Table S4,† the A_g states present the clear signature of the D_mix case (0.3 < Ω < 0.6, 2.5 < n_u,nl < 2.8). This analysis highlights that the D_mix case applies to a variety of molecules and is an important model to understand electronic excitations.

Fig. 8 Analysis of the lowest singly and doubly excited singlet states of various molecular systems: excitation energies (ΔE in eV), squared 1TDM norms (Ω) and numbers of unpaired electrons (n_u,nl).

Other interesting examples are the tetracene dimer and s-tetrazine. The tetracene dimer at 5 Å intermolecular separation is chosen as a more realistic illustration of dimer excited states relevant to singlet fission. Its lowest singlet excited state at 2.16 eV is a doubly excited ππ* state delocalized over both molecules. With values of Ω = 0.00 and n_u,nl = 4.05 it almost perfectly aligns with the D_OS case highlighting the relevance of this case for general dimer excited states. Finally, we investigate s-tetrazine. This molecule possesses a doubly excited nπ* state at 5.11 eV. With values of Ω = 0.00 and n_u,nl = 0.14 it is a close match to the D_CS case. Its energy is about twice the energy of S₁, the singly excited nπ* state. Thus, in agreement with the previous discussions, we find that D_CS type states are only found well above the first singly excited state.

5. Conclusions

This work presents a rigorous and transferable classification scheme for doubly excited states. We propose to define a doubly (or higher) excited state as a state that cannot be coupled to the ground state with any conceivable one-electron operator. This physically meaningful definition can be readily evaluated numerically using the 1TDM norm.

Within the manifold of doubly excited states, we define three cases: the closed-shell (D_CS) case, where two electrons are promoted together from one orbital to another; the open-shell case (D_OS), where the excitations occur between two independent orbital pairs; and the multiconfigurational mixed case (D_mix) possessing only partial doubly excited character and conforming with neither of the above definitions. The underlying energetics are presented, highlighting that the D_OS limiting case can occur as a low-lying excited state in realistic computations. Conversely, the pure D_CS case is expected at significantly higher energies, and only D_mix is a viable model for low-lying intramolecular doubly excited states. The differentiation between D_OS, D_CS, and D_mix in practical calculations using natural orbital occupations and other readily available density matrix descriptors is outlined. We find that the assignment of the D_OS and D_CS cases is generally clear and unambiguous. By contrast, D_mix type states are sometimes discussed quite controversially in the literature. Other authors^13,14 have refrained from using the term “doubly excited” for D_mix type states reserving it for the D_OS and D_CS cases. Ultimately, this is a question of terminology. Nonetheless, this study shows that D_mix states possess unique wavefunction properties that are clearly differentiated from traditional singly excited states and we, therefore, suggest labelling them as states with partial doubly excited character.

Several practical examples are presented to study the occurrence of different types of doubly excited states and their description with different computational methods. First, we highlight that D_OS-type states can be readily constructed in dimers where they occur as different combinations of the monomer singlet and triplet states, such as ¹(TT), ³(TT), ³(ST). Computations of the formaldehyde dimer at the ADC(3) level are presented, illustrating that even at this highly correlated level, one obtains 1(T)DM descriptors close to the idealized results.

We proceed to butadiene and larger polyenes to present results on their controversially discussed 2¹A_g excited states. Computations at the ADC(3), ab initio MRCI, and DFT/MRCI levels all agree that admixture of doubly excited character via the HOMO² → LUMO² transition plays an important role in the 2¹A_g excited states of these molecules. At the same time, the description is never close to a D_CS limiting case, but multiconfigurational character and admixture of singly excited configurations play a significant role in line with the D_mix case.

Unable to locate the D_CS case in polyenes, we proceed to a third model, the [2 + 2] cycloaddition of ethylene. At large intermolecular separations, a doubly excited D_OS (¹TT, ) state is found. Upon dimerization, this state converts into a near-perfect D_CS state that is strongly dominated by the HOMO² → LUMO² transition. However, this state lies at very high energies (above 18 eV), demonstrating, again, that a low-lying D_CS state cannot be achieved. Finally, we applied our scheme to an extended set of molecules to highlight that the three archetypes D_CS, D_OS, and D_mix are transferable models. Thus, we demonstrate the generality of our scheme for understanding doubly excited states in various systems.

In summary, we present a physically motivated definition of doubly excited character and a classification scheme able to distinguish between its limiting cases, providing a new approach to a long-standing problem. More specifically, we have highlighted challenges in the computational description of doubly excited states of various kinds, outlining the requirements for computational methods to describe them accurately. We emphasize that reproducing excitation energies is not enough to ensure the quality of a computational method to a determined system. A computation should only be deemed accurate if its wave functions and operator expectation values also comply with the reference. Nevertheless, the availability of well-defined and transferable descriptors can provide a solid basis for further discussions of the computational description of the doubly excited character. Furthermore, we hope the presented work can provide new ideas in the science surrounding doubly excited states and that the underlying physics discussed here can provide an improved language to discuss experimental results.

Data availability

The underlying computational research data (molecular coordinates and input/output files for Q-Chem, Columbus, DFT/MRCI and TheoDORE) are available at Loughborough University's data repository at https://doi.org/10.17028/rd.lboro.22303765.

Author contributions

Conceptualization, methodology, software: F. P.; data curation, visualization: M. T. do C.; funding acquisition and resources: M. T. do C. and M. B.; investigation, writing (original draft preparation): M. T. do C., F. P.; supervision: F. P., M. B., and J. M. T.; writing (review & editing): M. T. do C., F. P., J. M. T., M. B.

Conflicts of interest

The authors declare no competing financial interest.

Acknowledgements

F. P. thanks Ignacio Fernández Galván and Anna I. Krylov for discussions about wave function analysis methods. F. P. has received funding from the UK EPSRC (reference EP/V048686/1) and acknowledges access to Loughborough's high-performance computing system lovelace. M. T. do C., M. B., and J. M. T. have received funding from the European Union's Horizon 2020 Research and Innovation Programme under grant agreement No 828753 (BoostCrop). M. B. and J. M. T. thank the European Research Council (ERC) advanced grant SubNano (Grant agreement 832237). M. T. do C. is thankful for the support of the PhD mobility program funded by AMUtech. Centre de Calcul Intensif d'Aix-Marseille is acknowledged for granting access to its high performance computing resources.

References

1. M. E. Sandoval-Salinas and D. Casanova, ChemPhotoChem, 2021, 5, 282–293.
2. T. Zeng, R. Hoffmann and N. Ananth, J. Am. Chem. Soc., 2014, 136, 5755–5764.
3. M. B. Smith and J. Michl, Chem. Rev., 2010, 110, 6891–6936.
4. K. W. Stone, D. B. Turner, K. Gundogdu, S. T. Cundiff and K. A. Nelson, Acc. Chem. Res., 2009, 42, 1452.
5. B. P. Fingerhut, M. Richter, J. W. Luo, A. Zunger and S. Mukamel, Phys. Rev. B: Condens. Matter Mater. Phys., 2012, 86, 1–10.
6. C. Heshmatpour, P. Malevich, F. Plasser, M. Menger, C. Lambert, F. Šanda and J. Hauer, J. Phys. Chem. Lett., 2020, 11, 7776–7781.
7. P. De Silva, J. Phys. Chem. Lett., 2019, 10, 5674–5679.
8. J. H. Starcke, M. Wormit, J. Schirmer and A. Dreuw, Chem. Phys., 2006, 329, 39–49.
9. B. Hudson and B. Kohler, Annu. Rev. Phys. Chem., 1974, 25, 437–460.
10. K. Schulten, I. Ohmine and M. Karplus, J. Chem. Phys., 1976, 64, 4422–4441.
11. H. A. Frank and R. J. Cogdell, Photochem. Photobiol., 1996, 63, 257–264.
12. X. Hu, A. Damjanović, T. Ritz and K. Schulten, Proc. Natl. Acad. Sci. U. S. A., 1998, 95, 5935–5941.
13. Y. Shu and D. G. Truhlar, J. Am. Chem. Soc., 2017, 139, 13770–13778.
14. G. M. J. Barca, A. T. B. Gilbert and P. M. W. Gill, J. Chem. Theory Comput., 2018, 14, 9–13.
15. P. F. Loos, M. Boggio-Pasqua, A. Scemama, M. Caffarel and D. Jacquemin, J. Chem. Theory Comput., 2019, 15, 1939–1956.
16. M. E. Casida and M. Huix-Rotllant, Annu. Rev. Phys. Chem., 2012, 63, 287–323.
17. A. Dreuw and M. Wormit, WIREs Computational Molecular Science, 2015, 5, 82–95.
18. O. Christiansen, H. Koch and P. Jørgensen, J. Chem. Phys., 1995, 103, 7429–7441.
19. J. Leitner, A. L. Dempwolff and A. Dreuw, J. Chem. Phys., 2022, 157, 184101.
20. P. H. P. Harbach, M. Wormit and A. Dreuw, J. Chem. Phys., 2014, 141, 17.
21. P. F. Loos, A. Scemama, A. Blondel, Y. Garniron, M. Caffarel and D. Jacquemin, J. Chem. Theory Comput., 2018, 14, 4360–4379.
22. H. Lischka, D. Nachtigallová, A. J. A. Aquino, P. G. Szalay, F. Plasser, F. B. C. MacHado and M. Barbatti, Chem. Rev., 2018, 118, 7293–7361.
23. D. Casanova and A. I. Krylov, Phys. Chem. Chem. Phys., 2020, 22, 4326–4342.
24. Z. Rinkevicius, O. Vahtras and H. Ågren, J. Chem. Phys., 2010, 133, 114104.
25. G. M. J. Barca, A. T. B. Gilbert and P. M. W. Gill, J. Chem. Theory Comput., 2018, 14, 1501–1509.
26. D. Hait and M. Head-Gordon, J. Phys. Chem. Lett., 2021, 12, 4517–4529.
27. C. D. Lin, Nucl. Instrum. Methods Phys. Res., 1985, 240, 572–576.
28. C. D. Lin, Adv. Atom. Mol. Phys., 1986, 22, 77–142.
29. N. T. Maitra, Annu. Rev. Phys. Chem., 2022, 73, 117–140.
30. J. P. Coe and M. J. Paterson, Mol. Phys., 2014, 112, 733–739.
31. F. Plasser, M. Wormit and A. Dreuw, J. Chem. Phys., 2014, 141, 024106.
32. F. Plasser, A. I. Krylov and A. Dreuw, WIREs Computational Molecular Science, 2022, 1–15.
33. A. I. Krylov, J. Chem. Phys., 2020, 153, 080901.
34. S. Matsika, X. Feng, A. V. Luzanov and A. I. Krylov, J. Phys. Chem. A, 2014, 118, 11943–11955.
35. F. Plasser, S. A. Bäppler, M. Wormit and A. Dreuw, J. Chem. Phys., 2014, 141, 24107.
36. M. Head-Gordon, A. M. Grana, D. Maurice and C. A. White, J. Phys. Chem., 1995, 99, 14261–14270.
37. M. Head-Gordon, Chem. Phys. Lett., 2003, 372, 508–511.
38. T. Minami and M. Nakano, J. Phys. Chem. Lett., 2012, 3, 145–150.
39. F. Plasser, S. A. Mewes, A. Dreuw and L. González, J. Chem. Theory Comput., 2017, 13, 5343–5353.
40. A. V. Luzanov, A. A. Sukhorukov and V. É. Umanskii, Theor. Exp. Chem., 1976, 10, 354–361.
41. R. L. Martin, J. Chem. Phys., 2003, 118, 4775–4777.
42. F. Plasser and H. Lischka, J. Chem. Theory Comput., 2012, 8, 2777–2789.
43. F. Plasser, J. Chem. Phys., 2016, 144, 194107.
44. M. Nakano, H. Fukui, T. Minami, K. Yoneda, Y. Shigeta, R. Kishi, B. Champagne, E. Botek, T. Kubo, K. Ohta and K. Kamada, Theor. Chem. Acc., 2011, 130, 711–724.
45. I. Mayer, Chem. Phys. Lett., 2007, 443, 420–425.
46. F. Bettanin, L. F. A. Ferrão, M. Pinheiro, A. J. A. Aquino, H. Lischka, F. B. C. Machado and D. Nachtigallova, J. Chem. Theory Comput., 2017, 13, 4297–4306.
47. P. Kimber and F. Plasser, Classification and Analysis of Molecular Excited States, in Comprehensive Computational Chemistry, ed. R. Boyd and M. Yanez, 2023 (in print).
48. P. Kimber and F. Plasser, Phys. Chem. Chem. Phys., 2020, 22, 6058–6080.
49. A. T. B. Gilbert, N. A. Besley and P. M. W. Gill, J. Phys. Chem. A, 2008, 112, 13164–13171.
50. X. K. Chen, D. Kim and J. L. Brédas, Acc. Chem. Res., 2018, 51, 2215–2224.
51. W. Barford, Phys. Rev. B, 2022, 106, 035201.
52. S. Bai and M. Barbatti, J. Chem. Theory Comput., 2019, 15, 1503–1513.
53. T. Yanai, D. P. Tew and N. C. Handy, Chem. Phys. Lett., 2004, 393, 51–57.
54. T. H. Dunning, J. Chem. Phys., 1989, 90, 1007–1023.
55. S. Grimme, J. Comput. Chem., 2006, 27, 1787–1799.
56. M. J. Frisch, G. W. Trucks, H. B. Schlegel, G. E. Scuseria, M. A. Robb, J. R. Cheeseman, G. Scalmani, V. Barone, G. A. Petersson, H. Nakatsuji, X. Li, M. Caricato, A. V. Marenich, J. Bloino, B. G. Janesko, R. Gomperts, B. Mennucci, H. P. Hratchian, J. V. Ortiz, A. F. Izmaylov, J. L. Sonnenberg, D. Williams-Young, F. Ding, F. Lipparini, F. Egidi, J. Goings, B. Peng, A. Petrone, T. Henderson, D. Ranasinghe, V. G. Zakrzewski, J. Gao, N. Rega, G. Zheng, W. Liang, M. Hada, M. Ehara, K. Toyota, R. Fukuda, J. Hasegawa, M. Ishida, T. Nakajima, Y. Honda, O. Kitao, H. Nakai, T. Vreven, K. Throssell Jr, J. A. Montgomery, J. E. Peralta, F. Ogliaro, M. J. Bearpark, J. J. Heyd, E. N. Brothers, K. N. Kudin, V. N. Staroverov, T. A. Keith, R. Kobayashi, J. Normand, K. Raghavachari, A. P. Rendell, J. C. Burant, S. S. Iyengar, J. Tomasi, M. Cossi, J. M. Millam, M. Klene, C. Adamo, R. Cammi, J. W. Ochterski, R. L. Martin, K. Morokuma, O. Farkas, J. B. Foresman and D. J. Fox, Gaussian 16, Revision A. 03, Gaussian, Inc., Wallingford, CT, 2016.
57. S. Grimme and M. Waletzke, J. Chem. Phys., 1999, 111, 5645–5655.
58. F. Weigend and R. Ahlrichs, Phys. Chem. Chem. Phys., 2005, 7, 3297–3305.
59. A. D. Becke, J. Chem. Phys., 1993, 98, 1372–1377.
60. S. G. Balasubramani, G. P. Chen, S. Coriani, M. Diedenhofen, M. S. Frank, Y. J. Franzke, F. Furche, R. Grotjahn, M. E. Harding, C. Hättig, A. Hellweg, B. Helmich-Paris, C. Holzer, U. Huniar, M. Kaupp, A. Marefat Khah, S. Karbalaei Khani, T. Müller, F. Mack, B. D. Nguyen, S. M. Parker, E. Perlt, D. Rappoport, K. Reiter, S. Roy, M. Rückert, G. Schmitz, M. Sierka, E. Tapavicza, D. P. Tew, C. van Wüllen, V. K. Voora, F. Weigend, A. Wodyński and J. M. Yu, J. Chem. Phys., 2020, 152, 184107.
61. C. M. Marian, A. Heil and M. Kleinschmidt, Wiley Interdiscip. Rev.: Comput. Mol. Sci., 2019, 9, 1–31.
62. W. Fuß, Y. Haas and S. Zilberg, Chem. Phys., 2000, 259, 273–295.
63. M. R. Silva-Junior, M. Schreiber, S. P. A. Sauer and W. Thiel, J. Chem. Phys., 2008, 129, 104103.
64. A. D. Becke, Phys. Rev. A, 1988, 38, 3098.
65. C. Lee, W. Yang and R. G. Parr, Phys. Rev. B: Condens. Matter Mater. Phys., 1988, 37, 785.
66. R. Ditchfield, W. J. Hehre and J. A. Pople, J. Chem. Phys., 1971, 54, 720–723.
67. Y. Shao, Z. Gan, E. Epifanovsky, A. T. B. Gilbert, M. Wormit, J. Kussmann, A. W. Lange, A. Behn, J. Deng, X. Feng, D. Ghosh, M. Goldey, P. R. Horn, L. D. Jacobson, I. Kaliman, R. Z. Khaliullin, T. Kus, A. Landau, J. Liu, E. I. Proynov, Y. M. Rhee, R. M. Richard, M. A. Rohrdanz, R. P. Steele, E. J. Sundstrom, H. L. Woodcock, P. M. Zimmerman, D. Zuev, B. Albrecht, E. Alguire, B. Austin, G. J. O. Beran, Y. A. Bernard, E. Berquist, K. Brandhorst, K. B. Bravaya, S. T. Brown, D. Casanova, C. M. Chang, Y. Chen, S. H. Chien, K. D. Closser, D. L. Crittenden, M. Diedenhofen, R. A. Distasio, H. Do, A. D. Dutoi, R. G. Edgar, S. Fatehi, L. Fusti-Molnar, A. Ghysels, A. Golubeva-Zadorozhnaya, J. Gomes, M. W. D. Hanson-Heine, P. H. P. Harbach, A. W. Hauser, E. G. Hohenstein, Z. C. Holden, T. C. Jagau, H. Ji, B. Kaduk, K. Khistyaev, J. Kim, J. Kim, R. A. King, P. Klunzinger, D. Kosenkov, T. Kowalczyk, C. M. Krauter, K. U. Lao, A. D. Laurent, K. V. Lawler, S. V. Levchenko, C. Y. Lin, F. Liu, E. Livshits, R. C. Lochan, A. Luenser, P. Manohar, S. F. Manzer, S. P. Mao, N. Mardirossian, A. V. Marenich, S. A. Maurer, N. J. Mayhall, E. Neuscamman, C. M. Oana, R. Olivares-Amaya, D. P. Oneill, J. A. Parkhill, T. M. Perrine, R. Peverati, A. Prociuk, D. R. Rehn, E. Rosta, N. J. Russ, S. M. Sharada, S. Sharma, D. W. Small, A. Sodt, T. Stein, D. Stück, Y. C. Su, A. J. W. Thom, T. Tsuchimochi, V. Vanovschi, L. Vogt, O. Vydrov, T. Wang, M. A. Watson, J. Wenzel, A. White, C. F. Williams, J. Yang, S. Yeganeh, S. R. Yost, Z. Q. You, I. Y. Zhang, X. Zhang, Y. Zhao, B. R. Brooks, G. K. L. Chan, D. M. Chipman, C. J. Cramer, W. A. Goddard, M. S. Gordon, W. J. Hehre, A. Klamt, H. F. Schaefer, M. W. Schmidt, C. D. Sherrill, D. G. Truhlar, A. Warshel, X. Xu, A. Aspuru-Guzik, R. Baer, A. T. Bell, N. A. Besley, J. Da Chai, A. Dreuw, B. D. Dunietz, T. R. Furlani, S. R. Gwaltney, C. P. Hsu, Y. Jung, J. Kong, D. S. Lambrecht, W. Liang, C. Ochsenfeld, V. A. Rassolov, L. V. Slipchenko, J. E. Subotnik, T. Van Voorhis, J. M. Herbert, A. I. Krylov, P. M. W. Gill and M. Head-Gordon, Mol. Phys., 2015, 113, 184–215.
68. H. Lischka, T. Müller, P. G. Szalay, I. Shavitt, R. M. Pitzer and R. Shepard, Wiley Interdiscip. Rev.: Comput. Mol. Sci., 2011, 1, 191–199.
69. H. Lischka, R. Shepard, R. M. Pitzer, I. Shavitt, M. Dallos, T. Müller, P. G. Szalay, M. Seth, G. S. Kedziora, S. Yabushita and Z. Zhang, Phys. Chem. Chem. Phys., 2001, 3, 664–673.
70. H. Lischka, R. Shepard, I. Shavitt, R. M. Pitzer, M. Dallos, Th. Müller, P. G. Szalay, F. B. Brown, R. Ahlrichs, H. J. Böhm, A. Chang, D. C. Comeau, R. Gdanitz, H. Dachsel, C. Ehrhardt, M. Ernzerhof, P. Höchtl, S. Irle, G. Kedziora, T. Kovar, V. Parasuk, M. J. M. Pepper, P. Scharf, H. Schiffer, M. Schindler, M. Schüler, M. Seth, E. A. Stahlberg, J.-G. Zhao, S. Yabushita, Z. Zhang, M. Barbatti, S. Matsika, M. Schuurmann, D. R. Yarkony, S. R. Brozell, E. V. Beck, J.-P. Blaudeau, M. Ruckenbauer, B. Sellner, F. Plasser, J. J. Szymczak, R. F. K. Spada and A. Das, COLUMBUS, an ab initio electronic structure program, release 7.0, 2019.
71. J. A. Pople, R. Seeger and R. Krishnan, Int. J. Quantum Chem., 1977, 12, 149–163.
72. J. Schirmer, Phys. Rev. A, 1982, 26, 2395–2416.
73. F. Plasser, J. Chem. Phys., 2020, 152, 084108.
74. R. J. Cave and E. R. Davidson, Chem. Phys. Lett., 1988, 148, 190–196.
75. C. M. Marian and N. Gilka, J. Chem. Theory Comput., 2008, 4, 1501–1515.
76. M. T. do Casal, J. M. Toldo, F. Plasser and M. Barbatti, Phys. Chem. Chem. Phys., 2022, 23279–23288.
77. H. S. Yu, S. L. Li and D. G. Truhlar, J. Chem. Phys., 2016, 145, 130901.
78. F. Bernardi, M. Olivucci and M. A. Robb, Acc. Chem. Res., 1990, 23, 405–412.
79. N. D. Epiotis and S. Shaik, J. Am. Chem. Soc., 1978, 100, 9–17.
80. S. Ghosh, S. Shankar, D. S. Philips and A. Ajayaghosh, Mater. Today Chem., 2020, 16, 100242.
81. N. Luo, G. Zhang and Z. Liu, Org. Chem. Front., 2021, 8, 4560–4581.
82. A. Shukla, S. K. M. McGregor, R. Wawrzinek, S. Saggar, E. G. Moore, S. C. Lo and E. B. Namdas, Adv. Funct. Mater., 2021, 31, 1–11.
83. L. Shen, Z. Tang, X. Wang, H. Liu, Y. Chen and X. Li, Phys. Chem. Chem. Phys., 2018, 20, 22997–23006.
84. I. Papadopoulos, M. J. Álvaro-Martins, D. Molina, P. M. McCosker, P. A. Keller, T. Clark, Á. Sastre-Santos and D. M. Guldi, Adv. Energy Mater., 2020, 10, 1–10.
85. S. Klaiman and L. S. Cederbaum, J. Chem. Phys., 2014, 141, 194102.

---

Footnote

† Electronic supplementary information (ESI) available: An alternative discussion of Ω (S1); CO₂ example with 1TDM norm larger than 1 (S2); generalization of the excitation number (S3); DFT/MRCI parametrization comparison (S4, S5); formaldehyde dimer orbitals and discussion of doubly excited states within ADC(3) (S6); density descriptors of the polyenes' series (S7); density descriptors of ethylene's cycloaddition (S8); density descriptors of tetracene, diketopyrrolopyrrole derivatives and s-tetrazine (S9). See DOI: https://doi.org/10.1039/d2sc06990c

---