On the notion of strong correlation in electronic structure theory

Brad Ganoe and James Shee*
Department of Chemistry, Rice University, Houston, TX 77005, USA. E-mail: james.shee@rice.edu

First published on 10th April 2024

---

---

1 Introduction

There is no literature consensus about what exactly it means for a many-electron system to be strongly correlated.^1–4 Leading notions include the following:

1. The wavefunction cannot qualitatively be described by a single Slater determinant, and is “multiconfigurational”. A chemical example is when two states of different character – e.g., ionic and covalent configurations that can be non-orthogonal along much of a dissociation coordinate – are of competing significance (energetically nearly degenerate).⁵ Another example from physics of such an intermediate regime involves fluctuations in materials such as transition metal oxides close to a Mott metal–insulator transition⁶ or f-block compounds,⁷ in which electrons “hesitate” between delocalized (wave-like) and localized (particle-like) character.⁴

2. The two-body interaction is large relative to one-body terms in the Hamiltonian. This has its origin in the Hubbard model:

	(1)

specifically the parameter regime in which U ≫ t. The ground-state, when half-filled, involves strongly localized electrons coupled antiferromagnetically; doping gives rise to more exotic orders.⁸

3. Electrons are substantially “entangled”.⁹ From an information theory perspective, this can refer to entanglement entropy, e.g. Shannon or Renyi formulations.^10,11 Concepts such as 1- and 2-orbital entanglement entropies have also been proposed.^12,13

4. The overlap of the Hartree–Fock (HF) state or the largest-weighted configuration in a multi-determinant expansion with the true ground-state is small.

5. The eigenvalues of the one-electron reduced density matrix – the natural orbital occupation numbers (NOONs) – are non-integer, deviating substantially from 0 and 1 (or 2 for restricted models).¹⁴

6. The amplitudes of Coupled Cluster (CC) theory, especially T₁ and T₂ of CCSD,^15–17 become large. As the connected portions of triples amplitudes and higher become significant, the corresponding residuals remain large such that no finite truncation is well posed.¹⁸

7. Two-body spin or particle correlation functions, e.g. 〈S₀·S_i〉 or 〈_0i〉 in a real-space site basis¹⁹ that do not decay exponentially. This idea of long correlation lengths is familiar to statistical physics.

8. The energies from Møller–Plesset perturbation theory through low orders with the generalized Hartree–Fock (GHF) reference state qualitatively differ from that of Full Configuration Interaction (FCI).²⁰ The causal relationship between strong correlation and the breakdown of many-body perturbation theory is relatively well recognized in the physics literature.⁴

The failure of approximate density functional theory (DFT) is another common characterization. Alternatively, a colleague recently exclaimed of strong correlation that “you know it when you see it”!

While Small and Head-Gordon have defined strong correlation in terms of the breakdown of a specific perturbation theory, there is general agreement that a consequence of strong correlation is that feasibly-scaling electronic structure models such as low orders of single-reference perturbation or coupled cluster theories and approximate density functionals (especially those from the higher rungs of Jacob’s ladder) break down.^4,21,22 Thus it is of practical value to computational scientists to reliably diagnose when a system shows symptoms of strong correlation. For example, a recent work has delineated correlation regimes inside of which CCSD(T) is expected to produce ionization potentials of transition metal complexes to within kcal mol⁻¹ accuracy, and outside of which CCSD(T) is expected to qualitatively fail (no matter which orbitals are used in the reference determinant).²³ Such predictive classification relied on a composite diagnostic involving both spin-symmetry breaking in unrestricted CCSD wavefunctions and the largest CI coefficient in selected CI wavefunctions (in the basis of natural orbitals).

In contrast, the present work is motivated by a more fundamental question concerning the physical nature/characteristics of a strongly correlated electronic system: given a Hamiltonian and an exact eigenstate, what quantity or property can be associated with the regime of strong correlation? Furthermore, what is implied about physical reality when this quantity is large? We require such a quantity to be rigorously defined, physically interpretable, and generally meaningful for both molecules and quantum lattice models (other more technical properties are discussed in the final section). Note that this type of inquiry is distinct from classifying correlation regimes based on which approximate methods are broken. These two approaches to elucidating strong correlation are often conflated.

Another motivating goal of this work is to build a bridge between the chemistry and physics communities, with regard to perspectives and even language surrounding electron correlation. As frequently highlighted by R. Hoffmann,²⁴ a paradigmatic example of a concept that is fundamental, chemically-consequential, and beloved by chemists, but whose definition proves more evasive the more one tries to formalise it, is aromaticity and, more interestingly in our view, anti-aromaticity. According to one of the oldest and most widely-used definitions, aromaticity arises in cyclic compounds with 4N + 2 electrons in a conjugated π-orbital system, while antiaromatic compounds have 4N electrons. Baird’s rule further states that a system with an aromatic (antiaromatic) ground-state will have anti-aromatic (aromatic) low-lying vertical excited states.²⁵ Admittedly, such definitions and counting rules are not mathematically rigorous, yet empirically these concepts are seen to hold. For example, evidence of aromaticity includes literal aromaticity (smell), unusual stability/inertness to chemical reactivity, and magnetic shielding effects due to ring currents; while antiaromaticity implies an unusually unstable state that tends to rapidly react or distort geometrically.²⁶ A beautiful and unexpected result of this work is that antiaromatic molecules share a formal quantum-mechanical footprint with seemingly disparate phenomena ranging from bond-stretching and super exchange physics to antiferromagnetic order in some polynuclear transition metal catalysts and cuprate materials. Indeed, we venture to propose that the footprint is due to the elephant in the room – strong correlation.

2 Theory

2.1 Critical comments on metrics of strong correlation

To narrow the scope of our study, we briefly comment on the conceptions of strong correlation listed above that will not be examined below in our numerical tests.

Regarding 1: the expansion coefficients of a multiconfigurational wavefunction depend on the choice of single-particle basis; ref. 2 shows a particularly stark example in which delocalized vs. localized orbitals imply qualitatively different assessments of multiconfigurational character. It has also been shown that most mono-nuclear transition metal coordination complexes are qualitatively well described by a single determinant plus an appropriate level of dynamic correlation,³ and that exchange-coupled polynuclear transition metal compounds can often be described by a single spin-adapted configuration state function.^2,27

Regarding 3: in our view, the characterization of a system as strongly correlated should not depend on the partitioning of a Hilbert space into distinguishable subspaces, as is often required to compute metrics such as entanglement entropy. As in 1, orbital entanglement entropy metrics have been shown to be strongly dependent on the choice of single-particle basis employed.¹²

Regarding 4: while the overlap of a given wavefunction with FCI, 〈Φ|Ψ_FCI〉, is useful for finite systems such as molecules (as it reports on the percentage of the FCI Hilbert space that is relevant to that wavefunction), the focus on the case where |Φ〉 is the mean-field state is not relevant in the thermodynamic limit in light of the well-known “orthogonality catastrophe” initially described by Van Vleck:²⁸ given N infinitely separated He atoms, and assuming that 〈Φ_i|Ψ_FCI〉 = 1 − ε for each of the He atoms, then the overlap for the composite system (1 − ε)^N will decay exponentially with N, even for weakly correlated systems.

