one last thing

This commit is contained in:
Pierre-Francois Loos 2019-04-16 22:13:41 +02:00
parent a2ca1ee6c4
commit 408f25b528

View File

@ -51,6 +51,7 @@
\newcommand{\X}{\text{X}} \newcommand{\X}{\text{X}}
\newcommand{\UEG}{\text{UEG}} \newcommand{\UEG}{\text{UEG}}
\newcommand{\HF}{\text{HF}} \newcommand{\HF}{\text{HF}}
\newcommand{\ROHF}{\text{ROHF}}
\newcommand{\LDA}{\text{LDA}} \newcommand{\LDA}{\text{LDA}}
\newcommand{\PBE}{\text{PBE}} \newcommand{\PBE}{\text{PBE}}
\newcommand{\FCI}{\text{FCI}} \newcommand{\FCI}{\text{FCI}}
@ -144,7 +145,7 @@ As illustrative examples, we show how this density-based correction allows us to
\maketitle \maketitle
%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%
\section{Introduction} %\section{Introduction}
%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%
Contemporary quantum chemistry has developed in two directions --- wave function theory (WFT) \cite{Pop-RMP-99} and density-functional theory (DFT). \cite{Koh-RMP-99} Contemporary quantum chemistry has developed in two directions --- wave function theory (WFT) \cite{Pop-RMP-99} and density-functional theory (DFT). \cite{Koh-RMP-99}
Although both spring from the same Schr\"odinger equation, each of these philosophies has its own \textit{pros} and \textit{cons}. Although both spring from the same Schr\"odinger equation, each of these philosophies has its own \textit{pros} and \textit{cons}.
@ -156,7 +157,7 @@ One of the most fundamental drawbacks of conventional WFT methods is the slow co
This undesirable feature was put into light by Kutzelnigg more than thirty years ago. \cite{Kut-TCA-85} This undesirable feature was put into light by Kutzelnigg more than thirty years ago. \cite{Kut-TCA-85}
To palliate this, following Hylleraas' footsteps, \cite{Hyl-ZP-29} Kutzelnigg proposed to introduce explicitly the interelectronic distance $r_{12} = \abs{\br{1} - \br{2}}$ to properly describe the electronic wave function around the coalescence of two electrons. \cite{Kut-TCA-85, KutKlo-JCP-91, NogKut-JCP-94} To palliate this, following Hylleraas' footsteps, \cite{Hyl-ZP-29} Kutzelnigg proposed to introduce explicitly the interelectronic distance $r_{12} = \abs{\br{1} - \br{2}}$ to properly describe the electronic wave function around the coalescence of two electrons. \cite{Kut-TCA-85, KutKlo-JCP-91, NogKut-JCP-94}
The resulting F12 methods yields a prominent improvement of the energy convergence, and achieve chemical accuracy for small organic molecules with relatively small Gaussian basis sets. \cite{Ten-TCA-12, TenNog-WIREs-12, HatKloKohTew-CR-12, KonBisVal-CR-12} The resulting F12 methods yields a prominent improvement of the energy convergence, and achieve chemical accuracy for small organic molecules with relatively small Gaussian basis sets. \cite{Ten-TCA-12, TenNog-WIREs-12, HatKloKohTew-CR-12, KonBisVal-CR-12}
For example, at the CCSD(T) level, it is advertised that one can obtain quintuple-$\zeta$ quality correlation energies with a triple-$\zeta$ basis, \cite{TewKloNeiHat-PCCP-07} although computational overheads are introduced by the large auxiliary basis used to resolve three- and four-electron integrals. \cite{BarLoo-JCP-17} For example, at the CCSD(T) level, one can obtain quintuple-$\zeta$ quality correlation energies with a triple-$\zeta$ basis, \cite{TewKloNeiHat-PCCP-07} although computational overheads are introduced by the large auxiliary basis used to resolve three- and four-electron integrals. \cite{BarLoo-JCP-17}
To reduce further the computational cost and/or ease the transferability of the F12 correction, approximated and/or universal schemes have recently emerged. \cite{TorVal-JCP-09, KonVal-JCP-10, KonVal-JCP-11, BooCleAlaTew-JCP-2012, IrmHumGru-arXiv-2019, IrmGru-arXiv-2019} To reduce further the computational cost and/or ease the transferability of the F12 correction, approximated and/or universal schemes have recently emerged. \cite{TorVal-JCP-09, KonVal-JCP-10, KonVal-JCP-11, BooCleAlaTew-JCP-2012, IrmHumGru-arXiv-2019, IrmGru-arXiv-2019}
Present-day DFT calculations are almost exclusively done within the so-called Kohn-Sham (KS) formalism, which corresponds to an exact dressed one-electron theory. \cite{KohSha-PR-65} Present-day DFT calculations are almost exclusively done within the so-called Kohn-Sham (KS) formalism, which corresponds to an exact dressed one-electron theory. \cite{KohSha-PR-65}
@ -168,13 +169,13 @@ In the context of the present work, one of the interesting feature of density-ba
Progress toward unifying WFT and DFT are on-going. Progress toward unifying WFT and DFT are on-going.
In particular, range-separated DFT (RS-DFT) (see Ref.~\onlinecite{TouColSav-PRA-04} and references therein) rigorously combines these two approaches via a decomposition of the electron-electron (e-e) interaction into a smooth long-range part and a (complementary) short-range part treated with WFT and DFT, respectively. In particular, range-separated DFT (RS-DFT) (see Ref.~\onlinecite{TouColSav-PRA-04} and references therein) rigorously combines these two approaches via a decomposition of the electron-electron (e-e) interaction into a smooth long-range part and a (complementary) short-range part treated with WFT and DFT, respectively.
As the WFT method is relieved from describing the short-range part of the correlation hole around the e-e coalescence points, the convergence with respect to the one-electron basis set is greatly improved. \cite{FraMusLupTou-JCP-15} As the WFT method is relieved from describing the short-range part of the correlation hole around the e-e coalescence points, the convergence with respect to the one-electron basis set is greatly improved. \cite{FraMusLupTou-JCP-15}
Therefore, a number of approximate RS-DFT schemes have been developed using either single-reference \cite{AngGerSavTou-PRA-05, GolWerSto-PCCP-05, TouGerJanSavAng-PRL-09,JanHenScu-JCP-09} or multi-reference \cite{LeiStoWerSav-CPL-97, FroTouJen-JCP-07, FroCimJen-PRA-10, HedKneKieJenRei-JCP-15, FerGinTou-JCP-18} WFT approaches. Therefore, a number of approximate RS-DFT schemes have been developed within single-reference \cite{AngGerSavTou-PRA-05, GolWerSto-PCCP-05, TouGerJanSavAng-PRL-09,JanHenScu-JCP-09} or multi-reference \cite{LeiStoWerSav-CPL-97, FroTouJen-JCP-07, FroCimJen-PRA-10, HedKneKieJenRei-JCP-15, FerGinTou-JCP-18} WFT approaches.
Very recently, a major step forward has been taken by some of the present authors thanks to the development of a density-based basis set correction for WFT methods. \cite{GinPraFerAssSavTou-JCP-18} Very recently, a major step forward has been taken by some of the present authors thanks to the development of a density-based basis set correction for WFT methods. \cite{GinPraFerAssSavTou-JCP-18}
The present work proposes an extension of these new methodological development together with the first numerical tests on molecular systems. The present work proposes an extension of this new methodological development together alongside the first numerical tests on molecular systems.
%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%
\section{Theory} %\section{Theory}
%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%
The present basis set correction relies on the RS-DFT formalism to capture the missing part of the short-range correlation effects, a consequence of the incompleteness of the one-electron basis set. The present basis set correction relies on the RS-DFT formalism to capture the missing part of the short-range correlation effects, a consequence of the incompleteness of the one-electron basis set.
Here, we only provide the main working equations. Here, we only provide the main working equations.
@ -199,24 +200,24 @@ is the basis-dependent complementary density functional, $\hT$ is the kinetic op
In Eq.~\eqref{eq:E_funcbasis}, $\wf{}{\Bas}$ and $\wf{}{}$ are two general $\Ne$-electron wave functions belonging to the Hilbert space spanned by $\Bas$ and the complete basis, respectively. In Eq.~\eqref{eq:E_funcbasis}, $\wf{}{\Bas}$ and $\wf{}{}$ are two general $\Ne$-electron wave functions belonging to the Hilbert space spanned by $\Bas$ and the complete basis, respectively.
Both wave functions yield the same target density $\n{}{}$. Both wave functions yield the same target density $\n{}{}$.
Importantly, in the limit of a complete basis set (which we refer to as $\Bas \to \infty$), we have, for any density $\n{}{}$, $\lim_{\Bas \to \infty} \bE{}{\Bas}[\n{}{}] = 0$, which implies that Importantly, in the limit of a complete basis set (which we refer to as $\Bas \to \infty$), we have, for any density $\n{}{}$, $\lim_{\Bas \to \infty} \bE{}{\Bas}[\n{}{}] = 0$.
This implies that
\begin{equation} \begin{equation}
\label{eq:limitfunc} \label{eq:limitfunc}
\lim_{\Bas \to \infty} \qty( \E{\modY}{\Bas} + \bE{}{\Bas}[\n{\modZ}{\Bas}] ) = \E{\modY}{} \approx E, \lim_{\Bas \to \infty} \qty( \E{\modY}{\Bas} + \bE{}{\Bas}[\n{\modZ}{\Bas}] ) = \E{\modY}{} \approx E,
\end{equation} \end{equation}
where $\E{\modY}{}$ is the energy associated with the method $\modY$ in the complete basis set. where $\E{\modY}{}$ is the energy associated with the method $\modY$ in the complete basis set.
In the case where $\modY = \FCI$ in Eq.~\eqref{eq:limitfunc}, we have a strict equality as $\E{\FCI}{} = \E{}{}$. In the case where $\modY = \FCI$ in Eq.~\eqref{eq:limitfunc}, we have a strict equality as $\E{\FCI}{} = \E{}{}$.
Provided that the functional $\bE{}{\Bas}[\n{}{}]$ is known exactly, the only source of error at this stage lies in the potential approximate nature of the methods $\modY$ and $\modZ$.% for the FCI energy and density within $\Bas$, respectively. Provided that the functional $\bE{}{\Bas}[\n{}{}]$ is known exactly, the only source of error at this stage lies in the potential approximate nature of the methods $\modY$ and $\modZ$.
The functional $\bE{}{\Bas}[\n{}{}]$ is obviously \textit{not} universal as it depends on $\Bas$. The functional $\bE{}{\Bas}[\n{}{}]$ is obviously \textit{not} universal as it depends on $\Bas$.
Moreover, as $\bE{}{\Bas}[\n{}{}]$ aims at fixing the incompleteness of $\Bas$, its main role is to correct Moreover, as $\bE{}{\Bas}[\n{}{}]$ aims at fixing the incompleteness of $\Bas$, its main role is to correct
for the lack of cusp in $\wf{}{\Bas}$ (i.e.~a discontinuous derivative) at the e-e coalescence points, a universal condition of exact wave functions. for the lack of cusp in $\wf{}{\Bas}$ (i.e.~a discontinuous derivative) at the e-e coalescence points, a universal condition of exact wave functions.
Because the e-e cusp originates from the divergence of the Coulomb operator at $r_{12} = 0$, a cuspless wave function could equivalently originate from a Hamiltonian with a non-divergent two-electron interaction at $r_{12} = 0$. Because the e-e cusp originates from the divergence of the Coulomb operator at $r_{12} = 0$, a cuspless wave function could equivalently originate from a Hamiltonian with a non-divergent two-electron interaction at coalescence.
Therefore, as we shall do later on, it feels natural to approximate $\bE{}{\Bas}[\n{}{}]$ with short-range density functionals which deal with a smooth long-range electron interaction. Therefore, as we shall do later on, it feels natural to approximate $\bE{}{\Bas}[\n{}{}]$ with short-range density functionals where one deals with smooth long-range interactions.
Contrary to the conventional RS-DFT scheme which requires a range-separated \textit{parameter} $\rsmu{}{}$, here we use a range-separated \textit{function} $\rsmu{\Bas}{}(\br{})$ which automatically adapts to quantify the incompleteness of $\Bas$ in $\mathbb{R}^3$. Contrary to the conventional RS-DFT scheme which requires a range-separation \textit{parameter} $\rsmu{}{}$, here we use a range-separation \textit{function} $\rsmu{\Bas}{}(\br{})$ that automatically adapts to quantify the incompleteness of $\Bas$ in $\mathbb{R}^3$.
The first step of the present basis set correction consists of obtaining an effective two-electron interaction $\W{\Bas}{}(\br{1},\br{2})$ ``mimicking'' the Coulomb operator in a incomplete basis $\Bas$. The first step of the present basis set correction consists of obtaining an effective two-electron interaction $\W{\Bas}{}(\br{1},\br{2})$ ``mimicking'' the Coulomb operator in an incomplete basis $\Bas$.
%The present definition ensures that $\W{\Bas}{}(\br{1},\br{2})$ is finite at the e-e coalescence point as long as an incomplete basis set is used, and tends to the genuine, unbounded $r_{12}^{-1}$ Coulomb interaction as $\Bas \to \infty$.
In a second step, we shall link $\W{\Bas}{}(\br{1},\br{2})$ to $\rsmu{\Bas}{}(\br{})$. In a second step, we shall link $\W{\Bas}{}(\br{1},\br{2})$ to $\rsmu{\Bas}{}(\br{})$.
In the final step, we employ short-range density functionals \cite{TouGorSav-TCA-05} with $\rsmu{\Bas}{}(\br{})$ as range separation. In the final step, we employ short-range density functionals \cite{TouGorSav-TCA-05} with $\rsmu{\Bas}{}(\br{})$ as range separation.
@ -257,7 +258,7 @@ Because Eq.~\eqref{eq:int_eq_wee} can be rewritten as
\iint r_{12}^{-1} \n{2}{}(\br{1},\br{2}) \dbr{1} \dbr{2} = \iint \W{\Bas}{}(\br{1},\br{2}) \n{2}{}(\br{1},\br{2}) \dbr{1} \dbr{2}, \iint r_{12}^{-1} \n{2}{}(\br{1},\br{2}) \dbr{1} \dbr{2} = \iint \W{\Bas}{}(\br{1},\br{2}) \n{2}{}(\br{1},\br{2}) \dbr{1} \dbr{2},
\end{equation} \end{equation}
it intuitively motivates $\W{\Bas}{}(\br{1},\br{2})$ as a potential candidate for an effective interaction. it intuitively motivates $\W{\Bas}{}(\br{1},\br{2})$ as a potential candidate for an effective interaction.
Note that the divergence condition of $\W{\Bas}{}(\br{1},\br{2})$ in Eq.~\eqref{eq:def_weebasis} ensures that one-electron systems are free of correction as the present approach must only correct the basis set incompleteness error originates from the e-e cusp. Note that the divergence condition of $\W{\Bas}{}(\br{1},\br{2})$ in Eq.~\eqref{eq:def_weebasis} ensures that one-electron systems are free of correction as the present approach must only correct the basis set incompleteness error originating from the e-e cusp.
As already discussed in Ref.~\onlinecite{GinPraFerAssSavTou-JCP-18}, $\W{\Bas}{}(\br{1},\br{2})$ is symmetric, \textit{a priori} non translational, nor rotational invariant if $\Bas$ does not have such symmetries. As already discussed in Ref.~\onlinecite{GinPraFerAssSavTou-JCP-18}, $\W{\Bas}{}(\br{1},\br{2})$ is symmetric, \textit{a priori} non translational, nor rotational invariant if $\Bas$ does not have such symmetries.
A key quantity is the value of the effective interaction at coalescence of opposite-spin electrons A key quantity is the value of the effective interaction at coalescence of opposite-spin electrons
@ -266,14 +267,12 @@ A key quantity is the value of the effective interaction at coalescence of oppos
\W{\Bas}{}(\br{}) = \W{\Bas}{}(\br{},{\br{}}), \W{\Bas}{}(\br{}) = \W{\Bas}{}(\br{},{\br{}}),
\end{equation} \end{equation}
which is necessarily \textit{finite} for an incomplete basis set as long as the on-top pair density $\n{2}{}(\br{})$ is non vanishing. which is necessarily \textit{finite} for an incomplete basis set as long as the on-top pair density $\n{2}{}(\br{})$ is non vanishing.
%Of course, there exists \textit{a priori} an infinite set of functions in $\mathbb{R}^6$ satisfying \eqref{eq:int_eq_wee}, but
Thanks to its definition one can show that (see Appendix B of Ref.~\onlinecite{GinPraFerAssSavTou-JCP-18}) Thanks to its definition one can show that (see Appendix B of Ref.~\onlinecite{GinPraFerAssSavTou-JCP-18})
\begin{equation} \begin{equation}
\label{eq:lim_W} \label{eq:lim_W}
\lim_{\Bas \to \infty}\W{\Bas}{}(\br{1},\br{2}) = r_{12}^{-1}\ \lim_{\Bas \to \infty}\W{\Bas}{}(\br{1},\br{2}) = r_{12}^{-1}\
\end{equation} \end{equation}
for any $(\br{1},\br{2})$ such that $\n{2}{}(\br{1},\br{2}) \ne 0$.% and for any $\wf{}{\Bas}$, which guarantees a physically satisfying limit. for any $(\br{1},\br{2})$ such that $\n{2}{}(\br{1},\br{2}) \ne 0$.% and for any $\wf{}{\Bas}$, which guarantees a physically satisfying limit.
%An important point here is that, with the present definition of $\W{\Bas}{}(\br{1},\br{2})$, one can quantify the effect of the incompleteness of $\Bas$ on the Coulomb operator itself as a removal of the divergence of the two-electron interaction near the electron coalescence.
%================================================================= %=================================================================
%\subsection{Range-separation function} %\subsection{Range-separation function}
@ -290,13 +289,11 @@ such that the long-range interaction
\w{}{\lr,\rsmu{\Bas}{}}(\br{1},\br{2}) = \frac{1}{2} \qty{ \frac{\erf[ \rsmu{\Bas}{}(\br{1}) r_{12}]}{r_{12}} + \frac{\erf[ \rsmu{\Bas}{}(\br{2}) r_{12}]}{ r_{12}} } \w{}{\lr,\rsmu{\Bas}{}}(\br{1},\br{2}) = \frac{1}{2} \qty{ \frac{\erf[ \rsmu{\Bas}{}(\br{1}) r_{12}]}{r_{12}} + \frac{\erf[ \rsmu{\Bas}{}(\br{2}) r_{12}]}{ r_{12}} }
\end{equation} \end{equation}
coincides with the effective interaction at coalescence, i.e.~$\w{}{\lr,\rsmu{\Bas}{}}(\br{},\br{}) = \W{\Bas}{}(\br{})$. coincides with the effective interaction at coalescence, i.e.~$\w{}{\lr,\rsmu{\Bas}{}}(\br{},\br{}) = \W{\Bas}{}(\br{})$.
%\PFL{This expression looks like a cheap spherical average.
%What about $\rsmu{\Bas}{}(\br{1},\br{2}) = \sqrt{\rsmu{\Bas}{}(\br{1}) \rsmu{\Bas}{}(\br{2})}$ and a proper spherical average to get $\rsmu{\Bas}{}(r_{12})$?}
%================================================================= %=================================================================
%\subsection{Short-range correlation functionals} %\subsection{Short-range correlation functionals}
%================================================================= %=================================================================
Once defined a range separation function $\rsmu{\Bas}{}(\br{})$, we can use RS-DFT functionals to approximate $\bE{}{\Bas}[\n{}{}]$. Once $\rsmu{\Bas}{}(\br{})$ is defined, it can be used within RS-DFT functionals to approximate $\bE{}{\Bas}[\n{}{}]$.
As in Ref.~\onlinecite{GinPraFerAssSavTou-JCP-18}, we consider here a specific class of short-range correlation functionals known as ECMD whose general definition reads \cite{TouGorSav-TCA-05} As in Ref.~\onlinecite{GinPraFerAssSavTou-JCP-18}, we consider here a specific class of short-range correlation functionals known as ECMD whose general definition reads \cite{TouGorSav-TCA-05}
\begin{multline} \begin{multline}
\label{eq:ec_md_mu} \label{eq:ec_md_mu}
@ -329,7 +326,7 @@ with $\hWee{\lr,\rsmu{}{}} = \sum_{i<j} \w{}{\lr,\rsmu{}{}}(r_{ij})$.
% \w{}{\lr,\rsmu{}{}}(r_{12}) = \frac{\erf(\rsmu{}{} r_{12})}{r_{12}}. % \w{}{\lr,\rsmu{}{}}(r_{12}) = \frac{\erf(\rsmu{}{} r_{12})}{r_{12}}.
%\end{equation} %\end{equation}
%is the long-range part of the Coulomb operator. %is the long-range part of the Coulomb operator.
The ECMD functionals admit, for any density $\n{}{}(\br{})$, the following two limiting forms The ECMD functionals admit, for any $\n{}{}(\br{})$, the following two limiting forms
\begin{subequations} \begin{subequations}
\begin{align} \begin{align}
\label{eq:large_mu_ecmd} \label{eq:large_mu_ecmd}
@ -340,18 +337,18 @@ The ECMD functionals admit, for any density $\n{}{}(\br{})$, the following two l
\end{align} \end{align}
\end{subequations} \end{subequations}
where $\Ec[\n{}{}(\br{})]$ is the usual universal correlation functional defined in KS-DFT. where $\Ec[\n{}{}(\br{})]$ is the usual universal correlation functional defined in KS-DFT.
The choice of the ECMD in the present scheme is motivated by the analogy between the definition of $\bE{}{\Bas}[\n{}{}]$ [Eq.~\eqref{eq:E_funcbasis}] and the ECMD functionals [Eq.~\eqref{eq:ec_md_mu}]. The choice of the ECMD in the present scheme is motivated by the analogy between the definition of $\bE{}{\Bas}[\n{}{}]$ [Eq.~\eqref{eq:E_funcbasis}] and the ECMD functional [Eq.~\eqref{eq:ec_md_mu}].
Indeed, provided that $\w{}{\lr,\rsmu{\Bas}{}}(\br{1},\br{2}) = \W{\Bas}{}(\br{1},\br{2})$, then $\wf{}{\rsmu{\Bas}{}}$ and $\wf{}{\Bas}$ coincide. Indeed, provided that $\w{}{\lr,\rsmu{\Bas}{}}(\br{1},\br{2}) = \W{\Bas}{}(\br{1},\br{2})$, then $\wf{}{\rsmu{\Bas}{}}$ and $\wf{}{\Bas}$ coincide.
%The ECMD functionals differ from the standard RS-DFT correlation functional by the fact that the reference is not the KS Slater determinant but a multi-determinantal wave function. %The ECMD functionals differ from the standard RS-DFT correlation functional by the fact that the reference is not the KS Slater determinant but a multi-determinantal wave function.
%This makes them particularly well adapted to the present context where one aims at correcting a general WFT method. %This makes them particularly well adapted to the present context where one aims at correcting a general WFT method.
Therefore, we approximate $\bE{}{\Bas}[\n{}{}]$ by the ECMD functionals evaluated with the range-separation function $\rsmu{\Bas}{}(\br{})$. Therefore, we approximate $\bE{}{\Bas}[\n{}{}]$ by ECMD functionals evaluated with the range-separation function $\rsmu{\Bas}{}(\br{})$.
The LDA version of the ECMD complementary functional is defined as The LDA version of the ECMD complementary functional is defined as
\begin{equation} \begin{equation}
\label{eq:def_lda_tot} \label{eq:def_lda_tot}
\bE{\LDA}{\Bas}[\n{}{}(\br{}),\rsmu{\Bas}{}(\br{})] = \int \be{\LDA}{\sr}\qty(\n{}{}(\br{}),\rsmu{\Bas}{}(\br{})) \n{}{}(\br{}) \dbr{}, \bE{\LDA}{\Bas}[\n{}{}(\br{}),\rsmu{\Bas}{}(\br{})] = \int \be{\LDA}{\sr}\qty(\n{}{}(\br{}),\rsmu{\Bas}{}(\br{})) \n{}{}(\br{}) \dbr{},
\end{equation} \end{equation}
where $\be{\LDA}{\sr}(\n{}{},\rsmu{}{})$ is the reduced (i.e.~per electron) ECMD functional of the uniform electron gas (UEG) \cite{LooGil-WIRES-16} parametrized in Ref.~\onlinecite{PazMorGorBac-PRB-06}. where $\be{\LDA}{\sr}(\n{}{},\rsmu{}{})$ is the reduced (i.e.~per electron) ECMD functional of the uniform electron gas (UEG) \cite{LooGil-WIRES-16} parametrized in Ref.~\onlinecite{PazMorGorBac-PRB-06}.
The short-range LDA correlation functional relies on the transferability of the physics of the UEG which is certainly valid for large $\mu$ but is known to over correlate for small $\mu$. The short-range LDA correlation functional relies on the transferability of the physics of the UEG which is certainly valid for large $\mu$ but is known to over correlate for small $\mu$.
In order to correct such a defect, we propose here a new PBE ECMD functional In order to correct such a defect, we propose here a new PBE ECMD functional
\begin{equation} \begin{equation}
@ -374,11 +371,11 @@ This represents a major computational saving without loss of performance as we e
Depending on the functional choice, the complementary functional $\bE{}{\Bas}[\n{\modZ}{}]$ is then equal to $\bE{\LDA}{\Bas}[\n{\modZ}{}(\br{}),\rsmu{\Bas}{}(\br{})]$ or $\bE{\PBE}{\Bas}[\n{\modZ}{}(\br{}),\rsmu{\Bas}{}(\br{})]$ where $\rsmu{\Bas}{}(\br{})$ is given by Eq.~\eqref{eq:mu_of_r}. Depending on the functional choice, the complementary functional $\bE{}{\Bas}[\n{\modZ}{}]$ is then equal to $\bE{\LDA}{\Bas}[\n{\modZ}{}(\br{}),\rsmu{\Bas}{}(\br{})]$ or $\bE{\PBE}{\Bas}[\n{\modZ}{}(\br{}),\rsmu{\Bas}{}(\br{})]$ where $\rsmu{\Bas}{}(\br{})$ is given by Eq.~\eqref{eq:mu_of_r}.
%================================================================= %=================================================================
%\subsection{Valence approximation} %\subsection{Frozen-core approximation}
%================================================================= %=================================================================
As most WFT calculations are performed within the frozen-core (FC) approximation, it is important to define an effective interaction within a subset of MOs. As most WFT calculations are performed within the frozen-core (FC) approximation, it is important to define an effective interaction within a subset of MOs.
We then naturally split the basis set as $\Bas = \Cor \bigcup \BasFC$, where $\Cor$ is the set of core MOs, and define the FC version of the effective interaction as We then naturally split the basis set as $\Bas = \Cor \bigcup \BasFC$ (where $\Cor$ is the set of core MOs) and define the FC version of the effective interaction as
\begin{equation} \begin{equation}
\W{\Bas}{\FC}(\br{1},\br{2}) = \W{\Bas}{\FC}(\br{1},\br{2}) =
\begin{cases} \begin{cases}
@ -410,25 +407,21 @@ Defining $\n{\modZ}{\FC}$ as the FC (i.e.~valence-only) one-electron density obt
%================================================================= %=================================================================
%\subsection{Computational considerations} %\subsection{Computational considerations}
%================================================================= %=================================================================
One of the most computationally intensive task of the present approach is the evaluation of $\W{\Bas}{}(\br{})$ [see Eqs.~\eqref{eq:wcoal}] at each quadrature grid point. The most computationally intensive task of the present approach is the evaluation of $\W{\Bas}{}(\br{})$ [see Eq.~\eqref{eq:wcoal}] at each quadrature grid point.
Yet embarrassingly parallel, this step scales, in the general (multi-determinantal) case, as $\Ng \Nb^4$ (where $\Nb$ is the number of basis functions in $\Bas$) but is reduced to $\order*{ \Ng \Ne^2 \Nb^2}$ in the case of a single Slater determinant. Yet embarrassingly parallel, this step scales, in the general (multi-determinantal) case, as $\Ng \Nb^4$ (where $\Nb$ is the number of basis functions in $\Bas$) but is reduced to $\order*{ \Ng \Ne^2 \Nb^2}$ in the case of a single Slater determinant.
\titou{As shown in Ref.~\onlinecite{GinPraFerAssSavTou-JCP-18}, choosing a single Slater determinant wave function to compute $\W{\Bas}{}(\br{1},\br{2})$ already provides a quantitative representation of the incompleteness of $\Bas$ for weakly correlated systems, and therefore we use this framework all through this work. } As shown in Ref.~\onlinecite{GinPraFerAssSavTou-JCP-18}, this simple choice already provides, for weakly correlated systems, a quantitative representation of the basis set incompleteness.
%\begin{equation} Hence, unless otherwise stated, we will stick to this choice throughout the current study.
% \label{eq:fcoal} In our current implementation, the computational bottleneck is the four-index transformation to get the two-electron integrals in the MO basis which appear in Eqs.~\eqref{eq:n2basis} and \eqref{eq:fbasis}.
% \f{\Bas}{\HF}(\br{}) = \sum_{pq \in \Bas} \sum_{ij}^{\occ} \SO{p}{} \SO{q}{} \V{pq}{ij} \SO{i}{} \SO{j}{},
%\end{equation}
In our current implementation, the computational bottleneck of the present basis set correction is the four-index transformation to get the two-electron integrals in the molecular orbital basis which appear in Eqs.~\eqref{eq:n2basis} and \eqref{eq:fbasis}.
Nevertheless, this step usually has to be performed for most correlated WFT calculations. Nevertheless, this step usually has to be performed for most correlated WFT calculations.
Modern integral decomposition techniques (such as density fitting \cite{Whi-JCP-73}) or atomic-orbital-based algorithms could be employed to significantly speed up this step. Modern integral decomposition techniques (such as density fitting \cite{Whi-JCP-73}) or atomic-orbital-based algorithms could be employed to significantly speed up this step.
%When the four-index transformation become prohibitive, by performing successive matrix multiplications, one could rewrite the equations directly in the AO basis where it scales formally as $\order{\Ng \Nb^4}$ but where one can take advantage of the sparsity atomic-orbital-based algorithms to significantly speed up the calculations.
To conclude this section, we point out that, thanks of the definitions \eqref{eq:def_weebasis} and \eqref{eq:mu_of_r} as well as the properties \eqref{eq:lim_W} and \eqref{eq:large_mu_ecmd}, independently of the DFT functional, the present basis set correction To conclude this section, we point out that, thanks to the definitions \eqref{eq:def_weebasis} and \eqref{eq:mu_of_r} as well as the properties \eqref{eq:lim_W} and \eqref{eq:large_mu_ecmd}, independently of the DFT functional, the present basis set correction
i) can be applied to any WFT model that provides an energy and a density, i) can be applied to any WFT model that provides an energy and a density,
ii) does not correct one-electron systems, and ii) does not correct one-electron systems, and
iii) vanishes in the limit of a complete basis set, hence guaranteeing an unaltered CBS limit for a given WFT model. iii) vanishes in the limit of a complete basis set, hence guaranteeing an unaltered complete basis set (CBS) limit for a given WFT model.
%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%
\section{Results} %\section{Results}
%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%
%%% FIGURE 1 %%% %%% FIGURE 1 %%%
@ -490,48 +483,43 @@ iii) vanishes in the limit of a complete basis set, hence guaranteeing an unalte
We begin our investigation of the performance of the basis set correction by computing the atomization energies of \ce{C2}, \ce{N2}, \ce{O2} and \ce{F2} obtained with Dunning's cc-pVXZ basis sets (X $=$ D, T, Q and 5). We begin our investigation of the performance of the basis set correction by computing the atomization energies of \ce{C2}, \ce{N2}, \ce{O2} and \ce{F2} obtained with Dunning's cc-pVXZ basis sets (X $=$ D, T, Q and 5).
\ce{N2}, \ce{O2} and \ce{F2} are weakly correlated systems and belong to the G2 set \cite{CurRagTruPop-JCP-91} (see below), whereas \ce{C2} already contains a non-negligible amount of strong correlation. \cite{BooCleThoAla-JCP-11} \ce{N2}, \ce{O2} and \ce{F2} are weakly correlated systems and belong to the G2 set \cite{CurRagTruPop-JCP-91} (see below), whereas \ce{C2} already contains a non-negligible amount of strong correlation. \cite{BooCleThoAla-JCP-11}
In a second time, we compute the atomization energies of the entire G2 set \cite{CurRagTruPop-JCP-91} composed by 55 molecules with the cc-pVXZ family of basis sets. In a second time, we compute the atomization energies of the entire G2 set \cite{CurRagTruPop-JCP-91} composed by 55 molecules with the cc-pVXZ family of basis sets.
This molecular set has been exhausively studied in the last 20 years (see, for example, Refs.~\onlinecite{FelPetDix-JCP-08, Gro-JCP-09, FelPet-JCP-09, NemTowNee-JCP-10, FelPetHil-JCP-11, HauKlo-JCP-12, PetTouUmr-JCP-12, FelPet-JCP-13, KesSylKohTewMar-JCP-18}) \titou{and can be considered as a representative set for typical quantum chemical calculations on small organic molecules}. This molecular set has been exhausively studied in the last 20 years (see, for example, Refs.~\onlinecite{FelPetDix-JCP-08, Gro-JCP-09, FelPet-JCP-09, NemTowNee-JCP-10, FelPetHil-JCP-11, HauKlo-JCP-12, PetTouUmr-JCP-12, FelPet-JCP-13, KesSylKohTewMar-JCP-18}) and can be considered as a representative set of small organic and inorganic molecules.
As a method $\modY$ we employ either CCSD(T) or exFCI. As a method $\modY$ we employ either CCSD(T) or exFCI.
Here, exFCI stands for extrapolated FCI energies computed with the CIPSI algorithm. \cite{HurMalRan-JCP-73, GinSceCaf-CJC-13, GinSceCaf-JCP-15} Here, exFCI stands for extrapolated FCI energies computed with the CIPSI algorithm. \cite{HurMalRan-JCP-73, GinSceCaf-CJC-13, GinSceCaf-JCP-15}
We refer the interested reader to Refs.~\onlinecite{HolUmrSha-JCP-17, SceGarCafLoo-JCTC-18, LooSceBloGarCafJac-JCTC-18, SceBenJacCafLoo-JCP-18, LooBogSceCafJAc-JCTC-19} for more details. We refer the interested reader to Refs.~\onlinecite{HolUmrSha-JCP-17, SceGarCafLoo-JCTC-18, LooSceBloGarCafJac-JCTC-18, SceBenJacCafLoo-JCP-18, LooBogSceCafJAc-JCTC-19} for more details.
In the case of the CCSD(T) calculations, we have $\modZ = \HF$ as we use the restricted open-shell HF (ROHF) one-electron density to compute the complementary energy. In the case of the CCSD(T) calculations, we have $\modZ = \ROHF$ as we use the restricted open-shell HF (ROHF) one-electron density to compute the complementary energy.
For exFCI, the one-electron density is computed from a very large CIPSI expansion containing several million determinants. For exFCI, the one-electron density is computed from a very large CIPSI expansion containing several million determinants.
%For the definition of the interaction, we use a single Slater determinant built in the ROHF basis for the CCSD(T) calculation, and built with the natural orbitals of the converged variational wave function for the exFCI calculations. CCSD(T) energies are computed with Gaussian09 using standard threshold values, \cite{g09} while RS-DFT and exFCI calculations are performed with {\QP}. \cite{QP2}
CCSD(T) energies are computed with Gaussian09 using standard threshold values. \cite{g09}
RS-DFT and exFCI calculations are performed with {\QP}. \cite{QP2}
For the numerical quadratures, we employ the SG-2 grid. \cite{DasHer-JCC-17} For the numerical quadratures, we employ the SG-2 grid. \cite{DasHer-JCC-17}
Except for the carbon dimer where we have taken the experimental equilibrium bond length (\InAA{1.2425}), all geometries have been extracted from Ref.~\onlinecite{HauJanScu-JCP-09} and have been obtained at the B3LYP/6-31G(2df,p) level of theory. Except for the carbon dimer where we have taken the experimental equilibrium bond length (\InAA{1.2425}), all geometries have been extracted from Ref.~\onlinecite{HauJanScu-JCP-09} and have been obtained at the B3LYP/6-31G(2df,p) level of theory.
Frozen-core calculations are defined as such: an \ce{He} core is frozen from \ce{Li} to \ce{Ne}, while a \ce{Ne} core is frozen from \ce{Na} to \ce{Ar}. Frozen-core calculations are defined as such: a \ce{He} core is frozen from \ce{Li} to \ce{Ne}, while a \ce{Ne} core is frozen from \ce{Na} to \ce{Ar}.
In the context of the basis set correction, the set of MOs $\BasFC$ involved in the definition of the effective interaction refers to the non-frozen MOs. In the context of the basis set correction, the set of MOs $\BasFC$ involved in the definition of the effective interaction refers to the non-frozen MOs.
The FC density-based correction is set consistently when the FC approximation was applied in WFT methods. The FC density-based correction is set consistently when the FC approximation was applied in WFT methods.
\titou{In order to estimate the complete basis set (CBS) limit for each model, following the work of Ref.~\onlinecite{HalHelJorKloKocOlsWil-CPL-98}, To estimate the CBS limit of each method, following Ref.~\onlinecite{HalHelJorKloKocOlsWil-CPL-98}, we perform a two-point extrapolation of the correlation energies using the quadruple- and quintuple-$\zeta$ data that we add up to the HF energies obtained in the largest (i.e.~quintuple-$\zeta$) basis.
we employ the two-point extrapolation for the correlation energies for each model in quadruple- and quintuple-$\zeta$ basis sets, which is refered to as $\CBS$ correlation energies, and we add the HF energies in the largest basis sets (\textit{i.e.} quintuple-$\zeta$ quality basis sets) to the CBS correlation energies to estimate the CBS FCI energies.}
%\subsection{Convergence of the atomization energies with the WFT models } As the exFCI calculations are converged with a precision of about 0.1 {\kcal}, we can label these atomization energies as near-FCI.
As the exFCI calculations are converged with a precision of about 0.1 {\kcal}, we can consider these atomization energies as near-FCI values, and they will be our references for \ce{C2}, \ce{N2}, \ce{O2} and \ce{F2}. Hence, they will be our references for \ce{C2}, \ce{N2}, \ce{O2} and \ce{F2}.
The results for these diatomics are reported in Fig.~\ref{fig:diatomics}. The results for these diatomics are reported in Fig.~\ref{fig:diatomics}.
The corresponding numerical data can be found in the {\SI}. The corresponding numerical data can be found in the {\SI}.
As one can see, the convergence of the exFCI atomization energies is, as expected, slow with respect to the basis set: chemical accuracy (error below 1 {\kcal}) is barely reached for \ce{C2}, \ce{O2} and \ce{F2} even with a cc-pV5Z basis set, and the atomization energies are consistently underestimated. As one can see, the convergence of the exFCI atomization energies is, as expected, slow with respect to the basis set: chemical accuracy (error below 1 {\kcal}) is barely reached for \ce{C2}, \ce{O2} and \ce{F2} even with the cc-pV5Z basis set, and the atomization energies are consistently underestimated.
%Also, the atomization energies are consistently underestimated, reflecting that, in a given basis, the atom is always better described than the molecule due to the larger number of interacting electron pairs in the molecular system.
A similar trend holds for CCSD(T). A similar trend holds for CCSD(T).
%, and one can notice that the atomization energies of the CCSD(T) are always slightly underestimated with respect to the CIPSI ones, showing that the CCSD(T) ansatz is better suited for the atoms than for the molecule.
Regarding the effect of the basis set correction, several general observations can be made for both exFCI and CCSD(T). Regarding the effect of the basis set correction, several general observations can be made for both exFCI and CCSD(T).
First, in a given basis set, the basis set correction systematically improves the atomization energies (both at the LDA and PBE level). First, in a given basis set, the basis set correction systematically improves the atomization energies (both at the LDA and PBE levels).
A small overestimation can occur compared to the CBS values by a few tenths of a {\kcal} (the largest deviation being 0.6 {\kcal} for \ce{N2} at the CCSD(T)+PBE/cc-pV5Z level). A small overestimation can occur compared to the CBS value by a few tenths of a {\kcal} (the largest deviation being 0.6 {\kcal} for \ce{N2} at the CCSD(T)+PBE/cc-pV5Z level).
Nevertheless, the deviation observed for the largest basis set is typically within the extrapolation error of the CBS atomization energies, which is highly satisfactory knowing the marginal computation cost of the present correction. Nevertheless, the deviation observed for the largest basis set is typically within the CBS extrapolation error, which is highly satisfactory knowing the marginal computational cost of the present correction.
In most cases, the basis set corrected triple-$\zeta$ results are on par with the uncorrected quintuple-$\zeta$ ones. In most cases, the basis set corrected triple-$\zeta$ atomization energies are on par with the uncorrected quintuple-$\zeta$ ones.
%Also, the values obtained with the largest basis sets tends to converge toward a value close to the estimated CBS values.
Importantly, the sensitivity with respect to the SR-DFT functional is quite large for the double- and triple-$\zeta$ basis sets, where clearly the PBE functional performs better. Importantly, the sensitivity with respect to the SR-DFT functional is quite large for the double- and triple-$\zeta$ basis sets, where clearly the PBE functional performs better.
However, from the quadruple-$\zeta$ basis, the LDA and PBE functionals agree within a few tenths of a {\kcal}. However, from the quadruple-$\zeta$ basis, the LDA and PBE functionals agree within a few tenths of a {\kcal}.
Such weak sensitivity when reaching large basis sets shows the robustness of the approach. Such weak sensitivity when reaching large basis sets shows the robustness of the approach.
As a second set of numerical examples, we compute the error with respect to the CBS values of the atomization energies of the G2 test sets with $\modY=\CCSDT$, $\modZ=\HF$ and the cc-pVXZ basis sets. As a second set of numerical examples, we compute the error (with respect to the CBS values) of the atomization energies from the G2 test sets with $\modY=\CCSDT$, $\modZ=\ROHF$ and the cc-pVXZ basis sets.
\titou{Here we use HF/CBS values to compute the atomization energies which is equivalent to looking at the correlation energy contribution to the atomization energies. \titou{Here we use the near-CBS HF/cc-pV5Z energies to compute atomization energies.
This is equivalent to looking at the correlation energy contribution to the atomization energies.
Investigating the convergence of correlation energies or difference of such quantities is usually done to appreciate the performance of basis set corrections aiming at correcting two-electron effects\cite{Tenno-CPL-04,TewKloNeiHat-PCCP-07,IrmGru-arXiv-2019}, as these quantities do not contain the HF energy component whose rate of convergence is very different depending on the molecular system.} Investigating the convergence of correlation energies or difference of such quantities is usually done to appreciate the performance of basis set corrections aiming at correcting two-electron effects\cite{Tenno-CPL-04,TewKloNeiHat-PCCP-07,IrmGru-arXiv-2019}, as these quantities do not contain the HF energy component whose rate of convergence is very different depending on the molecular system.}
The ``plain'' CCSD(T) atomization energies as well as the corrected CCSD(T)+LDA and CCSD(T)+PBE values are depicted in Fig.~\ref{fig:G2_Ec}. The ``plain'' CCSD(T) atomization energies as well as the corrected CCSD(T)+LDA and CCSD(T)+PBE values are depicted in Fig.~\ref{fig:G2_Ec}.
The raw data can be found in the {\SI}. The raw data can be found in the {\SI}.
A statistical analysis of these data is also provided in Table \ref{tab:stats}, where we report the mean absolute deviation (MAD), root-mean-square deviation (RMSD), and maximum deviation (MAX) with respect to the CCSD(T)/CBS atomization energies. A statistical analysis of these data is also provided in Table \ref{tab:stats}, where we report the mean absolute deviation (MAD), root-mean-square deviation (RMSD), and maximum deviation (MAX) with respect to the CCSD(T)/CBS atomization energies.
\titou{Compared to the reference values extracted from Ref.~\onlinecite{HauKlo-JCP-12} which corresponds to frozen-core non-relativistic atomization energies obtained at the CCSD(T)(F12)/cc-pVQZ-F12 level of theory corrected for higher-excitation contributions ($E_\text{CCSDT(Q)/cc-pV(D+d)Z} - E_\text{CCSD(T)/cc-pV(D+d)Z})$, our CCSD(T)/CBS atomization energies differ by MAD = 0.37 {\kcal}.} Note that the MAD of our CCSD(T)/CBS atomization energies is only 0.37 {\kcal} compared to the values extracted from Ref.~\onlinecite{HauKlo-JCP-12} which corresponds to frozen-core non-relativistic atomization energies obtained at the CCSD(T)(F12)/cc-pVQZ-F12 level of theory corrected for higher-excitation contributions ($E_\text{CCSDT(Q)/cc-pV(D+d)Z} - E_\text{CCSD(T)/cc-pV(D+d)Z})$.
From double-$\zeta$ to quintuple-$\zeta$ basis, the MAD associated with the CCSD(T) atomization energies goes down slowly from 14.29 to 1.28 {\kcal}. From double-$\zeta$ to quintuple-$\zeta$ basis, the MAD associated with the CCSD(T) atomization energies goes down slowly from 14.29 to 1.28 {\kcal}.
For a commonly-used basis like cc-pVTZ, the MAD of CCSD(T) is still 6.06 {\kcal}. For a commonly-used basis like cc-pVTZ, the MAD of CCSD(T) is still 6.06 {\kcal}.
Applying the basis set correction drastically reduces the basis set incompleteness error. Applying the basis set correction drastically reduces the basis set incompleteness error.
@ -539,18 +527,19 @@ Already at the CCSD(T)+LDA/cc-pVDZ and CCSD(T)+PBE/cc-pVDZ levels, the MAD is re
With the triple-$\zeta$ basis, the MAD of CCSD(T)+PBE/cc-pVTZ is already below 1 {\kcal} with 36 cases (out of 55) where we achieve chemical accuracy. With the triple-$\zeta$ basis, the MAD of CCSD(T)+PBE/cc-pVTZ is already below 1 {\kcal} with 36 cases (out of 55) where we achieve chemical accuracy.
CCSD(T)+LDA/cc-pVQZ and CCSD(T)+PBE/cc-pVQZ return MAD of 0.33 and 0.31 kcal/mol (respectively) while CCSD(T)/cc-pVQZ still yields a fairly large MAD of 2.50 {\kcal}. CCSD(T)+LDA/cc-pVQZ and CCSD(T)+PBE/cc-pVQZ return MAD of 0.33 and 0.31 kcal/mol (respectively) while CCSD(T)/cc-pVQZ still yields a fairly large MAD of 2.50 {\kcal}.
\titou{Therefore, similar to F12 methods, \cite{TewKloNeiHat-PCCP-07} we can safely claim that the present basis set correction recovers quintuple-$\zeta$ quality atomization and correlation energies with triple-$\zeta$ basis sets for a much cheaper computational cost. Therefore, similar to F12 methods, \cite{TewKloNeiHat-PCCP-07} we can safely claim that the present basis set correction provides significant basis set reduction and recovers quintuple-$\zeta$ quality atomization and correlation energies with triple-$\zeta$ basis sets for a much cheaper computational cost.
Encouraged by these results obtained for weakly correlated systems, we are currently developing this theory towards the treatment of the basis set error for strongly correlated systems, excited states and the treatment of the one-electron error in the basis set incompleteness.} Encouraged by these promising results, we are currently pursuing various avenues toward basis set reduction for strongly correlated systems and electronically excited states.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section*{Supporting information} \section*{Supporting information}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
See {\SI} for raw data associated with the atomization energies of the four diatomics and the G2 atomization energies. See {\SI} for raw data associated with the atomization energies of the four diatomics and the G2 set.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\begin{acknowledgements} \begin{acknowledgements}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
The authors would like to thank the \textit{Centre National de la Recherche Scientifique} (CNRS) for funding. The authors would like to thank the \emph{Centre National de la Recherche Scientifique} (CNRS) for funding.
This work was performed using HPC resources from GENCI-TGCC (Grant No.~2018-A0040801738), and CALMIP (Toulouse) under allocation 2019-18005.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\end{acknowledgements} \end{acknowledgements}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%