Potential energy surfaces in atomic structure: The role of Coulomb correlation in the ground state of helium

TL;DR

This paper introduces a novel factorization of the two-electron atom's wavefunction, leading to a potential energy surface analogy that enhances understanding of Coulomb correlation effects in helium's ground state.

Contribution

It develops an exact factorization approach for two-electron atoms, deriving pseudoeigenvalue equations that reveal a molecular-like potential energy surface for atomic structure analysis.

Findings

Reveals a nonadiabatic potential energy surface for helium's ground state.

Shows correlation-induced symmetry breaking and quantum phase transition.

Provides a new interpretative framework for Coulomb correlation in atoms.

Abstract

For the $S$ states of two-electron atoms, we introduce an exact and unique factorization of the internal eigenfunction in terms of a marginal amplitude, which depends functionally on the electron-nucleus distances $r_1$ and $r_2$, and a conditional amplitude, which depends functionally on the interelectronic distance $r_{12}$ and parametrically on $r_1$ and $r_2$. Applying the variational principle, we derive pseudoeigenvalue equations for these two amplitudes, which cast the internal Schr\"odinger equation in a form akin to the Born-Oppenheimer separation of nuclear and electronic degrees of freedom in molecules. The marginal equation involves an effective radial Hamiltonian, which contains a nonadiabatic potential energy surface that takes into account all interparticle correlations in an averaged way, and whose unique eigenvalue is the internal energy. At each point $(r_1,r_2)$, such surface is, in turn, the unique eigenvalue in the conditional equation. Employing the ground state of He as prototype, we show that the nonadiabatic potential energy surface affords a molecularlike interpretation of the structure of the atom, and aids in the analysis of energetic and spatial aspects of the Coulomb correlation, in particular correlation-induced symmetry breaking and quantum phase transition.

Equations

$$\hat{H}\Phi (r_1,r_2,r_{12})=E\Phi (r_1,r_2,r_{12}),$$

$$\hat{H}=$$

$$D(r_1,r_2,r_{12})=r_1{r}_2{r}_{12}|\Phi (r_1,r_2,r_{12}){|}^2.$$

$$D(r_1,r_2,r_{12})=D_m(r_1,r_2)D_c(r_{12}|r_1,r_2),$$

$$D_m(r_1,r_2):={\int }_{|r_1-r_2|}^{r_1+r_2}d{r}_{12}D(r_1,r_2,r_{12})$$

$${\int }_0^{\infty }d{r}_1{\int }_0^{\infty }d{r}_2{\int }_{|r_1-r_2|}^{r_1+r_2}d{r}_{12}D(r_1,r_2,r_{12})=1,$$

$${\int }_0^{\infty }d{r}_1{\int }_0^{\infty }d{r}_2{D}_m(r_1,r_2)=1,$$

$${\int }_{|r_1-r_2|}^{r_1+r_2}d{r}_{12}{D}_c(r_{12}|r_1,r_2)=1.$$

$$\Phi (r_1,r_2,r_{12})=\psi (r_1,r_2)\chi (r_{12}|r_1,r_2),$$

$$\psi (r_1,r_2)$$

$$\chi (r_{12}|r_1,r_2):=e^{-i\alpha (r_1,r_2)}\frac{\Phi (r_1,r_2,r_{12})}{{⟨\Phi |\Phi ⟩}^{1/2}},$$

$$D_m(r_1,r_2)=r_1{r}_2|\psi (r_1,r_2){|}^2,$$

$$D_c(r_{12}|r_1,r_2)=r_{12}|\chi (r_{12}|r_1,r_2){|}^2.$$

$${\int }_0^{\infty }d{r}_1{r}_1{\int }_0^{\infty }d{r}_2{r}_2⟨\Phi |\Phi ⟩=1,$$

$${\int }_0^{\infty }d{r}_1{r}_1{\int }_0^{\infty }d{r}_2{r}_2|\psi {|}^2=1,$$

$$⟨\chi |\chi ⟩=1.$$

$$F[\Phi ]$$

$$\left({\hat{T}}_0+U(r_1,r_2)\right)\psi (r_1,r_2)=ϵ\psi (r_1,r_2),$$

$${\hat{\Omega }}_{\psi }(r_1,r_2)\chi (r_{12}|r_1,r_2)=U(r_1,r_2)\chi (r_{12}|r_1,r_2),$$

$${\hat{T}}_0\equiv -\sum _i\left(\frac{1}{2}\frac{{\partial }^2}{\partial r_i^2}+\frac{1}{r_i}\frac{\partial }{\partial r_i}\right),$$

