The coupling constant averaged exchange–correlation energy density

ABSTRACT The exchange–correlation energy, central to density-functional theory, may be represented in terms of the coupling constant averaged (CCA) exchange–correlation energy density. We present an approach to calculate the CCA energy density using accurate ab initio methods and its application to simple atomic systems. This function provides a link between intrinsically non-local, many-body electronic structure methods and simple local and semi-local density-functional approximations (DFAs). The CCA energy density is resolved into separate exchange and correlation terms and the features of each compared with those of quantities commonly used to construct DFAs. In particular, the more complex structure of the correlation energy density is found to exhibit features that align well with those present in the Laplacian of the density, suggesting its role as a key variable to be used in the construction of improved semi-local correlation functionals. The accurate results presented in this work are also compared with those provided by the Laplacian-dependent Becke–Roussel model for the exchange energy.


Introduction
A great deal of progress has been made in the development of increasingly accurate Kohn-Sham (KS) [1] density-functional approximations (DFAs) in recent decades. Over this time, there has been a general shift away from simple DFAs, constructed using the purely local information, such as the electron density and its derivatives, towards increasingly non-local models involving the KS orbitals directly [2][3][4][5]. In most cases, this is limited to the use of occupied orbitals only via orbitaldependent exchange or dependence on the kinetic energy density, however some models extend this to include virtual orbitals and their energies to allow a more accurate treatment of correlation.
Historically, the progression from local density approximation (LDA) models to generalised gradient CONTACT Andrew M. Teale andrew.teale@nottingham.ac.uk approximation (GGA) models delivered a step change in the accuracy of practical KS calculations. The development of hybrid functionals, that include a proportion of exact (orbital) exchange, has led to further significant improvements in the accuracy of thermochemical predictions. However, although several meta-GGA functionals have become well established [6][7][8][9], their performance in many other applications only amounts to an incremental refinement over that of standard GGAs. Questions remain as to how best local functionals may be constructed to give a higher level of accuracy than that typical of GGAs. The meta-GGAs go beyond the GGA level by including additional dependencies; most meta-GGAs utilise the kinetic energy density τ , however a small number also use the Laplacian of the density 2 ρ. In addition to the best choice of variables, it is not clear under which circumstances one might expect substantial improvements by using meta-GGA type constructions. In an attempt to shed light on the strengths and weaknesses of existing DFAs, a number of groups have developed techniques for the calculation of accurate densityfunctional quantities using systematically improvable ab initio methods. Such approaches include ab initio densityfunctional theory (DFT) [10] in which orbital-dependent exchange-correlation functionals are developed and the corresponding KS equations solved by the optimised effective potential method [11,12], inversion techniques that deliver KS potentials, orbitals and orbital energies [13][14][15], and techniques to compute the adiabatic connection (AC) [16][17][18][19] linking the non-interacting KS system to the physical system [15,20,21].
The success of the KS approach is largely due to the availability and broad applicability of simple, local DFAs. Whilst each of the techniques that employ a combination of ab initio and DFT methods have been useful in elucidating the shortcomings of simple DFAs, it is often difficult to use these results to gain direct insight into how to construct/modify simple local DFAs. In this work, we use the method of Lieb maximisation to calculate accurate energy densities for systems along the adiabatic connection, from which we can derive the coupling constant averaged (CCA) exchange-correlation energy density. This local function can be regarded as a target for modelling by simple DFAs.
We commence in Section 2.1 by defining energy densities in a convenient gauge and outlining their known properties. In Section 2.2, we outline how the energy densities can be determined for specific interaction strengths via Lieb maximisation [22]. We then outline how the CCA exchange-correlation energy density, the relevant target for DFAs, may be defined and calculated in Section 2.3. In Section 3, we describe the specific computational details for the calculation of these energy densities, before presenting results for atomic systems in Section 4.1. We then examine which aspects of these energy densities may be reflected in local quantities used in typical DFAs in Section 4.2. Properties of the local adiabatic connection integrand are examined in Section 4.3. In Section 4.1, we also compare our results with the Becke-Roussel (BR) exchange energy model [6]; a simple DFA in the correct gauge to allow direct comparison with our ab initio results. Finally in Section 5, we make some concluding remarks and discuss directions for future work.

