Kinetic Correlation Functionals from the Entropic Regularization of the Strictly Correlated Electrons Problem.

In this work, we study the entropic regularization of the strictly correlated electrons formalism, discussing the implications for density functional theory and establishing a link with earlier works on quantum kinetic energy and classical entropy. We carry out a very preliminary investigation (using simplified models) on the use of the solution of the entropic regularized problem to build approximations for the kinetic correlation functional at large coupling strengths. We also analyze lower and upper bounds to the Hohenberg-Kohn functional using the entropic regularized strictly correlated electrons problem.

Introduction

Despite all their successes, present approximations for the exchange–correlation (XC) functional of the Kohn–Sham (KS) density functional theory (DFT) are still plagued by the so-called strong correlation problem:[1] typically, the approximations fail when the physics of the system under study differs too much from the noninteracting one of the KS reference system.

The leading term of the strong-coupling limit of the DFT adiabatic connection (strictly correlated electrons (SCE) functional), equivalent to a semiclassical limit (ℏ → 0) at fixed one-electron density, gives access to the exact XC functional in the extreme case that the kinetic energy is neglected with respect to the electron–electron interactions.[2−5] This strictly correlated regime is complementary to the one described by the noninteracting KS system. By applying uniform-coordinate scaling, one sees that this limit captures the right physics for low-density systems, i.e., when the average electron–electron distance is much larger than the Bohr radius.[6,7] Indeed, when used as an approximation for the XC functional in the self-consistent KS scheme, SCE provides results that get closer and closer to the exact ones as the system is driven to lower and lower density.[8−11] However, with the exception of interesting models for electrons confined at the interface of semiconductor heterostructures,[9,10,12−14] chemical systems are never close to this extreme case. Yet, the SCE mathematical structure can be simplified and rescaled to design functionals for the electron–electron interaction at physical coupling strength[15−17] or can be used to build interpolations between the KS and the SCE limits.[18−29] While these strategies are both very promising, as, for example, they can describe accurately the H2 and H2+ dissociation curves in the KS spin-restricted formalism,[15] their main problem is that they do not capture the effects of the kinetic correlation energy, which is known to play a crucial role in the description of strongly correlated systems in the KS setting,[30−33] with its functional derivative displaying nonintuitive features such as “peaks” and “steps”.[30−32,34−36]

The next leading term in the strong-coupling expansion, corresponding to zero-point oscillations in a metric dictated by the density,[37] provides a “first-order” kinetic correlation energy correction,[38] but it is difficult to evaluate in the general case, with its functional derivative displaying features that are too extreme.[39] Moreover, this way to do the strong-coupling expansion is not the right one for problems such as bond-breaking excitations because in a molecular system, the density close to the atoms remains high: only when we drive the whole system to low density, the expansion is really able to capture the right physics.[40] The purpose of this work is to explore a different route, based on the entropic regularization of optimal transport,[41−44] which has been studied in mathematics and economics, as well as, more recently, has been applied in data sciences and statistical inference (see, for instance, ref (44) and references therein).

The OT formulation of the SCE functional[2,3] triggered cross-fertilization between two different research fields, which led to several formal proofs, setting the SCE limit on firm grounds,[4,5,45,46] as well as to new ideas and algorithms.[47−51] Here, we focus on the entropic regularization of the SCE problem[47,52,53] and explore whether this extension can be used to build approximations for the kinetic correlation energy functional and, more generally, to gain new insight into the problem of describing and understanding strong correlation within DFT. As we will explain, the entropic regularization of the SCE problem brings in a new link and perspective on the seminal work of Sears, Parr, and Dinur[54] on the relation between various definitions of entropy, information theory, and kinetic energy. Moreover, the formalism is quite general and could also be applied to other interactions and other kinds of particles, for example, if one wants to treat the nuclei in a quantum DFT setting.[55]

This paper is organized as follows: in Section , we introduce the theoretical aspects and describe the general form of the solution of the entropic regularization of the SCE functional. To illustrate its main properties, we present simple analytical and numerical examples in Section . We then compare, in Section , the entropic regularized SCE functional with the Hohenberg–Kohn functional, discussing inequalities and approximations, with the corresponding numerical and analytical studies in Section . Conclusions and future perspectives are discussed in Section .