$${\hat{\Omega }}_{\psi }(r_1,r_2)\equiv \hat{H}-\sum _{i\ne j}^2\left(\frac{1}{\psi }\frac{\partial \psi }{\partial r_i}\right)\left(\frac{r_i^2-r_j^2+r_{12}^2}{2r_i{r}_{12}}\frac{\partial }{\partial r_{12}}+\frac{\partial }{\partial r_i}\right).$$

$$U(r_1,r_2)=\frac{\lambda (r_1,r_2)}{{\rho }_m(r_1,r_2)}-\frac{{\hat{T}}_0\psi (r_1,r_2)}{\psi (r_1,r_2)}.$$

$$\Phi (r_1,r_2,r_{12})=N{e}^{-\zeta (r_1+r_2)},$$

$$\Phi (r_1,r_2,r_{12})=N{e}^{-\zeta (r_1+r_2)}(1+c_1{r}_{12}),$$

$$\Phi (r_1,r_2,r_{12})=$$

Full text

atomic shape, Coulomb correlation, helium atom, marginal-conditional factorization, nonadiabatic potential energy surface.

I Introduction

The notion that molecules and other atomic aggregates possess well-defined shapes, arguably the most basic paradigm of chemistry, is rooted in the classical structural theory Woolley1978 ; Claverie1980 . In the quantum-mechanical framework, this notion is put in on the basis of the topography of a Born-Oppenheimer (BO) potential energy surface (PES) Sutcliffe2012 .

In atomic physics, the success of the independent-particle central-field model led, for a long time, to the view that atoms are essentially spherical Berry1989-95 . However, the discovery that the intrashell supermultiplets in the double-excitation spectra of He can be empirically assigned to collective rotational-like and bendinglike motions of the electrons Kellman1980 , analogous to the rovibrational motions of a nonrigid linear $XYX$ molecule, challenged this view Berry1989-95 . Afterwards, it was found that more complex atoms can have nonspherical shapes even in their ground states Bao1994 . An analogous situation, referred to as “crystallization”, arises in “artificial atoms”, which are mesoscopic systems constituted by electrons confined in semiconductor quantum dots Puente2004 .

In nuclear physics, it was early recognized that nuclei can deviate from the spherical shape Nilsson1995 . Moreover, the concept of a BO PES was introduced, in analogy with the molecular case, albeit in the limited context of the empirical liquid drop model, with the roles of the adiabatic degrees of freedom played by the parameters that specify the surface deformation Bohr1939 ; Moller2009 .

Without recourse to the BO approximation, for any collection of particles, elements of structure can still be extracted from the patterns in suitably defined probability densities obtained from all-particle wavefunctions. For example, such densities reveal that in some states of few-valence-electron atoms Berry1989-95 ; Bao1994 ; Berry1981-87 and other few-body systems Matyus2011 ; Becerra2013 ; Birman2013 the correlations can induce such a collective behavior of the particles that the system acquires a moleculelike shape. Nevertheless, this probability-based notion of “quantum structure” is not as far reaching as the BO-PES-based notion of classical structure commonly employed in the molecular sciences Woolley1978 ; Claverie1980 .

Hence, it is highly desirable to generalize the concept of PES, starting from the quantum-mechanical Hamiltonian, to any assembly of particles, since this would allow the transfer of chemical-like notions, e.g., bond length, bond angle, and transition state, together with the conceptual and technical apparatus of quantum chemistry and molecular spectroscopy, e.g., the Franck-Condon principle, to other realms, thereby providing a unified treatment of atoms in molecules Berry1983 ; Kellman1997 , electrons in atoms Berry1983 ; Kellman1997 ; Essen1977 , quasiparticles in nanostructures Puente2004 ; Birman2013 , nucleons in nuclei Birman2013 ; Villars1987 , and even quarks in baryons Richard1992 . At first sight, this program might seem to be hampered by the fact that an adiabatic separation of degrees of freedom, the basic tenet of the BO approximation Sutcliffe2012 ; Essen1977 ; Villars1987 , is not always possible. However, in this paper we rigorously prove that such generalization is indeed feasible, by introducing an exact and unique marginal-conditional factorization (MCF) of the wavefunction Hunter1975 .

As a prototype, we show how to define a nonadiabatic PES (NAPES) for the radial motions of the electrons in $S$ states of two-electron atoms, which amounts to formulating an exact central-field model for these states. For the case of the ground state of He, the topography of this surface allows us to extract a classical molecularlike interpretation of the structure of the atom. In addition, we analyze the contributions of the different aspects of the electron correlation to this topography. Furthermore, we show that this NAPES provides a convenient conceptual framework for addressing issues of correlation-induced symmetry breaking and quantum phase transition.

