ORIGINAL ARTICLE: Sara Giarrusso and Aurora Pribram-Jones, Møller–Plesset and Density-Fixed Adiabatic Connections for a Model Diatomic System at Different Correlation Regimes, J. Chem. Theory Comput. 2023.
DOI: https://doi.org/10.1021/acs.jctc.3c00529
Article Abstract: In recent years, adiabatic connection (AC) interpolations developed within density functional theory (DFT) have been found to provide good performances in the calculation of interaction energies when used with Hartree–Fock (HF) ingredients. The physical and mathematical reasons for such unanticipated performance have been clarified, to some extent, by studying the strong-interaction limit of the Møller–Plesset (MP) AC. In this work, we calculate both the MP and the DFT AC integrand for the asymmetric Hubbard dimer, which allows for a systematic investigation of different correlation regimes by varying two simple parameters in the Hamiltonian: the external potential, Δν, and the interaction strength, U. Notably, we find that, while the DFT AC integrand appears to be convex in the full parameter space, the MP integrand may change curvature twice. Furthermore, we discuss different aspects of the second-order expansion of the correlation energy in each AC, and we demonstrate why the derivative of the λ-dependent density in the MP AC at
λ=0
(i.e., at the HF density) is zero in the model. Concerning the strong-interaction limit of both ACs in the Hubbard dimer setting, we show that the asymptotic value of the MP AC,
HF
W
∞
, is lower than (or equal to) its DFT analogue,
KS
W
∞
, if the two are compared at a given density, just like in real space. However, we also show that this is not always the case if the two quantities are compared at a given external potential.

Summary

This Mathematica notebook provides the analytical expressions and numerical calculations in support of the accompanying manuscript:
​“Møller-Plesset and density-fixed adiabatic connections for a model diatomic system at different correlation regimes” https://t.co/9KueSjCjdy
In this work, we provide a detailed analysis of the Møller−Plesset (MP) and Density-Fixed Adiabatic Connections (ACs) for a model diatomic system at different correlation regimes: the asymmetric Hubbard dimer. This model allows to control the correlation regime by changing the parameters U and Δv in the Hamiltonian.
AC methods rely on the idea of gradually switching from a formally non-interacting Hamiltonian, which is comparatively simple, to a “fully interacting” one, which is more complicated and describes the system of interest, by multiplying the interaction operator by an interaction-strength parameter.
The main motivation for this work is that in recent years, AC interpolations developed within density functional theory (DFT), that interpolate between the weak-interaction and the strong-interaction limits of the AC integrand, have been found to provide good performances in the calculation of interaction energies when used with Hartree–Fock (HF) ingredients. Such unanticipated performance has been clarified, to some extent, by studying the strong-interaction limit of the MP AC, which uses the HF Hamiltonian and the HF state as reference (i.e. when λ=0). However, many things still need to be better understood: for instance, the shape of the MP AC integrand is quite different than that of the DFT AC which is reckoned to be at least piecewise convex. The main result we report is that, while the DFT AC integrand appears to be convex for any value of the parameters
{U,Δν}
(in line with expectations), the MP integrand shows a double change of curvature for a continuous range of Δv, at any given
U>0
.
Authors: Sara Giarrusso and Aurora Pribram-Jones
​Keywords: electronic-structure theories, density-functional theory, adiabatic connection, asymmetric Hubbard dimer, strong-interaction limit.
Readers are invited to reproduce, and extend, the results contained in this notebook. To fully understand the notebook, reading of the manuscript is recommended.
Each cell must be evaluated sequentially to ensure the correct functionality.

MaTeX

This subroutine allows to use
L
A
T
E
X
style in plots’ labels.
In[]:=
Needs["MaTeX`"]​​SetOptions[MaTeX,"Preamble"->{"\\usepackage{amssymb,amsmath,latexsym,amsfonts,amsthm,xcolor,bm,mhchem,physics}"}];

Interacting Hamiltonian and site-occupation operators of the asymmetric Hubbard dimer

The asymmetric Hubbard dimer is a model system made up of only two sites, labelled 0 and 1. Its Hamiltonian in second quantization reads:

H
=-t
Σ
σ=,

†

a
0σ

a
1σ
+
†

a
1σ