Regarding 5, Mazziotti⁵³ has pointed out that the 1-RDM, and thus its eigenvalues, cannot capture spin-correlations which need not decay with the 1/R interaction (and therefore can be non-zero even when the correlation energy is zero). We note also that the NOON spectrum will be identical in, e.g., H₂ in the dissociation limit for both singlet and triplet states; this contradicts the qualitative difference in the wavefunctions: and |Ψ_T〉 = |σ,σ*〉. If two subspaces are assumed distinguishable, one can write: (L = left, R = right), which exhibits the famous “spooky action at a distance”, also known as spin correlations. Although the UHF wave function recovers the correct energy with a spin-symmetry broken single determinant, it does so for reasons other than having the exact many-particle wave function. Despite the ability of fractional NOONs to signal a multiconfigurational wave function, we are interested in quantities that can also capture spin-correlations/entanglements even in the dissociation limit. In what follows, we also highlight an example where the NOONs fail to distinguish two geometric configurations of H₄ with different electronic structures.

Regarding 6, it is possible that for a weakly correlated state, the T₁ norm from a CCSD wave function can be large; this would simply reflect a reference determinant with suboptimal orbitals. While a large T₂ norm from CCSD can possibly be a symptom of strong correlation, in this work we are interested in how strong correlation is manifested in the true (FCI) wave function.

Regarding 7, the definition of correlation functions with respect to real-space sites is natural for lattice models and specific systems such as well-localized μ-oxo-bridged transition metals. However, such a construction is more challenging to systematically define for, e.g., Be₂ at short interatomic separation, and (more obviously) atoms like cerium with complicated low-lying states.²⁹

Regarding 8: as mentioned in our response to point 5, in this work we are interested in quantities and properties of the exact wave function or quantities derived from the 2-RDM which can, in principle, be exactly mapped to the exact wave function. Although the size-consistency of GHF and subsequent low-order perturbation theories indeed enables qualitatively correct energies at that level of theory to signal correlation regimes in which energies can be predicted tractably, this metric is out of the scope of this work since we are interested in physical properties of exact eigenstates (which, like spin-entanglements, can persist even when the correlation energy is zero). Additionally, there are many ways to break low-order single-reference perturbation theory that are not (at least directly) physically meaningful, including the idea of “artificial” symmetry breaking in the HF state.^3,30 In an effort to report on physically-meaningful, i.e. “essential”, symmetry breaking, in this work we will report numerical values of 〈S²〉 from an unrestricted CCSD wavefunction (vide infra). Finally, we note that, GHF aside, redefining the 0th order Hamiltonian can often lead to qualitatively correct energies from low-order perturbation theories, as in CASPT2 or NEVPT2.

2.2 Reduced density matrices and cumulants

The N-body density matrix, Γ = |Φ〉〈Φ|, can generate the reduced-density matrices:

	(2)

The 1-RDM,

1Γpq = 〈Φ|c†pcq|Φ〉	(3)

whose eigenvectors and eigenvalues define the natural orbitals and their occupation numbers (NOONs), is sufficient to describe all one-particle properties of the system, including the mean-field energy. Similarly, the physical information of pairwise interacting systems (which includes all molecules and materials) is contained within the 2-RDM,

2Γpqrs = 〈Φ|c†pc†qcrcs|Φ〉	(4)

As guaranteed by Rosina’s theorem, there exists a bijection from the exact wavefunction to ²Γ of a non-degenerate ground state. Therefore, for a system with 2-body terms and fewer in the Hamiltonian, the 2-RDM contains the information required to obtain all higher-order RDMs and the exact wave function.³¹ This property is not true of the 1-RDM nor the electron density despite both being able to yield the exact energies, in principle. A recent extension of the theorem is also valid for degenerate states differentiable by a two-particle operator.³² These features have made 2-RDM-based methodologies very compelling, especially with regard to strong correlation. Seminal works^33,34 by Yang and Coleman investigating the properties of the p-RDMs of fermionic systems have also demonstrated that large eigenvalues of the 2-RDM are connected to “off-diagonal long range order” (ODLRO), which has been associated with phenomena such as superconductivity due to electronic degrees of freedom and exciton condensation,^35,36 while diagonal long range order (DLRO) makes an appearance in magnetism.³⁷ In analogy with the 1-RDM and the natural orbitals, the eigenvectors of the 2-RDM correspond to the natural geminals of the system.

The reduced density matrices can be related to the p-matrices:³⁸

	(5)

which are normalized to one, and thus, when contracted with the position or two positions describe the probability or joint probability densities of finding one or two particles at particular positions in space (occasionally called the distribution functions). From these quantities, the correlation coefficient, τ, can be obtained, which is the ratio between the covariance with respect to the position vectors of two particles, and , and the scalar variance with respect to the position of one particle. Given the indistinguishability and choice of origin at the center of mass, we can write the correlation coefficient as

	(6)

The correlation coefficient can reveal when variables are statistically uncorrelated. In fact, Kutzelnigg observed³⁹ that the initial impetus for Wigner and Seitz’s coining of “electron correlation” was the observation of a lack of statistical correlation between the positions of parallel spin electrons within the Hartree Product.^40–42 In the case of Hartree–Fock, the Slater determinant form of the wave function adds Fermi correlation (HF exchange), and separately, the fact that the joint-density function is equivalent to the product of the density functions – consistent with idempotency – implies statistical independence.³⁹

Beyond the mean-field approximation, the introduction of non-independence between the particles may be evaluated in a statistical sense by a moment-generating function of stochastic variables to evaluate the cumulants of the expansion;⁴³ a mapping by way of the antisymmeterized logarithm⁴⁴ yields the cumulants of the ^pΓ, ^pΔ, a quantity which is generated from the connected portions of the respective ^pΓ from the generalized Wick’s theorem. The properties of these terms, especially up through the ⁴Δ with their role in establishing N-representability criteria for the reconstruction of the 2-RDM,³¹ have been well reviewed in the literature and have been explicitly connected to other quantities such as CC amplitudes.^44–52

As in traditional statistical theory, in which the first-order cumulant is identical to the mean, the ¹Δ is identical to the ¹Γ and yields no additional unique information over the above quantity. As the lowest unique quantity, of particular interest among these is the cumulant of the 2-RDM, which is defined in terms of the 2-RDM and the 1-RDM (¹Γ):

2Δijkl = 2Γijkl − 1Γik1Γjl + 1Γil1Γjk	(7)

The ²Δ contains the first instance of statistically non-independent behavior and is both size-extensive and mappable uniquely to the many-body wavefunction.⁵³ By setting ²Δ to zero, one obtains the usual idempotency requirement. The cumulant is a measure of the statistical dependence of particles; large values of the trace or square norm of ²Δ is what we would like to associate with the term “strong correlation”.