Entropic Regularization of the SCE Functional

Let ρ(x), with , be a density such that . The SCE functional is defined asi.e., as the infimum, over all possible fermionic wave functions having the prescribed density ρ, of the expectation value of the electron–electron repulsion operatorWe have an infimum in eq because the minimum is attained not on the space of wave functions Ψ (with ) but on the larger space of probability measures (in physicists/chemists language, by allowing also Dirac delta distributions).[3,4] We denote probability measures as γ(x1, ..., x). In a loose way, we identifyeven if γ lives in a larger space (i.e., it is allowed to become a distribution). To illustrate what is meant, consider the simple case of N = 2 and D = 3. Then, the minimizer of eq has been proven[2,3] to be always of the SCE form[56,57]which is zero everywhere except on the three-dimensional (3D) manifold x2 = f(x1), parametrized by the co-motion function (or optimal map) , with the position of the first electron dictating the position of the second (strict correlation). Note that the SCE functional has been recently proven to yield the asymptotic low-density (or strong-coupling, or ℏ → 0) limit of the universal Hohenberg–Kohn (HK) functional.[3−5]

On the one hand, the fact that the infimum in eq is attained on a probability measure (i.e., γSCE is concentrated on a low-dimensional manifold of the full configuration space) is exactly what makes the SCE mathematical structure and its density dependence much more accessible than the HK functional.[8−10,15,21,57,58] On the other hand, the challenge of including the effects of kinetic correlation energy stems exactly from the fact that γSCE has infinite kinetic energy. We know that in the exact HK functional, even when very close to the SCE limit, kinetic energy will “spread out” a little bit the optimal γ, making it a true |Ψ|2. The zero-point energy (ZPE) expansion gives a recipe for this spreading out, but, as mentioned, in a rather complicated way.[37−39] Here, we consider a particular definition of entropy, used in the OT as a computational regularization, to realize this “spreading”.

Since it has been proven[4,5] that the fermionic statistics has no effect on the value of VeeSCE[ρ], we work directly in terms of γ(x1, ..., x), which has the loose sense of eq . We then consider the following minimization problemwhere the “entropic” functional Eτ[γ] is defined for τ > 0 aswithWe stress that the entropy term is defined on the set of signed measures such that ∫γ = 1, and it is defined as S[γ] = −∫γ log γ if γ is a probability density and S[γ] = +∞ otherwise. These conditions force the probability measures to be a probability density γτ in and not a Dirac delta on a manifold as, for example, γSCE of eq , since minus S[γSCE] would be equal to +∞. The constraint γ → ρ reads explicitlywhere the notation means that we do not integrate over the variable x.

We point out that the problem (eq ), typically with N = 2 and vee(x, y) in eq equal to the p-distance |x – y| (p ≥ 1), is being studied in different fields, including probability theory (e.g., refs (59, 60)), machine learning (e.g., refs (42, 44)), scientific computing,[61] statistical physics,[62,63] and economics.[43] In the following, we want to analyze the entropic regularization (eq ) in the framework of the DFT formalism.