II Formalism

For the $S$ states of two-electron atoms, the internal Schrödinger equation,

[TABLE]

determines the nonrelativistic energy, $E$ Hylleraas1929-64 . Considering the center of mass located at the nucleus, the internal coordinates $r_{1},r_{2},r_{12}$, with $r_{1},r_{2}$ being the electron-nucleus distances and $r_{12}$ being the interelectronic distance, determine the shape and size of the electron-nucleus-electron triangle. Such coordinates are independent, except that they are constrained by the triangle condition $|r_{1}-r_{2}|\leq r_{12}\leq r_{1}+r_{2}$. The volume element of this configuration space is $dV=r_{1}r_{2}r_{12}dr_{1}dr_{2}dr_{12}$. In atomic units, the fixed-nucleus Hamiltonian is given by

[TABLE]

Evidently, in this coordinate system there appear kinetic couplings, i.e., terms containing derivatives with respect to $r_{12}$, besides the original Coulomb potential coupling, $r_{12}^{-1}$.

We wish to construct a PES for the radial degrees of freedom, $r_{1},r_{2}$, which we will accomplish by averaging over the $r_{12}$ variable, in the following way. The internal distribution function is given by

[TABLE]

According to Bayes’ rule, this function can be factorized as

[TABLE]

where

[TABLE]

is the marginal distribution function for finding the electrons at $(r_{1},r_{2})$ irrespective of $r_{12}$, and $D_{c}(r_{12}|r_{1},r_{2})$ is the conditional distribution function for finding the electrons at $r_{12}$ provided that they are found at $(r_{1},r_{2})$. (Here and henceforth, when we speak of the probability of finding the electrons at a point, we actually mean in an infinitesimal neighborhood around that point.) The normalization of $D(r_{1},r_{2},r_{12})$,

[TABLE]

automatically implies the normalization of $D_{m}(r_{1},r_{2})$,

[TABLE]

and the local normalization of $D_{c}(r_{12}|r_{1},r_{2})$,

[TABLE]

Following Hunter Hunter1975 , we introduce the MCF of the internal eigenfunction

[TABLE]

by defining marginal and conditional amplitudes Gross2012-14 ; Arce2012

[TABLE]

[TABLE]

such that

[TABLE]

[TABLE]

[From Eq. (10) onwards, angular brackets express integrals over $r_{12}$ with the Jacobian $r_{12}$.] The normalization conditions (6)–(8) now read

[TABLE]

[TABLE]

[TABLE]

The following observations about the MCF [Eqs. (9)-(16)] are in order. First, it does not presuppose an adiabatic separation of the degrees of freedom, i.e., that $r_{1}$ and $r_{2}$ are “slow” in comparison with $r_{12}$. Second, even if the last summation in Eq. (2), which explicitly couples all the variables, were neglected, $\Phi(r_{1},r_{2},r_{12})$ still could not be exactly factorized as, say, $\psi(r_{1},r_{2})\xi(r_{12})$, because of the triangle condition. Third, the phase $e^{i\alpha(r_{1},r_{2})}$, with $\alpha(r_{1},r_{2})$ real, is arbitrary. Thus, given the exchange symmetry of $\Phi$, this phase can be chosen to set the exchange symmetries of $\psi$ and $\chi$, in the following way: Noting that $\left\langle\Phi|\Phi\right\rangle^{1/2}$ is always symmetric, for singlet states, where $\Phi$ must be symmetric, $\alpha$ can be chosen as symmetric (antisymmetric) to make both $\psi$ and $\chi$ symmetric (antisymmetric); for triplet states, where $\Phi$ must be antisymmetric, $\alpha$ can be chosen as symmetric (antisymmetric) to make $\psi$ symmetric (antisymmetric) and $\chi$ antisymmetric (symmetric). Fourth, the local normalization of $\chi$ (Eq. (16)) guarantees it to be nontrivial and unique, within the phase $\alpha(r_{1},r_{2})$ Gross2012-14 . Fifth, $\chi$ is not globally normalizable despite the fact that $\Phi$ is, since $\int_{0}^{\infty}dr_{1}r_{1}\int_{0}^{\infty}dr_{2}r_{2}\left\langle\chi|\chi\right\rangle=\mathcal{V}\rightarrow\infty$, with $\mathcal{V}$ the volume of the $\{r_{1},r_{2}\}$ subspace Arce2012 . Finally, if $\psi(r_{1},r_{2})$ had a node at, say, $r_{1}=a$ then $\Phi(a,r_{2},r_{12})$ would have to vanish at all $r_{12}$ [see Eq. (10)]. Moreover, for $\chi(r_{12}|r_{1},r_{2})$ to remain finite as $r_{1}\rightarrow a$, $\Phi(r_{1},r_{2},r_{12})$ would have to approach zero faster than $\left\langle\Phi|\Phi\right\rangle^{1/2}$ at all $r_{12}$ in this limit [see Eq. (11)], which is impossible. This constitutes a proof that marginal amplitudes must be nodeless, alternative to the one presented by Hunter Hunter1981 .