a
0σ
+U
Σ
i=0,1

n
i↑

n
i↓
+
Σ
i=0,1
v
i

n
i
with
†

a
,

a
the usual creation and annihilation operators, σ the index for the spin variable,

n
iσ
=
†

a
iσ

a
iσ
and

n
i
=
Σ
σ=,

n
iσ
. t, U and
{
v
i
}
are positive parameters. We expand the Hamiltonian in the site basis for two-particle singlet states: |0↑0↓〉, |1↑1↓〉, and
1
2
(|0↑1↓〉-|0↓1↑〉) to form the following matrix
In[]:=
H3D[U_,t_,v0_,v1_]=U+2v0,0,-
2
t,0,U+2v1,-
2
t,-
2
t,-
2
t,v0+v1;
We further define (v1-v0)=Δv=V, while the parameter t is set to 1/2 throughout the work.
The site-occupation operators

n
i
in the site basis read
In[]:=
n03D=DiagonalMatrix[{2,0,1}];​​n13D=DiagonalMatrix[{0,2,1}];

Ground State (GS) quantities

GS energy

In[]:=
ϵ[U_,V_]=
2
3
U-
3+
2
U
+3
2
V
Cos
1
3
ArcCos
U(9+2
2
U
-18
2
V
)
2
3/2
(3+
2
U
+3
2
V
)
;
The above energy expression assumes that the gauge of the potential is v1+v0=0. It corresponds to the lowest eigenvalue of matrix H3D[U,1/2,-V/2,V/2].

GS solution


Interacting site occupation difference


Møller-Plesset adiabatic connection

This section contains the operators and ingredients used to build the MP AC Hamiltonian (Eq. 36 of the paper)
HF

H
λ
=-t
Σ
σ=,

†

a
0σ

a
1σ
+
†

a
1σ

a
0σ
+
Σ
i=0,1

v
i

n
i
+λU
Σ
i=0,1

n
i↑

n
i↓
-
HF
n
i
2

n
i
​​=-t
Σ
σ=,

†

a
0σ

a
1σ
+
†

a
1σ

a
0σ
+λU
Σ
i=0,1

n
i↑

n
i↓
+
Σ
i=0,1
λ,HF
v
i

n
i
The λ-dependent GS associated to
HF

H
λ
is labelled
HF
Ψ
λ
.

RHF quantities


Møller-Plesset fluctuation potential and
HF
c
0

In[]:=
Umat3d[U_]=
U
0
0
0
U
0
0
0
0
;
In[]:=
MPFP3D[U_,x_]=Umat3d[U]-
U
2
1-
x
2
n03D+1+
x
2
n13D;
The matrix MPFP3D corresponds to the operator
U
Σ
i=0,1

n
i↑

n
i↓
-
HF
n
i
2

n
i
, which is the one multiplied by the interaction parameter λ in
HF

H
λ
and is therefore called MP fluctuation (or perturbation) potential operator.
In[]:=
c0hf[U_,x_]=
U
2
1+
2
x
4
;
The quantity c0hf is an energetic shift corresponding to

HF
Φ
U
Σ
i=0,1

n
i↑

n
i↓

HF
Φ

with
HF
Φ
=
HF
Ψ
0
.

AC integrand as a function of external potential

In[]:=
Wλhf[U_,V_,λ_]=Psi0[λU,Vλhfgen[U,V,λ]].MPFP3D[U,NhfvsV[U,V]].Psi0[λU,Vλhfgen[U,V,λ]]+c0hf[U,NhfvsV[U,V]];

AC integrand as a function of HF density (
"Nhf"=x=
HF
Δn
)

In[]:=
WλhfvsNhf[U_,x_,λ_]=-
U
3
+
1
8
U(4+
2
x
)+24U
3+
2
U
2
λ
+3
2
x
2
1
4-
2
x
+
Uλ
2
-64
4-
2
x
+
3
U
8
x
3
λ
+2U
6
x
λ3+U
4-
2
x
λ-6
2
U
2
λ
+16U
4
x
λ-6-2U
4-
2
x
λ+3
2
U
2
λ
+16
2
x
3
4-
2
x
+18Uλ+6
2
U
4-
2
x
2
λ
-4
3
U
3
λ
Cos
1
6
π+2ArcCos-
Uλ9+2
2
U
2
λ
-18
2
x
2
1
4-
2
x
+
Uλ
2
2
3/2
3+
2
U
2
λ
+3
2
x
2
1
4-
2
x
+
Uλ
2