First, we remark that the problem (eq ) admits a unique solution γτ since the functional Eτ[γ] is strictly convex in γ. Second, this unique solution can be fully characterized. In fact, as shown, for instance, in refs (59, 60) and[64], γτ is the solution of (eq ) if and only ifwhere is the so-called entropic weight and is fixed by the density constraintThe entropic weight aτ(x) can be written as an exponential of the entropic one-body potential uτ(x)with uτ(x) having the usual physical interpretation of DFT, as (minus) the potential that enforces the density constraint. The theorems behind eqs –2.12 are nontrivial, and we point to ref (60) for a rigorous proof in the case of bounded interactions vee and to Appendix for more details on how this potential appears as the dual variable with respect to the density, as in standard DFT. Here, to provide an intuitive idea of the role of the entropic weight, we consider the problem (eq ) in a box and we minimize Eτ[γ] with respect to γ without fixing the density constraint, obtaining the usual result, i.e., that γτ is a Gibbs stateThis clearly shows that the entropic weight is a Lagrange multiplier to enforce the constraint γ → ρ in eq . The solution γτ in eq can then be written asWe should remark at this point that the one-body potential uτ(x) is not gauged to approach zero when |x| → ∞, but it is shifted by a constant Cτ[ρ]When τ → 0, this constant ensures thatThis way, we see that γτ→0 of eq becomes increasingly more concentrated on the manifold where Vee(x1, ..., x) – ∑u0(x) is minimum and equal to 0. We can interpret Vee(x1, ..., x) – ∑u0(x) as a hamiltonian without kinetic energy, whose minimizing wave function is constrained to yield the given density ρ by the one-body potential u0(x). In fact, this is the hamiltonian that appears as a leading term in the strong-coupling limit of the usual density-fixed DFT adiabatic connection,[57,65] whose minimizing γ (if we relax the space in which we search for the minimum) will be zero everywhere except on the manifold where Vee(x1, ..., x) – ∑u0(x) has its global minimum. This is exactly the SCE manifold parameterized by the co-motion functions.

Note that the constant C0[ρ] = limτ→0Cτ[ρ] is precisely the same,[66] in the strong-coupling limit of DFT, as the one discussed by Levy and Zahariev in the context of KS DFT.[67] In fact, since the potential u0(x) is gauged at infinity to a constant that guarantees that the minimum of Vee(x1, ..., x) – ∑u0(x) is equal to zero, and since the optimal γτ→0 will be concentrated on the manifold where the minimum is attained, by simply taking the expectation value of Vee(x1, ..., x) – ∑u0(x) on γτ→0, we obtainMoreover, we also have that u0 is a functional derivative with respect to ρ (gauged to a constant at infinity) of VeeSCE[ρ].[8,9] If we use VeeSCE[ρ] as an approximation for the Hartree and exchange–correlation energy, as in the KS SCE approach,[8−10] then eq is exactly the condition imposed by Levy and Zahariev[67] to their constant shift.

Interpretation of the Parameter τ and of the Entropy S[γ]

One can simply regard τ > 0 as a parameter interpolating between two opposite regimes: the strictly correlated one (τ → 0) and the uncorrelated bosonic case (τ → ∞) with the prescribed density.

In fact, when τ → 0, the problem (eq ) becomes the one defined by the SCE functional of eq ,[53] and, as just discussed, γτ, given by eq , in this limit is increasingly more concentrated on the manifold on which Vee(x1, ..., x) – ∑u0(x) = 0. In the case N = 2, this is exactly the three-dimensional manifold {x1 = x, x2 = f(x)} parametrized by the co-motion function (or optimal map) f(x) of eq . To visualize this, in Figure , we show a simple example with N = 2 particles in one dimension (1D), having a Gaussian density, whose interaction is repulsive harmonic. In panel (a) of this figure, we show γτ→0(x1, x2), which is concentrated on the manifold x2 = f(x1), where for this special case f(x) = −x. For N > 2, we usually (but not always) also have a three-dimensional manifold parametrized by cyclical maps f(x).[52,57]

Figure 1

Optimal γ(x1, x2) for the interaction vee(x, y) = −(x – y)2, at different values of τ. Note that the marginals , remain the same at all τ, while γ evolves from the strictly correlated regime , with f(x) = −x, to the symmetric uncorrelated one . See Section for a fully analytical description of this example.

When τ → ∞, the problem (eq ) converges to the one of maximizing S[γ] alone under the constraint γ → ρ,This is equivalent to maximize the entropy of γ relative to the product state . In fact, with ρ̃(x) = ρ(x)/N,Since the density is held fixed, the second term in the last line is a constant during the maximization. Gibbs inequality applied to the relative entropy (first term in the last line) then gives , and the optimal γ that maximizes S[γ] is then the uncorrelated product state. Equation also shows that the entropy S[γ] with fixed one-electron density is a relative entropy (Kullback–Leibler divergence) with respect to the uncorrelated product, a.k.a. noninteracting bosonic state with the prescribed density. In other words, at fixed density, the uncorrelated product is the probability density whose support has the maximal volume. This is illustrated, again in the simple 1D case with repulsive harmonic interactions, in Figure c, where we also show, in panel (b), a case in between these two extremes.

