Molecular orbital framework is of central importance in chemistry. Often used by chemists and physicists to gain insight into molecular properties, Hartree–Fock or Kohn–Sham orbitals are obtained from rather crude treatments and, strictly speaking, are not observables. Yet, quantum mechanics offers a route for connecting general many-electron wavefunctions with reduced quantities—density matrices and orbitals—which give rise to observable properties. Such mapping makes possible, in principle, reconstruction of these objects from sufficiently detailed experimental data. This Perspective discusses Dyson orbitals and various types of natural transition orbitals and illustrates their role in modeling and interpreting different types of spectroscopic measurements.

## I. INTRODUCTION

Quantum mechanics explains chemical and physical phenomena by using abstract mathematical constructs, such as point charges and masses, spins, Hamiltonians, and wavefunctions. Today, it is rightfully considered to be the foundation of chemistry, and most general chemistry textbooks and courses begin with a brief foray into quantum mechanics. However, back in the 19th century, this trait of quantum mechanics, its reliance on mathematical language and abstractions to describe physics, was not unanimously embraced by the scientific community. Some of the best minds of the day denied quantum mechanics the legitimacy of being a scientific discipline because it operates with objects that cannot be observed.

Such was the opinion of one influential man from the 19th century, Ernst Mach,^{1} a founder of the Vienna Circle. Mach was not only an accomplished physicist (he even has a principle named after him) but also a philosopher of science. The school of thought that he helped to shape was called logical positivism, and one of its key principles was that science should be based entirely on directly observable phenomena.^{2} This sounds like something most of us subscribe to today, but there was a caveat in what Mach meant by observable phenomena. Mach defended a type of phenomenalism that recognized only sensations as real. That is, things you can see, hear, feel, or taste. Consequently, Mach considered atoms and molecules to be artificial constructs of the mind. He famously declared after an 1897 lecture by Ludwig Boltzmann at the Imperial Academy of Science in Vienna: “I don’t believe that atoms exist!”

These debates were eventually settled, and we now have a much broader view of what observables are. We do not need to rely on our unaugmented faculties to interact with the physical world. We accept that what we can see through microscopes, telescopes, or spectrometers is real. However, the main idea so ardently promoted by Mach is still valid—that science should only deal with phenomena that can be observed. That is, from the point of view of physical reality, a mathematical construct is only valid if one can design an experiment that can, in principle, probe it. With this in mind, let us talk about a very basic concept—orbitals. Are they real or just a mathematical construct?

When a common chemist hears “quantum mechanics,” she thinks “orbitals.” Orbitals are used to explain the structure of the Periodic Table of elements, the existence of chemical bonds, trends in reactivity, and the colors of compounds. The role of molecular orbitals (MOs) in providing rational explanation of chemical transformations was recognized by the 1981 Nobel Prize given to Kenichi Fukui and Roald Hoffmann.^{3}

Molecular orbital theory is now a truly native language of chemistry. Just as a small child, without understanding its grammar and rules, learns to speak his native language with great effectiveness, chemists learn how to speak orbitals before they learn quantum mechanics. Indeed, some chemists never learn quantum mechanics but use molecular orbital theory with great effectiveness.

In this Perspective, I wish to discuss molecular orbital theory with the experts, chemists who not only speak the language but also understand and appreciate its grammar. Hence, let us revisit the fundamentals and discuss the place of orbitals in quantitative and rigorous physical chemistry.

In one-electron systems, orbitals are wavefunctions. In a hydrogen-like atom, orbitals correspond to its electronic states. In molecules, atomic orbitals (AO) give rise to delocalized molecular orbitals, as shown in Fig. 1 for $H2+$. The delocalization lowers the kinetic energy of the electrons, gluing atoms into molecules. This is the quantum mechanical nature of what we call a chemical bond.^{4,5} Known as MO–LCAO (molecular orbital–linear combination of atomic orbitals) theory, a qualitative generalization of this exact result is used to explain trends in molecular properties at all levels from freshman textbooks to state-of-the-art research papers.

Being wavefunctions, orbitals of one-electron systems are connected to observables. Their squares give the probability of finding an electron at a certain position in space, known as electron density. Electron density can be probed directly by, for example, STM (scanning tunneling microscopy), classic and ultrafast electron diffraction, or x-ray scattering (see, for example, Refs. 6 and 7). Different shapes of electron densities give rise to concrete properties, such as multipole moments. Differences in the shapes of orbitals corresponding to different states encode optical properties, for example, one-photon transition between *σ* and *σ*^{*} orbitals in $H2+$ is induced only by light polarized along the molecular axis. Thus, experiments that probe these properties probe the wavefunctions and, consequently, orbital shapes. Given sufficiently detailed information, orbital shapes can be reconstructed from the experimental data. Hence, these orbitals are observables.

*N*-electron many-body wavefunction is approximated by the antisymmetrized product of

*N*orbitals (called the Slater determinant),

^{8}

The meaning of orbitals becomes even murkier when one considers strongly correlated wavefunctions, which contain many Slater determinants so that one can no longer rely on the pseudo-independent electron model to assign a physical meaning to the individual orbitals. Hence, one may wonder, do orbitals make any sense at all in many-electron systems?