The relation between the spin components of the two systems is shown in the following:

	(8)

and

	(9)

when the ¹Γ and ²Γ are traced to N and N(N − 1), respectively, and may be calculated explicitly as the non-terminating block components of the full 2-RDM or from the spin-traced, sometimes called charge, density matrix.^54,55 The cumulants of the density have several properties that make them useful in the construction of a metric for strong correlation. For example, NOONs have been proposed as a way to identify active space orbitals, since partial occupation of orbitals is considered an identifier of multiconfigurational character. This strategy first appeared with regard to the broken spin-symmetry mean-field level of calculation.^14,56,57 However, as first noted by Juhász and Mazziotti,^53,58 the ¹Δ contains insufficient data to describe the entire spectrum of possible strong correlation effects and offers no unique signature for infinitely separated hydrogen atoms with 〈S²〉 = 0. Interestingly, ¹Δ is able to yield the diagonal elements of ²Δ, and hence its trace and statistical dependence of pairwise orbitals, although problems were encountered for mixed/ensemble states.⁵⁹ Alternatively, the squared Frobenius norm of the ²Δ was proposed by Mazziotti as a size-extensive quantity invariant to unitary transformation of the orbitals which can capture the spin-entanglements of atom centers at separation, while the trace and eigenvalues of ²Δ may capture the rise of ODLRO.⁶⁰ Further studies quantifying strong correlation utilizing the ²Δ have been proposed, including metrics based on a two-electron two-fragment model in a minimal basis⁶¹ and other terms derivable from the ²Δ such as the cross-trace and infinity-norm.^62,63

2.3 Additional metrics

The physical basis of spin-symmetry breaking in approximately (weakly) correlated single-reference wave functions has been described elsewhere.^3,23,30,64 At present, we will highlight a few salient points.

• For the specific case of H₂ with a single Slater determinant wave function,⁶⁵

	(10)

where n_LUNO is the NOON corresponding to the lowest unoccupied molecular orbital (LUMO) from HF theory. 〈S²〉 can, at least in some cases, encode information related to relevant eigenvalues of the 1-RDM.

• For finite molecules, spin is a well-defined quantum number, and yet spin-symmetry breaking in a single-determinant wave function can recover more accurate energies in certain cases than the spin-pure counterpart. Let us consider again the case of stretched H₂ (which is isoelectronic to chemically relevant organic diradicaloids as well as many di-Cu(II) active sites). Here, 〈S²〉_UHF goes from zero before the Coulson–Fischer point to one at complete dissociation; in the latter regime the UHF wave function can be expanded as an equal mixture of pure singlet and triplet states. The triplet state in fact has a more accurate electron density distribution than that of the closed-shell spin-pure determinant (although the spin-density distribution is also qualitatively wrong). Thus, the spin-contamination enables a single-determinant to capture (to some extent) an open-shell singlet electronic structure that is required of biradicals, and can be viewed as physically motivated, or “essential” according to Head-Gordon and coworkers.

• Strong spin correlation implies spin-symmetry breaking, though the converse is not true.²⁰ In what follows, we will compute 〈S²〉_UCCSD, as UCCSD does not contain free parameters and captures to a large extent what is frequently referred to as weak correlation along with some parts of strong correlation. Ref. 23 corroborates our hypothesis that the spin-symmetry breaking in UHF that remains after UCCSD can be attributed to strong spin correlations that are known to break low-order coupled cluster theory.

• Strong spin correlation is one of possibly many components of the as-yet mysterious quantity of interest: strong correlation. Thus, 〈S²〉_UCCSD is, by construction, incomplete. In addition, there are a small number of known cases (e.g., two-electron systems including stretched H₂) where 〈S²〉_UCCSD will always be spin-pure, despite qualitatively different correlation regimes being traversed.

In addition to quantities derived from the 2-RDM and its cumulant, and 〈S²〉_UCCSD, we would like to explore an additional metric defined here. The Born–Oppenheimer Hamiltonian can be written as Ĥ = + _en + _ee + _nn = Ĥ_1e + _ee + E_nn, where is the kinetic energy operator, (V_pqrs are two-electron integrals), and E_nn is a constant energy shift. The first candidate metric is motivated by the large U/t limit of the Hubbard model, which is widely accepted by physicists to exhibit strong correlation, and also the Virial Theorem.² One possible analog for molecules is the ratio:

	(11)

The denominator can also reasonably be replaced by 〈Ψ|Ĥ_1e|Ψ〉.

3 Numerical studies

We start by considering the half-filled Hubbard model at various U/t ratios. To test the generality of the behavior of various computed quantities, we then switch to a number of molecular systems including ground-state dissociations and selected low-lying excited states. All calculations were done with PySCF^66,67 employing double-ζ basis sets (cc-pVDZ for H₂, HHeH, and H₄; 6-31G otherwise). All quantities except for 〈S²〉_UCCSD reflect an exact (FCI-level) underlying wavefunction. Given there are only three unique spin blocks in the systems for this work, we shall use ²Γ_αβ, ²Γ_ββ, ²Γ_αβ to denote the subblocks of the ²Γ where appropriate.

3.1 Hubbard model

The energy ratio, α (eqn (11)), the maximum eigenvalues (γ) of ²Γ_αα and ²Γ_αβ, the maximum eigenvalues (δ) of ²Δ_αα and ²Δ_αβ, the trace of ²Δ, and the square norm of ²Δ with respect to the exact ground-state wave function of the 6-site half-filled Hubbard model are shown in Fig. 1 as a function of U with t = 1. As the ratio U/t increases, we expect the electronic structure to change from the non-interacting limit with freely delocalized electrons to a strongly interacting regime in which electrons are localized on each site, exhibiting antiferromagnetic order. All of the ²Δ quantities are seen to monotonically increase in magnitude, as expected. The ratio α is always negative due to the denominator, and (almost by construction) increases in magnitude as the U/t ratio increases. The largest eigenvalues of ²Γ are found to decrease, rather surprisingly, though this will be discussed in detail in the final section.

3.2 Molecules

Hydrogen lattices are prototypical models of correlated electronic systems, pushing strongly interacting Hubbard-like models one step further into the chemical realm. Although these simple s-orbital systems are not themselves representative of real-world chemistry, such systems can often be solved exactly, and a rather surprising amount of correlated-electron physics can be uncovered⁶⁸ when, e.g., the lattice parameters are varied.