The problem (eq ) has been already used as an auxiliary functional to compute numerically the solutions of eq . In fact, the entropy term reinforces the uniqueness of the minimizer in eq . The parameter τ in this case regularizes the problem of eq (“spreading out” the support of γ, as in Figure ), which can be solved via the Sinkhorn algorithm.[42,61]

We should emphasize that, as eq clearly shows, the entropy S[γ] used here is different from the quantum mechanical entropy of finite-temperature DFT (see refs (68−71) and references therein), which is defined in terms of density matrices and favors mixed states. Here, S[γ] can be interpreted in terms of mutual information (or discrimination information), measuring how a probability γ differs from a reference distribution, in this case, the uncorrelated product. A related definition and interpretation in terms of the Kullback–Leibler divergence, including its link to kinetic energy, was considered by Sears, Parr, and Dinur[54] in the context of DFT. The link between various definitions of entropy and kinetic energy is also present in several works in the literature; in particular, the link with the kinetic correlation energy is conjectured in ref (72).

Before comparing the functional Fentrτ[ρ] with the Hohenberg–Kohn functional close to the strong-coupling regime, we find it important to illustrate the formalism just introduced with simple examples.

Analytic and Numerical Examples of the Entropic Regularization Problem

Harmonic Interactions Case

We start by considering the repulsive and attractive harmonic interaction vee(x, y) = ξ(x – y)2, with ξ = ±1. This interaction is interesting not only because it allows for analytic solutions with which one can fully illustrate the formalism, but also because it arises as a leading term in the effective interaction between electrons bound on two different distant neutral fragments (dispersion). In fact, if we keep the densities of the two fragments frozen at their isolated ground-state values (a variational constraint that has several computational advantages and can lead to very accurate or even exact results[73]), minimizing the dipolar interaction, which contains terms like x1x2 orthogonal to the bond axis and −z1z2 parallel to it, is equivalent to minimizing the repulsive and attractive harmonic interaction, respectively. This is simply because ±x1x2 differs from ∓1/2(x1 – x2)2 only by one-body terms, which do not affect the minimizer when the density is held fixed. Another case in which harmonic interactions could be interesting is if we want to treat (some) nuclei quantum mechanically.

(a) N = 2

To allow for a completely analytic solution, we fix the one-body density to be a Gaussian. This is exactly the Drude quantum oscillator model for the coarse-grained dispersion between two fragments[74,75] when we forbid the oscillator density to change with respect to its isolated value (a constraint that gives the exact result for the dispersion coefficient C6 between two oscillators, exactly like in the case of the H atom[73]). Since the dipolar interaction separates in the three spatial directions, we can consider the one-dimensional case withIn the following, we use the notation x = x1 and y = x2 for the coordinates of the two particles in 1D. By writing and dividing both sides of eq by aτ(x), we see that eq becomes, after writing As previously discussed, if we find the explicit form for uτ(x) that satisfies eq , then we have found the optimal one. We then first assume that the solution uτ can be restricted to a class of second-degree polynomialsand verify that indeed it is possible to obtain a solution of this kind, which amounts to solving the system of equationswhich yields, choosing the negative solutionDefining , the corresponding minimizing γτ(x, y) readsand it is shown at different values of τ (with σ = 1 and ξ = −1) in Figure , where, as anticipated in Section , we see the transition from the SCE-like state at small τ to the uncorrelated product state at large τ

(b) N > 2, D = 1

In the case when N > 2, the first equation of the system in eq readswithBy arguing similarly as in the previous paragraph, one can obtain that the solution of the equationis given by

Regularized Coulomb Interaction Case