In fact, orbitals do make sense, and we can connect them to experimental observables. We can also do this without invoking approximations—molecular orbital concepts are valid within rigorous many-body formalisms. Molecular orbital theory is fully compatible with the exact quantum mechanical description of many-electron systems, and orbitals are related to physical observables so that they are themselves observables. This Perspective reviews the formalism behind this rigorous connection and presents examples highlighting the utility of quantitative molecular orbital theory in modern quantum chemistry.

Our strategy is to follow a specific observable and to define orbitals from the perspective of how they relate to the chosen observable. This means that different observables give rise to different orbitals. For example, Dyson orbitals are related to photodetachment and photoionization cross sections probed by photoelectron spectroscopy or to electron-impact ionization cross sections, natural transition orbitals are related to linear spectroscopies (e.g., ultraviolet-visible or x-ray absorption spectroscopies), natural (density) difference orbitals (giving rise to attachment–detachment densities) are related to spectroscopies probing differences in states’ properties (e.g., Stark-effect spectroscopies), and so on. Thus, we abandon the idea of a single definition of the best molecular orbitals—in the author’s opinion, they do not exist. I also do not consider orbital definitions that cannot be related to potentially measurable quantities—in this respect, I share the sentiment that science should concern itself with observables. Different experiments probe different types of orbitals—hence, what you see depends on the tool you use to interrogate your system. This lack of a single definition of the true orbitals is not a flaw of the theory; rather, it is fully in accord with the quantum mechanical world in which the probe and the system that is probed are interconnected.

## II. TRUE ORBITALS FOR EVERY OCCASION

### A. Dyson orbitals

^{9}

^{,}

**is the dipole-moment operator,**

*μ***is a unit vector along the polarization of ionizing radiation, $\Psi IN$ and $\Psi FN\u22121$ denote the initial and the ionized target state, respectively, and $\Psi kel$ is the wavefunction of the ejected electron. Both $\Psi IN$ and $\Psi FN\u22121$ are many-body wavefunctions, but their explicit knowledge is not required for computing matrix elements of one-electron operators, such as, for example, the dipole moment**

*u*^{10}

*ϕ*

^{d}(

**) is a one-electron function, which contains all the necessary information about the molecular system before and after it is ionized,**

*r**p*

^{†}is the creation operator corresponding to orbital

*ϕ*

_{p}from the set of orthonormal spin-orbitals used to construct the many-electron basis functions (Slater determinants used to expand Ψ). In the derivation of this result,

^{11–14}we relied on the indistinguishability of the electrons, which allowed us to carry out the integration over

*N*electrons in two steps, thus reducing the many-body integral in Eq. (2) to integration over the coordinates of just one electron. The only approximation here is the assumption that $\Psi kel$ is orthogonal to the molecular core, which is usually justified by the large size of the free-electron state (this assumption can be lifted, giving rise to a slightly bulkier expression, with overlap corrections to the main term).

This function, *ϕ*^{d}(** r**), defined by Eq. (5), is called a Dyson orbital (or a generalized overlap amplitude).

^{11,14–20}Dyson orbitals can be defined for any type of wavefunction, including the exact one; thus, they do not rely on the pseudo-non-interacting electron model.

^{21}They simply tell you what the difference between the

*N*and

*N*− 1 electron system is. Thus, they can be interpreted as the initial state of the ejected electron.

*e*, 2

*e*) cross sections.

^{22,23}Electron-impact ionization is an example of electron-momentum spectroscopy, formulated in the momentum space. The (

*e*, 2

*e*) cross section is proportional to the energy-momentum spectral function,

*∑*

_{a$v$}denotes averaging over the initial and final state degeneracies and $\varphi id(p)$ is the Dyson orbital corresponding to ionized state

*i*(with ionization energy

*E*

_{i}) in the momentum space,

Because Dyson orbitals are directly related to the photoionization/photodetachment/electron-impact ionization cross sections, they can be reconstructed from the experimental data and are, therefore, observable quantities. Although the process of the reconstruction is not at all trivial (more on this below), here we merely point out that experimental imaging of Dyson orbitals is, in principle, possible.

Alternatively, we can compute experimental observables from a given Dyson orbital and judge the quality of the computed wavefunctions based on how well the theoretical and experimental results agree. An example in Fig. 2 illustrates this point. Absolute photoionization cross sections are very sensitive to the quality of Dyson orbitals. In a detailed benchmark study,^{14} we compared experimental absolute cross sections of the photoionization of small molecules with the theoretical values based on the Dyson orbitals computed from coupled-cluster (CC) and equation-of-motion coupled-cluster (EOM-CC) wavefunctions.^{24–26} Figure 2 shows such a comparison for formaldehyde. As one can see, theory and experiment are in perfect agreement. Importantly, using lower-quality orbitals (e.g., Kohn–Sham or Hartree–Fock) leads to noticeable discrepancies. Thus, the ability to reproduce the absolute photoionization cross sections provides a measure of the closeness of the shape of the Dyson orbital to the exact one [of course, the calculation also requires the photoelectron wavefunction, as per Eq. (4); hence, the quality of the computed cross section depends also on the quality of $\Psi kel$].