We derive the equations that govern $\psi$ and $\chi$ from the variational principle, as follows. First, we set up the constrained functional Gross2012-14

[TABLE]

where the first term is the expectation value of the energy, the second term ensures the local normalization of $\chi(r_{12}|r_{1},r_{2})$ at every point of the $\{r_{1},r_{2}\}$ subspace, and the third term ensures the normalization of $\psi(r_{1},r_{2})$, with $\lambda(r_{1},r_{2})$ and $\epsilon$ being Lagrange multipliers. Then, we impose the extremization condition $\delta F=0$, which yields, after a convenient rearrangement,

[TABLE]

[TABLE]

where

[TABLE]

[TABLE]

Here, the notation $\hat{\Omega}_{\psi}(r_{1},r_{2})$ indicates that this operator contains $\psi(r_{1},r_{2})$ and acts also on the $r_{1},r_{2}$ variables. The fact that in Eq. (19) $r_{1}$ and $r_{2}$ play the roles of parameters can be made more explicit by rewriting its left-hand side as $\hat{\Omega}_{\psi}(r_{1}^{\prime},r_{2}^{\prime})\chi(r_{12}|r_{1}^{\prime},r_{2}^{\prime})\vline_{r_{1}^{\prime}=r_{1},r_{2}^{\prime}=r_{2}}$, which means that $r_{1}$ and $r_{2}$ are fixed only after the operator has acted on $\chi$.

Equations (18) and (19) constitute a pair of exact coupled pseudoeigenvalue equations with the following characteristics. First, since the operators on their left-hand sides are Hermitian the eigenvalues $\epsilon$ and $U(r_{1},r_{2})$ are real. Second, since $\psi$ and $\chi$ are unique up to a phase for a given $\Phi$, each one of them possesses only one acceptable solution. Third, the presence of the logarithmic derivatives of $\psi$ in the second term at the right-hand side of Eq. (21) makes them nonlinear, which implies that their solution requires an iterative self-consistent scheme. (For mathematical caveats about this kind of problem, see Ref. Jecko2015 .) Finally, the eigenvalue $U(r_{1},r_{2})$ and the Lagrange multiplier $\lambda(r_{1},r_{2})$ are related by

[TABLE]

Consequently, by substituting Eq. (22) into Eq. (18) and taking into account Eq. (15) we see that $\epsilon=\int_{0}^{\infty}dr_{1}r_{1}\int_{0}^{\infty}dr_{2}r_{2}\lambda(r_{1},r_{2})$. Furthermore, with some additional manipulation we obtain that $\epsilon=\int_{0}^{\infty}dr_{1}r_{1}\int_{0}^{\infty}dr_{2}r_{2}\left\langle\Phi|\hat{H}|\Phi\right\rangle=E$. Therefore, $\lambda(r_{1},r_{2})=\left\langle\Phi|\hat{H}|\Phi\right\rangle$ can be interpreted as a local energy, i.e., the energy of the system when the electrons are positioned at $(r_{1},r_{2})$.

In Eq. (18) $\hat{T}_{0}$, as given by Eq. (20), represents the kinetic energy of two electrons in a central field. On the other hand, according to Eq. (19) $U(r_{1},r_{2})=\left\langle\chi|\hat{\Omega}_{\psi}(r_{1},r_{2})|\chi\right\rangle$, so this function carries all the information about the electron-nucleus attractions, electron-electron repulsion, and kinetic couplings, averaged over $r_{12}$. Hence, in Eq. (18) $U(r_{1},r_{2})$ plays the role of an effective radial potential that, nonetheless, correlates the electrons fully. Consequently, Eq. (18) constitutes an exact central-field model for $S$ states. In addition, since this is a Schrödinger equation for $\sqrt{D_{m}(r_{1},r_{2})/r_{1}r_{2}}$ (see Eq. (12)), we observe that the energy is a functional of the marginal distribution function, which from now on we will call the radial distribution.