For illustrative purposes, we now consider a 1D problem with N = 2 electrons interacting via the effective Coulomb repulsion, , which has been shown in ref (76) to yield results that closely mimic the 3D electronic structure. In Section , we will also consider another 1D interaction, with a long-range Coulomb tail, finding results qualitatively very similar. We fix the density to beThe reason to choose this particular density is that it has an exponential decay at large x (similar to an atomic density) and allows for an analytic solution in the SCE case.[38] For the entropic regularization case, however, the solution of the system of eq cannot be obtained analytically, and therefore, we have computed it numerically via the Sinkhorn algorithm[42] (POT library[77]). In Figure , we report our results for the support of γτ, as τ increases: in panel (a), corresponding to a small value of τ, we see that γτ(x1, x2) is different from zero only very close to the manifold parametrized by the co-motion function, x2 = f(x1), thus becoming a very good approximation for . We also show, as a tiny red line, the co-motion function f(x) computed analytically[38] from the SCE theory. Panel (c) corresponds to a relatively high value of τ, and we see that γτ is approaching the uncorrelated bosonic product , losing any resemblance with the SCE state. The central panel (b) is for us the most interesting: the system is still close to the SCE state, but it has a significant spreading, which could be used to approximate the quantum system close to (but not at) the SCE limit, mimicking the role of kinetic energy. We will explore this possibility in the next two sections.

Figure 2

Support of the optimal γτ(x1, x2) for the interaction , at different values of τ for the density (eq ). The red line represents the co-motion function x2 = f(x1).

Comparison with the Hohenberg–Kohn Functional

In this section, we compare the entropic functional Fentrτ[ρ] of eq with the Hohenberg–Kohn[78] functional (HK) in its extension of Levy and Lieb[79,80] as a constrained minimization problem, generalized to arbitrary coupling strength λ ≥ 0withand Vee[γ] defined in eq . Note that F0[ρ] is the Kohn–Sham functional and F1[ρ] is the Hohenberg–Kohn functional at physical coupling strength. In particular, we are interested in exploring how the large-λ expansion of λ–1Fλ[ρ] compares with the entropic functional Fentrτ[ρ] at small τ. We already know that the two limits are equalbut we want to compare how they behave when approaching the SCE limit, slightly spreading out the optimal γ into a |Ψ|2 around the SCE manifold as in Figure b. To begin with, we briefly recall how Fλ[ρ] behaves at large λ, namely[37,38,57]where FZPE[ρ] is the zero-point energy functional. Similarly to the functional S[γ], the zero-point oscillations performed by the N particles around the manifold parametrized by the co-motion functions (optimal maps) f(x) allow for the corresponding probability density γZPE to provide a finite kinetic energy. Calling the Hessian matrix of the function Vee(x1, ..., x) – ∑u0(x) evaluated on the manifold {x1 = x, x2 = f2(x), ..., x = f(x)}, the two functionals in eq can be written explicitly as[37]In particular, due to the virial theorem, we can write the λ-dependent expectation value of the electron–electron interaction and of the kinetic energy operator at large λwhere Ψλ[ρ] is the minimizer of eq . We should stress that, while for the leading term in eq there are rigorous mathematical proofs,[4,5] the term of order is a very plausible conjecture,[37] which has been confirmed numerically in some simple cases.[38]

Inequalities and Approximations

First of all, as shown in ref (52), as a simple consequence of the logarithmic Sobolev inequality for the Lebesgue measure,[81] it holdsHowever, this entropic lower bound to the HK functional can be very loose, as we will show in Figures and 4 with some numerical examples. We also havesimply because this way we have added a positive quantity to the value of Vee[γ] obtained with the γτ that minimizes eq .

Figure 3

Exact expansion of the solution to eq for the repulsive harmonic interaction and a Gaussian density as a function of ϵ = λ–1. A = ϵG1/ϵτ=ϵ[ρ], , C = ϵF1/ϵ[ρ], , E = Fentr(π/2)ϵ[ρ], F = VeeSCE. See text in Section for a detailed explanation.

Figure 4

Exact expansion of the solution to eq for different interactions and density ρC. A = ϵG1/ϵτ=ϵ[ρ], , C = VeeSCE[ρ] + ϵFZPE[ρ], , E = Fentr(π/2)ϵ[ρ]. The horizontal line represents VeeSCE[ρ]. The numerical value of a in eq reads, respectively, a = 0.27 (top) and a = 0.32 (bottom). See text in Section for a detailed explanation.