However, although total cross sections are highly sensitive to the shape of the Dyson orbital, we cannot really reconstruct the orbital from the total cross section, which is a too-deeply integrated quantity. Indirect information about orbital shapes can be deduced from vibrational progressions and Franck–Condon factors, but to really see the orbitals, one needs more detailed experimental observables. Differential cross sections, which provide angular-resolved photoelectron spectra, contain such information. For randomly oriented molecules, the photoelectron angular distribution (PAD) is completely characterized by the dipole anisotropy parameter *β*: *β* ∼ 0 corresponds to the isotropic distribution of photoelectrons, *β* ∼ 2 corresponds to the photoelectrons aligned parallel to the polarization of the ionizing radiation, and so on. Patterns in the PADs can be rationalized from dipole selection rules. As illustrated in Fig. 3, ionization from an *s*-orbital gives rise to a *p*-wave (with *β* = 2), ionization of a *p*-orbital gives rise to interfering *s*- and *d*-waves, and so on.

Thus, PADs (or *β*’s) encode the information about the shapes of the Dyson orbitals from which the photoelectrons originate.^{27–29} PADs can be used to quantify the extent of orbital deformation in the molecular environment due to hybridization and polarization. For example, PADs of gas-phase water molecules reveal that the molecular Dyson orbitals can be described as slightly deformed 2*p*_{O} orbitals, giving rise to *β* ≈ 1.5 at high photoelectron energies (above 100 eV). PADs are also used to interrogate the extent to which a solute’s wavefunction is perturbed by the solvent.^{29} A recent study^{30} used angular-resolved photoelectron spectroscopy and the microjets’ technique to understand the electronic structure of bulk water, i.e., to assess whether water molecules in bulk retain their electronic identities or give rise to qualitatively different delocalized band-like states. As shown in Fig. 4, the delocalized states are expected to give rise to an isotropic distribution of photoelectrons, in contrast to *β* ≈ 1.5 of localized states of individual water molecules. So which picture applies to bulk water? It turns out, the answer depends on the energy of the ionizing radiation: in the low-energy regime, one sees delocalized states, whereas in the high-energy regime, one observes the states of the individual waters. The two regimes and the transition between them can be understood in terms of interference between photoelectrons coming off individual molecules. When the de Broglie wavelength of photoelectrons is shorter than the distance between the individual molecules [20 eV, corresponding to 3.5 Å, first peak in *g*(*r*) of the liquid water], then there is no interference, and one observes an incoherent superposition of photoelectrons from the individual centers. Such an interplay between the de Broglie wavelength of photoelectrons and the distance between the emitting centers has been illustrated by Sanov’s molecular interferometer experiment, showing the transition from the delocalized and entangled molecular state of $I2\u2212$ to the incoherent atomic states in the course of dissociation.^{31} Surprisingly, the microjet measurements of neat water^{30} show persistent anisotropy of the photoelectrons, with only moderate reduction of *β* relative to the gas phase, due to transport scattering. In stark contrast to the energetic consequences of hydrogen bonding and electrostatic interactions, manifested by large (1 eV–2 eV) solvent-induced shifts and inhomogeneous broadening of the valence bands, these interactions do not result in reduction of anisotropy, meaning that the valence states of water in bulk can be thought of as the energy-shifted valence states of the isolated water molecules.

The utility of Dyson orbitals goes far beyond the interpretation of photoelectron or electron momentum (*e*, 2*e*) spectra. Because they show the difference between the *N* and *N* − 1 electron states, the shape of Dyson orbitals can be related to the differences in electron density and, consequently, bonding patters. In the formaldehyde example (Fig. 2), the Dyson orbital corresponding to the lowest ionization energy can be described as a lone pair on oxygen (with smaller contributions from *σ*_{OH}). Thus, the ionization of this mostly non-bonding orbital leads to relatively small changes in the geometry and short vibrational progression (a small feature at the onset is due to the respective Franck–Condon factors).

In open-shell species (e.g., radicals), Dyson orbitals can be used to describe the states of the unpaired electrons. Figure 5 shows Dyson orbitals of the four lowest electronic states in SrOH. The orbitals have clear atomic-like characters, revealing that these species can be described as a cationic core with the unpaired electron residing on the Sr. Because of this feature, the electronic transitions in SrOH are atomic-like and the change in the state of the unpaired electron does not perturb the core. Consequently, the geometry does not change upon excitation, leading to diagonal Franck–Condon factors. This type of electronic structure enables multiple optical transitions, which was instrumental in extending laser-cooling techniques to polyatomic molecules.^{33} The design of laser-coolable molecules has been guided by the analysis of orbitals involved in optical transitions,^{32,34–38} highlighting the utility of molecular orbital theory in the emerging field of quantum information science.

### B. Orbitals for excited states

It is difficult to imagine a conversation about excited states and photoinduced processes without invoking molecular orbitals.^{39,40} Yet, here again we face conceptual challenges: in many-electron systems, electrons interact with each other, and changing the state of one electron affects the states of the remaining ones. Hence, a naïve Koopmans’ picture of electronic excitations as removing an electron from an occupied orbital and attaching it to a virtual orbital is wrong quantitatively and often even qualitatively. For example, in a typical molecule, estimates of the lowest excited-state energy as orbital energy gaps are off by many electron-volts due to a large electron–hole interaction energy neglected within the pseudo-independent electron picture and a poor description of virtual orbitals. Qualitatively, the character of excited states also often differs from orbital-energy based predictions, especially when realistic basis sets with diffuse orbitals are used.

_{0}) and the target wavefunction is described as a linear combination of singly excited Slater determinants, $\Phi ia$,

*i*,

*j*,

*k*, … denote orbitals occupied in the reference determinant Φ

_{0}, and

*a*,

*b*,

