loos/Ec
2
0
Fork 0
Code Issues Pull Requests Releases Wiki Activity
11c9c7bb71
Ec/Manuscript/Ec.tex

634 lines
46 KiB
TeX
Raw Blame History

\documentclass[aps,prb,reprint,noshowkeys,superscriptaddress]{revtex4-1}
\usepackage{graphicx,dcolumn,bm,xcolor,microtype,multirow,amscd,amsmath,amssymb,amsfonts,physics,wrapfig,txfonts,siunitx}
\usepackage[version=4]{mhchem}
\newcommand{\ie}{\textit{i.e.}}
\newcommand{\eg}{\textit{e.g.}}
\newcommand{\alert}[1]{\textcolor{black}{#1}}
\usepackage[normalem]{ulem}
\newcommand{\titou}[1]{\textcolor{red}{#1}}
\newcommand{\trashPFL}[1]{\textcolor{red}{\sout{#1}}}
\newcommand{\PFL}[1]{\titou{(\underline{\bf PFL}: #1)}}
\newcommand{\toto}[1]{\textcolor{green}{#1}}
\newcommand{\trashAS}[1]{\textcolor{green}{\sout{#1}}}
\newcommand{\AS}[1]{\toto{(\underline{\bf AS}: #1)}}
% useful shortcuts
\newcommand{\mc}{\multicolumn}
\newcommand{\fnm}{\footnotemark}
\newcommand{\fnt}{\footnotetext}
\newcommand{\tabc}[1]{\multicolumn{1}{c}{#1}}
% key quantities
\newcommand{\Nel}{n}
\newcommand{\Norb}{N}
\newcommand{\Ndet}{N_\text{det}}
% operators
\newcommand{\hT}{\Hat{T}}
\newcommand{\hH}{\Hat{H}}
\newcommand{\hh}{\Hat{h}}
\newcommand{\hX}{\Hat{X}}
\newcommand{\hk}{\Hat{\kappa}}
\newcommand{\cre}[1]{\Hat{a}_{#1}^{\dagger}}
\newcommand{\ani}[1]{\Hat{a}_{#1}^{}}
% vectors and matrices
\newcommand{\bO}{\boldsymbol{0}}
\newcommand{\bI}{\boldsymbol{1}}
\newcommand{\bX}{\boldsymbol{X}}
\newcommand{\bH}{\boldsymbol{H}}
\newcommand{\bg}{\boldsymbol{g}}
\newcommand{\bk}{\boldsymbol{\kappa}}
\newcommand{\bc}{\boldsymbol{c}}
\newcommand{\br}{\boldsymbol{r}}
% curly letters
\newcommand{\cA}{\mathcal{A}}
\newcommand{\cI}{\mathcal{I}}
\newcommand{\cP}{\mathcal{P}}
% energies
\newcommand{\EHF}{E_\text{HF}}
\newcommand{\Ec}{E_\text{c}}
\newcommand{\Evar}{E_\text{var}}
\newcommand{\EPT}{E_\text{PT2}}
\newcommand{\ECIPSI}{E_\text{CIPSI}}
% wave functions and orbitals
\newcommand{\PsiO}{\Psi_0}
\newcommand{\Psivar}{\Psi_\text{var}}
\newcommand{\MO}[1]{\phi_{#1}}
\usepackage[
colorlinks=true,
citecolor=blue,
breaklinks=true
]{hyperref}
\urlstyle{same}
\begin{document}
\newcommand{\LCPQ}{Laboratoire de Chimie et Physique Quantiques (UMR 5626), Universit\'e de Toulouse, CNRS, UPS, France}
\newcommand{\CEISAM}{Universit\'e de Nantes, CNRS, CEISAM UMR 6230, F-44000 Nantes, France}
\title{Reference correlation energies in finite Hilbert spaces: five- and six-membered rings}
\author{Yann Damour}
\affiliation{\LCPQ}
\author{Micka\"el V\'eril}
\affiliation{\LCPQ}
\author{Michel Caffarel}
\affiliation{\LCPQ}
\author{Denis Jacquemin}
\affiliation{\CEISAM}
\author{Anthony Scemama}
\email{scemama@irsamc.ups-tlse.fr}
\affiliation{\LCPQ}
\author{Pierre-Fran\c{c}ois Loos}
\email{loos@irsamc.ups-tlse.fr}
\affiliation{\LCPQ}
% Abstract
\begin{abstract}
In the continuity of our recent work on the benzene molecule [\href{https://doi.org/10.1063/5.0027617}{J.~Chem.~Phys.~\textbf{153}, 176101 (2020)}], itself motivated by the blind challenge of Eriksen \textit{et al.} [\href{https://doi.org/10.1021/acs.jpclett.0c02621}{J.~Phys.~Chem.~Lett.~\textbf{11}, 8922 (2020)}] on the same system, we report reference frozen-core correlation energies for twelve five- and six-membered ring molecules (cyclopentadiene, furan, imidazole, pyrrole, thiophene, benzene, pyrazine, pyridazine, pyridine, pyrimidine, tetrazine, and triazine) in the standard correlation-consistent double-$\zeta$ Dunning basis set (cc-pVDZ).
This corresponds to Hilbert spaces with sizes ranging from $10^{28}$ (for thiophene) to $10^{36}$ (for benzene).
Our estimates are based on energetically optimized-orbital selected configuration interaction (SCI) calculations performed with the \textit{Configuration Interaction using a Perturbative Selection made Iteratively} (CIPSI) algorithm.
The performance and convergence properties of several series of methods are investigated.
In particular, we study the convergence properties of i) the M{\o}ller-Plesset perturbation series up to fifth-order (MP2, MP3, MP4, and MP5), ii) the iterative approximate single-reference coupled-cluster series CC2, CC3, and CC4, and iii) the single-reference coupled-cluster series CCSD, CCSDT, and CCSDTQ.
The performance of the ground-state gold standard CCSD(T) is also investigated.
\end{abstract}
% Title
\maketitle
%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Introduction}
\label{sec:intro}
%%%%%%%%%%%%%%%%%%%%%%%%%
Electronic structure theory relies heavily on approximations. \cite{Szabo_1996,Helgaker_2013,Jensen_2017}
Loosely speaking, to make any theory useful, three main approximations must be enforced.
The first fundamental approximation, known as the Born-Oppenheimer approximation, usually consists in assuming that the motion of nuclei and electrons are decoupled. \cite{Born_1927}
The nuclei coordinates can then be treated as parameters in the electronic Hamiltonian.
The second central approximation which makes calculations feasible by a computer is the basis set approximation where one introduces a set of pre-defined basis functions to represent the many-electron wave function of the system.
In most molecular calculations, a set of one-electron, atom-centered gaussian basis functions are introduced to expand the so-called one-electron molecular orbitals which are then used to build the many-electron Slater determinants.
The third and most relevant approximation in the present context is the ansatz (or form) of the electronic wave function $\Psi$.
For example, in configuration interaction (CI) methods, the wave function is expanded as a linear combination of Slater determinants, while in (single-reference) coupled-cluster (CC) theory, \cite{Cizek_1966,Paldus_1972,Crawford_2000,Piecuch_2002,Bartlett_2007,Shavitt_2009} a reference Slater determinant $\PsiO$ [usually taken as the Hartree-Fock (HF) wave function] is multiplied by a wave operator defined as the exponentiated excitation operator $\hT = \sum_{k=1}^\Nel \hT_k$ (where $\Nel$ is the number of electrons).
The truncation of $\hT$ allows to define a hierarchy of non-variational and size-extensive methods with improved accuracy:
CC with singles and doubles (CCSD), \cite{Cizek_1966,Purvis_1982} CC with singles, doubles, and triples (CCSDT), \cite{Noga_1987a,Scuseria_1988} CC with singles, doubles, triples, and quadruples (CCSDTQ), \cite{Oliphant_1991,Kucharski_1992} with corresponding computational scalings of $\order*{\Norb^{6}}$, $\order*{\Norb^{8}}$, and $\order*{\Norb^{10}}$, respectively (where $\Norb$ denotes the number of orbitals).
Parallel to the ``complete'' CC series presented above, an alternative series of approximate iterative CC models have been developed by the Aarhus group in the context of CC response theory \cite{Christiansen_1998} where one skips the most expensive terms and avoids the storage of the higher-excitation amplitudes: CC2, \cite{Christiansen_1995a} CC3, \cite{Christiansen_1995b,Koch_1997} and CC4. \cite{Kallay_2005,Matthews_2021,Loos_2021}
These iterative methods scale as $\order*{\Norb^{5}}$, $\order*{\Norb^{7}}$, and $\order*{\Norb^{9}}$, respectively, and can be seen as cheaper approximations of CCSD, CCSDT, and CCSDTQ.
Coupled-cluster methods have been particularly successful at computing accurately various properties for small- and medium-sized molecules.
\cite{Kallay_2003,Kallay_2004a,Gauss_2006,Kallay_2006,Gauss_2009}
A similar systematic truncation strategy can be applied to CI methods leading to the well-established family of methods known as CISD, CISDT, CISDTQ, \ldots~where one systematically increases the maximum excitation degree of the determinants taken into account.
Except for full CI (FCI) where all determinants from the Hilbert space (\ie, with excitation degree up to $\Nel$) are considered, truncated CI methods are variational but lack size-consistency.
The non-variationality of truncated CC methods being less of an issue than the size-inconsistency of the truncated CI methods, the formers have naturally overshadowed the latters in the electronic structure landscape.
However, a different strategy has recently made a come back in the context of CI methods. \cite{Bender_1969,Whitten_1969,Huron_1973}
Indeed, selected CI (SCI) methods, \cite{Booth_2009,Giner_2013,Evangelista_2014,Giner_2015,Holmes_2016,Tubman_2016,Liu_2016,Ohtsuka_2017,Zimmerman_2017,Coe_2018,Garniron_2018} where one iteratively selects the important determinants from the FCI space (usually) based on a perturbative criterion, has been recently shown to be highly successful in order to produce reference energies for ground and excited states in small- and medium-sized molecules \cite{Holmes_2017,Li_2018,Li_2020,Loos_2018a,Chien_2018,Loos_2019,Loos_2020b,Loos_2020c,Loos_2020e,Garniron_2019,Eriksen_2020,Yao_2020,Veril_2021,Loos_2021} thanks to efficient deterministic, stochastic or hybrid algorithms well suited for massive parallelization.
We refer the interested reader to Refs.~\onlinecite{Loos_2020a,Eriksen_2021} for recent reviews.
SCI methods are based on a well-known fact: amongst the very large number of determinants contained in the FCI space, only a relative small fraction of them significantly contributes to the energy.
Accordingly, the SCI+PT2 family of methods performs a sparse exploration of the FCI space by selecting iteratively only the most energetically relevant determinants of the variational space and supplementing it with a second-order perturbative correction (PT2). \cite{Huron_1973,Garniron_2017,Sharma_2017,Garniron_2018,Garniron_2019}
Although the formal scaling of such algorithms remain exponential, the prefactor is greatly reduced which explains their current attractiveness in the electronic structure community and much wider applicability than their standard FCI parent.
Note that, very recently, several groups \cite{Aroeira_2021,Lee_2021,Magoulas_2021} have coupled CC and SCI methods via the externally-corrected CC methodology, \cite{Paldus_2017} showing promising performances for weakly and strongly correlated systems.
A rather different strategy in order to reach the holy grail FCI limit is to resort to M{\o}ller-Plesset (MP) perturbation theory, \cite{Moller_1934}
which popularity originates from its black-box nature, size-extensivity, and relatively low computational scaling, making it easily applied to a broad range of molecular systems.
Again, at least in theory, one can obtain the exact energy of the system by ramping up the degree of the perturbative series. \cite{Marie_2021}
The second-order M{\o}ller-Plesset (MP2) method \cite{Moller_1934} [which scales as $\order*{\Norb^{5}}$] has been broadly adopted in quantum chemistry for several decades, and is now included in the increasingly popular double-hybrid functionals \cite{Grimme_2006} alongside exact HF exchange.
Its higher-order variants [MP3, \cite{Pople_1976}
MP4, \cite{Krishnan_1980} MP5, \cite{Kucharski_1989} and MP6 \cite{He_1996a,He_1996b} which scales as $\order*{\Norb^{6}}$, $\order*{\Norb^{7}}$, $\order*{\Norb^{8}}$, and $\order*{\Norb^{9}}$ respectively] have been investigated much more scarcely.
However, it is now widely recognised that the series of MP approximations might show erratic, slow, or divergent behavior that limit its applicability and systematic improvability. \cite{Laidig_1985,Knowles_1985,Handy_1985,Gill_1986,Laidig_1987,Nobes_1987,Gill_1988,Gill_1988a,Lepetit_1988,Malrieu_2003}
Again, MP perturbation theory and CC methods can be coupled.
The CCSD(T) method, \cite{Raghavachari_1989} known as the gold-standard of quantum chemistry for weakly correlated systems, is probably the most iconic example of such coupling.
Motivated by the recent blind test of Eriksen \textit{et al.}\cite{Eriksen_2020}~reporting the performance of a large panel of emerging electronic structure methods [the many-body expansion FCI (MBE-FCI), \cite{Eriksen_2017,Eriksen_2018,Eriksen_2019a,Eriksen_2019b} adaptive sampling CI (ASCI), \cite{Tubman_2016,Tubman_2018,Tubman_2020} iterative CI (iCI), \cite{Liu_2014,Liu_2016,Lei_2017,Zhang_2020} semistochastic heat-bath CI (SHCI), \cite{Holmes_2016,Holmes_2017,Sharma_2017} the full coupled-cluster reduction (FCCR), \cite{Xu_2018,Xu_2020} density-matrix renormalization group (DMRG), \cite{White_1992,White_1993,Chan_2011} adaptive-shift FCI quantum Monte Carlo (AS-FCIQMC), \cite{Booth_2009,Cleland_2010,Ghanem_2019} and cluster-analysis-driven FCIQMC (CAD-FCIQMC) \cite{Deustua_2017,Deustua_2018}] on the non-relativistic frozen-core correlation energy of the benzene molecule in the standard correlation-consistent double-$\zeta$ Dunning basis set (cc-pVDZ), some of us have recently investigated the performance of the \textit{Configuration Interaction using a Perturbative Selection made Iteratively} (CIPSI) method \cite{Huron_1973,Giner_2013,Giner_2015,Garniron_2018,Garniron_2019} on the very same system \cite{Loos_2020e} [see also Ref.~\onlinecite{Lee_2020} for a study of the performance of phaseless auxiliary-field quantum Monte Carlo (ph-AFQMC) \cite{Motta_2018}].
In the continuity of this recent work, we report here a significant extension by estimating the (frozen-core) FCI/cc-pVDZ correlation energy of twelve cyclic molecules (cyclopentadiene, furan, imidazole, pyrrole, thiophene, benzene, pyrazine, pyridazine, pyridine, pyrimidine, tetrazine, and triazine) with the help of CIPSI employing energetically-optimized orbitals at the same level of theory. \cite{Yao_2020,Yao_2021}
These systems are depicted in Fig.~\ref{fig:mol}.
This set of molecular systems corresponds to Hilbert spaces with sizes ranging from $10^{28}$ (for thiophene) to $10^{36}$ (for benzene).
In addition to CIPSI, the performance and convergence properties of several series of methods are investigated.
In particular, we study i) the MP perturbation series up to fifth-order (MP2, MP3, MP4, and MP5), ii) the CC2, CC3, and CC4 approximate series, and ii) the ``complete'' CC series up to quadruples (\ie, CCSD, CCSDT, and CCSDTQ).
The performance of the ground-state gold standard CCSD(T) is also investigated.
%%% FIG 1 %%%
\begin{figure*}
\includegraphics[width=\linewidth]{mol}
\caption{
Five-membered rings (top) and six-membered rings (bottom) considered in this study.
\label{fig:mol}}
\end{figure*}
%%% %%% %%%
The present manuscript is organized as follows.
In Sec.~\ref{sec:compdet}, computational details concerning geometries, basis sets, and methods are reported.
Section \ref{sec:OO-CIPSI} provides theoretical details about the CIPSI algorithm and the orbital optimization procedure that we have employed here.
In Sec.~\ref{sec:res}, we report our reference FCI correlation energies for the five-membered and six-membered cyclic molecules obtained thanks to extrapolated orbital-optimized CIPSI calculations (Sec.~\ref{sec:cipsi_res}).
These reference correlation energies are then used to benchmark and study the convergence properties of various perturbative and CC methods (Sec.~\ref{sec:mpcc_res}).
Finally, we draw our conclusions in Sec.~\ref{sec:ccl}.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Computational details}
\label{sec:compdet}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
The geometries of the twelve systems considered in the present study have been all obtained at the CC3/aug-cc-pVTZ level of theory and have been extracted from a previous study. \cite{Loos_2020a}
Note that, for the sake of consistency, the geometry of benzene considered here is different from one of Ref.~\onlinecite{Loos_2020e} which has been computed at a lower level of theory [MP2/6-31G(d)]. \cite{Schreiber_2008}
The MP2, MP3, MP4, CC2, CC3, CC4, CCSD, CCSDT, and CCSDTQ calculations have been performed with CFOUR, \cite{Matthews_2020} while the CCSD(T) and MP5 calculations have been computed with GAUSSIAN 09. \cite{g09}
The CIPSI calculations have been performed with QUANTUM PACKAGE. \cite{Garniron_2019}
In the current implementation, the selection step and the PT2 correction are computed simultaneously via a hybrid semistochastic algorithm \cite{Garniron_2017,Garniron_2019} (which explains the statistical error associated with the PT2 correction in the following).
Here, we employ the renormalized version of the PT2 correction which has been recently implemented and tested for a more efficient extrapolation to the FCI limit thanks to a partial resummation of the higher-order of perturbation. \cite{Garniron_2019}
We refer the interested reader to Ref.~\onlinecite{Garniron_2019} where one can find all the details regarding the implementation of the PT2 correction and the CIPSI algorithm.
For all these calculations, we consider Dunning's correlation-consistent double-$\zeta$ basis (cc-pVDZ).
Although the FCI energy has the enjoyable property of being independent of the set of one-electron orbitals used to construct the many-electron Slater determinants, as a truncated CI method, the convergence properties of CIPSI strongly dependent on this orbital choice.
In the present study, we investigate, in particular, the convergence behavior of the CIPSI energy for two sets of orbitals: natural orbitals (NOs) and optimized orbitals (OOs).
Following our usual procedure, \cite{Scemama_2018,Scemama_2018b,Scemama_2019,Loos_2018a,Loos_2019,Loos_2020a,Loos_2020b,Loos_2020c,Loos_2020e} we perform first a preliminary SCI calculation using HF orbitals in order to generate a SCI wave function with at least $10^7$ determinants.
Natural orbitals are computed based on this wave function and they are used to perform a new CIPSI run.
Successive orbital optimizations are then performed, which consist in minimizing the variational CIPSI energy at each iteration up to approximately $2 \times 10^5$ determinants.
When convergence is achieved in terms of orbital optimization, as our ``production'' run, we perform a new CIPSI calculation from scratch using this set of optimized orbitals.
Using optimized orbitals has the undeniable advantage to produce, for a given variational energy, more compact CI expansions.
For the benzene molecule, we also explore the use of localized orbitals (LOs) which are produced with the Boys-Foster localization procedure \cite{Boys_1960} that we apply to the natural orbitals in several orbital windows in order to preserve a strict $\sigma$-$\pi$ separation in the planar systems considered here. \cite{Loos_2020e}
Because they take advantage of the local character of electron correlation, localized orbitals have been shown to provide faster convergence towards the FCI limit compared to natural orbitals. \cite{Angeli_2003,Angeli_2009,BenAmor_2011,Suaud_2017,Chien_2018,Eriksen_2020,Loos_2020e}
As we shall see below, employing optimized orbitals has the advantage to produce an even smoother and faster convergence of the SCI energy toward the FCI limit.
Note both localized and optimized orbitals do break the spatial symmetry.
Unlike excited-state calculations where it is important to enforce that the wave functions are eigenfunctions of the $\Hat{S}^2$ spin operator, \cite{Chilkuri_2021} the present wave functions do not fulfill this property as we aim for the lowest possible energy of a closed-shell singlet state.
We have found that $\expval*{\Hat{S}^2}$ is, nonetheless, very close to zero for each system.
%%%%%%%%%%%%%%%%%%%%%%%%%
\section{CIPSI with optimized orbitals}
\label{sec:OO-CIPSI}
%%%%%%%%%%%%%%%%%%%%%%%%%
Here, we provide key details about the CIPSI method \cite{Huron_1973,Garniron_2019} as well as the orbital optimization procedure which has been shown to be highly effective in the context of SHCI by Umrigar and coworkers. \cite{Eriksen_2020,Yao_2020,Yao_2021}
Note that we focus on the ground state but the present discussion can be easily extended to excited states. \cite{Scemama_2019,Veril_2021}
At the $k$th iteration, the total CIPSI energy $\ECIPSI^{(k)}$ is defined as the sum of the variational energy
\begin{equation}
\Evar^{(k)} = \frac{\mel*{\Psivar^{(k)}}{\hH}{\Psivar^{(k)}}}{\braket*{\Psivar^{(k)}}{\Psivar^{(k)}}}
\end{equation}
and a second-order perturbative correction
\begin{equation}
\EPT^{(k)}
= \sum_{\alpha \in \cA_k} e_{\alpha}^{(k)}
= \sum_{\alpha \in \cA_k} \frac{\abs*{\mel*{\Psivar^{(k)}}{\hH}{\alpha}}^2}{\Evar^{(k)} - \mel*{\alpha}{\hH}{\alpha}}
\end{equation}
where $\hH$ is the (non-relativistic) electronic Hamiltonian,
\begin{equation}
\label{eq:Psivar}
\Psivar^{(k)} = \sum_{I \in \cI_k} c_I^{(k)} \ket*{I}
\end{equation}
is the variational wave function, $\cI_k$ is the set of internal determinants $\ket*{I}$ and $\cA_k$ is the set of external determinants (or perturbers) $\ket*{\alpha}$ which do not belong to the variational space at the $k$th iteration but are linked to it via a nonzero matrix element, \ie, $\mel*{\Psivar^{(k)}}{\hH}{\alpha} \neq 0$.
The sets $\cI_k$ and $\cA_k$ define, at the $k$th iteration, the internal and external spaces, respectively.
In the selection step, the perturbers corresponding to the largest $\abs*{e_{\alpha}^{(k)}}$ values are then added to the variational space at iteration $k+1$.
In our implementation, the size of the variational space is roughly doubled at each iteration.
Hereafter, we label these iterations over $k$ which consist in enlarging the variational space as \textit{macroiterations}.
In practice, $\Evar^{(k)}$ is computed by diagonalizing the $\Ndet^{(k)} \times \Ndet^{(k)}$ CI matrix with elements $\mel{I}{\hH}{J}$ via Davidson's algorithm. \cite{Davidson_1975}
The magnitude of $\EPT^{(k)}$ provides, at iteration $k$, a qualitative idea of the ``distance'' to the FCI limit. \cite{Garniron_2018}
We then linearly extrapolate, using large variational wave functions, the CIPSI energy to $\EPT = 0$ (which effectively corresponds to the FCI limit).
Further details concerning the extrapolation procedure are provided below (see Sec.~\ref{sec:res}).
Orbital optimization techniques at the SCI level are theoretically straightforward, but practically challenging.
Here, we detail our orbital optimization procedure within the CIPSI algorithm and we assume that the variational wave function is normalized, \ie, $\braket*{\Psivar}{\Psivar} = 1$.
As stated in Sec.~\ref{sec:intro}, $\Evar$ depends on both the CI coefficients $\{ c_I \}_{1 \le I \le \Ndet}$ [see Eq.~\eqref{eq:Psivar}] but also on the orbital rotation parameters $\{\kappa_{pq}\}_{1 \le p,q \le \Norb}$.
Here, motivated by cost saving arguments, we have chosen to optimize separately the CI and orbital coefficients by alternatively diagonalizing the CI matrix after each selection step and then rotating the orbitals until the variational energy for a given number of determinants is minimal.
To do so, we conveniently rewrite the variational energy as
\begin{equation}
\label{eq:Evar_c_k}
\Evar(\bc,\bk) = \mel{\Psivar}{e^{\hk} \hH e^{-\hk}}{\Psivar},
\end{equation}
where $\bc$ gathers the CI coefficients, $\bk$ the orbital rotation parameters, and
\begin{equation}
\hk = \sum_{p < q} \sum_{\sigma} \kappa_{pq} \qty(\cre{p\sigma} \ani{q\sigma} - \cre{q\sigma} \ani{p\sigma})
\end{equation}
is a real-valued one-electron antisymmetric operator, which creates an orthogonal transformation of the orbital coefficients when exponentiated, $\ani{p\sigma}$ ($\cre{p\sigma}$) being the second quantization annihilation (creation) operator which annihilates (creates) a spin-$\sigma$ electron in the (real-valued) spatial orbital $\MO{p}(\br)$. \cite{Helgaker_2013}
Applying the Newton-Raphson method by Taylor-expanding the variational energy to second order around $\bk = \bO$, \ie,
\begin{equation}
\label{eq:EvarTaylor}
\Evar(\bc,\bk) \approx \Evar(\bc,\bO) + \bg \cdot \bk + \frac{1}{2} \bk^{\dag} \cdot \bH \cdot \bk,
\end{equation}
one can iteratively minimize the variational energy with respect to the parameters $\kappa_{pq}$ by setting
\begin{equation}
\bk = - \bH^{-1} \cdot \bg,
\end{equation}
where $\bg$ and $\bH$ are the orbital gradient and Hessian, respectively, both evaluated at $\bk = \bO$.
Their elements are explicitly given by the following expressions: \cite{Bozkaya_2011,Henderson_2014a}
\begin{equation}
\begin{split}
g_{pq}
&= \left. \pdv{\Evar(\bc,\bk)}{\kappa_{pq}}\right|_{\bk=\bO}
\\
&= \sum_{\sigma} \mel{\Psivar}{\comm*{\cre{p\sigma} \ani{q\sigma} - \cre{q\sigma} \ani{p\sigma}}{\hH}}{\Psivar}
\\
&= \cP_{pq} \qty[ \sum_r \left( h_p^r \ \gamma_r^q - h_r^q \ \gamma_p^r \right) + \sum_{rst} \qty( v_{pt}^{rs} \Gamma_{rs}^{qt} - v_{rs}^{qt} \Gamma_{pt}^{rs} ) ],
\end{split}
\end{equation}
and
\begin{equation}
\begin{split}
H_{pq,rs}
& = \left. \pdv{\Evar(\bc,\bk)}{\kappa_{pq}}{\kappa_{rs}}\right|_{\bk=\bO}
\\
& = \cP_{pq} \cP_{rs} \Bigg\{
\frac{1}{2} \sum_{\sigma\sigma'} \mel*{\Psivar}{\comm*{\cre{r \sigma'} \ani{s \sigma'}}{\comm*{\cre{p \sigma} \ani{q \sigma}}{\hH}}}{\Psivar}
\\
& \phantom{\cP_{pq} \cP_{rs} \Bigg\{} + \frac{1}{2} \sum_{\sigma\sigma'} \mel*{\Psivar}{\comm*{\cre{p \sigma} \ani{q \sigma}}{\comm*{\cre{r \sigma'} \ani{s \sigma'}}{\hH}}}{\Psivar}
\Bigg\}
\\
& = \cP_{pq} \cP_{rs} \Bigg\{
\frac{1}{2} \sum_u \qty[ \delta_{qr}(h_p^u \gamma_u^s + h_u^s \gamma_p^u) + \delta_{ps}(h_r^u \gamma_u^q + h_u^q \gamma_u^r)]
\\
& \phantom{\cP_{pq} \cP_{rs} \Bigg\{} - (h_p^s \gamma_r^q + h_r^q \gamma_p^s)
\\
& \phantom{\cP_{pq} \cP_{rs} \Bigg\{} + \frac{1}{2} \sum_{tuv} \delta_{qr}(v_{pt}^{uv} \Gamma_{uv}^{st} + v_{uv}^{st} \Gamma_{pt}^{uv})
\\
& \phantom{\cP_{pq} \cP_{rs} \Bigg\{} + \frac{1}{2} \sum_{tuv} \delta_{ps}(v_{uv}^{qt} \Gamma_{rt}^{uv} + v_{rt}^{uv}\Gamma_{uv}^{qt})]
\\
& \phantom{\cP_{pq} \cP_{rs} \Bigg\{} + \sum_{uv} (v_{pr}^{uv} \Gamma_{uv}^{qs} + v_{uv}^{qs} \Gamma_{ps}^{uv})
\\
& \phantom{\cP_{pq} \cP_{rs} \Bigg\{} - \sum_{tu} (v_{pu}^{st} \Gamma_{rt}^{qu}+v_{pu}^{tr} \Gamma_{tr}^{qu}+v_{rt}^{qu}\Gamma_{pu}^{st} + v_{tr}^{qu}\Gamma_{pu}^{ts})]
\Bigg\},
\end{split}
\end{equation}
where $\delta_{pq}$ is the Kronecker delta, $\cP_{pq} = 1 - (p \leftrightarrow q)$ is a permutation operator,
\begin{subequations}
\begin{gather}
\gamma_p^q = \sum_{\sigma} \mel{\Psivar}{\hat{a}_{p \sigma}^{\dagger} \hat{a}_{q \sigma}^{}}{\Psivar},
\\
\Gamma_{pq}^{rs} = \sum_{\sigma \sigma'} \mel{\Psivar}{\cre{p\sigma} \cre{r\sigma'} \ani{s\sigma'} \ani{q\sigma}}{\Psivar}
\end{gather}
\end{subequations}
are the elements of the one- and two-electron density matrices, and
\begin{subequations}
\begin{gather}
\label{eq:one}
h_p^q = \int \MO{p}(\br) \, \hh(\br) \, \MO{q}(\br) d\br,
\\
\label{eq:two}
v_{pq}^{rs} = \iint \MO{p}(\br_1) \MO{q}(\br_2) \frac{1}{\abs*{\br_1 - \br_2}} \MO{r}(\br_1) \MO{s}(\br_2) d\br_1 d\br_2
\end{gather}
\end{subequations}
are the one- and two-electron integrals, respectively.
Because the size of the CI space is much larger than the orbital space, for each macroiteration, we perform multiple \textit{microiterations} which consist in iteratively minimizing the variational energy \eqref{eq:Evar_c_k} with respect to the $\Norb(\Norb-1)/2$ independent orbital rotation parameters.
Microiterations are stopped when a stationary point is found, \ie, $\norm{\bg}_\infty < \tau$, where $\tau$ is a user-defined threshold which has been set to $10^{-3}$ a.u.~in the present study, and a new CIPSI selection step is performed.
Note that a tight convergence is not critical here as a new set of microiterations is performed at each macroiteration and a new production CIPSI run is performed from scratch using the final set of orbitals.
It is also worth pointing out that, after each orbital rotation, the one- and two-electron integrals defined in Eqs.~\eqref{eq:one} and \eqref{eq:two} have to be updated for the next iteration.
%\begin{equation}
% \Evar = \sum_{pq} h_p^q \gamma_p^q + \frac{1}{2} \sum_{pqrs} v_{pq}^{rs} \Gamma_{pq}^{rs},
%\end{equation}
To enhance the convergence of the microiteration process, we employ a variant of the Newton-Raphson method known as ``trust region''. \cite{Nocedal_1999}
This popular variant defines a region where the quadratic approximation \eqref{eq:EvarTaylor} is an adequate representation of the objective energy function \eqref{eq:Evar_c_k} and it evolves during the optimization process in order to preserve the adequacy via a constraint on the step size preventing it from overstepping, \ie, $\norm{\bk^{(\ell+1)}} \leq \Delta^{(\ell)}$, where $\Delta^{(\ell)}$ is the trust radius at the $\ell$th microiteration.
By introduction a Lagrange multiplier $\lambda$ to control the trust-region size, one obtains $\bk^{(\ell+1)} = - (\bH^{(\ell)} + \lambda \bI)^{-1} \cdot \bg^{(\ell)}$.
The addition of the level shift $\lambda \geq 0$ removes the negative eigenvalues and ensure the positive definiteness of the Hessian matrix by reducing the step size.
By choosing the right value of $\lambda$, the step size is constraint into a hypersphere of radius $\Delta^{(\ell)}$ and is able to evolve from the Newton direction at $\lambda = 0$ to the steepest descent direction as $\lambda$ grows.
The evolution of the trust radius during the optimization and the use of a condition to cancel the step where the energy rises ensure the convergence of the algorithm.
More details can be found in Ref.~\onlinecite{Nocedal_1999}.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Results and discussion}
\label{sec:res}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\begin{figure*}
\includegraphics[width=0.24\textwidth]{Cyclopentadiene_EvsNdet}
\includegraphics[width=0.24\textwidth]{Furan_EvsNdet}
\includegraphics[width=0.24\textwidth]{Imidazole_EvsNdet}
\includegraphics[width=0.24\textwidth]{Pyrrole_EvsNdet}
\\
\includegraphics[width=0.24\textwidth]{Thiophene_EvsNdet}
\includegraphics[width=0.24\textwidth]{Benzene_EvsNdet}
\includegraphics[width=0.24\textwidth]{Pyrazine_EvsNdet}
\includegraphics[width=0.24\textwidth]{Pyridazine_EvsNdet}
\\
\includegraphics[width=0.24\textwidth]{Pyridine_EvsNdet}
\includegraphics[width=0.24\textwidth]{Pyrimidine_EvsNdet}
\includegraphics[width=0.24\textwidth]{Tetrazine_EvsNdet}
\includegraphics[width=0.24\textwidth]{Triazine_EvsNdet}
\caption{$\Delta \Evar$ (solid) and $\Delta \Evar + \EPT$ (dashed) as functions of the number of determinants $\Ndet$ in the variational space for the twelve cyclic molecules represented in Fig.~\ref{fig:mol}.
Two sets of orbitals are considered: natural orbitals (NOs, in red) and optimized orbitals (OOs, in blue).
The CCSDTQ correlation energy is represented as a thick black line.
\label{fig:vsNdet}}
\end{figure*}
\begin{figure*}
\includegraphics[width=0.24\textwidth]{Cyclopentadiene_EvsPT2}
\includegraphics[width=0.24\textwidth]{Furan_EvsPT2}
\includegraphics[width=0.24\textwidth]{Imidazole_EvsPT2}
\includegraphics[width=0.24\textwidth]{Pyrrole_EvsPT2}
\\
\includegraphics[width=0.24\textwidth]{Thiophene_EvsPT2}
\includegraphics[width=0.24\textwidth]{Benzene_EvsPT2}
\includegraphics[width=0.24\textwidth]{Pyrazine_EvsPT2}
\includegraphics[width=0.24\textwidth]{Pyridazine_EvsPT2}
\\
\includegraphics[width=0.24\textwidth]{Pyridine_EvsPT2}
\includegraphics[width=0.24\textwidth]{Pyrimidine_EvsPT2}
\includegraphics[width=0.24\textwidth]{Tetrazine_EvsPT2}
\includegraphics[width=0.24\textwidth]{Triazine_EvsPT2}
\caption{$\Delta \Evar$ as a function of $\EPT$ for the twelve cyclic molecules represented in Fig.~\ref{fig:mol}.
Two sets of orbitals are considered: natural orbitals (NOs, in red) and optimized orbitals (OOs, in blue).
The four-point weighted linear fit using the four largest variational wave functions for each set is depicted as a dashed black line.
The weights are taken as the inverse square of the perturbative corrections.
The CCSDTQ correlation energy is also represented as a thick black line.
\label{fig:vsEPT2}}
\end{figure*}
\begin{squeezetable}
\begin{table*}
\caption{Total energy $E$ (in \SI{}{\hartree}) and correlation energy $\Delta E$ (in \SI{}{\milli\hartree}) for the frozen-core ground state of five-membered rings in the cc-pVDZ basis set.
For the CIPSI estimates of the correlation energy, the fitting error associated with the four-point linear fit is reported in parenthesis.
\label{tab:Tab5-VDZ}}
\begin{ruledtabular}
\begin{tabular}{lcccccccccc}
& \mc{2}{c}{Cyclopentadiene} & \mc{2}{c}{Furan} & \mc{2}{c}{Imidazole} & \mc{2}{c}{Pyrrole} & \mc{2}{c}{Thiophene} \\
\cline{2-3} \cline{4-5} \cline{6-7} \cline{8-9} \cline{10-11}
Method & $E$& $\Delta E$ & $E$ & $\Delta E$ & $E$ & $\Delta E$ & $E$ & $\Delta E$ & $E$ & $\Delta E$ \\
\hline
HF & $-192.8083$ & & $-228.6433$ & & $-224.8354$ & & $-208.8286$ & &$-551.3210$ & \\
\hline
MP2 & $-193.4717$ & $-663.4$ & $-229.3508$ & $-707.5$ & $-225.5558$ & $-720.4$ & $-209.5243$ & $-695.7$ & $-551.9825$ & $-661.5$ \\
MP3 & $-193.5094$ & $-701.0$ & $-229.3711$ & $-727.8$ & $-225.5732$ & $-737.8$ & $-209.5492$ & $-720.6$ & $-552.0104$ & $-689.4$ \\
MP4 & $-193.5428$ & $-734.5$ & $-229.4099$ & $-766.6$ & $-225.6126$ & $-777.2$ & $-209.5851$ & $-756.5$ & $-552.0476$ & $-726.6$ \\
MP5 & $-193.5418$ & $-733.4$ & $-229.4032$ & $-759.9$ & $-225.6061$ & $-770.8$ & $-209.5809$ & $-752.3$ & $-552.0426$ & $-721.6$\\
\hline
CC2 & $-193.4782$ & $-669.9$ & $-229.3605$ & $-717.2$ & $-225.5644$ & $-729.0$ & $-209.5311$ & $-702.5$ & $-551.9905$ & $-669.5$ \\
CC3 & $-193.5449$ & $-736.6$ & $-229.4090$ & $-765.7$ & $-225.6115$ & $-776.1$ & $-209.5849$ & $-756.3$ & $-552.0473$ & $-726.3$ \\
CC4 & $-193.5467$ & $-738.4$ & $-229.4102$ & $-766.9$ & $-225.6126$ & $-777.2$ & $-209.5862$ & $-757.6$ & $-552.0487$ & $-727.7$ \\
\hline
CCSD & $-193.5156$ & $-707.2$ & $-229.3783$ & $-735.0$ & $-225.5796$ & $-744.2$ & $-209.5543$ & $-725.7$ & $-552.0155$ & $-694.5$ \\
CCSDT & $-193.5446$ & $-736.2$ & $-229.4076$ & $-764.3$ & $-225.6099$ & $-774.6$ & $-209.5838$ & $-755.2$ & $-552.0461$ & $-725.1$ \\
CCSDTQ & $-193.5465$ & $-738.2$ & $-229.4100$ & $-766.7$ & $-225.6123$ & $-776.9$ & $-209.5860$ & $-757.4$ & $-552.0485$ & $-727.5$ \\
\hline
CCSD(T) & $-193.5439$ & $-735.6$ & $-229.4073$ & $-764.0$ & $-225.6099$ & $-774.5$ & $-209.5836$ & $-754.9$ & $-552.0458$ & $-724.8$
\\
\hline
CIPSI & & $-739.3(1)$ & & $-768.1(2)$ & & $-778.3(1)$ & & $-758.4(2)$ & & $-729.1(3)$\\
\end{tabular}
\end{ruledtabular}
\end{table*}
\end{squeezetable}
\begin{squeezetable}
\begin{table*}
\caption{Total energy $E$ (in \SI{}{\hartree}) and correlation energy $\Delta E$ (in \SI{}{\milli\hartree}) for the frozen-core ground state of six-membered rings in the cc-pVDZ basis set.
For the CIPSI estimates of the correlation energy, the fitting error associated with the four-point linear fit is reported in parenthesis.
\label{tab:Tab6-VDZ}}
\begin{ruledtabular}
\begin{tabular}{lcccccccccccccc}
& \mc{2}{c}{Benzene} & \mc{2}{c}{Pyrazine} & \mc{2}{c}{Pyridazine} & \mc{2}{c}{Pyridine} & \mc{2}{c}{Pyrimidine} & \mc{2}{c}{Tetrazine} & \mc{2}{c}{Triazine} \\
\cline{2-3} \cline{4-5} \cline{6-7} \cline{8-9} \cline{10-11} \cline{12-13} \cline{14-15}
Method & $E$ & $\Delta E$ & $E$ & $\Delta E$ & $E$ & $\Delta E$ & $E$ & $\Delta E$
& $E$ & $\Delta E$ & $E$ & $\Delta E$ & $E$ & $\Delta E$ \\
\hline
HF & $-230.7222$ & & $-262.7030$ & & $-262.6699$ & & $-246.7152$ & & $-262.7137$ & & $-294.6157$ & & $-278.7173$ \\
\hline
MP2 & $-231.5046$ & $-782.3$ & $-263.5376$ & $-834.6$ & $-263.5086$ & $-838.7$ & $-247.5227$ & $-807.5$ & $-263.5437$ & $-830.1$ & $-295.5117$ & $-895.9$ & $-279.5678$ & $-850.5$\\
MP3 & $-231.5386$ & $-816.4$ & $-263.5567$ & $-853.7$ & $-263.5271$ & $-857.3$ & $-247.5492$ & $-834.0$ & $-263.5633$ & $-849.6$ & $-295.5152$ & $-899.5$ & $-279.5809$ & $-863.6$ \\
MP4 & $-231.5808$ & $-858.5$ & $-263.6059$ & $-902.9$ & $-263.5778$ & $-907.9$ & $-247.5951$ & $-879.9$ & $-263.6129$ & $-899.3$ & $-295.5743$ & $-958.6$ & $-279.6340$ & $-916.7$ \\
MP5 & $-231.5760$ & $-853.8$ & $-263.5968$ & $-893.8$ & $-263.5681$ & $-898.3$ & $-247.5881$ & $-872.9$ & $-263.6036$ & $-890.0$ & $-295.5600$ & $-944.3$ & $-279.6228$ & $-905.4$ \\
\hline
CC2 & $-231.5117$ & $-789.4$ & $-263.5475$ & $-844.5$ & $-263.5188$ & $-848.9$ & $-247.5315$ & $-816.3$ & $-263.5550$ & $-841.3$ & $-295.5247$ & $-909.0$ & $-279.5817$ & $-864.4$ \\
CC3 & $-231.5814$ & $-859.1$ & $-263.6045$ & $-901.5$ & $-263.5761$ & $-906.2$ & $-247.5948$ & $-879.6$ & $-263.6120$ & $-898.4$ & $-295.5706$ & $-954.9$ & $-279.6329$ & $-915.6$ \\
CC4 & $-231.5828$ & $-860.6$ & $-263.6056$ & $-902.6$ & $-263.5773$ & $-907.5$ & $-247.5960$ & $-880.8$ & $-263.6129$ & $-899.3$ & $-295.5716$ & $-955.9$ & $-279.6334$ & $-916.1$ \\
\hline
CCSD & $-231.5440$ & $-821.8$ & $-263.5640$ & $-861.0$ & $-263.5347$ & $-864.9$ & $-247.5559$ & $-840.7$ & $-263.5716$ & $-858.0$ & $-295.5248$ & $-909.1$ & $-279.5911$ & $-873.8$ \\
CCSDT & $-231.5802$ & $-857.9$ & $-263.6024$ & $-899.4$ & $-263.5739$ & $-904.0$ & $-247.5931$ & $-877.9$ & $-263.6097$ & $-896.1$ & $-295.5673$ & $-951.6$ & $-279.6300$ & $-912.7$ \\
CCSDTQ & $-231.5826$ & $-860.4$ & $-263.6053$ & $-902.3$ & $-263.5770$ & $-907.1$ & $-247.5960$ & $-880.8$ & $-263.6126$ & $-899.0$ & $-295.5712$ & $-955.4$ & $-279.6331$ & $-915.8$ \\
\hline
CCSD(T) & $-231.5798$ & $-857.5$ & $-263.6024$ & $-899.4$ & $-263.5740$ & $-904.1$ & $-247.5929$ & $-877.7$ & $-263.6099$ & $-896.2$ & $-295.5680$ & $-952.2$ & $-279.6305$ & $-913.1$ \\
\hline
CIPSI & & $-863.0(4)$ & & $-904.6(6)$ & & $-908.8(2)$ & & $-883.4(0)$ & & $-900.4(4)$ & & $-957.3(2)$ & & $-918.5(5)$\\
\end{tabular}
\end{ruledtabular}
\end{table*}
\end{squeezetable}
\begin{squeezetable}
\begin{table}
\caption{
Extrapolated correlation energies $\Delta \Evar$ (in \SI{}{\milli\hartree}) for the twelve cyclic molecules represented in Fig.~\ref{fig:mol} and their associated fitting errors (in \SI{}{\milli\hartree}) obtained via weighted linear fits with a varying number of points.
The weights are taken as the inverse square of the perturbative corrections.
For a $n$-point fit, the $n$ largest the largest variational wave functions are used.
\label{tab:fit}}
\begin{ruledtabular}
\begin{tabular}{lccc}
Molecule & Number of & \mc{2}{c}{Fitting parameters} \\
\cline{3-4}
& fitting points & $\Delta \Evar$ & Fitting error \\
\hline
Cyclopentadiene & 3 & $-739.295$ & $0.199$ \\
& 4 & $-739.309$ & $0.088$ \\
& 5 & $-739.230$ & $0.074$ \\
& 6 & $-739.304$ & $0.072$ \\
& 7 & $-739.292$ & $0.055$ \\
& 8 & $-739.119$ & $0.131$ \\
\hline
Furan & 3 & $-767.790$ & $0.064$ \\
& 4 & $-768.104$ & $0.196$ \\
& 5 & $-768.194$ & $0.135$ \\
& 6 & $-768.060$ & $0.131$ \\
& 7 & $-768.086$ & $0.101$ \\
& 8 & $-767.904$ & $0.154$ \\
\hline
Imidazole & 3 & $-778.295$ & $0.356$ \\
& 4 & $-778.270$ & $0.150$ \\
& 5 & $-778.178$ & $0.105$ \\
& 6 & $-778.174$ & $0.072$ \\
& 7 & $-778.051$ & $0.099$ \\
& 8 & $-777.992$ & $0.089$ \\
\hline
Pyrrole & 3 & $-758.650$ & $0.321$ \\
& 4 & $-758.389$ & $0.174$ \\
& 5 & $-758.460$ & $0.110$ \\
& 6 & $-758.352$ & $0.100$ \\
& 7 & $-758.347$ & $0.075$ \\
& 8 & $-758.357$ & $0.059$ \\
\hline
Thiophene & 3 & $-728.744$ & $0.691$ \\
& 4 & $-729.052$ & $0.331$ \\
& 5 & $-728.948$ & $0.203$ \\
& 6 & $-728.987$ & $0.140$ \\
& 7 & $-729.067$ & $0.117$ \\
& 8 & $-728.876$ & $0.162$ \\
\hline
Benzene & 3 & $-862.325$ & $0.279$ \\
& 4 & $-863.024$ & $0.424$ \\
& 5 & $-862.890$ & $0.266$ \\
& 6 & $-862.360$ & $0.383$ \\
& 7 & $-862.083$ & $0.339$ \\
& 8 & $-861.711$ & $0.370$ \\
\hline
Pyrazine & 3 & $-904.867$ & $1.420$ \\
& 4 & $-904.588$ & $0.650$ \\
& 5 & $-904.550$ & $0.385$ \\
& 6 & $-903.982$ & $0.439$ \\
& 7 & $-903.746$ & $0.359$ \\
& 8 & $-903.549$ & $0.311$ \\
\hline
Pyridazine & 3 & $-909.292$ & $0.024$ \\
& 4 & $-908.808$ & $0.230$ \\
& 5 & $-908.820$ & $0.133$ \\
& 6 & $-908.342$ & $0.303$ \\
& 7 & $-908.368$ & $0.224$ \\
& 8 & $-908.229$ & $0.198$ \\
\hline
Pyridine & 3 & $-883.363$ & $0.047$ \\
& 4 & $-883.413$ & $0.029$ \\
& 5 & $-882.700$ & $0.405$ \\
& 6 & $-882.361$ & $0.341$ \\
& 7 & $-882.023$ & $0.330$ \\
& 8 & $-881.732$ & $0.322$ \\
\hline
Pyrimidine & 3 & $-900.817$ & $0.726$ \\
& 4 & $-900.383$ & $0.356$ \\
& 5 & $-900.496$ & $0.214$ \\
& 6 & $-900.698$ & $0.190$ \\
& 7 & $-900.464$ & $0.206$ \\
& 8 & $-900.226$ & $0.227$ \\
\hline
Tetrazine & 3 & $-957.559$ & $0.246$ \\
& 4 & $-957.299$ & $0.160$ \\
& 5 & $-957.869$ & $0.349$ \\
& 6 & $-957.744$ & $0.247$ \\
& 7 & $-957.709$ & $0.183$ \\
& 8 & $-957.558$ & $0.176$ \\
\hline
Triazine & 3 & $-919.596$ & $0.105$ \\
& 4 & $-918.457$ & $0.538$ \\
& 5 & $-918.355$ & $0.312$ \\
& 6 & $-918.206$ & $0.226$ \\
& 7 & $-917.876$ & $0.267$ \\
& 8 & $-917.533$ & $0.308$ \\
\end{tabular}
\end{ruledtabular}
\end{table}
\end{squeezetable}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsection{CIPSI estimates}
\label{sec:cipsi_res}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
We first study the convergence of the CIPSI correlation energy $\Delta \Evar = \Evar - \EHF$ as a function of the number of determinants.
Our motivation here is to generate FCI-quality reference correlation energies for the twelve cyclic molecules represented in Fig.~\ref{fig:mol} in order to benchmark, in a second time, the performance and convergence properties of various mainstream MP and CC methods (see Sec.~\ref{sec:mpcc_res}).
For the natural and optimized orbital sets, we report, in Fig.~\ref{fig:vsNdet}, the evolution of the variational correlation energy $\Delta \Evar$ and its perturbatively corrected value $\Delta \Evar + \EPT$ with respect to the number of determinants for each cyclic molecule.
As compared to natural orbitals (solid red lines), one can see that, for a given number of determinants, the use of optimized orbitals greatly lowers $\Delta \Evar$ (solid blue lines).
Adding the perturbative correction $\EPT$ yields similar curves for both sets of orbitals (dashed lines).
This indicates that, for a given number of determinants, $\EPT$ (which provides a qualitative idea to the distance to the FCI limit) is much smaller for optimized orbitals than for natural orbitals.
This is further evidenced in Fig.~\ref{fig:vsEPT2} where we show the behavior of $\Delta \Evar$ as a function of $\EPT$ for both sets of orbitals.
The four-point weighted linear fit using the four largest variational wave functions are also represented (dashed black lines), while the CCSDTQ correlation energy (solid black line) is also reported for comparison purposes.
From Fig.~\ref{fig:vsEPT2}, it is clear that the behavior of $\Delta \Evar$ is much more linear and produces smaller $\EPT$ values, hence facilitating the extrapolation procedure to the FCI limit.
\begin{figure}
\includegraphics[width=\linewidth]{Benzene_EvsNdetLO}
\caption{$\Delta \Evar$ (solid) and $\Delta \Evar + \EPT$ (dashed) as functions of the number of determinants $\Ndet$ in the variational space for the benzene molecule.
Three sets of orbitals are considered: natural orbitals (NOs, in red), localized orbitals (LOs, in red), and optimized orbitals (OOs, in blue).
The CCSDTQ correlation energy is represented as a thick black line.
\label{fig:BenzenevsNdet}}
\end{figure}
Figure~\ref{fig:BenzenevsNdet} compares the convergence of $\Delta \Evar$ for the natural, localized, and optimized sets of orbitals.
As mentioned in Sec.~\ref{sec:compdet}, although both sets break the spatial symmetry to take advantage of the local nature of electron correlation, optimized orbitals further improve on the use of localized orbitals.
More quantitatively, optimized orbitals produce the same variational energy as localized orbitals with, roughly, a ten-fold reduction in the number of determinants.
A similar improvement is observed going from natural to localized orbitals.
\titou{T2: comment on PT2 for localized orbitals?}
Accordingly to these observations, all our FCI correlation energy estimates have been produced with the set of optmized orbitals.
To do so, we have then extrapolated to the orbital-optimized CIPSI calculations to $\EPT = 0$ via a weighted linear extrapolation using the four largest variational wave functions.
These estimates are reported in Tables \ref{tab:Tab5-VDZ} and \ref{tab:Tab6-VDZ} for the five- and six-membered rings, respectively.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsection{Benchmark of CC and MP methods}
\label{sec:mpcc_res}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\begin{figure*}
\includegraphics[width=0.32\textwidth]{Cyclopentadiene_MPCC}
\includegraphics[width=0.32\textwidth]{Furan_MPCC}
\includegraphics[width=0.32\textwidth]{Imidazole_MPCC}
\\
\includegraphics[width=0.32\textwidth]{Pyrrole_MPCC}
\includegraphics[width=0.32\textwidth]{Thiophene_MPCC}
\includegraphics[width=0.32\textwidth]{Benzene_MPCC}
\\
\includegraphics[width=0.32\textwidth]{Pyrazine_MPCC}
\includegraphics[width=0.32\textwidth]{Pyridazine_MPCC}
\includegraphics[width=0.32\textwidth]{Pyridine_MPCC}
\\
\includegraphics[width=0.32\textwidth]{Pyrimidine_MPCC}
\includegraphics[width=0.32\textwidth]{Tetrazine_MPCC}
\includegraphics[width=0.32\textwidth]{Triazine_MPCC}
\caption{Convergence of the correlation energy (in \SI{}{\milli\hartree}) as a function of the computational cost for the twelve cyclic molecules represented in Fig.~\ref{fig:mol}.
Three series of methods are considered: i) MP2, MP3, MP4, and MP5 (blue), ii) CC2, CC3, and CC4 (green), and iii) CCSD, CCSDT, CCSDTQ (red)
The CIPSI estimate of the correlation energy is represented as a black line.
\label{fig:MPCC}}
\end{figure*}
%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Conclusion}
\label{sec:ccl}
%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%
\begin{acknowledgements}
This work was performed using HPC resources from GENCI-TGCC (2020-gen1738) and from CALMIP (Toulouse) under allocation 2021-18005.
This project has received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme (Grant agreement No.~863481).
\end{acknowledgements}
%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%
\section*{Data availability statement}
%%%%%%%%%%%%%%%%%%%%%%%%%
The data that support the findings of this study are openly available in Zenodo at \href{http://doi.org/XX.XXXX/zenodo.XXXXXXX}{http://doi.org/XX.XXXX/zenodo.XXXXXXX}.
%%%%%%%%%%%%%%%%%%%%%%%%%
\bibliography{Ec}
%%%%%%%%%%%%%%%%%%%%%%%%%
\end{document}
Reference in New Issue View Git Blame Copy Permalink