Exchange-correlation energy densities
In this work, we consider the adiabatic connection [16][17][18][19] between the KS and physical systems as a function of the interaction strength λ. The ground state wave function for each λ, λ can be determined as the solution to the Schrödinger equation with the Hamiltonian, Here,T is the kinetic energy operator,Ŵ is the electronelectron interaction operator, scaled by the coupling constant λ and v λ is the external potential that binds the density, ρ λ yielded by λ , at a given λ ࢠ [0, 1], such that ρ λ = ρ 1 . As such at λ = 0, we have the KS system and at λ = 1 the physical interacting system. For a review of the adiabatic connection formalism, see Ref. [23].
A key result from the adiabatic connection formalism is that the exchange-correlation energy E xc [ρ] may be expressed as where the exchange-correlation integrand W xc,λ = λ |Ŵ | λ − J[ρ] is derived using the Hellmann-Feynman theorem and J[ρ] is the classical Coulomb energy. In the second equality, we have introduced the exchange-correlation energy density at coupling constant λ, w xc,λ (r), defined such that W xc,λ = ρ(r)w xc,λ (r)dr. Clearly, this energy density cannot be uniquely defined -it is only defined within a gauge, any function that integrates to zero over all space may be added to this integrand. For further discussion of this point, see Refs [24,25]. In the context of the adiabatic connection however, a natural definition for this energy density is [24,26] w xc,λ (r) = 1 2 where P 2,λ (r, r ) is the pair density which, in contrast to the one-electron density, varies with interaction strength over the adiabatic connection (as discussed in Ref. [23]) and h xc,λ (r, r ) is the exchange-correlation hole, which is also dependent on interaction strength. The energy density may then be expressed in the form w xc,λ (r) = 1 2ρ(r) This definition of the energy density has been widely used in the literature, for example, Refs [27][28][29][30], and has recently been computed for λ = 0, 1 and Ý by Mirtschink et al. [24]. The first term of Equation (5) for λ = 1 may also be identified as the potential v cond (r) discussed in Ref. [31]. The choice of gauge implied by the definition of Equation (3) is particularly appealing for several reasons. First, it is the electrostatic potential of the exchangecorrelation hole for a given interaction strength -this gives a somewhat intuitive picture of exchange and correlation effects. Second, many of the properties of the exchange-correlation energy density in this gauge are well understood -in particular it is known that this w xc (r) → −1/2r as r → Ý for finite systems [32,33]. Third, for λ = 0, the KS system is described by a single Slater determinant, hence the corresponding energy density reduces to the exchange-only energy density (which is independent of λ), where h x (r, r ) is the exchange-hole. As a result, the exchange-correlation energy density may be resolved into exchange and correlation components at any given point along the adiabatic connection, leading to the following definition for correlation energy density: The final and most compelling reason for the choice of gauge in Equation (3) is that the only quantity required to compute the energy density is the λ-dependent pair density P 2,λ (r, r ), which is readily accessible from methods for calculating the adiabatic connection integrand based on ab initio theories.