*c*, … denote orbitals from the virtual subspace.) The analysis of the leading amplitudes $cia$ in this expansion enables the interpretation of electronic transitions in terms of hole and particle orbitals. For example, a CIS/6-31G

^{*}calculation of ethylene shows that the lowest excited state can indeed be described as a

*π*→

*π*

^{*}transition (the weight of the respective $\Phi ia\u2248$ 0.94), as expected from Koopmans’ theorem (HOMO → LUMO transition); this is illustrated in Fig. 6. Unfortunately, this type of analysis (based on the wavefunction amplitudes) is prone to several problems. First, it is not unique. The CIS energy and excited-state properties are invariant with respect to any unitary transformation of orbitals within the occupied or virtual orbital spaces, but the wavefunction amplitudes are not. Thus, choosing different orbitals changes the apparent character of the wavefunction. One consequence of this is that the weights of the leading amplitudes are strongly basis-set dependent. Often, in large bases, the wavefunction features several amplitudes with similar weights, which may or may not reflect the true character of the transition. Second, this analysis quickly becomes intractable when both states acquire multi-configurational character. For example, already at the CIS level, the analysis of the transitions between two

*excited*states is much more challenging than between the ground and excited states. In correlated methods, both the ground-state and excited-state wavefunctions become multiconfigurational and include multiply excited determinants.

^{41}The transition density is an object that connects the two states,

*x*

_{e}and

*x*

_{h}denote electron (particle) and hole coordinates, respectively, and $\gamma pqFI$ is a one-particle transition density matrix (1PTDM) between the initial and final states,

*p*

^{†}and

*q*are creation and annihilation operators corresponding to the

*ϕ*

_{p}and

*ϕ*

_{q}molecular orbitals, respectively. The transition density

*ρ*(

*x*

_{e},

*x*

_{h}) is related to observable properties because it defines the probability of transition. For example, the transition dipole moment between two many-body states is

*ρ*(

*x*

_{e},

*x*

_{h}) is mapped into the observable, it can, in principle, be reconstructed from the experimental data. When plotted in spatial coordinates,

*ρ*(

*x*

_{e},

*x*

_{h}) provides a visual map of the changes in the electron density. By virtue of Eq. (11), one can interpret

*ρ*(

*x*

_{e},

*x*

_{h}) as an exciton’s wavefunction and use it to rationalize the computed properties.

^{42}For example, the shape of

*ρ*(

*x*

_{e},

*x*

_{h}) for different excited states in a bacteriochlorophyll dimer shown in Fig. 7 explains the relative intensities and polarization of the main spectral features and clearly distinguishes local, entangled, and charge-transfer transitions.

^{43}

*γ*

^{FI}, the one-particle transition density matrix (1PTDM), contains the coefficients of the expansion of

*ρ*over molecular orbitals. It is a much simpler object than two many-body wavefunctions. One can think of it as of a map connecting the two states by moving one electron at a time, thus extending the CIS-like picture to correlated states,

^{F}from Ψ

^{I},

*γ*contains all the information needed to compute the transition dipole moment and, hence, relates to an experimentally observable property. For the CIS transitions from the ground state, i.e., when Ψ

^{F}= Ψ

^{CIS}and Ψ

^{I}= Φ

_{0}, only the occupied-virtual block of

*γ*is non-zero and $\gamma ia=cia$.

^{44}

*γ*and Δ are closely related and contain identical information;

^{45}however, for correlated wavefunctions (or transitions between the CIS states) this is no longer the case. For example, Δ can be used to characterize the differences between the states that are multiply excited with respect to each other, whereas

*γ*only captures the differences in terms of one-electron excitations.

While *ρ*(*x*_{e},*x*_{h}) and Δ(*r*) are invariant with respect to different orbital choices—any set of orthonormal orbitals can be used in Eqs. (11)–(14)—the values of the individual elements of *γ*_{pq} and Δ_{pq} depend on the specific orbital choice in the same way the wavefunction amplitudes do. Rather than obscuring the intrinsic nature of electronic transitions, this arbitrariness is the key to deriving a rigorous and clear orbital picture. One can find the best set of orbitals, which represents Δ_{pq} and *γ*^{FI} in the most compact way.

An important difference between the two objects is that the state density matrix is Hermitian, whereas the transition density matrix is not.^{46} Hence, the most compact representation of Δ_{pq} is achieved by simple diagonalization, yielding natural density-difference orbitals (NDOs) and their occupations, which add up to zero. By interpreting negative occupations as electron detachment from state Ψ^{I} and positive occupations as electron attachment to state Ψ^{F}, one can compute promotion numbers, which quantify how many electrons need to be simultaneously excited to describe the transition. NDOs and their occupation numbers yield attachment and detachment densities, which compactly represent the difference between the two states.^{44,47}

*γ*

^{FI}, the most compact representation is accomplished by the singular value decomposition,

^{91,92}

**Σ**is the diagonal matrix of singular values,

*σ*

_{K}, and matrices

**V**and

**U**contain the hole and particle (electron) NTOs according to

*ρ*has the most compact form,

*σ*

_{K}are non-vanishing. The squares of the

*σ*

_{K}s can be interpreted as the weights of the respective NTO pair when divided by the square of the Frobenius norm (Ω) of

*γ*,

^{48}

^{44,45,47,49–57}which provide a way to describe electronic excitations in terms of molecular orbitals for any type of wavefunction. NTOs are also directly related to the observables,