5/2
48+12U
2
x
4-
2
x
λ+
2
U
(16+8
2
x
-3
4
x
)
2
λ

1-
2
U
2
λ
2
9+2
2
U
2
λ
-18
2
x
2
1
4-
2
x
+
Uλ
2
4
3
3+
2
U
2
λ
+3
2
x
2
1
4-
2
x
+
Uλ
2
-
U2Uλ+3
2
x
1
4-
2
x
+
Uλ
2
Sin
1
6
π+2ArcCos-
Uλ9+2
2
U
2
λ
-18
2
x
2
1
4-
2
x
+
Uλ
2
2
3/2
3+
2
U
2
λ
+3
2
x
2
1
4-
2
x
+
Uλ
2

3
3+
2
U
2
λ
+3
2
x
2
1
4-
2
x
+
Uλ
2
;

Density-fixed adiabatic connection

This section contains the operators and ingredients used to build the DFT AC Hamiltonian (Eq.40 of the paper)
KS

H
λ
=-t
Σ
σ=,

†

a
0σ

a
1σ
+
†

a
1σ

a
0σ
+λU
Σ
i=0,1

n
i↑

n
i↓
+
λ,KS
Δv
2


n
1
-

n
0

The λ-dependent GS associated to
KS

H
λ
is labelled
KS
Ψ
λ
.

Lieb Maximisation

In[]:=
VλLM[U_,Nn_,λ_]:=ReturnV/.NMaximizeϵ[λU,V]-
V
2
Nn,V[[2]]
The function VλLM corresponds to
λ,KS
Δv
: the λ-dependent external potential that delivers the prescribed density Δn (=”Nn”) at all λ, i.e. such that
Δn=
KS
Ψ
λ


n
1
-

n
0

KS
Ψ
λ

.

AC integrand as a function of external potential

In[]:=
PsiλKS[U_,V_,λ_]:=Return[Psi0[λU,VλLM[U,nn[U,V],λ]]]
In[]:=
c0KSvsV[U_,V_]=-
U
2
1+
2
nn[U,V]
2
;
In[]:=
WλKS[U_,V_,λ_]:=PsiλKS[U,V,λ].Umat3d[U].PsiλKS[U,V,λ]+c0KSvsV[U,V]

AC integrand as a function of interacting density (
"N"=x=
HF
Δn
)

In[]:=
PsiλKSvsN[U_,x_,λ_]:=Return[Psi0[λU,VλLM[U,x,λ]]]
In[]:=
c0KSvsN[U_,x_]=-
U
2
1+
2
x
2
;
In[]:=
WλKSvsN[U_,x_,λ_]:=PsiλKSvsN[U,x,λ].Umat3d[U].PsiλKSvsN[U,x,λ]+c0KSvsN[U,x]

Correlation energies, their separate contributions and their weak- and strong-interaction expansions [see also J. Chem. Phys. 157, 054102 (2022)]

In this section we report the total and separate contributions to the correlation energy and we analyse their weak-(small) and strong-(large) interaction expansions. Note that, for each of the 3 following ingredients,
HF
I
c
(λU,x)=
HF
I
c,λ
(U,x)
with I=U,T,V

HF
U
c
[eq (53a)]


HF
T
c
[eq (53b)]


HF
V
c
[eq (53c)]


Figure 3

Figure 6

Figure 7

Comparing MP and DFT AC integrands [Figure 1]

Relative error vs accuracy predictor [Figure 4]

Calculation of the tables is not immediate (takes more than few minutes)

Appendix

AC integrand [Figure 9]

Global Correlation energy [Figure 10]

CITE THIS NOTEBOOK

Møller–Plesset and Density-Fixed Adiabatic Connections for a Model Diatomic System at different correlation regimes​
by Sara Giarrusso and Aurora Pribram-Jones
Wolfram Community, STAFF PICKS, September 5, 2023
​https://community.wolfram.com/groups/-/m/t/3005864