Computing energy densities from ab initio theory
We now turn our attention to the determination of the quantities, as a function of the coupling constant λ, required to compute the energy densities of Equations (5)- (7). The ground state energy at a given λ in the external potential v may be expressed in the form where the minimisation is over all ensemble density matricesγ containing N electrons. As shown by Lieb [22], For all v-representable input densities, the supremum in this equation becomes a maximum. Recently, it has become possible to perform the maximisation of Equation (10) practically using accurate ab initio methods [15,20,21,34] by choosing an ansatz for the computation of E [v] and expanding the external potential v in a finite basis.
To study the adiabatic connection relevant to the KS DFT, one may then proceed as follows. First, choose an accurate ansatz for the computation of E 1 [v] and determine the corresponding electron density ρ 1 . This provides a good approximation to the physical density at λ = 1. For a series of interaction strengths along the adiabatic connection, the maximisation of Equation (10) is then performed using the fixed input density ρ 1 . Fixing the input density at the physical value naturally leads to the relevant KS reference system at λ = 0 and an optimising potential v 0 = v KS . In general, the optimising potential, v λ , is the Lagrange multiplier function that imposes the density constraint.
Throughout this work, we employ the coupled-cluster singles and doubles with perturbative triples [CCSD(T)] method [35,36] for the determination of E [v]. The implementation of Ref. [21] is used to calculate F λ [ρ] with the finite basis expansion proposed by Wu and Yang [15,37] for where v ext (r) is the standard external potential due to the nuclei, v ref (r) is a fixed reference potential chosen to ensure correct appropriate asymptotic behaviour of v (r) (in these calculations, the Fermi-Amaldi potential [38] was used for this purpose) and {g t (r)} are a set of Gaussian basis functions with expansion coefficients {b t } to be determined through optimisation. Since F λ [ρ] is convex in v, this maximisation can be performed practically by Newton methods as shown by Wu and Yang [15,37]. To employ this approach for general λ and a variety of non-variational wave functions, one can determine the gradient and Hessian of this functional using the relaxed Lagrangian formulation of Helgaker and Jørgensen [39] and second-order response theory as described in Refs [21,34]. Once the maximisation of Equation (10) is complete, the optimising potential v λ is such that ρ λ = ρ 1 . In addition, E λ [v] and the relaxed λ-interacting one and two particle Lagrangian density matrices are readily obtained. As a result, the quantities required for the construction of w xc,λ (r) in Equation (3) can be easily extracted. In this work, we have implemented the construction of this function on a numerical grid. As discussed in Section 4.1, the behaviour of w xc,λ (r) is similar to that of typical exchange correlation potentials, making it well suited to evaluation and integration on standard girds used in densityfunctional calculations. We use the standard Lebedev angular and Lindh-Malmqvist-Gagliardi radial quadrature [40] implemented in the Dalton program [41,42]. For plotting and analysis, we also compute the same functions on evenly spaced Cartesian grids.