We begin our analysis with the well-studied dihydrogen molecule. At equilibrium, a single large-weighted configuration corresponds to the HF solution, with electron fluctuations into the doubly-excited space accounting for the bulk of the correlation energy. Upon stretching, this molecular orbital picture becomes inadequate as two equally-weighted configurations come to dominate the ground-state wave function. Despite the correlation energy tending towards zero at dissociation (in fact, the singlet and triplet states are isoenergetic, and the UHF energy is exact), the entanglement characteristic of the singlet state does not follow such a trend as mentioned above. As shown in Fig. 2, this is well reflected by many of the metrics presently under consideration: ²γ_αβ only has its connected portions survive, ²δ_αβ, captured in both the trace and squared Frobenius norm of ²Δ. As must follow from the exactness of UCCSD, 〈S²〉_UCCSD here does not reflect any spin entanglement. ²Γ_αα takes on no values while the elements of ¹Γ_αα contribute to the increasing squared Frobenius norm of ²Δ as it approaches at dissociation, the previously derived limit for this system.⁵⁴ There is a notable non-monotonicity in the ratio α, which declines as the bond length is stretched beyond 1 Å. This ratio does not appear to be applicable beyond models like the Hubbard model, possibly due to the delicate balance between phenomena such as node-induced confinement, orbital contractions, and polarization effects that can modulate the kinetic energy in molecular systems.⁶⁹

H₄ is a challenging molecule that has been the subject of many studies in a variety of geometric configurations.^70,71 A previous report has shown that the NOONs are identical for square and rhombic geometric configurations, yet notable differences occur in the underlying electronic structures.⁷² Indeed, that report motivated us to look beyond quantities related to ¹Γ. We investigate two reaction coordinates – the P4 rectangular D_2h dissociation into hydrogen dimers and the S4 square D_4h dissociation into hydrogen atoms. In the former, the HF determinant and the lowest energy doubly-excited configuration contribute significantly as the system crosses into the square geometry. There is a region of both spin and spatial frustration as we approach a crossing of the S₀ and S₁ surfaces, alongside a near-degenerate T₁.⁷⁰ Fig. 3 (left) shows that in the singlet ground-state, the ²Δ quantities correctly flag the electronically challenging region surrounding the square geometry with lattice parameter approaching 1.06 Å (2 Bohr). Spin-symmetry breaking in UCCSD shows similar behavior (though with a slightly wider peak), and we note that at 1.05 Å the UCCSD wave function is in equal parts of singlet and triplet character (averaging out to 〈S²〉 = 1). This is consistent with the open-shell singlet (two electrons in two degenerate orbitals) picture implied by Huckel/MO theory. We again find that the largest eigenvalue of ²Γ decreases in the square regime. In contrast, these metrics corresponding to the quintet state (shown in the right panel of Fig. 3) are all flat. As was the case for H₂, the energy ratio α does not seem to contain interpretable information here.

In the case of the S4 dissociation to separate hydrogen atoms (Fig. 4), the ²Δ quantities reveal that, in contrast to the P4 coordinate, S4 does not involve a closed-shell region aptly described by mean-field theory (invoking the statistical sense with which we interpret correlation strength, indeed we can describe such a state as weakly correlated). As was the case in P4, the S4 quintet state (Fig. 4 right panel) is always weakly correlated, with all of our investigated metrics – except α – flat-lining. Interestingly, in the singlet case the squared Frobenius norm and the magnitude of the trace of ²Δ are the only metrics that increase with lattice parameter; 〈S²〉_UCCSD decreases (and appears to plateau at a value less than the maximum value of unity from the P4 coordinate), as do the maximum eigenvalues of the cumulant, albeit ever so slightly.

As a brief aside, it has been noted that the S4 model of equivalent side lengths, R_S4, has an identical spectrum of NOONs as some trapezoidal configurations in the C_2v geometry with one side length t₂ and the other parallel side consisting of infinitely separated hydrogen atoms.⁷³ In the cases of R_S4 = 2.0 Å and t₂ = 1.793 Å, although neither can be considered weakly correlated, the S4 example is known to exhibit additional correlation effects. Metrics based on the NOON spectrum of ¹Γ are unable to detect any difference (in other words, the ¹γ take on identical values). Despite this, differences in the 2-RDM yield differing values derived from the ²Δ in these examples: the trace of the ²Δ evaluates to −0.844 and −0.344 in the S4 and trapezoidal geometries, respectively, while the squared Frobenius norms are 1.37 and 0.52, respectively. Thus, in contrast to the NOONs these ²Δ-based quantities are able to distinguish between two very different regimes of electron correlation.

The six-membered hydrogen ring has been studied extensively, in part because it is isoelectronic to the π system of an aromatic molecule such as benzene (although of course H₆ involves s-orbitals while the π system of benzene involves out-of-plane p orbitals). While the always-square S4 state of the 4-membered hydrogen ring is generally expected to have a strongly-correlated character across the board, the expected weakly-correlated regime that corresponds to aromaticity in the six-membered ring is well reflected by the ²Δ and 〈S²〉_UCCSD metrics in the singlet ground-state before and in the neighborhood of the equilibrium geometry of 1 Å (Fig. 5 left panel). When the lattice parameter is stretched beyond 1 Å the largest eigenvalues, magnitude of the trace, and squared Frobenius norm of ²Δ increase, as does the spin-symmetry breaking from UCCSD. As before, these metrics are largely flat for the highest spin state – in this case a septet (Fig. 5 right panel) – confirming weak correlation. At the end of this section we will return to the low-lying excited states of this model.

The linear symmetric dissociation of HHeH functions as a simple model of superexchange with tunable biradical character, in which the paramagnetic hydrogens are coupled through a diamagnetic helium intermediate. The parallel or anti-parallel coupling of the spins localized on the H atoms results in two prominent low-lying states, ³A_σu⁺ and ¹X_σg⁺. While even UHF may recapture the singlet–triplet energy splitting (and thus the exchange coupling constant, J) at infinite distance, the wave function itself must include a large weight on the lowest-lying doubly-excited configuration. The moderate to large biradicaloid character after 1.25 Å results in the breakdown of CCSD and CCSD(T).⁷⁴ Fig. 6 shows that the coupling of the hydrogens mediated by the helium center leads to a similar trend as found in stretched H₂, though in the singlet HHeH case ²γ_αα and ²γ_ββ arise from the additional electrons. In addition, the ²γ_αβ and ²δ_αβ in HHeH follow a similar trend as in stretched H₂, with quantitative differences introduced by the intermediate coupling of a third center. Likewise, the squared Frobenius norm of ²Δ approaches a value incrementally higher. While 〈S²〉_UCCSD is now able to register a signal, it is slower to reflect the emergence of biradicaloid character in the region of 1.25–1.5 Å, compared to ²Δ quantities.