*π*→

*π*

^{*}has the transition dipole along the molecular axis, whereas

*π*→

*Ry*(3

*s*) has the transition dipole perpendicular to the molecular plane.

Figure 8 illustrates the ability of NTOs to unambiguously determine the character of the excited states. The CIS/aug-cc-pVDZ wavefunction of the lowest singlet state uracil features five leading amplitudes (which collectively account for only 75% of the wavefunction), suggesting a naïve interpretation of this transition as having multi-configurational Rydberg-valence character. However, the NTOs computed for this wavefunction reveal a pure (99.1%) valence *lp*(*O*) → *π*^{*} transition. Thus, the mixed multi-configurational appearance in the basis of canonical orbitals is just an artifact of the representation.

The ability of NTOs to provide an essential picture of the transition becomes even more important when using correlated wavefunctions. For example, in EOM-CCSD (EOM-CC with single and double excitations), both the ground-state and excited-state wavefunctions include single and double excitations, making amplitude-based analysis even more ambiguous. Yet, NTOs distill these complex wavefunctions into a simple orbital picture. Figure 9 shows the NTOs for the bright XA_{g} → 1B_{u} transition and dark XA_{g} → 4A_{g} transition in *trans*-stilbene computed using EOM-CCSD/daug-cc-pVDZ wavefunctions. The analysis shows that 1B_{u} can be described as a simple *π* → *π*^{*} excitation, whereas 4A_{g} has more complex character and requires two pairs of NTOs. Being independent of the correlation treatment and the basis, NTOs enable comparisons of wavefunctions computed by different methods,^{55} including comparisons between wavefunction treatments and time-dependent density functional theory.

*PR*

_{NTO}) is defined as

*γ*∥

^{2}. For a transition between a reference and a CIS target state,

*PR*

_{NTO}equals 1 for a CIS target state with exactly one configuration, 2 for a CIS state represented by two configurations with equal weights, and so on. Thus, pure excitations [such as

*ππ*

^{*}or

*π*→

*Ry*(

*s*) transitions in ethylene] are revealed by

*PR*

_{NTO}≈ 1, whereas states of mixed character (e.g., Rydberg-valence) would have larger values. Plasmonic excitations are identified by very large

*PR*

_{NTO}.

*PR*

_{NTO}is closely related to the number of entangled states, an alternative metric

^{58}also defined on the basis of the 1PTDMs,

*S*

_{HE}is the hole–electron entanglement entropy,

*S*

_{HE}is the sum of the

*αα*and

*ββ*parts of the transition). Thus, for singlet–singlet transitions,

*Z*

_{HE}accounting for spin-entanglement is the square of

*Z*

_{HE}computed from the spinless transition density. A recent study

^{36}of novel molecular frameworks for quantum information science used these metrics to investigate the entanglement between two optical cycling centers attached to the same molecular scaffold. The calculations on model systems (Me–O–CC–OMe′, Me = Mg, Ca, Ba, Sr) showed that the extent of entanglement and electronic communication can be controlled by the diradical character.

NTOs obtained from 1PTDM can always be used to describe the difference between the two states; however, their connection to spectroscopic observables requires further discussion. As per Eq. (20), NTOs can be used to understand the trends in one-particle transition properties. For example, from Fig. 9, one can easily see why the XA_{g} → 1B_{u} transition is bright and the XA_{g} → 4A_{g} is dark in linear spectroscopy. However, it is not at all obvious why the XA_{g} → 4A_{g} has a very large two-photon absorption (2PA) cross section simply because the NTOs shown in Fig. 9 do not directly relate to the 2PA cross section.

^{59–61}

_{ni}=

*E*

_{n}−

*E*

_{i}and

*ω*

_{1}and

*ω*

_{2}are the frequencies of the two absorbed photons (polarized along

*x*- and

*y*-directions) such that

*ω*

_{1}+

*ω*

_{2}= Ω

_{fi}.

Equation (24) reflects important physics of the 2PA and other nonlinear optical processes—the properties of the transition can no longer be understood by considering only the initial and final states. Rather, one needs to consider all states of the system. Except for resonant regimes (such as a 2PA transition in which one of the photons is resonant with a bright excited state Ψ_{n}), the SOS expressions cannot be truncated to a few terms without significant loss of accuracy, presenting an obstacle to deriving molecular orbital pictures of non-linear phenomena, such as 2PA or RIXS (resonant inelastic x-ray scattering).

*X*and $X\u0303$, are computed by solving auxiliary response equations.

^{60–63}Within the formalism outlined in Refs. 61 and 63, the EOM-CC response states are

_{I}s denote Slater determinants from the target-state manifold. The physical meaning of response states becomes clear once we recognize that they describe the first-order response

^{64}of Ψ

_{k}to the perturbation with frequency

*ω*

_{x},

*X*, which assimilate the contributions from all real excited states of the system, provide concrete meaning to the so-called “virtual states,” which are commonly invoked in non-linear spectroscopies.

Following Eq. (25), one can define and compute the respective 1PTDMs (called perturbed 1PTDMs) and extract the NTOs. Such analysis delivers simple molecular orbital pictures of nonlinear processes. Figure 10 shows NTOs for the bright 2PA transition in ethylene (X*A*_{g} → 2A_{g}). As one can see from the regular (linear or one-photon) NTOs, the 2A_{g} state has *πRy*(*p*_{z}) character. The NTO analysis of the M_{zz} and M_{xx} 2PA moments shows that the former moment is dominated by the pathway via the *Ry*(*s*) virtual state, whereas the latter moment involves the *π*^{*} virtual state.