Equations (18) and (19) are arranged in a form analogous to Hunter’s nonadiabatic electronic-nuclear separation Hunter1975 ; Gross2012-14 ; Arce2012 ; Cederbaum2013-2014 , which, in turn, is a sort of exact version of the BO approximation. Thus, $r_{12}$ and $r_{1},r_{2}$ are the analogs of the electronic and the nuclear coordinates, respectively, and $\hat{\Omega}$, $\hat{T}_{0}$ and $U$ are the analogs of the clamped-nuclei Hamiltonian (with the electron-nucleus attractions included), the nuclear kinetic energy operator, and the molecular NAPES. Consequently, we will refer to our $U(r_{1},r_{2})$ as the atomic NAPES. By the same token, $\chi$ and $\psi$ are analogous to the electronic and the vibrational eigenfunctions. Nevertheless, it must be kept in mind that, due to the uniqueness of the MCF, there can be only one marginal (“vibrational”) eigenfunction, with energy $E$, associated with this atomic NAPES, in contrast with the molecular BO case, where the PES can sustain several vibrational states with different energies.

III Illustrative Calculations and Discussion

We applied the foregoing formalism to the ground state of He. The solution to the nonlinear system of Eqs. (18) and (19) is a difficult problem Jecko2015 that lies beyond the scope of this work. Since our main goal is to learn about the topography of the atomic NAPES and its dependence upon the interparticle correlations, we employed these equations for the purposes of analysis only. Therefore, we extracted approximate marginal and conditional amplitudes, using Eqs. (10) and (11), choosing $\alpha\equiv 0$, from variationally optimized trial functions, and then evaluated approximate atomic NAPES’s by $U\approx\left\langle\chi|\hat{\Omega}_{\psi}|\chi\right\rangle$, in accordance with Eq. (19).

In order of increasing account of the Coulomb correlation, the optimized trial functions selected were

[TABLE]

with $N=13.59$, $\zeta=27/16$;

[TABLE]

with $N=12.26$, $\zeta=1.857270$, $c_{1}=0.390807$; and

[TABLE]

with $N=12.27$, $\zeta=1.755656$, $c_{1}=0.337294$, $c_{2}=-0.037024$, $c_{3}=-0.145874$, $c_{4}=0.023634$, $c_{5}=0.112519$. Expression (23) is the familiar Kellner uncorrelated function, whereas expressions (24) and (25) are of the Hylleraas type, with linear and quadratic polynomial correlation factors, respectively Hylleraas1929-64 . For these trial functions the expectation values of the energy turn out to be $\langle E\rangle=-2.848,-2.891,-2.903$ hartree, respectively Hylleraas1929-64 , which increasingly approach the reported “exact” energy Nakashima2007 . Of course, much better trial functions than (25) can be devised Nakashima2007 , but this level of approximation suffices for our largely qualitative analysis.

Before presenting the results, we must qualitatively compare the spatial behavior of the functions (23)-(25). First, we recall that, strictly, the correlation factor in a trial function should go to a constant as $r_{12}\rightarrow\infty$ (which can occur only when $r_{1}\rightarrow\infty$ and/or $r_{2}\rightarrow\infty$), because in this limit the electrons are uncorrelated. Hence, the only function that fulfills this condition is (23). In fact, the correlation factors of functions (24) and (25) diverge in this limit. (However, the full trial functions remain well behaved because the exponential factor damps the divergence.) Thus, in regions of large $r_{1}$ and/or large $r_{2}$ the uncorrelated function (23) actually provides a much better approximation to the exact wavefunction. On the other hand, close to the nucleus the correlated functions (24) and (25) provide better approximations. Very close to the nucleus, function (25) should provide a better approximation than function (24), but at sufficiently long distance from the nucleus the situation should reverse, since the divergence of (25) is stronger. Hence, we will confine our comparison of the NAPES’s extracted from these functions to a region relatively close to the nucleus, which is where the electron correlation plays an important role, anyway. Second, the correlation factor of function (25) also has corrections in $r_{1}$ and $r_{2}$, which (24) does not have. Therefore, due to the interplay of these characteristics, at a particular point $(r_{1},r_{2})$ it is difficult to predict which NAPES will be more accurate, although, evidently, the one evaluated with function (25) should be more accurate overall, since it provides the best energy.

Trial functions with polynomial correlation factors, in spite of their wrong asymptotic behaviors, are advantageous to us because they permit analytical evaluations of $\psi$, $\chi$ and $U$. For this task, we employed Mathematica 10.3.1 Mathematica . However, we do not provide the expressions here, since some of them are too formidable. As a check, we also reevaluated the energies by means of the expression $\langle E\rangle=\int_{0}^{\infty}dr_{1}r_{1}\int_{0}^{\infty}dr_{2}r_{2}\psi^{*}(r_{1},r_{2})\left(\hat{T}_{0}+U(r_{1},r_{2})\right)\psi(r_{1},r_{2})$ (see Eq. (18)), which had to be performed numerically, employing the same software; the values obtained turned out to be the same as the ones reported above.