The beryllium dimer has been the focus of a barrage of experimental and theoretical studies, due in part to the incredibly small measured bond energy of 2.67 kcal mol⁻¹. Another interesting feature is the quasi-degeneracy of the Be 2s and 2p atomic orbitals;⁷⁵ the mixing of these orbitals as the atom centers approach one another gives rise to a notoriously pathological electronic regime. Predictions on Be₂ were found to be sensitive to choices of active space – RHF and RCAS(4e,4o) yield qualitatively incorrect repulsive curves while RCAS(2e,2o) recovers the correct behavior at the cost of quantitatively incorrect energies. The CCSD curve is also repulsive, and CCSD(T) yields both less consistent and less accurate correlation energies recovering less than 70% of the total.^76,77 Furthermore, in small basis sets FCI does not bind Be₂. In sum, mean-field and even low-order coupled cluster theory are inadequate in describing the beryllium dimer, and the most challenging regime is not at dissociation (though this is still non-trivial due to the near-degeneracy of the 2s and 2p atomic orbitals) but at the equilibrium bond length. As shown in Fig. 7, the pathological behavior at shorter bond lengths is correctly reflected by the ²Δ and 〈S²〉_UCCSD; however, the latter decays to zero in the limit of separated Be atoms, which in our view is unphysical: despite being essentially a two-electron problem (assuming large energy separation between the 1s and 2s orbitals), there is non-negligible correlation (“correlation energy” and statistical two-electron correlation) that arises due to the presence of the competing 2p orbitals. ‖²Δ‖² and Tr[²Δ] clearly reflect this near the dissociation limit (as does ²δ_αβ to a lesser extent).

4 Discussion and conclusions

Having focused our analysis mostly on the cumulant-based quantities, we briefly circle back to the RDMs themselves to discuss prominent features and their relation to the cumulants. For the sake of clarity, we take minimal basis (STO-6G) H₂ as an illustrative example, since the RDM eigenvalues correspond to the absolute values of the density matrix elements. Being a two-electron system, the block structure is simplified in the singlet such that ¹Γ_αα = ¹Γ_ββ and only the sub-blocks related to the ²Γ_αβ portions contribute, while for the triplet state only the ¹Γ_αα and ²Γ_αα blocks are non-zero.

As seen in Table 3, the squared Frobenius norm of the ²Γ is constant despite its changing eigenvalues ²γ. The change of the cumulant is dictated by the decrease in the ¹γ values, given their inverse relationship. Meanwhile, the triplet, here constrained by the minimal basis into a simple form for all geometries, formally never registers any signal across its spectrum because of the constrained nature of ¹γ. Beyond lacking size-extensivity, the ²Γ squared Frobenius norm is incapable of capturing all the correlation reflected in the eigenvalues, while the ²Δ norm and trace incorporate this in a straightforward way. We remark that larger triple- and quadruple-ζ basis sets show similar trends, with only minor numerical differences in the metrics.

---

---

We summarize the properties of the metrics investigated herein in Table 4. A few points are in order to clarify the listed properties: although UCCSD is of course size-extensive, 〈S²〉_UCCSD does not scale with particle number, in contrast to the rest of the metrics. This is not necessarily a bad thing, given that comparisons among molecules of different sizes may require something along the lines of division by N for size-extensive metrics. By orbital-invariant, we mean that a fixed wave function will yield a fixed metric no matter which single-particle orbitals are used to represent it. In addition, an ideal metric that tracks the correlation strength will flag strong correlations that arise due to the need to restore symmetries of the system – not just spin-symmetry but spatial and others as well – and should also be meaningful in the thermodynamic limit (〈S²〉 is not a good quantum number in this limit).

---

A nuanced point about the orbital invariance property of ²Δ quantities: as pointed out recently,⁵² the hydrogen dimer offers a counter-example to when the diagonal elements of the ²Δ alone may serve as a signal of strong correlation. The minimal wave function that incorporates mixing of the configurations, θ, and orbital rotations, ϕ

|Φ〉 = (cosθcos2ϕ + sinθsin2ϕ)|0α0β〉 + (sinθcos2ϕ + cosθsin2ϕ)|1α1β〉 + (cosθ − sinθ)cosϕsinϕ(|0α1β〉 + |1α0β〉)	(12)

is able to describe the entirety of the space. When ϕ = 0 the 1-RDM is always diagonal at every geometry when expressed in the molecular orbital basis. Yet at the maximally-entangled stretched geometry, , the 1-RDM is diagonal for all orbital optimization parameter values, ϕ. When such non-uniqueness occurs, the elements of the RDMs and their cumulants may then likewise take on non-unique values, even as their norms remain invariant. In this example, the diagonal elements themselves become , evaluating to zero when the four determinants have equivalent contributions at . In these rare situations where non-uniqueness of the 1-RDM may occur at a single geometry, the trace does not have physical content, and the (more generally) orbital-invariant square norm of the ²Δ is to be preferred.

In light of our numerical results and our assessment of the satisfaction of the properties listed above for most chemically-relevant electronic structure problems, both the trace and square norm of the ²Δ are apt metrics that reflect the statistical two-electron interdependence that we are proposing here as a meaningful interpretation of the term “strong correlation”. Given that Tr[²Δ] = (1 − ¹Γ)¹Γ, i.e. is computable from the 1-RDM alone, this quantity is arguably the most practical for electronic structure calculations of pure states which do not suffer from a non-uniqueness of the natural orbitals at a particular geometry as in the case of stretched H₂. In the natural spin orbital basis, , where it is evident that the trace takes on its maximal value when all NOONs are equal (i.e. when the wave function is a superposition of equi-probable many-body basis states). As pointed out previously, the squared Frobenius norm of ²Δ avoids this issue entirely and is generalizable to open/mixed systems while the trace is not. Importantly, both Tr[²Δ] and ‖²Δ‖ remain robust metrics of strong correlation even when the NOON spectrum alone fails to distinguish cases such as square and trapezoidal H₄.^52,59

In conclusion, after giving a critical overview of different conceptions of strong correlation, we find that the statistical interpretation offered by the trace or square norm of the cumulant of the 2-RDM is rigorously defined, physically meaningful, and computable. We have shown that for a strongly interacting 6-site Hubbard model (half-filled, at large U/t) and for a number of molecular reaction coordinates involving challenging electronic structure regimes, these metrics based on ²Δ indicate what can be justifiably viewed as strong correlation. In addition, a salient result of our work is that we build a bridge between the rather mysterious chemical concept of antiaromaticity and our view of strong correlation.

It is of interest to us whether the physical interpretation of strong correlation, as substantiated by the statistical meaning of the cumulants, can be bridged with chemical notions (in addition to aromaticity and antiaromaticity) or other observables of interest. Another important challenge that we are currently exploring involves obtaining approximate but qualitatively informative 1- and 2-RDMs from scalable models such as quantum Monte Carlo or selected configuration interaction. This will enable reliable assessments of the nature of electron correlations in larger systems of interest to the chemistry and physics communities.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

We thank Gustavo Scuseria for insightful discussions regarding 2-RDMs and comments on the manuscript. We also thank Alistair Sterling for preliminary UHF-based investigations of 〈V〉/〈T〉 ratios. This work was motivated by many interesting conversations in Znojmo (Summer 2023) at the “Strong Correlation in Molecules” conference.

References