A route we explore in this work is the use of the γτ[ρ] from an entropic calculation at finite τ to compute an approximate many-body kinetic energy in the λ → ∞ limitwhere γτ[ρ] is the optimum in the problem (eq ) with the given ρ. Since γτ has the explicit form eq (in terms of the entropic potential uτ(x) that needs to be computed numerically), we obtainObviously, γτ will not have the right nodal surface and will miss the fermionic character. However, the fermionic statistics is expected[6,37] to appear in Fλ[ρ] at large λ only through orders , a conjecture that was supported by numerical evidence.[38] The idea is to use the large-λ functional as an approximation for the Hartree exchange–correlation functional so that the fermionic character will be captured by the KS kinetic energy, similarly to the KS SCE scheme.[8−10,50] More generally, we will analyze the functional Gλτ[ρ] defined aswith γτ[ρ] the minimizer of Fentrτ[ρ]. As a consequence of the variational principle, we have for the special case of a N = 2 closed-shell systemHowever, for N > 2, the inequality will not be valid in general, as does not have the right fermionic antisymmetry. We still expect it to hold for large λ with τ ∝ λ–1/2, where the energetic difference between fermionic and bosonic minimizers should become exponentially small,[38] of orders . In Section , we provide a first explorative study into different ways to find an optimal relation between τ and λ, to make Gλτ[ρ] as close as possible to Fλ[ρ]. Note that by looking at eq , one may expect that Fentrτ[ρ] diverges as 1/τ2 for small τ. However, the divergence is milder because when τ → 0, the integrand in eq tends to zero, as γτ→0 is increasingly more concentrated on the manifold where Vee(x1, ..., x) – ∑u0(x) is minimum (and stationary, i.e., where its gradient, contained in the modulus square inside the integrand, is zero). We believe that Fentrτ[ρ] diverges only as 1/τ for small τ, implying that τ should be proportional to λ–1/2 to match the large-λ expansion of the HK functional, a conjecture that seems to be confirmed by our analytical and numerical results in Section . However, we have no rigorous proof for this statement.

Analytical and Numerical Investigation

In Section , a specific relation between τ and λ was used to establish a rigorous inequality, namely, eq , which holds ∀λ ≥ 0 when . The question we want to address here is whether for a given λ (and in particular for large λ), the inequality (eq ) can be sharpened into an equality by tuning τ according to a general dependence τ(λ). We thus look for τ that solvesAlthough this equation can probably be always solved, at least for large λ, the real question is whether we can find a reasonably accurate general approximation for the relation between τ and λ, as, obviously, we do not want to compute the exact HK functional each time to determine the proper τ(λ). Here, we make a very preliminary numerical and analytic exploration, which supports the already conjectured relation τ(λ) ∼ λ–1/2 at large λ. Finding an approximate τ(λ) that is generally valid, however, remains for now an open challenge, which requires further investigations.

Repulsive Harmonic Interaction

Equation can be solved explicitly for the example discussed in Section , where N = 2, the density is a Gaussian and the electron–electron interaction is repulsive harmonic. In fact, we start by noting that the exact wave function minimizing Fλ[ρ] with repulsive harmonic electron–electron interaction and a Gaussian density has the form (see, e.g., the appendix of ref (57))whileimplying that γτ can always be mapped to γexactλ by settingThis implies that, being γτ essentially of the exact form for this specific case, we can just evaluate the functional and minimize it with respect to the coefficients Aτ, Bτ. The constraint γτ → ρ impliesEquation reads thenand we obtain the optimal Bτ as a function of λ by settingThe only positive solution, Bτ(λ), provides the answer. In fact, direct comparison of eq with eq shows thatorEquation ,5.9b has the following asymptotic expansionsconfirming τ(λ) ∼ λ–1/2 for λ → ∞, as discussed in Section . Both series at small and large λ will have a finite radius of convergence since the function τ(λ), eq ,5.9b, has several branch cuts. The exact τ(λ) of eq ,5.9b can be very accurately represented as λ–1 plus a correction in the form of a simple Padé approximant that interpolates between the two limits of eq In Figure , we compare, as a function of ϵ = λ–1, the exact HK functional ϵF1/ϵ[ρ] (curve labeled “C”) with the results obtained from the functional Gλτ[ρ] of eq by using for τ(λ) different approximations. In the curve labeled “A”, we have used the λ → 0 leading term of eq , τ(λ) = λ–1, and in the curve labeled “B”, we have used the λ → ∞ leading term, τ(λ) = λ–1/2. We see that, this way, we approximate Fλ[ρ] at different correlation regimes. We also show in the same figure the left-hand side of the inequality (eq ) when we set , which was found in the inequality (eq ), curve labeled “D”. As it should, this curve stays above the value of VeeSCE[ρ] (horizontal line, labeled “F”), but, in this case, it also stays below the HK functional, which is a nice feature, although probably peculiar to the harmonic interaction (see Section ). We also show the right-hand side of the inequality (eq ) (curve labeled “E”), which, as anticipated, is a very loose lower bound. The result obtained by using the Padé approximant τPad(λ) of eq in Gλτ[ρ] is, on the scale of Figure , indistinguishable from the exact curve.