Figure 1 displays the approximate atomic NAPES evaluated with the trial function (25). Several features stand out. First, $U$ has negative values everywhere, except for the steep repulsive walls at small $r_{1}$ or $r_{2}$. Second, it contains a ridge along $r_{1}=r_{2}$, which asymptotically goes to zero very slowly. (Strictly speaking, exactly at $r_{1}=r_{2}$ the value of $U$ is undetermined, an artifact of the approximate $\Phi$ employed. We obtained these values by a simple interpolation. The same observation applies to the NAPES that will be obtained later from the trial function (24)). Third, this ridge separates two basins with minima of $U_{e}=-7.73$ hartree positioned at $r_{e,i}=0.26,r_{e,j}=0.43$ bohr, and asymptotic values of $-5.57$ hartree. Fourth, there is a saddle point of $U^{\ddagger}=-7.42$ hartree positioned at $r_{1}^{\ddagger}=r_{2}^{\ddagger}=0.32$ bohr, which, along the minimum-energy path, lies at the top of a barrier of height $U_{b}\equiv U^{\ddagger}-U_{e}=0.31$ hartree located in between the two minima.

A classical molecularlike interpretation of this topography is as follows. The atomic NAPES is associative, as it should since we are dealing with a bound state. The minima correspond to two versions of the same equilibrium structure, with the two electrons on respective “Bohrian” orbits of radii 0.26 and 0.43 bohr. These versions are interconvertible by a degenerate rearrangement through the saddle point Bunker1998 , the latter corresponding to a transition structure with the two electrons on the same (unstable) orbit of radius 0.32 bohr. The stable orbits can perform highly anharmonic “breathing” motions, which are analogous to molecular bond stretching modes. Since the spatial part of the eigenfunction is symmetric under permutation of $r_{1}$ and $r_{2}$ (see Eq. (25)), these breathings occur in phase, analogously to a molecular symmetric stretching mode.

To examine the sensitivity of the topography of the NAPES to the degree of accuracy with which the electron correlation is taken into account, we also evaluated NAPES’s with the two other trial functions. Figure 2 shows the results. (We chose the size of the region exhibited in accordance with the effective range of the radial distribution, as discussed below.) We see that the three surfaces are very similar, which means that the topography of the NAPES is not very sensitive to the degree of correlation in the wavefunction. As the electron correlation in the trial function increases, the only trend we can discern is that the basins become more anharmonic, which implies that the amplitudes of the breathings of the orbits get larger.

To assess the roles of the Coulomb interactions in shaping the topography of the atomic NAPES, first we calculated surfaces with the electron-nucleus attractions, $-2r_{i}^{-1}$, removed from $\hat{\Omega}$ [see Eqs. (21) and (2)], for the three trial functions. As observed in Fig. 3, the wells disappear, i.e., the basins become valleys, and, consequently, $U$ becomes dissociative. Moreover, as the degree of correlation carried by the trial function increases, the maximum on the ridge slightly shifts towards shorter distance from the nucleus. Then, we calculated surfaces with both $-2r_{i}^{-1}$ and $r_{12}^{-1}$ removed from $\hat{\Omega}$, which means that the only couplings remaining are the kinetic ones. Now, in Fig. 4 we see that the ridge disappears and a single shallow well remains, whose bottom corresponds to an equilibrium configuration with both electrons on the same orbit of radius \sim\ 0.9 bohr. Furthermore, as the degree of correlation increases, such well becomes more anharmonic. Thus, the kinetic couplings play a small stabilizing role in the atom, since this surface is associative, and contribute to the anharmonicity of the basins present in the NAPES (see Fig. 2). Naturally, the electron-nucleus attractions play the major stabilizing role, manifested by the presence of the deep wells in the NAPES (see Fig. 2), and the electron-electron repulsion plays a destabilizing role, manifested by the presence of the ridge (see Figs. 2 and 3). Therefore, the Coulomb repulsion is responsible for inducing the symmetry breaking Birman2013 that gives rise to the split-orbit equilibrium configuration. Interestingly, the concept of a correlation-induced orbit splitting has been introduced in an ad hoc manner by defining inner and outer radial probability densities Koga2006 .