1. J. W. Hollett and P. M. W. Gill, The two faces of static correlation, J. Chem. Phys., 2011, 134, 114111 CrossRef PubMed .
2. R. Izsák, A. V. Ivanov, N. S. Blunt, N. Holzmann and F. Neese, Measuring Electron Correlation: The Impact of Symmetry and Orbital Transformations, J. Chem. Theory Comput., 2023, 19, 2703–2720 CrossRef PubMed .
3. J. Shee, M. Loipersberger, D. Hait, J. Lee and M. Head-Gordon, Revealing the nature of electron correlation in transition metal complexes with symmetry breaking and chemical intuition, J. Chem. Phys., 2021, 154, 194109 CrossRef CAS PubMed .
4. C. Ahn, A. Cavalleri, A. Georges, S. Ismail-Beigi, A. J. Millis and J.-M. Triscone, Designing and controlling the properties of transition metal oxide quantum materials, Nat. Mater., 2021, 20, 1462–1468 CrossRef CAS PubMed .
5. N. J. Mayhall, P. R. Horn, E. J. Sundstrom and M. Head-Gordon, Spin–flip non-orthogonal configuration interaction: a variational and almost black-box method for describing strongly correlated molecules, Phys. Chem. Chem. Phys., 2014, 16, 22694–22705 RSC .
6. M. Imada, A. Fujimori and Y. Tokura, Metal-insulator transitions, Rev. Mod. Phys., 1998, 70, 1039 CrossRef CAS .
7. S. G. Stewart, Heavy-fermion systems, Rev. Mod. Phys., 1984, 56, 755 CrossRef .
8. B.-X. Zheng, C.-M. Chung, P. Corboz, G. Ehlers, M.-P. Qin, R. M. Noack, H. Shi, S. R. White, S. Zhang and G. K.-L. Chan, Stripe order in the underdoped region of the two-dimensional Hubbard model, Science, 2017, 358, 1155–1160 CrossRef CAS PubMed .
9. L. Ding, S. Mardazad, S. Das, S. Szalay, U. Schollwock, Z. Zimborás and C. Schilling, Concept of orbital entanglement and correlation in quantum chemistry, J. Chem. Theory Comput., 2021, 17, 79–95 CrossRef CAS PubMed .
10. T. Shen, H. Barghathi, A. Del Maestro and B. Rubenstein, Disentangling the Physics of the Attractive Hubbard Model via the Accessible and Symmetry-Resolved Entanglement Entropies. arXiv, 2023, preprint, arXiv:2312.11746,  DOI:10.48550/arXiv.2312.11746.
11. L. Ding, S. Knecht and C. Schilling, Quantum Information-Assisted Complete Active Space Optimization (QICAS), J. Phys. Chem. Lett., 2023, 14, 11022–11029 CrossRef CAS PubMed .
12. C. J. Stein and M. Reiher, Measuring multi-configurational character by orbital entanglement, Mol. Phys., 2017, 115, 2110–2119 CrossRef CAS .
13. K. Boguslawski, P. Tecmer, G. Barcza, O. Legeza and M. Reiher, Orbital entanglement in bond-formation processes, J. Chem. Theory Comput., 2013, 9, 2959–2973 CrossRef CAS PubMed .
14. P. Pulay and T. P. Hamilton, UHF natural orbitals for defining and starting MC-SCF calculations, J. Chem. Phys., 1988, 88, 4926–4933 CrossRef CAS .
15. T. J. Lee and P. R. Taylor, A diagnostic for determining the quality of single-reference electron correlation methods, Int. J. Quantum Chem., 1989, 36, 199–207 CrossRef .
16. I. M. B. Nielsen and C. L. Janssen, Double-substitution-based diagnostics for coupled-cluster and Møller–Plesset perturbation theory, Chem. Phys. Lett., 1999, 310, 568–576 CrossRef CAS .
17. I. W. Bulik, T. M. Henderson and G. E. Scuseria, Can single-reference coupled cluster theory describe static correlation?, J. Chem. Theory Comput., 2015, 11, 3171–3179 CrossRef CAS PubMed .
18. M. Degroote, T. M. Henderson, J. Zhao, J. Dukelsky and G. E. Scuseria, Polynomial similarity transformation theory: A smooth interpolation between coupled cluster doubles and projected BCS applied to the reduced BCS Hamiltonian, Phys. Rev. B, 2016, 93, 125124 CrossRef .
19. J. Hachmann, J. J. Dorando, M. Avilés and G. K. Chan, The radical character of the acenes: A density matrix renormalization group study, J. Chem. Phys., 2007, 127, 134309 CrossRef PubMed .
20. D. W. Small and M. Head-Gordon, Post-modern valence bond theory for strongly correlated electron spins, Phys. Chem. Chem. Phys., 2011, 13, 19285–19297 RSC .
21. Y. A. Aoto, A. P. de Lima Batista, A. Kohn and A. G. de Oliveira-Filho, How to arrive at accurate benchmark values for transition metal compounds: Computation or experiment?, J. Chem. Theory Comput., 2017, 13, 5291–5316 CrossRef CAS PubMed .
22. D. Hait, N. M. Tubman, D. S. Levine, K. B. Whaley and M. Head-Gordon, What levels of coupled cluster theory are appropriate for transition metal systems? A study using near-exact quantum chemical values for 3d transition metal binary compounds, J. Chem. Theory Comput., 2019, 15, 5370–5385 CrossRef CAS PubMed .
23. H. Neugebauer, H. T. Vuong, J. L. Weber, R. A. Friesner, J. Shee and A. Hansen, Toward Benchmark-Quality Ab Initio Predictions for 3d Transition Metal Electrocatalysts: A Comparison of CCSD(T) and ph-AFQMC, J. Chem. Theory Comput., 2023, 19, 6208–6225 CrossRef CAS PubMed .
24. R. Hoffmann, The many guises of aromaticity, Am. Sci., 2015, 103, 18 CrossRef .
25. M. Rosenberg, C. Dahlstrand, K. Kilsa and H. Ottosson, Excited state aromaticity and antiaromaticity: opportunities for photophysical and photochemical rationalizations, Chem. Rev., 2014, 114, 5379–5425 CrossRef CAS PubMed .
26. A. I. Boldyrev and L.-S. Wang, All-metal aromaticity and antiaromaticity, Chem. Rev., 2005, 105, 3716–3757 CrossRef CAS PubMed .
27. G. Li Manni, W. Dobrautz and A. Alavi, Compression of spin-adapted multiconfigurational wave functions in exchange-coupled polynuclear spin systems, J. Chem. Theory Comput., 2020, 16, 2202–2215 CrossRef CAS PubMed .
28. W. Kohn, Nobel Lecture: Electronic structure of matter—wave functions and density functionals, Rev. Mod. Phys., 1999, 71, 1253 CrossRef CAS .
29. J. Ning and D. G. Truhlar, Spin–Orbit Coupling Changes the Identity of the Hyper-Open-Shell Ground State of Ce⁺, and the Bond Dissociation Energy of CeH⁺ Proves to Be Challenging for Theory, J. Chem. Theory Comput., 2021, 17, 1421–1434 CrossRef CAS PubMed .
30. J. Lee and M. Head-Gordon, Distinguishing artificial and essential symmetry breaking in a single determinant: approach and application to the C₆₀, C₃₆, and C₂₀ fullerenes, Phys. Chem. Chem. Phys., 2019, 21, 4763–4778 RSC .
31. D. A. Mazziotti, Complete reconstruction of reduced density matrices, Chem. Phys. Lett., 2000, 326, 212–218 CrossRef CAS .
32. D. A. Mazziotti, Contracted Schrödinger equation: Determining quantum energies and two-particle density matrices without wave functions, Phys. Rev. A, 1998, 57, 4219 CrossRef CAS .
33. C. N. Yang, Concept of off-diagonal long-range order and the quantum phases of liquid He and of superconductors, Rev. Mod. Phys., 1962, 34, 694 CrossRef CAS .
34. A. Coleman, Structure of fermion density matrices. II. Antisymmetrized geminal powers, J. Math. Phys., 1965, 6, 1425–1431 CrossRef .
35. D. Jérome, T. Rice and W. Kohn, Excitonic insulator, Phys. Rev., 1967, 158, 462 CrossRef .
36. S. Safaei and D. A. Mazziotti, Quantum signature of exciton condensation, Phys. Rev. B, 2018, 98, 045122 CrossRef CAS .
37. V. Tsurkan, S. Zherlitsyn, L. Prodan, V. Felea, P. T. Cong, Y. Skourski, Z. Wang, J. Deisenhofer, H.-A. K. von Nidda and J. Wosnitza, et al., Ultra-robust high-field magnetization plateau and supersolidity in bond-frustrated MnCr₂S₄, Sci. Adv., 2017, 3, e1601982 CrossRef PubMed .
38. A. J. Coleman and V. I. Yukalov, Reduced Density Matrices: Coulson’s Challenge, Springer Science & Business Media, 2000, vol. 72 Search PubMed .
39. W. Kutzelnigg, G. Del Re and G. Berthier, Correlation coefficients for electronic wave functions, Phys. Rev., 1968, 172, 49 CrossRef .
40. E. Wigner, On the interaction of electrons in metals, Phys. Rev., 1934, 46, 1002 CrossRef CAS .
41. E. Wigner and F. Seitz, On the constitution of metallic sodium, Phys. Rev., 1933, 43, 804 CrossRef CAS .
42. E. Wigner and F. Seitz, On the constitution of metallic sodium. II, Phys. Rev., 1934, 46, 509 CrossRef CAS .
43. R. Kubo, Generalized cumulant expansion method, J. Phys. Soc. Jpn., 1962, 17, 1100–1120 CrossRef .
44. W. Kutzelnigg and D. Mukherjee, Cumulant expansion of the reduced density matrices, J. Chem. Phys., 1999, 110, 2800–2809 CrossRef CAS .
45. F. E. Harris, Cumulant-based approximations to reduced density matrices, Int. J. Quantum Chem., 2002, 90, 105–113 CrossRef CAS .
46. M. Hanauer and A. Köhn, Meaning and magnitude of the reduced density matrix cumulants, Chem. Phys., 2012, 401, 50–61 CrossRef CAS .
47. J. P. Misiewicz, J. M. Turney and H. F. Schaefer III, Reduced density matrix cumulants: The combinatorics of size-consistency and generalized normal ordering, J. Chem. Theory Comput., 2020, 16, 6150–6164 CrossRef CAS PubMed .
48. J. E. Harriman, Grassmann products, cumulants, and two-electron reduced density matrices, Phys. Rev. A, 2002, 65, 052507 CrossRef .
49. J. M. Herbert, Magnitude and significance of the higher-order reduced density matrix cumulants, Int. J. Quantum Chem., 2007, 107, 703–711 CrossRef CAS .
50. E. Pérez-Romero, L. M. Tel and C. Valdemoro, Traces of spin-adapted reduced density matrices, Int. J. Quantum Chem., 1997, 61, 55–61 CrossRef .
51. K. Shamasundar, Cumulant decomposition of reduced density matrices, multireference normal ordering, and Wicks theorem: A spin-free approach, J. Chem. Phys., 2009, 131, 174109 CrossRef CAS PubMed .
52. L. Kong and E. F. Valeev, A novel interpretation of reduced density matrix and cumulant for electronic structure theories, J. Chem. Phys., 2011, 134, 214109 CrossRef PubMed .
53. T. Juhász and D. A. Mazziotti, The cumulant two-particle reduced density matrix as a measure of electron correlation and entanglement, J. Chem. Phys., 2006, 125, 174105 CrossRef .
54. A. Luzanov and O. Prezhdo, High-order entropy measures and spin-free quantum entanglement for molecular problems, Mol. Phys., 2007, 105, 2879–2891 CrossRef CAS .
55. M. Piris, J. Matxain, X. Lopez and J. M. Ugalde, Spin conserving natural orbital functional theory, J. Chem. Phys., 2009, 131, 021102 CrossRef CAS PubMed .
56. O. J. Fajen and K. R. Brorsen, Multicomponent CASSCF revisited: Large active spaces are needed for qualitatively accurate protonic densities, J. Chem. Theory Comput., 2021, 17, 965–974 CrossRef CAS PubMed .
57. M. L. Abrams and C. D. Sherrill, Natural orbitals as substitutes for optimized orbitals in complete active space wavefunctions, Chem. Phys. Lett., 2004, 395, 227–232 CrossRef CAS .
58. M. D. Benayoun, A. Lu and D. A. Mazziotti, Invariance of the cumulant expansion under 1-particle unitary transformations in reduced density matrix theory, Chem. Phys. Lett., 2004, 387, 485–489 CrossRef CAS .
59. J. T. Skolnik and D. A. Mazziotti, Cumulant reduced density matrices as measures of statistical dependence and entanglement between electronic quantum domains with application to photosynthetic light harvesting, Phys. Rev. A, 2013, 88, 032517 CrossRef .
60. A. Raeber and D. A. Mazziotti, Large eigenvalue of the cumulant part of the two-electron reduced density matrix as a measure of off-diagonal long-range order, Phys. Rev. A, 2015, 92, 052502 CrossRef .
61. E. Ramos-Cordoba, P. Salvador and E. Matito, Separation of dynamic and nondynamic correlation, Phys. Chem. Chem. Phys., 2016, 18, 24015–24023 RSC .
62. Y. Pavlyukh and J. Berakdar, Accessing electronic correlations by half-cycle pulses and time-resolved spectroscopy, Phys. Rev. A, 2014, 90, 053417 CrossRef .
63. D. R. Alcoba, R. C. Bochicchio, L. Lain and A. Torre, On the measure of electron correlation and entanglement in quantum chemistry based on the cumulant of the second-order reduced density matrix, J. Chem. Phys., 2010, 133, 144104 CrossRef CAS PubMed .
64. A. Rettig, J. Shee, J. Lee and M. Head-Gordon, Revisiting the orbital energy-dependent regularization of orbital-optimized second-order Møller–Plesset theory, J. Chem. Theory Comput., 2022, 18, 5382–5392 CrossRef CAS PubMed .
65. D. Hait, A. Rettig and M. Head-Gordon, Beyond the Coulson–Fischer point: characterizing single excitation CI and TDDFT for excited states insingle bond dissociations, Phys. Chem. Chem. Phys., 2019, 21, 21761–21775 RSC .
66. Q. Sun, T. C. Berkelbach, N. S. Blunt, G. H. Booth, S. Guo, Z. Li, J. Liu, J. D. McClain, E. R. Sayfutyarova and S. Sharma, et al., PySCF: the Python-based simulations of chemistry framework, Wiley Interdiscip. Rev.: Comput. Mol. Sci., 2018, 8, e1340 Search PubMed .
67. Q. Sun, X. Zhang, S. Banerjee, P. Bao, M. Barbry, N. S. Blunt, N. A. Bogdanov, G. H. Booth, J. Chen and Z.-H. Cui, et al., Recent developments in the PySCF program package, J. Chem. Phys., 2020, 153, 024109 CrossRef CAS PubMed .
68. M. Motta, C. Genovese, F. Ma, Z.-H. Cui, R. Sawaya, G. K.-L. Chan, N. Chepiga, P. Helms, C. Jiménez-Hoyos and A. J. Millis, et al., Ground-state properties of the hydrogen chain: Dimerization, insulator-to-metal transition, and magnetic phases, Phys. Rev. X, 2020, 10, 031058 CAS .
69. A. J. Sterling, D. S. Levine, A. Aldossary and M. Head-Gordon, Chemical Bonding and the Role of Node-Induced Electron Confinement, J. Am. Chem. Soc., 2024, 146(14), 9532–9543 CrossRef CAS PubMed .
70. K. Jankowski and J. Paldus, Applicability of coupled-pair theories to quasidegenerate electronic states: A model study, Int. J. Quantum Chem., 1980, 18, 1243–1269 CrossRef CAS .
71. U. Baek, D. Hait, J. Shee, O. Leimkuhler, W. J. Huggins, T. F. Stetina, M. Head-Gordon and K. B. Whaley, Say no to optimization: A nonorthogonal quantum eigensolver, PRX Quantum, 2023, 4, 030307 CrossRef .
72. O. V. Gritsenko, J. Wang and P. J. Knowles, Symmetry dependence and universality of practical algebraic functionals in density-matrix-functional theory, Phys. Rev. A, 2019, 99, 042516 CrossRef CAS .
73. J. Wang and P. J. Knowles, Nonuniqueness of algebraic first-order density-matrix functionals, Phys. Rev. A, 2015, 92, 012520 CrossRef .
74. M. Włoch, J. R. Gour and P. Piecuch, Extension of the renormalized coupled-cluster methods exploiting left eigenstates of the similarity-transformed Hamiltonian to open-shell systems: A benchmark study, J. Phys. Chem. A, 2007, 111, 11359–11382 CrossRef PubMed .
75. M. El Khatib, G. L. Bendazzoli, S. Evangelisti, W. Helal, T. Leininger, L. Tenti and C. Angeli, Beryllium dimer: A bond based on non-dynamical correlation, J. Phys. Chem. A, 2014, 118, 6664–6673 CrossRef CAS PubMed .
76. W. Helal, S. Evangelisti, T. Leininger and A. Monari, A FCI benchmark on beryllium dimer: The lowest singlet and triplet states, Chem. Phys. Lett., 2013, 568, 49–54 CrossRef .
77. I. Magoulas, N. P. Bauman, J. Shen and P. Piecuch, Application of the CC(P;Q) hierarchy of coupled-cluster methods to the beryllium dimer, J. Phys. Chem. A, 2018, 122, 1350–1368 CrossRef CAS PubMed .
78. T. M. Krygowski and M. K. Cyrański, Structural aspects of aromaticity, Chem. Rev., 2001, 101, 1385–1420 CrossRef CAS PubMed .
79. I. Fernández, Aromaticity: Modern Computational Methods and Applications, Elsevier, 2021 Search PubMed .
80. R. F. Bader, Atoms in molecules, Acc. Chem. Res., 1985, 18, 9–15 CrossRef CAS .
81. Z. Chen, C. S. Wannere, C. Corminboeuf, R. Puchta and P. V. R. Schleyer, Nucleus-independent chemical shifts (NICS) as an aromaticity criterion, Chem. Rev., 2005, 105, 3842–3888 CrossRef CAS .
82. J. Kruszewski and T. M. Krygowski, Harmonic Oscillator Approach to the Definition of Aromaticity, 1972 Search PubMed.
83. E. Matito, M. Solà, P. Salvador and M. Duran, Electron sharing indexes at the correlated level. Application to aromaticity calculations, Faraday Discuss., 2007, 135, 325–345 RSC .
84. A. Kekulé, Sur la constitution des substances aromatiques, Bulletin mensuel de la Société Chimique de Paris, 1865, 3, 98 Search PubMed .
85. O. Werba, A. Raeber, K. Head-Marsden and D. A. Mazziotti, Signature of van der Waals interactions in the cumulant density matrix, Phys. Chem. Chem. Phys., 2019, 21, 23900–23905 RSC .
86. M. Motta, K. J. Sung, K. B. Whaley, M. Head-Gordon and J. Shee, Bridging physical intuition and hardware efficiency for correlated electronic states: the local unitary cluster Jastrow ansatz for electronic structure, Chem. Sci., 2023, 14, 11213–11227 RSC .
87. P. B. Karadakov, Ground-and excited-state aromaticity and antiaromaticity in benzene and cyclobutadiene, J. Phys. Chem. A, 2008, 112, 7303–7309 CrossRef CAS PubMed .
88. P. B. Karadakov, P. Hearnshaw and K. E. Horner, Magnetic shielding, aromaticity, antiaromaticity, and bonding in the low-lying electronic states of benzene and cyclobutadiene, J. Org. Chem., 2016, 81, 11346–11352 CrossRef CAS PubMed .
89. R. Papadakis and H. Ottosson, The excited state antiaromatic benzene ring: a molecular Mr Hyde?, Chem. Soc. Rev., 2015, 44, 6472–6493 RSC .
90. Y. Jung, T. Heine, P. v. R. Schleyer and M. Head-Gordon, Aromaticity of four-membered-ring 6π-electron systems: N₂S₂ and Li₂C₄H₄, J. Am. Chem. Soc., 2004, 126, 3132–3138 CrossRef CAS PubMed .
91. J. Pedersen and K. V. Mikkelsen, A Benchmark Study of Aromaticity Indexes for Benzene, Pyridine, and the Diazines - II. Excited State Aromaticity, J. Phys. Chem. A, 2023, 127, 122–130 CrossRef CAS PubMed .

This journal is © The Royal Society of Chemistry 2024