Recently, following the same steps, we extended the NTO analysis of perturbed density matrices to the RIXS process.^{65} Figure 11 shows the orbital picture of the RIXS transitions in benzene derived from the NTOs computed from perturbed density matrices. As in the case of 2PA, these NTOs are directly mapped onto the cross sections. Without making simplifying assumptions, they reveal the character of the virtual state responsible for the transitions. In the benzene example, the character of the virtual state changes from valence *π*^{*}-like to Rydberg *s*-like depending on the frequency of the excitation pulse.

### C. Orbital picture of spin-forbidden transitions and magnetic properties

Molecular orbital theory has been used to explain relativistic phenomena, ranging from orbital contraction in heavy elements due to scalar relativistic effects to inter-system crossing (ISC), intensity borrowing, and magnetic anisotropies facilitated by spin–orbit interactions. The last effect plays an important role even in molecules composed of light atoms, as it couples otherwise non-interacting states (e.g., singlets and triplets), opening up new relaxation pathways.

In photochemistry, the rates of ISC are often explained by El-Sayed’s rules,^{66,67} which state that a large SOC (spin–orbit coupling) can be expected for transitions involving a flip of an orbital, such as in $[(\pi )2(n)1(\pi *)1]\u2192[(\pi )2(n)2(\pi *)0]$, whereas transitions that do not involve orbital flip, such as in $[(\pi )1(n)2(\pi *)1]\u2192[(\pi )2(n)2(\pi *)0]$, would have much smaller SOC and, consequently, slower ISC. El-Sayed’s rules are based on the analysis of the one-electron part of the Breit–Pauli Hamiltonian, which is responsible for the SOC. As illustrated in Fig. 12, the SOC term has the form of a local angular momentum operator whose matrix representation reveals the need for orbital torque, for example, in the basis of three *p*-orbitals, the only non-zero matrix elements of $L^z$ are between *p*_{x} and *p*_{y}.

Can we use the concept of NTOs to develop a quantitative version of El-Sayed’s rules applicable to the correlated many-body wavefunctions? The conceptual difficulty is that the SOC is a tensorial property that depends on the matrix elements between all multiplet components. For example, the SOC constant (SOCC), which enters the expressions for ISC rates and spectral line splittings, involves the sum over all spin-projections within the two coupled multiplets (e.g., the three M_{s} components of a triplet state, in the case of singlet–triplet couplings). Thus, the calculation of NTOs between one multiplet component (such as the M_{s} = 0 component, most often available from electronic structure calculations) is not sufficient for quantitative analysis of the SOCC. This obstacle can be overcome by using Wigner–Eckart’s theorem and introducing a new type of NTOs appropriate for the analysis of spin-forbidden phenomena.^{68} The definition of spin–orbit NTOs^{68} is based on the spinless density matrices that are used for the calculation of the SOCs for the entire multiplet from just one transition.^{69,70} This approach allows one to treat the transitions between states with arbitrary spin projections in a uniform way and to quantitatively describe the contributions of specific orbital pairs to the overall SOC; hence, these NTOs are linked to experimentally observable quantities (e.g., ISC and oscillator strengths). Figure 12 shows NTOs computed for the SOC within the quintet manifold in an iron-contained single-molecule magnet described by the EOM-EA-CC ansatz.^{68,70} A large spin-reversal barrier^{71} of this single-molecule magnet can be explained by the shapes of these NTOs. In the spirit of El-Sayed’s rule, these two orbitals, which can be described as *d*_{zy} and $dz2$, give rise to non-zero matrix elements of the angular momentum operator.

### D. Orbital concepts in non-Hermitian quantum chemistry

Standard quantum chemistry deals with electronically bound states (i.e., states stable with respect to electron ejection). These states are discrete, and their wavefunctions are finite (*L*^{2}-integrable). Ground states and low-lying excited states of many molecules are of this type. Yet, many important physical phenomena take us outside the bound part of the spectrum into the electronic continuum.^{72–75} Examples include highly excited states (any excited state above the ionization threshold), core-level states, excited and sometimes even ground states of anions (i.e., transient anions), and molecules on metal surfaces. In the continuum, there are no discrete states, and all wavefunctions are decaying asymptotically. Yet, the continuum has a structure, revealed spectroscopically. It contains states, called resonances, which behave as bound states, but with a finite lifetime. In Hermitian quantum mechanics, these states correspond to an increased density of states in certain parts of the continuum. Hence, standard quantum chemistry is not applicable, neither for a quantitative treatment of these states nor for their conceptual analysis. The increased density of states does not have a potential energy surface governing nuclear motion. The increased density of states cannot be easily interpreted using molecular orbital theory. Yet, spectral signatures of resonances, such as vibrational progressions and angular distributions of photoelectrons, suggest that the concepts of Born–Oppenheimer separation of nuclear and electronic motion and a molecular orbital picture of the transitions do apply to these states hidden in the continuum.