Because of the averaging over $r_{12}$, this atomic NAPES provides explicit information about the radial correlations only. To obtain explicit information about the angular correlations, we turn to the conditional distribution function, $D_{c}(r_{12}|r_{1},r_{2})$, which, for convenience of interpretation, we transformed into $D_{c}(\theta_{12}|r_{1},r_{2})$ using the relation $r_{12}=(r_{1}^{2}+r_{2}^{2}-2r_{1}r_{2}\cos\theta_{12})^{1/2}$. Figure 5 displays this function, extracted from the trial function (25), evaluated at selected points $(r_{1},r_{2})$ of the NAPES. Its behavior reveals the presence of the Coulomb hole: First, the probabilities of finding the system in the collinear electron-nucleus-electron ($\theta_{12}=\pi$) and nucleus-electron-electron ($\theta_{12}=0$) configurations are maximal and minimal, respectively. [Note that the probability of finding the electrons at zero separation vanishes because of the Jacobian $r_{12}$ present in Eq. (13).] Second, the values and curvature of this function decrease with the distance from the nucleus, which must be due to the weakening of the electron-electron correlation as the size of the electron-nucleus-electron triangle grows.

The corresponding radial distribution, $D_{m}(r_{1},r_{2})$, is shown in Fig. 6. It turns out to be unimodal, with its maximum located on the ridge, at $r_{1}=r_{2}=0.60$ bohr. (By r_{i}\sim 3\ bohr, this function has practically decayed to zero. This is why above we showed all the surfaces within the region 0\leq r_{i}\leq 3.0\ bohr.) Hence, quantum mechanics has frustrated the symmetry breaking “latent” in the NAPES, and placed the most probable configuration with both electrons on the same orbit of radius 0.60 bohr. Taking into account that the profile of the NAPES along the minimum-energy path looks like a double well, we suspect this has happened because the “zero-point energy”, $E_{0}\equiv\langle E\rangle-U_{e}=4.83$ hartree, is very much above the barrier, $U_{b}=0.31$ hartree. Consequently, the atom behaves analogously to a fluxional molecule. Had $E_{0}$ been below the barrier, $D_{m}$ would have displayed two humps, each associated mainly with the breathing of one orbit. It appears that something like this occurs in some of the doubly-excited states of He, for instance the nominal 2s3s\ ${}^{1}S^{e}$ state [see Fig. 2 in Ref. Berry1981-87 (d)], where the electron-electron correlation is more effective and, consequently, the ridge must be higher.

To support this suspicion, we considered fictitious atoms with the electron-electron repulsion increased, i.e., with $r_{12}^{-1}$ replaced by $\kappa r_{12}^{-1}$ ($\kappa>1$) in the Hamiltonian (2). For each value of $\kappa$ we used a trial function of the form (25), with all the parameters reoptimized. In the NAPES, we observed that, as $\kappa$ increases, hills continuously develop on each side of the ridge and the basins become more anharmonic. Concomitantly, at $\kappa\sim 2.5$ the radial distribution begins to develop two humps along the basins. Figures 7 and 8 display the NAPES and $D_{m}$, respectively, for $\kappa=4$. In this case the energy turns out to be $\langle E\rangle=-1.302$ hartree, which indicates that this fictitious atom is much less stable than the real one, as expected. The minima of the basins are now located at $r_{e,i}=0.27,r_{e,j}=2.28$ bohr with potential values of $U_{e}=-4.53$ hartree, and the saddle point is now positioned at $r_{1}^{\ddagger}=r_{2}^{\ddagger}=1.60$ bohr with a potential value of $U^{\ddagger}=-0.97$ hartree. The barrier height becomes $U_{b}=3.56$ hartree, which is much higher than the one for the $\kappa=1$ case. The two maxima of $D_{m}$ are positioned at $r_{i}=0.65,r_{j}=3.64$, very displaced from the potential minima due to the high anharmonicity of the basins. Along the ridge this function practically vanishes. Thus, the most probable configurations now have the electrons on respective orbits of radii 0.65 and 3.64 bohr, with the probability of interconverting by tunneling practically vanishing. These orbits can be envisioned undergoing essentially independent in-phase breathing motions. Hence, the symmetry breaking latent in the NAPES has become actual, due to the stronger electron-electron correlation Birman2013 . With confidence, we can associate this phenomenon with the fact that the zero-point energy, $E_{0}=3.23$ hartree, is now 0.33 hartree below the barrier.

The smooth evolution of the radial distribution from unimodal to bimodal as the Coulomb strength parameter, $\kappa$, increases is the hallmark of a continuous quantum phase transition, where $\kappa$ and the double-well profile of the minimum-energy path are analogous to the temperature and the free energy in thermodynamics, respectively. Phase transitions of this type in Coulomb three-body systems have been characterized Kais2000 ; the atomlike (unimodal) to moleculelike (bimodal) evolution of the particle density in the sequence H${}^{-}\rightarrow$He$\rightarrow$Ps${}^{-}\rightarrow$H${}_{2}^{+}$ Matyus2011 ; Becerra2013 is a case in point.