Effective Coulomb Interaction

For an interaction that mimics the electron–electron repulsion in quasi-1D systems, there is no analytical computation available. As anticipated in Section , we resort to the Sinkhorn algorithm to obtain the quantities of interest and repeat the computation just done for the harmonic cost. We tested two different interaction forms for vee, namely, a regularized Coulomb interaction and the exponential interaction already used in Section to compute γτ at various regimeswith the same density of eq . In Figure , we compare Gλτ(λ)[ρ], using different approximations for τ(λ), with the expansion (curve labeled C), which for N = 2 electrons in 1D has been shown[38] to approximate very accurately the exact HK functional at large λ. The analogous of eq cannot be derived analytically, but we use for the asymptotics of τ(λ) at high couplings the dependence discussed in Section and confirmed in eq , namelyand we optimize a to match the expansion of the HK functional. We get a ≈ 0.27 for veereg(x) and a very similar value, a ≈ 0.32, for veeexp(x). The curve labeled B shows the corresponding Gλτ(λ)[ρ] when we set τ(λ) equal to eq . In the curve labeled A, we have simply set τ = λ–1, which was the small-λ expansion found for the harmonic interaction case. We also show in the same figure the left-hand side of the inequality (eq ) when we set , which was found in the inequality (eq ), curve labeled D. As it should, this curve stays above the value of VeeSCE[ρ] (horizontal line, labeled F), but, contrary to the harmonic case of Figure , this time, this curve does not stay below the λ-dependent HK functional. We also show the right-hand side of the inequality (eq ) (curve labeled E), which again is found to be a very loose lower bound.

Conclusions and Outlook

In this paper, we introduced and studied structural properties of a new class of density functionals based on the entropic regularization of the SCE functional. Although the entropic regularization of the OT-SCE problem has been previously used as a numerical tool to compute the SCE energy via the Sinkhorn algorithm, here we have investigated whether it could also provide a route to build and study approximations of the Hohenberg–Kohn functional at large coupling constant λ. We have first focused on the link between the (classical) entropy with fixed marginals used here, the quantum kinetic energy, and the Kullback–Leibler divergence, with links to the seminal work of Sears, Parr, and Dinur,[54] and with other recent works in the same spirit.[72,82−89]

We have performed a very preliminary investigation on whether the minimizing wave function of the regularized SCE entropic problem, which has an explicit form, could be used to estimate the kinetic energy. A more extensive investigation is needed, to assess whether it is possible to find an approximate general relation between τ and λ, at least for large λ. We conjectured here, and we have numerical evidence in very simple cases, that when λ → ∞, it holds τ ∼ aλ–1/2, with a probably a density-dependent constant.

We should remark that, from a computational viewpoint, a challenging problem is to face the very unfavorable scaling with respect to the number of electrons (marginals) N of the Sinkhorn algorithm when solving the entropic-SCE problem.[90] This implies that to provide functionals for routine applications, we might need to construct approximations inspired to the mathematical form of eq , similar to what has been done for the leading SCE term.[15−17,21,91] To this purpose, it will be essential to further study properties of uτ at small τ, also comparing and testing it as a candidate for the Hartree exchange–correlation potential.