Non-Hermitian quantum mechanics^{72,73,75} offers an elegant framework for tackling resonances. By recasting the Schrödinger equation into a non-Hermitian form, one can project out the continuum and describe the resonances as isolated states with *L*^{2}-integrable wavefunctions and complex energies. The real part of the resonance energy corresponds to the resonance position, and the imaginary part corresponds to the resonance width (inverse lifetime). By using complex potential energy surfaces,^{76–79} one can describe nuclear dynamics associated with resonances. Similar to bound states, the observable properties of resonances are encoded in the respective state and transition density matrices. Within the non-Hermitian formalism, these densities become complex-valued.^{80,81} For example, Dyson orbitals and NTOs acquire real and imaginary components.

What can we learn from these complex-valued orbitals? In order to understand the physical meaning of real and imaginary densities, one needs to relate them to experimental observables. The lifetime of a resonance is the key observable. It is related to the imaginary density, opening a route to interpretation. Further insight into the meaning of imaginary densities and orbitals can be obtained by following the Feshbach–Fano formalism.^{82} Within this framework,^{81} wavefunctions of metastable excited states (such as excited states in many anions^{83}) can be analyzed in terms of real and imaginary excitons via direct extension of Eq. (12). Figure 13 shows such real and imaginary NTOs for two excited states of CN^{−}. The real NTO pairs reveal that the two states can be described as *π* → *π*^{*} and *σ* → *π*^{*} excitations. The imaginary NTOs characterize the decay channels due to the coupling to the continuum—the hole orbitals are very similar to the real hole orbitals, but the particle orbitals have different nodal structure and are more diffuse. Figure 14 shows the real and imaginary NTOs for C_{3}N^{−} and C_{7}N^{−}. The NTOs in the ^{1}Δ state of C_{3}N^{−} are similar to CN^{−}, revealing the *π* → *π*^{*} character of the bound part of the exciton and a *p*-like decay channel. In the ^{1}Σ^{+} state of C_{7}N^{−}, the NTOs have more complex character: the bound part of the exciton is a mixture of *π* → *π*^{*} and *σ* → *σ*(*Ry*) excitation and the imaginary part shows two different decay channels. Note that the relative weights of the configurations involving the *π* and *σ* holes are flipped in the real and imaginary parts—in the real exciton, the contribution of the *σ* hole is minor, whereas in the imaginary part, it becomes the leading component. Recall that, while the singular values for the real part have similar meaning as in the bound-state calculations (they describe the character of the transition), the singular values of the imaginary part are related to the partial widths of the resonance.^{81} Hence, this flip of leading *σ*_{K} reveals that this state of *ππ*^{*} character decays predominantly via the *σσ*(*Ry*) channel.

## III. CAN WE REALLY SEE THEM?

The theory provides a solid foundation for experimental orbital imaging, but design and implementation of such experiments are far from trivial. The interpretation of orbital imaging experiments is still controversial, highlighting the need for more robust theoretical modeling and data processing frameworks. Despite these challenges, there are many exciting examples illustrating our progress toward direct experimental orbital imaging.

STM enables atomic-scale imaging of molecules adsorbed on surfaces. By measuring changes in the local molecular density of states under different applied voltages, STM can provide real-space images of the electron density of the Dyson orbitals that connect the neutral molecules with their respective negative and positive ions. Figure 15 shows the images for cobalt phthalocyanine (CoPC). The experimental differential conductance maps reveal a clear nodal structure of the Dyson orbitals, in good agreement with the computed maps. However, the phase information is not available, and at present, STM imaging is limited to the two-dimensional picture. Additional complications involve the perturbation of the electronic structure by the STM tip and the surface.

Advanced light sources offer new possibilities for orbital imaging, such as the use of high-order harmonic generation.^{85} In addition, strong fields can be used to align molecules in space, which eliminates orientational averaging and, therefore, increases the amount of information that can be used for reconstruction. Angular-resolved photoelectron experiments, especially when executed for oriented molecules, can deliver three-dimensional images of orbitals;^{86–89} however, the exact details of the reconstruction process remain a key challenge. Given that the interpretation of these experiments often involves certain simplifying assumptions (such as a plane-wave treatment of photoelectrons) and multiple steps of converting the raw data into the image, one may argue that we are not yet truly able to see the orbitals, at least not in a way that would convince Mach and his fellow members of the Vienna Circle. Challenges of the full tomographic orbital imaging aside, the mapping between orbitals and observables is indisputable, and this connection provides a sound guiding principle for responsible use of orbital concepts in physics and chemistry.

## IV. CONCLUSION

In this Perspective, I have discussed molecular orbital concepts from the point of view of many-body electronic structure and observable molecular properties. The examples illustrate how simple molecular orbital pictures can be rigorously distilled from correlated wavefunctions by using reduced quantities such as one-particle transition density matrices. By using this formalism, molecular orbitals can be related to molecular properties and spectroscopic observables. Thus, orbitals are not a mere mathematical abstraction; they are real and can be observed. However, the type of orbital you see depends on which tool you use to look at them. This is an intrinsic feature of quantum mechanics—the result of the observation depends not only on the state of the system but also on the act of observation. Orbital concepts are instrumental in deriving physical insight from spectroscopy. Here again, different spectroscopies relate to different types of orbitals. I conclude with a general guiding principle for using molecular orbital concepts: Follow the observable or design a new one.

## DATA AVAILABILITY

Data sharing is not applicable to this article as no new data were created or analyzed in this work.

## ACKNOWLEDGMENTS