Coupling constant averaged energy densities
Whilst the structure of the energy densities at each interaction strength λ is of interest, integration of these functions with the density ρ over space yields only the global adiabatic integrand W xc,λ of Equation (2). This quantity can then be integrated over λ to obtain the total exchange-correlation energy of the system. Here, we instead focus on the CCA energy densityw λ (r). This reduces the calculation of the exchange-correlation energy to a single spatial integration. Using the fact that (ρ λ = ρ 1 ) λ, the CCA energy density can be defined by this CCA energy density may be easily resolved into exchange and correlation parts, since w 0 (r) = w x (r). The functionw xc (r) and its components may be considered as targets for modelling by practical DFAs.
The λ-dependent w xc,λ (r) defines a local adiabatic connection integrand, as opposed to the global quantity W xc,λ in Equation (5). Although this integrand must then be evaluated for each grid point, the local approach offers some important advantages. First, as outlined in Ref. [43], this local quantity can be more amenable to producing size-consistent AC models in the absence of degeneracy, something that is not easily achieved using non-linear interpolation models of the global AC [34,[44][45][46]. Second, the local nature of this energy density is much more akin to energy densities often calculated by conventional DFAs and, although care must be taken to ensure a consistent gauge in comparisons, this information should provide valuable insight to guide to construction of new DFAs. The key quantities provided by the correlation component of this local AC to be considered in this study are illustrated in Figure 1. The integral of w c,λ (r), shown as the shaded area, defines the value of CCA energy densityw c (r). However, further information is encoded in the shape of this local integrand, given that the exchange contribution at each grid point is constant in λ. As a result, the shape of the local integrand contains information about how the energy density departs from exchange-only behaviour at each point in space. In an attempt to quantify this behaviour, we will also consider the initial slope of this plot at λ = 0.

Computational details
The additional elements required for the evaluation of the energy densities described by Equations (3) and (5) and the CCA energy density described by Equation (12) have been implemented in a development version of the Dalton quantum chemistry program [41,42]. In particular, we have modified the DFT module to allow the evaluation of these energy densities on standard quadrature grids and also on evenly spaced Cartesian grids. Here, we consider the application of this method for calculating energy densities to the closed shell atoms He, Be, Ne and Ar. Accurate electron densities ρ 1 for these systems have been computed at the CCSD(T) level of theory using the u-aug-cc-pCVQZ basis set for He, Be and Ne and the u-aug-cc-pCVTZ basis set for Ar, where the prefix u signifies that the basis set has been fully uncontracted. These densities are used as the input for Lieb maximisation calculations at all λ. The Lieb maximisation uses the implementation described in Ref. [34], in which the basis set used to expand the external potential v in Equation (11) is the same as that used to model the molecular orbitals.
In these calculations, the Hessian of F λ [ρ] is numerically regularised via a singular value decomposition with threshold 10 −7 . The calculations were considered converged if the 2-norm of the gradient of F λ [ρ] fell below 10 −6 . To obtain an accurate estimate of the CCA energy densitiesw xc (r), calculations were carried out for the following values of λ: These values were selected to provide sufficient samples over the range λ = 0 → 1 for integration over the coupling constant. In addition, a number of densely spaced values of λ close to 0 were included so that the behaviour of the local adiabatic connection in the noninteracting limit can be examined more closely (see Section 4.3).

... The exchange energy density
The exchange energy density, w x , is independent of λ at each point in space. This function is examined for the He, Be, Ne and Ar atoms in Figure 2. It is well documented [47,48] that the exchange-energy density in atomic systems varies radially from the nucleus in a relatively smooth and uncomplicated way; the data presented in Figure 2 mirror those observations. In the individual cases, the exchange-energy density for He shows a smooth monotonic decay with increasing distance from the nucleus, tending towards the asymptotic limit. Those for Be, Ne and Ar exhibit more varied structures, with features in regions typically associated with the atomic shell structures. However, it is important to note that the magnitude of the exchange-energy density at each point in the system exceeds that of the correlation-energy density by often several orders of magnitude. Therefore, it is essential that DFAs model this component accurately independent of the correlation component if error cancellations are to be avoided. The vast majority of DFAs do not attempt to model the exchange-energy density in the gauge of the exchangehole. Instead typical functionals include GGA contributions in which Laplacian terms have been eliminated using integration by parts and as a result their energy densities are in a different gauge. As a result, they cannot be easily compared with the accurate energy densities in this work, for further discussion of this point, see, for example, Ref. [25]. An alternative approach is taken in the construction of the BR model for the exchange energy. In this approach, the Taylor expansion of the spherically averaged exchange-hole to second-order [32] is employed, which refers to a shell of radius s centred on reference point r. This exact form is utilised to parameterise a model exchange hole derived from a consideration of hydrogenic systems. The exchange-hole curvature for spin σ is given by where τ σ is the everywhere positive-definite local kinetic energy density, The BR exchange-hole model takes the form where the parameters a and b are determined by equating the values of Equations (13) and (16) at s = 0 and equating the curvature Q σ (r) with that of the model in Equation (16). These parameters can be obtained by solving the 1D non-linear equation for x σ = a σ b σ , at each grid point. Setting then the exchange-energy density at a given reference point r is given by and the total exchange energy is given by For closed shell systems as considered throughout this work, w BR x (r) = w BR x,α (r) = w BR x,β (r). The exchange-energy densities w x (r) and w BR x (r), evaluated on the CCSD(T) electron density, are shown in Figure 2 for the atoms He, Be, Ne and Ar. On the scale of w x (r), the BR exchange-energy density is reasonably accurate. By construction, the model has the correct −1/2r asymptotic behaviour (shown by the dashed lines) of the energy density and is close to w x (r) in most regions is shown in more detail below the plots for each system. In general, the BR model gives a too negative energy density near to the nucleus in these systems. Some oscillatory behaviour is also present in this region, arising from the Gaussian basis functions used in this work. For the systems Be, Ne and Ar, the largest differences between the BR model and the accurate energy density occur in regions where the energy density is most rapidly changing, such as features close to inter-shell transitions.
To quantify the quality of the BR model, we consider the exchange energies obtained by integration with the CCSD(T) electronic density in Table 1. The errors in the exchange energy E x are also presented. In general, the BR model is accurate to within 1%-2% for E x . For the He, Be and Ne atoms, the too negative behaviour of w BR x (r) close to the nucleus dominates and leads to overly negative values of E x . For Ar, the behaviour is different, with positive errors in the valence region close to the inter-shell features in the electron density dominating to give an overall too positive value for E x .

... The correlation energy density
The λ-dependent correlation energy densities w c,λ (r) and the CCA correlation energy densityw c (r) are plotted for each of the four atoms He, Be, Ne and Ar in Figures 3-6, respectively. The CCA correlation energy density is computed by numerical integration over λ. It is clear that whilst the magnitude of the correlation energy density changes with λ, its structure is not sensitive to the interaction strength. The structure of the correlation energy densities is clearly more complex than that of the exchange-energy density andw c (r) is the target for DFAs to model. In Section 4.2, we discuss to what extent these features may be recognised using quantities typically used in simple local DFAs. The variation of w c,λ (r) with λ defines a local adiabatic connection at each point in space and this is discussed in Section 4.3.
The KS energy components corresponding to the CCSD(T) electronic density are presented in Table 2. The non-interacting kinetic (T s ), Coulomb (J) exchange (K)   and electron-nuclear (E ne ) components are computed by performing the Lieb maximisation at λ = 0 and using the resulting molecular orbitals. Subtracting these from the total CCSD(T) energy (the value of E 1 [v]) gives the value of the KS correlation energy E c . To confirm the accuracy of the calculated CCA correlation energy densities, we have also evaluated the correlation energy by numerical integration of this quantity with the electronic density, denoted in Table 2 as E w c c . Both estimates of E c agree to within the accuracy expected, confirming the quality of the calculated CCA correlation energy densities.

Comparison with local quantities
Given that the correlation energy densities show a much more complex structure than the corresponding exchange-only energy densities, it is interesting to examine local quantities, such as the density and its derivatives in comparison withw c to investigate how it may be reconstructed by simple DFAs. In this work, we inspected ρ, ρ, τ and 2 ρ and have looked for any correspondences between features in these local quantities and those in w c (r). We have found that, whilst the information contained in the density and its gradient may be sufficient to determine the main features ofw x , point-wise, these alone do not provide enough information for an accurate reconstruction of the more detailed structure ofw c . For the kinetic energy density τ , a number of forms may be considered; the most commonly employed forms are those based on the Laplacian of the orbitals and the (everywhere positive) form based on the gradient of the orbitals. Neither of these forms show features that align well with those in the correlation energy density.
However, comparing 2 ρ with features ofw c (r) is more informative, as the nodes of the Laplacian delimit fluctuations in the energy density. We note that for the full exchange-correlation energy density, it has previously been observed that deviations from the LDA energy density are well correlated with 2 ρ, see, for example, Refs [49,50]. The utility of 2 ρ as an indicator for features of the correlation energy density is illustrated in Figure 7 for the helium, beryllium, neon and argon atoms. Here, the CCA correlation energy density and 2 ρ are plotted together. The vertical dashed lines indicate the positions of nodes in 2 ρ and can be seen to delimit the core and valence regions of the atom. The correlation energy density clearly undergoes transitions moving between these regions.  For the helium atom, the only node in 2 ρ separates the valence region into two parts, this node (and the outer nodes in Be, Ne and Ar) does not correspond directly to shell structure in the density, but does help in determining the point beyond which the CCA energy density starts to smoothly approach its asymptotic limit. For the Be, Ne and Ar atoms, the inner nodes of 2 ρ reflect the shell structure of the atom. In general, the correlation of features in 2 ρ with those inw c (r) is only qualitative, but does show that the CCA energy density has a structure Table . Kohn-Sham electronic energy contributions in Hartree for He, Be, Ne and Ar. The correlation energy is computed by differences between the λ =  and λ =  energies and by integration of the CCA correlation energy density.

Quantity
He Be Ne Ar where v  and v  are the optimising potentials of Equation () at λ =  and λ = , respectively.
that can be sensibly divided according to subtle features of the topology of the density, as captured by 2 ρ. In general, this correlation is more pronounced for higher Z atoms and in core regions, see, for example, the first and second turning points inw c (r) for Ne and Ar. This indicates that some of the seemingly complex structure of the correlation energy density can be understood from a relatively simple quantity depending only on the electronic density, namely 2 ρ. Other measures may be considered for understanding the structure ofw c (r), but are left to future work. In the construction of most standard GGA, functionals dependence on the Laplacian is eliminated by integration by parts, which also changes the gauge of the energy density -making it difficult to compare them with energy densities in the gauge of the exchange-correlation hole. Cancio et al. [51] has suggested that the construction of functionals including the Laplacian should be revisited -noting that the slowly varying limit can still be correctly recovered when the Laplacian is included, but that its inclusion may substantially improve the description of strongly inhomogeneous systems such as atoms, molecules and materials. The results presented here support this idea and suggest that in particular an adequate description of the correlation energy density may be easiest to achieve by introducing a dependence on the Laplacian. The investigation of energy densities for Laplaciandependent correlation models and their performance will be presented in future work.

The local adiabatic connection
The challenge for practical DFAs is to model the CCA energy density. We have seen that the Laplacian of the density shows promise as a local quantity that can be used in DFAs, with a clear correspondence to the structure of w c . To gain further insight into the nature of the correlation at each point in space, we now consider the local AC.
As representative examples, we consider the Ne atom, where the correlation is expected to be mainly dynamic in nature, and Be where static correlation effects play a role. For the global integrand W xc,λ , it is understood that as the nature of the correlation shifts from predominantly dynamic to more static in character, the curvature of this integrand becomes more and more pronounced. We now consider the behaviour of the local integrands in these systems.
In Figures 8 and 9, we presentw c (r) and the local adiabatic connection at several points in space for Be and Ne, respectively. These points are marked on thew c (r) plots by a coloured symbol, with the corresponding local AC curve plotted with the same marker in the right-hand panel. The first feature observed for both species is that the corresponding local ACs do not vary in magnitude (i.e. from bottom to top of the right-hand panels of Figures 8 and 9) sequentially with variation in r. This is because of the oscillatory nature ofw c (r) and the fact that the area between each local AC curve and the x-axis in the right-hand panel represents the value ofw c (r) in  the left-hand panel. More interestingly, we see that the local shape of the integrands varies in an intuitive way with r. For both systems, the local AC in the core region is strongly linear in λ, indicating that the nature of the correlation is relatively simple. For the Ne atom, the local AC is reasonably linear for all r values, indicating that the correlation is relatively simple throughout the system. As a result, the main challenge for a DFA is then to capture the structure ofw c (r) as discussed in the previous section. For the Be atom, the behaviour ofw c (r) is more complex; the first two points in the core and inner valence regions show simple, relatively linear local ACs.
However, the remaining two points in the valence and outer valence regions show much more pronounced curvature -indicating a more complex higher order dependence on λ, in line with the expected global behaviour for a system exhibiting static correlation effects. However, in the present example, it can be seen that these more complex effects are essentially localised to the valence and outer valence regions.
This connection between the nature of the correlation and the location in the system should be of significant utility in the construction and testing of DFAs. In particular, it suggests the possibility of a topological-DFT in which DFA energy densities are constructed by taking into account both the local λ-variation and the position in space in the system considered.
As indicated in the schematic of Figure 1, one could also consider the initial slope of the local AC as a measure of how rapidly the local model must depart from exchange-only behaviour at each point in space. In Figure 10, we present this derivative for each of the atomic systems considered. In general, the shape of these functions closely resembles the energy density itself. This arises because a largerw c is obtained when the area enclosed by w c is larger, generally involving a more negative initial-slope. For a more insightful understanding of the nature of the correlation using derivatives at λ = 0, it may then be necessary to proceed to higher orders. These data should be of particular interest for constructing DFAs via interpolation models based on the local AC and will be explored in future work.

Conclusions
We have presented an approach for calculating CCA energy densities in the gauge of the exchange-correlation hole. These quantities may be considered as a target for modelling by DFAs, defined in this gauge. The exchangeenergy density is relatively smooth with minor features close to the inter-shell structure in atomic densities. Furthermore, we have confirmed that the BR model, which defines an exchange-energy density in the gauge of the exchange-correlation hole, delivers a reasonable description of this component of the energy density.
The correlation component of the energy density shows much more complex structure. Interestingly, this structure aligns with the nodes of 2 ρ, suggesting that this quantity should be of greater use in the construction of DFAs. This is in agreement with recommendations by Cancio et al. [51] in studies of the exchangeonly DFAs, based observations of the full exchangecorrelation energy density in Silicon [50]. Further insight from this work is offered by the ability to resolve the energy density into its exchange and correlation components separately, revealing the more elaborate structure of the latter and its alignment with features of 2 ρ, suggesting its possible role as a key variable in the formulation of new DFAs for the more accurate treatment of dynamical correlation.
A local adiabatic connection at each point in space is defined by the λ-dependence of the energy densities, and the observations for Be and Ne suggest that this local connection captures the nature of the correlation effects. The area enclosed by the local AC gives the value of the energy density at each point in space and, of course, varies substantially across the system. However, the shape of the local AC still retains extra information about the nature of the correlation; higher order dependence on the interaction strength introduces more rapid curvature in the local AC. Observations for the Be and Ne atoms show that such effects can be localised to different regions of space, with less pronounced curvature in the core regions and more pronounced curvature in the valence regions. This feature is more pronounced in Be than the other atoms, consistent with the near degeneracy effects present in this system. Knowledge of the nature of the local AC should be useful for constructing DFAs via local interpolation models, which offer size-consistent models for the exchange-correlation energy in the absence of degeneracy [52,53].
A particularly interesting prospect arising from this work is the exploration of this connection between spatial resolution of the correlation energy and its dependence on the interaction strength at each point in space. The results here suggest that the development of a topological-DFT, in which the λ-dependence of the functional is explicitly connected to the topology of the local density, may be fruitful. Such an approach would allow more careful adaptation of the λ-dependence of DFAs to the regions in which they are applied. Furthermore, our results suggest that 2 ρ should be a key variable for identifying appropriate regions.
A number of areas for future work are now being explored. These include the generalisation of this procedure to deal with generalised range-dependent ACs [54][55][56], for example, with the widely used error-function attenuated interactions [57]; the application of this method to a wide range of molecular systems to provide benchmark data for the development of DFAs; the testing of existing models for exchange-correlation energy densities defined in the gauge of the exchange-correlation hole and the development of new models in this context. We expect that the methodology outlined here will be an essential tool for these developments.