IV Conclusions and Outlook

The particular MCF (9) defined in this work produced an exact central-field model of the two-electron atom, embodied in the atomic NAPES $U(r_{1},r_{2})$ [see Eq. (18)]. For the ground state of He, we found that the topography of this surface (see Figs. 1 and 7), together with the radial distribution (see Figs. 6 and 8), provide a convenient framework for discussing issues of electron correlation, symmetry breaking, and quantum phase transition.

Three alternative MCF’s can be defined, which allow one to focus on complementary aspects. First, the MCF $\Phi(r_{1},r_{2},r_{12})=\psi(r_{2},r_{12})\chi(r_{1}|r_{2},r_{12})$ yields the NAPES $U(r_{2},\theta_{12})$, which is useful for analyzing the spatial behavior of the Coulomb hole Boyd1982 . In addition, such NAPES, together with the marginal distribution $D_{m}(r_{2},\theta_{12})$, should also provide a more complete picture of the correlation-induced symmetry breaking underlying the quantum phase transitions that have been observed in Coulomb three-body systems Matyus2011 ; Becerra2013 ; Birman2013 ; Kais2000 . Second, from the MCF $\Phi(r_{1},r_{2},r_{12})=\psi(r_{1})\chi(r_{2},r_{12}|r_{1})$ the density $\rho(r_{2},\theta_{12}|r_{1})$ used extensively in an ad hoc way by Berry and coworkers Berry1989-95 ; Berry1981-87 ; Berry1983 can be obtained simply as $|\chi(r_{2},r_{12}|r_{1})|^{2}$ followed by the transformation $r_{12}\rightarrow\theta_{12}$. In this case, an atomic nonadiabatic potential energy curve (NAPEC) $U(r_{1})$ can also be extracted Hunter1986 . Finally, the MCF $\Phi(r_{1},r_{2},r_{12})=\psi(r_{12})\chi(r_{1},r_{2}|r_{12})$ generates the NAPEC $U(r_{12})$, which can be used to generalize the molecular-orbital-like method of Feagin and Briggs, initially devised for the classification of doubly excited states of He Feagin1986-88 , to cases where the $r_{12}$ coordinate cannot be treated as adiabatic, as in the ground and singly excited states of this atom and other three-particle systems Hunter1967 .

Further insight can be gained by working with the hyperspherical coordinates $R\equiv(r_{1}^{2}+r_{2}^{2})^{1/2}$, $\alpha\equiv\tan^{-1}(r_{2}/r_{1})$ Lin1995 , instead of the Hylleraas coordinates $r_{1},r_{2}$. In particular, the MCF $\Phi(R,\alpha,\theta_{12})=\psi(R)\chi(\alpha,\theta_{12}|R)$ produces the NAPEC $U(R)$, whose adiabatic counterpart has been extensively used for the classification of states of Coulomb three-particle systems Lin1995 . (We would like to clarify that what is often called the “potential surface” in this connection is not a PES in the sense used in this paper, but rather the sum of all the Coulomb potentials appearing in the Hamiltonian expressed in hyperspherical or Jacobi coordinates.)

Ludeña and coworkers Becerra2013 have discovered that the topologies of non-BO one-particle nuclear and electron densities for few-body systems depend upon the reference points selected for the definitions of these quantities. This non-uniqueness weakens the notion that elements of molecular structure can be extracted from the patterns exhibited by the densities Matyus2011 . Here, the MCF can come to the rescue since one can define a density for any particle by marginalizing the variables associated with the other particles, and, as we have seen, marginal distributions are uniquely defined. All these constitute topics for future research.

We must emphasize that the implementation of this methodology does not rely on finding solutions of coupled marginal and conditional equations, such as (18) and (19), which appears to be a very difficult nonlinear problem Jecko2015 . Instead, these equations are used for extracting and analyzing information contained in wavefunctions that can be generated by other, more practical means, such as the variational trial functions (23)–(25). Thus, the extension of this methodology to more complex systems seems quite feasible. One of the issues we plan to address in the near future is whether atoms with more than two valence electrons in their ground states possess well-defined geometrical shapes or behave analogously to fluxional molecules Bao1994 .

In conclusion, we believe we have taken a first step in developing a research program aimed at unifying the treatment of atoms, molecules, and other collections of quantum particles, where the concept of PES, with the associated notion of geometrical shape, in the sense employed in the molecular sciences, can play a key role.

Acknowledgements.

L.D. Salas acknowledges financial support from Colciencias through a doctoral fellowship.