I wish to dedicate this work to Dr. Vladimir Ivanovich Pupyshev (August 06, 1949 to June 16, 2020), who passed away^{90} as this manuscript was being finalized. Dr. Pupyshev was an inspiring teacher and a gentle mentor, who ignited my interest in quantum chemistry and encouraged my initial clumsy attempts in theoretical work. Dr. Pupyshev was able to convey the conceptual beauty of quantum mechanics while maintaining the highest level of mathematical rigor, always spiced with a good joke, a literary quote, or a thought-provoking sharp remark. His scientific style has provided a constant inspiration for my work.

This work was supported by the U.S. National Science Foundation (Grant No. CHE-1856342).

A.I.K. is the president and a part-owner of Q-Chem, Inc.

## REFERENCES

People of my generation may find amusing the connection between Mach and Lenin. In his essay “Materialism and Empirio-criticism,” Lenin scorned Mach and his followers, along with several other “reactionist philosophers.” Despite having extensive schooling in dialectical materialism and Lenin’s work—we had to write extensive conspectuses of this particular treatise to pass exams—I am completely blank on which part of Mach’s teaching Lenin deemed to be at odds with Marxist’s philosophy. Its low intellectual content and poor scholarly quality notwithstanding, “Materialism and Empirio-criticism” impacted science in a profound way. The arrogance of a poorly educated revolutionary zealot attacking competent scientists heralded the beginning of a new era in which ideology trumps the science and politicians censor it. In the spirit of “Materialism and Empirio-criticism,” Soviet commissars declared entire scientific disciplines as ideologically impure (reactionary, renegade, bourgeoisie, imperialistic, colonial, nationalistic, etc.) and antagonistic to the cause of the superior class (proletariat, in Soviet Russia), shutting down research and sanctioning (up to physically eliminating) the researchers.

Royal Swedish Academy of Sciences, Press Release, 1981, “The Royal Swedish Academy of Sciences has decided to award the 1981 Nobel Prize in chemistry jointly to Professor Kenichi Fukui, Kyoto University, Kyoto, Japan, and the other half to Professor Roald Hoffmann, Cornell University, Ithaca, NY, USA, for their theories, developed independently, concerning the course of chemical reactions.”

*Molecular Electronic Structure Theory*

*Quantum Mechanics of One and Two Electron Atoms*

*Ab initio*molecular dynamics and time-resolved photoelectron spectroscopy of electronically excited uracil and thymine

When computed using the Hartree–Fock wave function of the neutral and neglecting orbital relaxation in the cation, Dyson orbitals are canonical Hartree–Fock orbitals, as stipulated by Koopmans’s theorem. If orbitals are allowed to relax, as in the ΔSCF treatment of ionization, Dyson orbitals [as defined by Eq. (5)] can be computed by singular value decomposition of the overlap matrix between the orbitals of the neutral and the ionized systems, as was done by Martin and Davidson^{17} and described formally by Goscinski and Lindner.^{15} Reference 17, which represents one of the earlier appearances of Dyson orbitals in mainstream quantum chemistry,^{17} used the term “corresponding ionization orbitals” and did not make an explicit connection with the Dyson orbitals (also called generalized overlap amplitudes) appearing in Green’s function theory.^{15,16} Rather, the theory was described in terms of the “corresponding orbital transformation of Amos and Hall.”^{91}

*Electron Momentum Spectroscopy*

*Electronic Aspects of Organic Photochemistry*

In the same sense and with the same caveats as discussed in the case of Dyson orbitals—it can, in principle, be reconstructed from sufficiently detailed experimental information, i.e., angular resolved absorption cross sections for molecules fixed in space.

_{2}antenna complex of purple bacteria

The state density matrix is Hermitian in the exact theory and can become slightly non-Hermitian within approximate treatments, such as coupled-cluster and equation-of-motion coupled-cluster theory. However, since only the Hermitian part of the state density matrix contributes to state properties, such as dipole moments, one can work with symmetrized state density matrices for the purpose of orbital analyses and calculations of molecular properties.

Promotion numbers computed as the sum of positive (or negative) occupations of the density difference matrix Δ provide a more general metric. For CIS wave functions, both quantities yield the same answer. A large doubly excited character of a transition would manifest itself in a reduced value of Ω, whereas density-difference eigenvalues would deliver a more direct answer (i.e., ∼2) (see examples in Ref. 53).

The analysis of excited states in terms of transition densities and density differences was pioneered by Luzanov and co-workers^{49,50} in the 1970s. The term “Natural Transition Orbitals” was coined by Martin,^{51} who independently arrived at this concept in 2003, following Amos and Hall’s work.^{91} He illustrated NTOs by using the CIS ansatz, overlooking the fact that in this case NTOs are formally and numerically equivalent to the natural orbitals of density difference whose squares yield attachment-detachment densities, which were described by Head-Gordon and co-workers in 1995.^{44} NTOs (in the context of CIS) were rediscovered yet again in 2007 by Mayer.^{52} Responding to Mayer’s paper, Surjan has shown^{45} that CIS NTOs are equivalent to the NOs, again overlooking all prior work on the subject.^{44,49–51} Russian-speaking readers can enjoy a humorous account of the history of NTOs and singular value decomposition in quantum chemistry in the essay by Luzanov,^{92} presented as an illustration of the (extended) Arnold’s principle: “If a notion bears a personal name, then this name is not the name of the discoverer.”

*Non-Hermitian Quantum Mechanics*

*About V. Arnold’s Principle and Its Application to Quantum Chemistry and Related Fields*