In this paper, we discuss the extension of the recently introduced subsystem embedding subalgebra coupled cluster (SES-CC) formalism to unitary CC formalisms. In analogy to the standard single-reference SES-CC formalism, its unitary CC extension allows one to include the dynamical (outside the active space) correlation effects in an SES induced complete active space (CAS) effective Hamiltonian. In contrast to the standard single-reference SES-CC theory, the unitary CC approach results in a Hermitian form of the effective Hamiltonian. Additionally, for the double unitary CC (DUCC) formalism, the corresponding CAS eigenvalue problem provides a rigorous separation of external cluster amplitudes that describe dynamical correlation effects—used to define the effective Hamiltonian—from those corresponding to the internal (inside the active space) excitations that define the components of eigenvectors associated with the energy of the entire system. The proposed formalism can be viewed as an efficient way of downfolding many-electron Hamiltonian to the low-energy model represented by a particular choice of CAS. In principle, this technique can be extended to any type of CAS representing an arbitrary energy window of a quantum system. The Hermitian character of low-dimensional effective Hamiltonians makes them an ideal target for several types of full configuration interaction type eigensolvers. As an example, we also discuss the algebraic form of the perturbative expansions of the effective DUCC Hamiltonians corresponding to composite unitary CC theories and discuss possible algorithms for hybrid classical and quantum computing. Given growing interest in quantum computing, we provide energies for H_{2} and Be systems obtained with the quantum phase estimator algorithm available in the Quantum Development Kit for the approximate DUCC Hamiltonians.

## I. INTRODUCTION

Even though quantum chemistry and materials science communities have expended a great deal of effort designing numerous methods to describe collective behavior of electrons in correlated systems, the fundamental understanding of these processes in many systems is still inhibited by an exponential growth in computation associated with representing the many-body electronic wave function. These problems usually occur for systems characterized by small energy gaps between occupied and unoccupied one-particle states, where the use of advanced tools to describe electron correlation effects is a prerequisite. Several classes of many-body approaches based on the inclusion of higher-rank collective excitations,^{1,2} multireference concepts,^{3–13} and symmetry-breaking mechanisms^{14–20} have resulted in various approximations trying to describe static and dynamic correlation effects. Coupled cluster,^{21–25} density matrix renormalization group (DMRG),^{26–28} and density matrix methods^{29,30} have already demonstrated their efficiency in coping with complicated electron correlation effects in modest size molecular and materials systems. One should also emphasize the emergence of several new classes of methodologies for treating strong correlation effects (both static and dynamical) in molecular systems that include the full CI Quantum Monte Carlo (FCIQMC) method,^{31,32} adaptive configuration interaction formulations,^{33,34} and the novel idea of mixing deterministic CC and stochastic CIQMC calculations to obtain highly accurate or even numerically exact energetics in single- and multireference systems with largely polynomial costs.^{35,36} However, the applicability of these methods for larger systems is still defined by a trade-off between accuracy and computational costs. In many cases, this situation evolves into a scientific stalemate. While strongly correlated systems successfully elude the mainstream theoretical modeling exemplified by low-rank methods, mainly density functional theory (DFT), the applicability of very accurate yet very expensive many-body formalisms is significantly limited by computational resources offered by conventional computers. In this regard, enabling mathematically rigorous models where correlation effects are downfolded into a low-dimensionality space offers a unique chance to permanently eliminate the inherent bias/biases of currently employed many-body theories. The diagonalization of the resulting low-dimensionality effective Hamiltonians is also an ideal target for various algorithms including algorithms for classical computers as well as novel algorithms capable of taking advantage of emerging quantum information systems.^{37–45}

Of the various quantum chemistry methods, CC theory has become the *de facto* standard high accuracy calculations for nuclear, atomic, and molecular systems. The framework of CC theory, when combined with recent developments of our own, makes it well suited to address systems in which correlation effects are downfolded onto smaller spaces. We have recently shown in Refs. 46 and 47 that we can go beyond solving the ground-state CC equations in the conventional iterative manner, by decoupling the excitations into two disjoint sets as shown in Fig. 1. The set of excitations *A* is obtained from predefined classes of excitations [or subalgebra(s)], and the *B* set contains all the remaining wave function parameters needed to describe the whole system. This technique, known as subsystem embedding subalgebra coupled cluster (SES-CC), introduces the concept of active spaces in a natural way and provides a mathematically rigorous procedure for downfolding many-body effects for a subset of excitations into an effective Hamiltonian for tractable eigenvalue problems that provide the ground state energy for the full conventional CC calculation (see Fig. 1). The authors of this paper view the SES formalism, which has been shown to work well compared to standard CC approximations, as potential generalized extension of previously developed embedding methods (see Ref. 48) and Löwdin partitioning techniques^{49,50} which will be discussed in future work. However, as a consequence of the single-reference CC ansatz, the standard SES-CC effective Hamiltonians are not Hermitian, which precludes them from efficient utilization with full configuration interaction (FCI) type diagonalizers, such as those used in quantum algorithms and the DMRG method. So there is a natural need for extending SES-CC approach to the unitary CC formalisms which assure the Hermitian character of resulting effective Hamiltonians. Several key problems associated with this extension are related to the following questions: (1) What is the optimal way of utilizing diagonalizers for enabling accurate many-body formulations given current and near-term time perspectives? (2) How can many-body formulations be tuned (using mathematically rigorous procedures and operating on the parameters defining a given many-body approach) to existing and future computational systems to assure applications to more realistic and challenging problems than currently possible?

We address the above questions by developing a new downfolding strategy pertaining to the SES-CC formalism which enables one to represent the effective (or downfolded) many-body Hamiltonian in a much smaller orbital space (or active space), where the correlation effects from outside the active space are included in the form of a similarity unitary transformation involving only parameters corresponding to high-energy components of the corresponding wave function. The proposed computational scheme can give rise to an efficient hybrid classical/quantum computational approach: while the high-energy components of the wave function and the second quantized form of the effective Hamiltonian are calculated using classical computers, the diagonalization of the effective Hamiltonian is achieved by employing quantum algorithms. As a specific example of the above formalism, in this paper, we will consider a formalism where similarity unitary transformation is defined to downfolding an essential part of the virtual orbital space. We will introduce simple and numerically affordable models based on the inclusion of single and double excitations in relevant unitary operators and derive the algebraic form of downfolded Hamiltonian by combining particle-hole and physical vacuum representations of second quantized operators. We will also provide a numerical illustration of these procedures by applying quantum phase algorithm as implemented in the Quantum Development Kit (QDK)^{51} for H_{2} and Be systems. We will also illustrate the efficiency of the double unitary CC (DUCC) algorithms in compressing out-of-active-space correlation effects for several basis sets of increasing size for simple approximations of DUCC Hamiltonians.

## II. STANDARD SINGLE-REFERENCE FORMULATION

The standard single-reference CC formulations are predicated on the assumption that there exists a reasonable choice of a Slater determinant |Φ⟩ that can provide an adequate starting point for the exponential (CC) expansion of the ground-state electronic wave function |Ψ⟩,

where the so-called cluster operator *T* can be expanded in terms of its many-body components *T*_{k},

In the second quantized form, each *T*_{k} component can be expressed as

where indices *i*_{1}, *i*_{2}, … (*a*_{1}, *a*_{2}, …) refer to occupied (unoccupied) spin orbitals in the reference function |Φ⟩, $ta1\u2026aki1\u2026ik$ are cluster amplitudes, and $Ei1\u2026ika1\u2026ak$ are excitation operators. The excitation operators are defined through strings of standard creation ($ap\u2020$) and annihilation operators (*a*_{p}),

In the particle-hole formalism, the excitation operators are expressed in terms of particle/hole creation operators only. When acting on the reference function |Φ⟩, the $Ei1\u2026ika1\u2026ak$ operators produce the so-called excited configurations $|\Phi i1\u2026ika1\u2026ak\u2009$ defined as

Upon the substitution of the ansatz (1) into the Schrödinger equation, one gets the energy-dependent form of the CC equations,

where *P* and *Q* are projection operators onto the reference function (*P* = |Φ⟩⟨Φ|) and onto the excited configurations generated by the *T* operator when acting on the reference function,

A careful diagrammatic analysis^{1} leads to an equivalent (*at the solution*), energy-independent form of the CC equations for the cluster amplitudes,

where $H\xaf=e\u2212THeT$ is referred to as the similarity transformed Hamiltonian. It can also be shown that $H\xaf$ is expressible in terms of connected diagrams, i.e., $H\xaf=(HeT)C$, where *C* designates a connected form of a given operator expression. Once the *T* operator is determined by solving Eq. (8), the corresponding energy is given by the standard formula

The second-quantized form of the many-body Hamiltonian defined by up to pairwise interactions is given by the formula

where *p*, *q*, *r*, and *s* indices run over a spin-orbital basis involved in a given algebraic approximation (associated with the use of finite-dimensional one-particle basis set to discretize the Schrödinger equation), and $hqp$ and $vrspq$ represent one- and two-electron integrals (in the above representation, the $vrspq$ tensor is antisymmetric with respect to permutations among the sets of lower and upper spin-orbital indices). In typical molecular applications, based, for example, on the delocalized Hartree-Fock molecular orbitals expanded in the Gaussian basis set, the number of terms defining second-quantized Hamiltonian is proportional to *N*^{4}, where *N* stands for the number of basis set functions. It has recently been shown that using a different kind of basis set, namely, the plane wave dual basis set, the number of terms can be reduced from *N*^{4} to *N*^{2}.^{42} Similar reduction can also be achieved for Gaussian basis sets when combined Cholesky and singular value decompositions are employed to represent two-electron integrals (see Refs. 43 and 52 for details).

From the point of quantum computing applications, the net effect of the number of basis set functions and the number of nonvanishing terms in Hamiltonian define the circuit depth that determines the efficiency of quantum algorithms. The reduction in the number of non-negligible terms may also be achieved by employing localization techniques for Gaussian basis sets.^{53–55} An interesting alternative to the localized basis sets that may especially impact the choice of the initial state is the use of the Brückner orbitals^{56–61} that maximize the overlap |⟨Φ_{B}|Ψ⟩| between the normalized lowest energy Slater determinant |Φ_{B}⟩ and the correlated wave function |Ψ⟩. Given the form of this condition, one should also expect more efficient utilization of phase estimation techniques when Brückner orbitals are employed in the context of various quantum algorithms such as Trotterization.

In the exact wave function limit, the excitation level *m* is equal to the number of correlated electrons (*N*_{e}), while in the approximate CC formulations, *m* ≪ *N*_{e}. In this way, one can define standard approximations such as CC approach with singles and doubles (CCSD) (*m* = 2),^{25} CC approach with singles, doubles, and triples (CCSDT) (*m* = 3),^{62–64} and CC approach with singles, doubles, triples, and quadruples (CCSDTQ) (*m* = 4).^{65,66} Various standard CC approximations have been successfully applied to describe various many-body systems across energy and spatial scales ranging from nuclear matter to molecular and extended/periodic systems.^{67–76} The success of the CC methods in capturing correlation effects can be attributed to two main factors: (1) its size-extensivity, i.e., proper scaling of the energy with number of the particles, which is a direct consequence of connected character of diagrams contributing to the CC equations, and (2) possibility of approximating higher-order excitations by products of low-rank cluster operators.

## III. PROPERTIES OF CC SUBSYSTEM EMBEDDING SUBALGEBRAS

Certain properties of CC equations are inextricably associated with the possibility of partitioning of cluster operators in the CC wave function into components corresponding to various subalgebras of excitation Lie algebra denoted here as $g(N)$, which is generated by all excitation operators $Eilal=aal\u2020ail$. In a recent paper,^{46} we have analyzed properties of the CC equations stemming from the presence of CC-approximation-specific subalgebras of excitations that can naturally be identified with the so-called active spaces that are frequently used in many areas of quantum chemistry and physics. Algebraic properties of these subalgebras provide a means to recast the CC equations in the form of a set of eigenvalue problems and a set of equations that couple these eigenvalue problems. Although there were several attempts to re-express the CC equations as a nonlinear eigenvalue problem (either in the context of dressed configuration interaction Hamiltonian,^{77} inclusion of high-order excitations,^{78} or the analysis of multiple solutions of CC equations^{79}), in contrast to earlier efforts, all parameters (a subset of the cluster amplitudes) defining the matrices to be diagonalized in the eigenvalue subproblems are entirely decoupled from those parameters (also a subset of the cluster amplitudes) that define components of the corresponding eigenvectors. In particular, it was shown that through the similarity transformation of the electronic Hamiltonian, it is possible to downfold it to the form that acts in the active space and provides corresponding CC energy as its eigenvalue value. In contrast to the full electronic Hamiltonian, its effective active-space normal-product form representation involves only creation/annihilation operators carrying active-space indices. In typical applications, the number of active spin orbitals (*N*_{act}) is significantly smaller compared to the total number of spin orbitals *N*_{S}, i.e., *N*_{act} ≪ *N*_{S}.

Let us start the discussion by introducing basic notions defining subalgebras. An important class of subalgebras of the $g(N)$ excitation algebra is closely related to ideas underlying the active-space concepts in quantum chemistry, where one can define subalgebras $h$ through all possible excitations $Ei1\u2026ima1\u2026am$ that excite electrons from a subset of active occupied orbitals (denoted as R) to a subset of active virtual orbitals (denoted as S). These subalgebras will be denoted as $g(N)(R,S)$. The $g(N)(R,S)$ subalgebras can also be viewed as generators of various complete active spaces [CAS(R,S)] spanned by the reference function |Φ⟩ and all excited configurations obtained by acting with elements of $g(N)(R,S)$ onto |Φ⟩ (see Fig. 2). Specific examples of these subalgebras contain subalgebras involving all occupied and selected (S) virtual orbitals or selected set of occupied orbitals (R) and all virtual orbitals, which will be denoted for short as $g(N)(S)$ and $g(N)(R)$, respectively. In an alternative notation, we will denote $g(N)(R,S)$ as $g(N)(xR,yS)$, where numbers *x* and *y* refer to the number of orbitals included in sets R and S, respectively. Special subalgebras $g(N)(R)$ and $g(N)(S)$ will be denoted as $g(N)(xR)$ and $g(N)(yS)$. In this paper, we will entirely focus on active-space excitation subalgebras for the closed-shell single-reference CC formulations.

An important property of the excitation subalgebras is the fact that, in general, an arbitrary cluster operator *T* can be decomposed in a part that belongs to subalgebra of interest $h$ [denoted as an internal part, $Tint(h)$] and a part that belongs to $g(N)\u2212h$ [denoted as an external part, $Text(h)$], i.e.,

The *T*_{int} operator is a generator of the CAS(R,S), and the corresponding amplitudes of the *T*_{k,int} components are labeled exclusively by active-space spin-orbital indices, while the indices for the *T*_{k,ext} components of *T*_{ext} consist of one or more spin-orbitals outside of the active space. This decomposition entails decomposition of the corresponding CC wave function,

where the CAS-type wave function $|\Psi (h)\u2009$ is defined as

In (12), we use the fact that [*T*_{ext}, *T*_{int}] = 0, which can be seen immediately from the fact that the terms in the cluster operator in particle hole representation can be expressed as $Ei1,\u2026,ika1,\u2026,ak=ba1\u2020\cdots bak\u2020bik\u2020\cdots bi1\u2020$, for anticommuting $b\u2113\u2020$, where *b*_{ℓ}/$b\u2113\u2020$ operators are defined as

and *O* and *V* represent sets of occupied and unoccupied spin orbitals in |Φ⟩, respectively.

The subalgebras $h$ that satisfy the following two important requirements:

$|\Psi (h)\u2009=eTint(h)|\Phi \u2009$ is characterized by the same symmetry properties as |Ψ⟩ and |Φ⟩ vectors (for example, spin and spatial symmetries).

The $eTint(h)|\Phi \u2009$ ansatz generates the FCI expansion for the subsystem defined by the CAS corresponding to the $h$ subalgebra,

play an important role in further analysis of the structure of the CC equations, and we will refer to all subalgebras satisfying requirements (1) and (2) as subsystem embedding subalgebras (SESs).

In Ref. 46, it has been shown that $g(N)(1R,yS)$ along with $g(N)(xR,1S)$ form SESs for CCSD and $g(N)(2R,yS)$ along with $g(N)(xR,2S)$ are SESs for CCSDTQ. The largest subalgebras in these classes are $g(N)(1R)$ and $g(N)(2R)$ subalgebras. In Ref. 46, we have also demonstrated that for any SES $h$ (corresponding to some CC approximation), the standard,

and hybrid,

representations of CC equations are equivalent at the solution. In Eqs. (15)–(19), we used the following notations: (1) the projection operators $Qint(h)$ and $Qext(h)$ ($Q=Qint(h)+Qext(h)$) project onto subspaces spanned by all excited configurations generated by acting $Tint(h)$ and $Text(h)$ onto the reference function, respectively, (2) $H\xafN=H\xaf\u2212\u2009\Phi |H\xaf|\Phi \u2009$ is the normal product form of the similarity transformed Hamiltonian $H\xaf$, (3) the $H\xafext$ operator is defined as $H\xafext\u2261H\xafext(h)=e\u2212Text(h)HeText(h)$, and (4) for notational simplicity, we used the following notational convention:

The above mentioned equivalence means that cluster amplitudes corresponding to excitations included in $h$ can be obtained in a diagonalization procedure. Moreover, the standard form of the CC energy expression [given by Eq. (9)] is a special case of Eq. (18) corresponding to $h=g(N)(0)$ [where $g(N)(0)$ contains no excitations]—in this case, $Tint(h)=0$.

An immediate consequence of the above equivalence is the fact that the energy of the entire systems can be obtained at the solution as an eigenvalue of the effective Hamiltonian operator $H\xafexteff(h)$,

in the corresponding complete active space. By construction, the cluster amplitudes *T*_{ext}, used to define $H\xafext(h)$, are decoupled from cluster amplitudes *T*_{int} that define the components of the corresponding eigenvector. One should also notice that the many-body expansion of $H\xafext(h)$ may contain effective interactions involving higher-than-pairwise interactions.

Properties of SESs induced eigenvalue problems can also be utilized to define alternative ways of forming CC approximations and corresponding CC equations. For example, the CCSD equations can be recast in the form shown in Fig. 3. The form of the decomposition shown in Fig. 3 can also be viewed as an “echo” of the fact that CCSD theory is an exact theory for subsystem decomposed into noninteracting two-electron systems (an example is shown in Fig. 4). It is also interesting to notice that for the exact CC theory for closed-shell systems discussed here, there exist a chain of various types of SESs that meet requirements (1) and (2). For instance, for the exact formulation, one can consider a chain of SESs defined as

which results in separate eigenvalue problems corresponding to subsystems of various sizes. This observation can be used to define a new CC approximation where instead of referring to adding higher and higher ranks of excitations as a design principle (used to define standard approximations such as CCD, CCSD, CCSDT, and CCSDTQ), one can envision a strategy based on the inclusion excitations in the cluster operator that belong to a specific class (or classes) of SESs. For example, in Ref. 46, we discussed an approximation (subalgebra flow CCSD(2) approximation [SAF-CCSD(2)] that employs all amplitudes contained in all $g(N)(2R)$ SESs, which leads to a CC model that contains all singly and doubly excited cluster amplitudes and selected subsets of triply and quadruply excited ones and where the CC equations can be represented as a conglomerate of coupled eigenvalue problems shown in Fig. 5 (for more detailed discussion, see Ref. 46). It can also be shown that selecting cluster amplitudes based on subsystem embedding algebras is a natural way of introducing a notion of seniority number intensively studied in the context of CC applications with strongly correlated systems.^{80–83}

The character of the expansion (12) is identical with the one discussed in the original formulation of the active-space coupled cluster formalism^{84,85} (see also Refs. 86 and 87), which also utilizes the decomposition of the cluster operator into internal and external parts. The active-space CC methods were broadly utilized as an efficient tool for the selection of the most important high-rank excitations in the ground-state CC formulations^{84,85} and their excited-state equation-of-motion extension.^{88,89} In this context, the SES-CC approach may be viewed as an extension of the active-space CC ideas. Among new elements of theoretical design that differ SES-CC and active-space CC methods, one should mention that: (1) the SES-CC formalism is formulated in terms of the effective Hamiltonian formalism and defines conditions that need to be met to furnish the energy of the entire system as an eigenvalue of the effective Hamiltonian, (2) the effective Hamiltonian provides a rigorous decoupling between internal and external excitations, (3) amplitude selection mechanism in the SES-CC formulations is different from the one in the active-space CC methods and is based on subalgebra flow-equations, and (4) in the subalgebra flow equations, one uses multiple active space characteristic for various standard CC formulations. Since SES-CC is an effective Hamiltonian theory, one should also stress that within this formulation, even standard CC formulations can be viewed as a renormalization technique, where the dimensionality of the Hamiltonian describing the system is replaced by reduced-size Hamiltonian in the active space and is reminiscent of the renormalization-group techniques for interacting fermions.^{90}

To summarize this section, techniques based on the utilization of the subsystem embedding subalgebras can be used to downfold the full electronic Hamiltonian to an arbitrary active space corresponding to some subsystem embedding subalgebra. For example, using SES, one could downfold electronic Hamiltonian to the active space (usually containing highest occupied and lowest unoccupied orbitals) that contains the most important contributions to the electronic wave function of interest. This fact may be especially appealing in the context of quantum algorithms such as quantum phase estimation (QPE)^{39,41,91–95} or variational quantum eigensolver (VQE).^{44,45,96–100} Unfortunately, the main caveat is related to the fact that effective Hamiltonians $H\xafext(h)$ are non-Hermitian, which is a major obstacle in using these algorithms. The main goal of Sec. IV is to redefine unitary CC formalism to provide efficient downfolding algorithm, which yields Hermitian effective Hamiltonian(s) and at the same time assures the same properties as single-reference CC effective Hamiltonians $H\xafext(h)$ discussed in this section.

## IV. UNITARY CC METHOD

The unitary CC (UCC) formulations have been introduced in Refs. 101–109 and have been intensively studied in the recent years in the context of various quantum algorithms.^{43,44,96,98,110}

The standard UCC ansatz (the generalization of standard UCC ansatz has recently been discussed in the literature; for details, see Ref. 111)

is reminiscent of standard single-reference expansion (1) with the difference that the *σ* operator is anti-Hermitian, i.e.,

This property of the *σ* operator leads to a nonterminating character of many-body expansion for the wave function (26). The anti-Hermitian character of *σ* can be assured by the specific form of the *σ* operator, which in most formulations is represented as

where *T* has exactly the same structure as discussed in Eq. (2). In this and following Sec. V, we will focus on exact (i.e., including all possible excitations) formulation of the UCC theory.

In the analysis of the UCC formalisms, we will refer to two standard formulas for operator exponentials:

The Baker-Campbell-Hausdorff (BCH) formula:

where commutators $CBCH(k)$ are defined as $CBCH(1)(X,Y)=12[X,Y]$, $CBCH(2)(X,Y)=112([X,[X,Y]]+[Y,[Y,X]])$, etc.

The transposed variant of the Zassenhaus formula:

^{112,113}

where $CZ(2)=12[X,Y]$, $CZ(3)=13[Y,[[X,Y]]+16[X,[X,Y]]$, etc. The *R*_{Z}(*X*, *Y*) function is defined as

where its inverse is given by the formula

The equations for *σ* operator in the exact limit can be obtained by substituting ansatz (26) into the Schrödinger equation,

The equivalent representation, which provides the explicitly connected form of the equations, can be obtained by multiplying both sides by *e*^{−σ} operator and decoupling equations for amplitudes from the equation for energy (see Refs. 105 and 107 for details), i.e.,

In analogy to the standard single-reference approach, let us explore if partitioning of *σ* into part belonging to some SES $h$ ($\sigma int\u2261\sigma int(h)$) and its complement ($\sigma ext\u2261\sigma ext(h)$) leads to a downfolding of Hamiltonian in the form discussed in Sec. III. For this purpose, we will (1) apply Zassenhaus formula (30) to factorize $e\sigma int+\sigma ext$, (2) premultiply (33) from the left by $e\u2212\sigma extRZ\u22121(\sigma int,\sigma ext)$, and (3) project resulting equations onto *P* + *Q*_{int} space. This procedure leads to the following equations:

For simplicity, the above equation can be rewritten as

where the $H\xafextUCC$ operator is defined as

or equivalently

where the effective Hamiltonian $H\xafexteff(UCC)$ in Eq. (39) is defined as

In the above equation, we use the fact that

It can be shown that the form of the active-space eigenvalue-type problem (39) is equivalent to the (*P* + *Q*_{int}) projections of the connected form of the UCC equations described by Eq. (34). This can be shown by introducing the resolution of identity $e\u2212\sigma inte\sigma int$ to the left of $H\xafextUCC$ in Eq. (37) and noticing that $e\u2212\sigma intH\xafextUCCe\sigma int=e\u2212\sigma He\sigma $. Then, the resulting equation takes the form

Using the matrix representation of the *σ*_{int} operator in the CAS space denoted as *σ*_{int}, this equation can be rewritten as

where the first component of vector ** x** corresponds to ⟨Φ|

*e*

^{−σ}

*He*

^{σ}|Φ⟩ −

*E*, while the remaining components correspond to projections of

*e*

^{−σ}

*He*

^{σ}|Φ⟩ onto excited (with respect to the reference determinant |Φ⟩) configurations $\u2009\Phi \Delta (CAS)|$ belonging to $h$-induced CAS. The $[e\sigma int]$ matrix is also nonsingular, which is a consequence of the algebraic structure of

*σ*

_{int}rather than a particular values of cluster amplitudes defining

*T*

_{int}and $Tint\u2020$. To see it, let us calculate its determinant

where *T*_{int} and $Tint\u2020$ are matrix representations of *T*_{int} and $Tint\u2020$ in CAS. Since *T*_{int} and $Tint\u2020$ contain either excitations or de-excitations, we have

which means that $det(e\sigma int)=1$ and $(e\sigma int)$ is a nonsingular operator. Therefore, the only solution of Eq. (43) corresponds to ** x** = 0, which proves the equivalence of Eq. (37) with

*P*and

*Q*

_{int}projections of Eqs. (34) and (35).

Although $H\xafextUCC$ ($H\xafexteff(UCC)$) is Hermitian, in contrast to $H\xafext$ from Eq. (18), the $H\xafextUCC$ ($H\xafexteff$) operator does not decouple *σ*_{ext} amplitudes from the *σ*_{int} ones. If $h$ is chosen to contain highest/lowest lying occupied/virtual orbitals, we can view this property of $H\xafextUCC$ as mixing low- and high-energy components of the wave function (or using quantum chemical lingua, as a mixing of static and dynamic correlation effects). In Sec. V, we will discuss how to re-instate the separation of effective Hamiltonian while maintaining its Hermitian character. This will have important consequences on how the approximate formulations can be constructed.

## V. ALTERNATIVE UCC EXPANSIONS

In this section, we will discuss properties of a UCC ansatz given by the product of two unitary transformations [double unitary CC (DUCC) expansion or, for the reasons explained later, tailored unitary CC formulation] which explicitly employs the partitioning induced by some SES $h$,

The above expansion is driven by similar ideas as a class of CC methods that utilize double similarity transformations and have been extensively discussed in the literature.^{11,114–117} Moreover, ansatz (46) can also be viewed as a unitary generalization of the active-space CC methods^{85} and tailored CC formalism (TCC),^{118–121} where equations for *σ*_{int} are represented in the form of eigenvalue problem discussed in Sec. IV.

In analogy to Sec. IV, we will focus on the exact formulation of the DUCC approach that includes all possible excitations in the *σ*_{ext} and *σ*_{int} operators. Since using BCH expansion (29), the DUCC expansion can be transformed to the alternative single-reference ansatz,

where *D* is anti-Hermitian (*D*^{†} = −*D*), the DUCC formalism can also be viewed as a special case of a unitary CC ansatz.

In analogy to the UCC formulation, when the double UCC ansatz (46) is introduced into Schrödinger equation,

it can be rewritten in the equivalent form which decouples equations for cluster amplitudes from the equation for energy,

One can show that Eq. (49) corresponding to projections onto (*P* + *Q*_{int}) subspace can be written in equivalent form as an eigenvalue problem,

or equivalently, using effective Hamiltonian language,

where

and

To show this fact, it suffices to introduce the resolution of identity $e\sigma inte\u2212\sigma int$ to the left of the $H\xafextDUCC$ operator in Eq. (51) and notice that $e\u2212\sigma intH\xafextDUCCe\sigma int=e\u2212\sigma inte\u2212\sigma extHe\sigma inte\sigma ext$. Next, in analogy to Eqs. (42) and (43), Eq. (53) can be represented as

where the first component of the [** y**] vector is equivalent to $\u2009\Phi |e\u2212\sigma inte\u2212\sigma extHe\sigma inte\sigma ext|\Phi \u2009\u2212E$, while the remaining components correspond to projections of $e\u2212\sigma inte\u2212\sigma extHe\sigma inte\sigma ext|\Phi \u2009$ onto excited configurations belonging to

*Q*

_{int}. Given the nonsingular character of the $[e\sigma int]$ matrix, this proves the equivalence of these two representations.

By construction, the DUCC effective Hamiltonian is Hermitian and in contrast to the UCC case, it is expressible in terms of the external *σ*_{ext} amplitudes only, providing in this way in Eq. (52) a rigorous decoupling of degrees of freedom corresponding to *σ*_{ext} and *σ*_{int} in the sense of discussion of Sec. III. Since amplitudes defining *σ*_{ext} are characterized by larger perturbative denominator compared to *σ*_{int}, where small denominators may occur, it is much safer to determine *σ*_{ext} using perturbative techniques. Equation (52) also offers a possibility of downfolding the Hamiltonian to the (*P* + *Q*_{int}) space where the correlation effects from *Q*_{ext} can be included through the *σ*_{ext} operator, which makes calculations with the downfolded Hamiltonian amenable for quantum computing even for larger systems.

In the following part of the paper, we will discuss an approximate form of the second quantized representation of the $H\xafexteff(DUCC)$ operator (denoted as the Γ operator) for a specific choice of SES containing all occupied and lowest-lying virtual spin orbitals to define the SES active space (i.e., all occupied indices *i*, *j*, … and some small subset of virtual spin orbitals *a*, *b*, … are deemed active). In this case, amplitudes defining the *σ*_{ext} operator must carry at least one inactive virtual orbital. One can view this procedure as a downfolding of an essential part of the virtual spin-orbital space. This process will consist of two steps: (1) expansion of $H\xafexteff(DUCC)$ in powers of *σ*_{ext} operator and (2) approximation of *σ*_{ext} operator. The Γ operator is defined by strings of creation/annihilation operators that carry only active spin orbitals and can be written as

where the subscript act designates terms of a given operator expression that contains creation/annihilation operators carrying active-space spin-orbital labels. Using the Baker-Campbell-Hausdorff formula, one can recast the above equation in the form of infinite expansion,

Using the particle-hole (ph) formalism, one can also expand the Γ operator into the sum of its many-body components Γ_{i},

where *N*_{e,act} designates the number of active electrons and

where *N*[…]^{ph} designates the particle-hole normal product form of a given operator expression, *P*, *Q*, *R*, *S*, … designate general active spin-orbital indices, and Γ_{scalar} designates full contracted (scalar) part of the Γ operator with respect to the particle-hole vacuum |Φ⟩. In the above expansion, we also assume that all multidimensional tensors $\gamma Q1\u2026QiP1\u2026Pi$ are antisymmetric with respect to the permutations among the sets of lower and upper spin orbital indices.

In practical realizations of the DUCC formalism, one has to truncate both the many-body Γ expansion given by Eq. (61) and the excitation level included in the *σ*_{ext} operator. Below, as a specific example, we describe a variant of Γ-operator approximations based on the inclusion of one- and two-body interactions/excitations. This will also illustrate the benefits of using a hybrid particle-hole and physical vacuum second-quantized representations of all operators involved in the approximation. Following similar ideas as discussed in Ref. 103, where energy functionals were constructed based on the order of the energy contributions, one can select terms in expansion (57) based on the perturbative analysis of contributing terms. For example, using particle-hole formalism, the normal product form of the electronic Hamiltonian *H*_{N} (*H*_{N} = *H* − ⟨Φ|*H*|Φ⟩) can be split into its one-particle part *F*_{N} and two-particle component *V*_{N} and one can retain elements in (57) that are correct through the second order, i.e.,

An important aspect related to the approximate form of the Γ operator is related to the excitation order included in the *σ*_{ext} operator, which can include single, double, triple, etc., many-body components. Since *σ*_{ext} is mainly responsible for dynamic correlation effects, various iterative and noniterative approximations can be used to evaluate these terms. In this paper, we will consider an approximation where *σ*_{ext} is represented by external components of singly and doubly excited cluster operators (*T*_{ext,1} and *T*_{ext,2}), i.e.,

For the simplest approximation of *T*_{ext,1} and *T*_{ext,2} operators, one can employ external parts of the CCSD *T*_{1} and *T*_{2} operators. More sophisticated approximations may involve external cluster amplitudes corresponding to single, double, triple, quadruple, etc., excitations obtained from genuine UCC models.

Since we are interested in making DUCC formalism amenable for classical/quantum computing, in the remaining part of the paper, we will derive the algebraic form for the matrix elements defining Γ operator in the second-quantized form in the particle-hole and physical vacuum representations. The interest in latter representation is mainly caused by the fact that most of the classical (DMRG)/quantum diagonalizers utilize physical vacuum representation of the second-quantized forms of electronic Hamiltonian. In our analysis, we will focus our attention on one- and two-body interactions. To derive these formulas, we use a combined approach:

### A. Determination of the algebraic form of $\Gamma $ Eq. (62) using particle-hole formalism

Applying the particle-hole variant of Wick’s theorem to the operator expressions [defining expansion (62) for *T*_{ext} and $Text\u2020$ given by Eq. (63)],

and retaining terms only through two-body interactions, one obtains

where again *P*, *Q*, *R*, *S* designate general active spin-orbital indices (in the forthcoming analysis, we will designate occupied and virtual active orbitals by *I*, *J*, *K*, … and *A*, *B*, *C*, …, respectively), and (…)_{C,open} designates connected and open (i.e., having external lines in diagrammatic representation) parts of a given operator expression. In Eq. (65), (Γ)_{scalar} corresponds to the scalar part of Γ in particle-hole representation, i.e., Γ_{scalar} = ⟨ Φ|Γ|Φ⟩. One should also notice that first and second, third and sixth, and fourth and fifth operators in expression (64) are pairs of Hermitian conjugate operators [for example, $(HNText)C,open=(Text\u2020HN)C,open\u2020$]. The utilization of the particle-hole formalism helps in keeping track and including a broader class of correlation effects compared to the physical vacuum representation. The explicit expressions defining $\gamma QP$ and $\gamma RSPQ$ amplitudes are given in Tables I–III. For these tables, we assume that a restricted or unrestricted Hartree-Fock (RHF/UHF) reference is employed (i.e., all nondiagonal elements of the Fock matrix are zero), and as a result, the third and sixth terms in Eq. (64) vanish.

Amplitude . | Expression . | Amplitude . | Expression . | Amplitude . | Expression . | Amplitude . | Expression . |
---|---|---|---|---|---|---|---|

H_{N} term | |||||||

$\gamma AB=fAB$ | $\gamma IJ=fIJ$ | $\gamma AI=fAI$ | $\gamma IA=fIA$ | ||||

$\gamma IABC=vIABC$ | $\gamma IJKA=vIJKA$ | $\gamma ABCI=vABCI$ | $\gamma KAIJ=vKAIJ$ | ||||

$\gamma IAJB=vIAJB$ | $\gamma ABCD=vABCD$ | $\gamma IJKL=vIJKL$ | $\gamma ABIJ=vABIJ$ | ||||

$\gamma IJAB=vIJAB$ |

Amplitude . | Expression . | Amplitude . | Expression . | Amplitude . | Expression . | Amplitude . | Expression . |
---|---|---|---|---|---|---|---|

H_{N} term | |||||||

$\gamma AB=fAB$ | $\gamma IJ=fIJ$ | $\gamma AI=fAI$ | $\gamma IA=fIA$ | ||||

$\gamma IABC=vIABC$ | $\gamma IJKA=vIJKA$ | $\gamma ABCI=vABCI$ | $\gamma KAIJ=vKAIJ$ | ||||

$\gamma IAJB=vIAJB$ | $\gamma ABCD=vABCD$ | $\gamma IJKL=vIJKL$ | $\gamma ABIJ=vABIJ$ | ||||

$\gamma IJAB=vIJAB$ |

Amplitude . | Expression . |
---|---|

$(HNText)C,open$ term | |

$\gamma AB+=\u2212fAMsMB+veAMBsMe\u221212veAMNsMNeB$ | |

$\gamma IJ+=feJsIe+veIMJsMe+12vefMJsMIef$ | |

$\gamma AI+=veAMIsMe$ | |

$\gamma IA+=feAsIe\u2212fIMsMA+veIMAsMe+feMtMIeA\u221212veIMNsMNeA+12vefMAsMIef$ | |

$\gamma IABC+=veABCsIe\u2212vAIMBsMC+vAIMCsMB+fAMsMIBC\u2212veAMBsMIeC+veAMCsMIeB+12vIAMNsMNBC$ | |

$\gamma IJKA+=veJKAsIe\u2212veIKAsJe+vIJMKsMA+feKsIJeA+12vefKAsIJef\u2212veJMKsMIeA+veIMKsMJeA$ | |

$\gamma ABCI+=\u2212vABMIsMC$ | |

$\gamma KAIJ+=veAIJsKe$ | |

$\gamma IAJB+=veAJBsIe+vIAMJsMB\u2212veAMJsMIeB$ | |

$\gamma ABCD+=vABMCsMD\u2212vABMDsMC+12vABMNsMNCD$ | |

$\gamma IJKL+=veJKLsIe\u2212veIKLsJe+12vefKLsIJef$ | |

$\gamma IJAB+=veJABsIe\u2212veIABsJe+vIJMAsMB\u2212vIJMBsMA+feAsIJeB\u2212feBsIJeA+fJMsMIAB\u2212fIMsMJAB$ | |

$+12vefABsIJef+12vIJMNsMNAB\u2212veJMAsMIeB+veIMAsMJeB+veJMBsMIeA\u2212veIMBsMJeA$ | |

$(Text\u2020HN)C,open$ term | |

$\gamma BA+=\u2212fMAsBM+vMBeAseM\u221212vMNeAseBMN$ | |

$\gamma JI+=fJeseI+vMJeIsMe+12vMJefsefMI$ | |

$\gamma IA+=vMIeAsMe$ | |

$\gamma AI+=fAeseI\u2212fMIsAM+vMAeIseM+fMeseAMI\u221212vMNeIseAMN+12vMAefsefMI$ | |

$\gamma BCIA+=vBCeAseI\u2212vMBAIsCM+vMCAIsBM+fMAsBCMI\u2212vMBeAseCMI+vMCeAseBMI+12vMNIAsBCMN$ | |

$\gamma KAIJ+=vKAeJseI\u2212vKAeIseJ+vMKIJsAM+fKeseAIJ+12vKAefsefIJ\u2212vMKeJseAMI+vMKeIseAMJ$ | |

$\gamma CIAB+=\u2212vMIABsCM$ | |

$\gamma IJKA+=vIJeAseK$ | |

$\gamma JBIA+=vJBeAseI+vMJIAsBM\u2212vMJeAseBMI$ | |

$\gamma CDAB+=vMCABsDM\u2212vMDABsCM+12vMNABsCDMN$ | |

$\gamma KLIJ+=vKLeJseI\u2212vKLeIseJ+12vKLefsefIJ$ | |

$\gamma ABIJ+=vABeJseI\u2212vABeIseJ+vMAIJsBM\u2212vMBIJsAM+fAeseBIJ\u2212fBeseAIJ+fMJsABMI\u2212fMIsABMJ$ | |

$+12vABefsefIJ+12vMNIJsABMN\u2212vMAeJseBMI+vMAeIseBMJ+vMBeJseAMI\u2212vMBeIseAMJ$ |

Amplitude . | Expression . |
---|---|

$(HNText)C,open$ term | |

$\gamma AB+=\u2212fAMsMB+veAMBsMe\u221212veAMNsMNeB$ | |

$\gamma IJ+=feJsIe+veIMJsMe+12vefMJsMIef$ | |

$\gamma AI+=veAMIsMe$ | |

$\gamma IA+=feAsIe\u2212fIMsMA+veIMAsMe+feMtMIeA\u221212veIMNsMNeA+12vefMAsMIef$ | |

$\gamma IABC+=veABCsIe\u2212vAIMBsMC+vAIMCsMB+fAMsMIBC\u2212veAMBsMIeC+veAMCsMIeB+12vIAMNsMNBC$ | |

$\gamma IJKA+=veJKAsIe\u2212veIKAsJe+vIJMKsMA+feKsIJeA+12vefKAsIJef\u2212veJMKsMIeA+veIMKsMJeA$ | |

$\gamma ABCI+=\u2212vABMIsMC$ | |

$\gamma KAIJ+=veAIJsKe$ | |

$\gamma IAJB+=veAJBsIe+vIAMJsMB\u2212veAMJsMIeB$ | |

$\gamma ABCD+=vABMCsMD\u2212vABMDsMC+12vABMNsMNCD$ | |

$\gamma IJKL+=veJKLsIe\u2212veIKLsJe+12vefKLsIJef$ | |

$\gamma IJAB+=veJABsIe\u2212veIABsJe+vIJMAsMB\u2212vIJMBsMA+feAsIJeB\u2212feBsIJeA+fJMsMIAB\u2212fIMsMJAB$ | |

$+12vefABsIJef+12vIJMNsMNAB\u2212veJMAsMIeB+veIMAsMJeB+veJMBsMIeA\u2212veIMBsMJeA$ | |

$(Text\u2020HN)C,open$ term | |

$\gamma BA+=\u2212fMAsBM+vMBeAseM\u221212vMNeAseBMN$ | |

$\gamma JI+=fJeseI+vMJeIsMe+12vMJefsefMI$ | |

$\gamma IA+=vMIeAsMe$ | |

$\gamma AI+=fAeseI\u2212fMIsAM+vMAeIseM+fMeseAMI\u221212vMNeIseAMN+12vMAefsefMI$ | |

$\gamma BCIA+=vBCeAseI\u2212vMBAIsCM+vMCAIsBM+fMAsBCMI\u2212vMBeAseCMI+vMCeAseBMI+12vMNIAsBCMN$ | |

$\gamma KAIJ+=vKAeJseI\u2212vKAeIseJ+vMKIJsAM+fKeseAIJ+12vKAefsefIJ\u2212vMKeJseAMI+vMKeIseAMJ$ | |

$\gamma CIAB+=\u2212vMIABsCM$ | |

$\gamma IJKA+=vIJeAseK$ | |

$\gamma JBIA+=vJBeAseI+vMJIAsBM\u2212vMJeAseBMI$ | |

$\gamma CDAB+=vMCABsDM\u2212vMDABsCM+12vMNABsCDMN$ | |

$\gamma KLIJ+=vKLeJseI\u2212vKLeIseJ+12vKLefsefIJ$ | |

$\gamma ABIJ+=vABeJseI\u2212vABeIseJ+vMAIJsBM\u2212vMBIJsAM+fAeseBIJ\u2212fBeseAIJ+fMJsABMI\u2212fMIsABMJ$ | |

$+12vABefsefIJ+12vMNIJsABMN\u2212vMAeJseBMI+vMAeIseBMJ+vMBeJseAMI\u2212vMBeIseAMJ$ |

Amplitude . | Expression . |
---|---|

$12(Text\u2020(FNText)C,open)C,open$ term | |

$\gamma IJ+=12seJffes\u2009If\u221212seMfMJsIe\u221214sefMJfINsMNef+12segMJffes\u2009MIfg+14sefIMfMNsNJef$ | |

$\gamma BA+=12sBMfMNsNA\u221212seMfBesMA+14sfBMNfeAsMNef+14seBMNffesMNAf\u221212seBMKfKNsNMeA$ | |

$\gamma IA+=12seMffAsMIef\u221212seMfINsMNeA\u221212seMfMNsINAe+12seMffesIMAf$ | |

$\gamma KLIJ+=14sefIJfLMsMKef\u221214sefIJfKMsMLef+12segIJffes\u2009KLfg$ | |

$\gamma JBIA+=12seBMIfJNsMNeA\u221212seBMIffAsMJef\u221212seBMIffes\u2009MJfA+12seBMIfMNsNJeA$ | |

$\gamma CDAB+=14sCDMNfeAsMNeB\u221214sCDMNfeBsMNeA\u221212sCDMKfMNsNKAB$ | |

$\gamma KAIJ+=\u221212seAIJfKMsMe$ | |

$\gamma ABCI+=\u221212sABMIfeCsMe$ | |

$\gamma IJKA+=12seKffes\u2009IJfA\u221212seKfJMsIMeA+12seKfIMsJMeA+12seKffAsIJef$ | |

$\gamma CIAB+=12sCMfMNsNIAB\u221212sCMfeBsIMeA+12sCMfeAsIMeB+12sCMfINsMNAB$ | |

$12((Text\u2020FN)C,openText)C,open$ term | |

$\gamma JI+=12sJefefsfI\u221212sMe\u2009fJMseI\u221214sMJeffNIsefMN+12sMJegfefsfgMI+14sIMeffNMsefNJ$ | |

$\gamma AB+=12sMBfNMsAN\u221212sMe\u2009feBsAM+14sMNfBfAesefMN+14sMNeBfefsAfMN\u221212sMKeBfNKseANM$ | |

$\gamma AI+=12sMefAfsefMI\u221212sMefNIseAMN\u221212sMefNMsAeIN+12sMefefsAfIM$ | |

$\gamma IJKL+=14sIJeffMLsefMK\u221214sIJeffMKsefML+12sIJegfefsfgKL$ | |

$\gamma IAJB+=12sMIeB\u2009fNJseAMN\u221212sMIeBfAfsefMJ\u221212sMIeBfefsfAMJ+12sMIeBfNMseANJ$ | |

$\gamma ABCD+=14sMNCDfAeseBMN\u221214sMNCDfBeseAMN\u221212sMKCDfNMsABNK$ | |

$\gamma IJKA+=\u221212sIJeAfMKseM$ | |

$\gamma CIAB+=\u221212sMIAB\u2009fCeseM$ | |

$\gamma KAIJ+=12sKefefsfAIJ\u221212sKefMJseAIM+12sKefMIseAJM+12sKefAfsefIJ$ | |

$\gamma ABCI+=12sMCfNMsABNI\u221212sMCfBeseAIM+12sMCfAeseBIM+12sMCfNIsABMN$ |

Amplitude . | Expression . |
---|---|

$12(Text\u2020(FNText)C,open)C,open$ term | |

$\gamma IJ+=12seJffes\u2009If\u221212seMfMJsIe\u221214sefMJfINsMNef+12segMJffes\u2009MIfg+14sefIMfMNsNJef$ | |

$\gamma BA+=12sBMfMNsNA\u221212seMfBesMA+14sfBMNfeAsMNef+14seBMNffesMNAf\u221212seBMKfKNsNMeA$ | |

$\gamma IA+=12seMffAsMIef\u221212seMfINsMNeA\u221212seMfMNsINAe+12seMffesIMAf$ | |

$\gamma KLIJ+=14sefIJfLMsMKef\u221214sefIJfKMsMLef+12segIJffes\u2009KLfg$ | |

$\gamma JBIA+=12seBMIfJNsMNeA\u221212seBMIffAsMJef\u221212seBMIffes\u2009MJfA+12seBMIfMNsNJeA$ | |

$\gamma CDAB+=14sCDMNfeAsMNeB\u221214sCDMNfeBsMNeA\u221212sCDMKfMNsNKAB$ | |

$\gamma KAIJ+=\u221212seAIJfKMsMe$ | |

$\gamma ABCI+=\u221212sABMIfeCsMe$ | |

$\gamma IJKA+=12seKffes\u2009IJfA\u221212seKfJMsIMeA+12seKfIMsJMeA+12seKffAsIJef$ | |

$\gamma CIAB+=12sCMfMNsNIAB\u221212sCMfeBsIMeA+12sCMfeAsIMeB+12sCMfINsMNAB$ | |

$12((Text\u2020FN)C,openText)C,open$ term | |

$\gamma JI+=12sJefefsfI\u221212sMe\u2009fJMseI\u221214sMJeffNIsefMN+12sMJegfefsfgMI+14sIMeffNMsefNJ$ | |

$\gamma AB+=12sMBfNMsAN\u221212sMe\u2009feBsAM+14sMNfBfAesefMN+14sMNeBfefsAfMN\u221212sMKeBfNKseANM$ | |

$\gamma AI+=12sMefAfsefMI\u221212sMefNIseAMN\u221212sMefNMsAeIN+12sMefefsAfIM$ | |

$\gamma IJKL+=14sIJeffMLsefMK\u221214sIJeffMKsefML+12sIJegfefsfgKL$ | |

$\gamma IAJB+=12sMIeB\u2009fNJseAMN\u221212sMIeBfAfsefMJ\u221212sMIeBfefsfAMJ+12sMIeBfNMseANJ$ | |

$\gamma ABCD+=14sMNCDfAeseBMN\u221214sMNCDfBeseAMN\u221212sMKCDfNMsABNK$ | |

$\gamma IJKA+=\u221212sIJeAfMKseM$ | |

$\gamma CIAB+=\u221212sMIAB\u2009fCeseM$ | |

$\gamma KAIJ+=12sKefefsfAIJ\u221212sKefMJseAIM+12sKefMIseAJM+12sKefAfsefIJ$ | |

$\gamma ABCI+=12sMCfNMsABNI\u221212sMCfBeseAIM+12sMCfAeseBIM+12sMCfNIsABMN$ |

### B. Determination of the physical-vacuum representation of the $\Gamma $ operator

In order to find an equivalent characterization of the Γ operator given by Eq. (65) using the physical vacuum parameterization, we will employ the set of identities from Table IV that translate particle-hole normal product forms for typical strings of creation/annihilation operators to physical vacuum expressions where all creation operators are placed to the left of the annihilation operators. Using these identities, the physical vacuum Hamiltonian Γ, in one- and two-body interaction approximation, takes the form

where all $\chi QP$ and $\chi RSPQ$ coefficients are listed in Table V. These matrix elements can be implemented and used as an input for full configuration interaction type diagonalizers (in this case, limited to the diagonalization in the corresponding active space) including various “full” CC approaches (CC approaches involving all possible excitations within the active space), density matrix renormalization group, and quantum simulators (employing either QPE or VQE).

Normal product form: . | Normal product form: . |
---|---|

particle-hole formalism . | physical vacuum formalism . |

$N[aB\u2020aA]ph=aB\u2020aA$ | |

$N[aJ\u2020aI]ph=aJ\u2020aI\u2212\delta IJ$ | |

$N[aI\u2020aA]ph=aI\u2020aA$ | |

$N[aA\u2020aI]ph=aA\u2020aI$ | |

$N[aB\u2020aC\u2020aAaI]ph=\u2212aB\u2020aC\u2020aIaA$ | |

$N[aK\u2020aA\u2020aJaI]ph=\u2212aA\u2020aK\u2020aJaI\u2212\delta IKaA\u2020aJ+\delta JKaA\u2020aI$ | |

$N[aC\u2020aI\u2020aBaA]ph=aC\u2020aI\u2020aBaA$ | |

$N[aI\u2020aJ\u2020aAaK]ph=\u2212aI\u2020aJ\u2020aKaA\u2212\delta IKaJ\u2020aA+\delta JKaI\u2020aA$ | |

$N[aJ\u2020aB\u2020aAaI]ph=aB\u2020aJ\u2020aIaA\u2212\delta IJaB\u2020aA$ | |

$N[aC\u2020aD\u2020aBaA]ph=aC\u2020aD\u2020aBaA$ | |

$N[aK\u2020aL\u2020aJaI]ph=aK\u2020aL\u2020aJaI\u2212\delta JLaK\u2020aI+\delta ILaK\u2020aJ+\delta JKaL\u2020aI\u2212\delta IKaL\u2020aJ+\delta JL\delta IK\u2212\delta JK\delta IL$ | |

$N[aI\u2020aJ\u2020aBaA]ph=aI\u2020aJ\u2020aBaA$ | |

$N[aA\u2020aB\u2020aJaI]ph=aA\u2020aB\u2020aJaI$ |

Normal product form: . | Normal product form: . |
---|---|

particle-hole formalism . | physical vacuum formalism . |

$N[aB\u2020aA]ph=aB\u2020aA$ | |

$N[aJ\u2020aI]ph=aJ\u2020aI\u2212\delta IJ$ | |

$N[aI\u2020aA]ph=aI\u2020aA$ | |

$N[aA\u2020aI]ph=aA\u2020aI$ | |

$N[aB\u2020aC\u2020aAaI]ph=\u2212aB\u2020aC\u2020aIaA$ | |

$N[aK\u2020aA\u2020aJaI]ph=\u2212aA\u2020aK\u2020aJaI\u2212\delta IKaA\u2020aJ+\delta JKaA\u2020aI$ | |

$N[aC\u2020aI\u2020aBaA]ph=aC\u2020aI\u2020aBaA$ | |

$N[aI\u2020aJ\u2020aAaK]ph=\u2212aI\u2020aJ\u2020aKaA\u2212\delta IKaJ\u2020aA+\delta JKaI\u2020aA$ | |

$N[aJ\u2020aB\u2020aAaI]ph=aB\u2020aJ\u2020aIaA\u2212\delta IJaB\u2020aA$ | |

$N[aC\u2020aD\u2020aBaA]ph=aC\u2020aD\u2020aBaA$ | |

$N[aK\u2020aL\u2020aJaI]ph=aK\u2020aL\u2020aJaI\u2212\delta JLaK\u2020aI+\delta ILaK\u2020aJ+\delta JKaL\u2020aI\u2212\delta IKaL\u2020aJ+\delta JL\delta IK\u2212\delta JK\delta IL$ | |

$N[aI\u2020aJ\u2020aBaA]ph=aI\u2020aJ\u2020aBaA$ | |

$N[aA\u2020aB\u2020aJaI]ph=aA\u2020aB\u2020aJaI$ |

Amplitude . | Expression . | Amplitude . | Expression . |
---|---|---|---|

$\chi AB=\gamma AB\u2212\u2211M\gamma MAMB$ | $\chi IJ=\gamma IJ\u2212\u2211M\gamma MIMJ$ | ||

$\chi AI=\gamma AI\u2212\u2211M\gamma MAMI$ | $\chi IA=\gamma IA\u2212\u2211M\gamma MIMA$ | ||

$\chi IABC=\gamma IABC$ | $\chi IJAK=\gamma IJAK$ | ||

$\chi ABCI=\gamma ABCI$ | $\chi AKIJ=\gamma AKIJ$ | ||

$\chi AIBJ=\gamma AIBJ$ | $\chi ABCD=\gamma ABCD$ | ||

$\chi IJKL=\gamma IJKL$ | $\chi ABIJ=\gamma ABIJ$ | ||

$\chi IJAB=\gamma IJAB$ |

Amplitude . | Expression . | Amplitude . | Expression . |
---|---|---|---|

$\chi AB=\gamma AB\u2212\u2211M\gamma MAMB$ | $\chi IJ=\gamma IJ\u2212\u2211M\gamma MIMJ$ | ||

$\chi AI=\gamma AI\u2212\u2211M\gamma MAMI$ | $\chi IA=\gamma IA\u2212\u2211M\gamma MIMA$ | ||

$\chi IABC=\gamma IABC$ | $\chi IJAK=\gamma IJAK$ | ||

$\chi ABCI=\gamma ABCI$ | $\chi AKIJ=\gamma AKIJ$ | ||

$\chi AIBJ=\gamma AIBJ$ | $\chi ABCD=\gamma ABCD$ | ||

$\chi IJKL=\gamma IJKL$ | $\chi ABIJ=\gamma ABIJ$ | ||

$\chi IJAB=\gamma IJAB$ |

## VI. DISCUSSION

The accuracy of truncations behind the DUCC(2) formalism is contingent upon several factors:

the size of the active space (or equivalently the number of active virtual orbitals included in $h$),

the accuracy of

*σ*_{ext}estimates represented by single and double excitations,the role of missing higher-rank many-body effects in the

*σ*_{ext}and Γ operators.

As stated earlier, the size of the active space may affect the accuracy of *σ*_{ext} amplitudes, which is a consequence of the fact that utilization of larger spaces, whose choice is driven by the value of orbital energies, prevents perturbative denominators from being near singular. For situations where the energy separation between virtual active and inactive spin orbitals is sufficiently large, one should also expect that the role of higher-rank excitations in *σ*_{ext} is proportionally smaller. Otherwise, one needs to include higher many-body components of *σ*_{ext} and Γ operators, for example, three- and/or four-body components, i.e., *σ*_{ext,3}, *σ*_{ext,4}, …, and Γ_{3}, Γ_{4}, …. In such cases, one should also expect that standard single reference CC formulations including *T*_{3} and/or *T*_{4} cluster may not be a viable source of the information about exact *σ*_{ext,3} and/or *σ*_{ext,4} operators. Instead, one should resort to using genuine UCC formulations. For example, the *T*_{ext,3} amplitudes can be extracted from the UCC(4) model discussed in Ref. 103, where the sufficiency conditions for cluster amplitudes,

where *Q*_{1}, *Q*_{2}, and *Q*_{3} are projection operators onto singly, doubly, and triply excited configurations, allow to generate *T*_{3} amplitudes in on-the-fly manner, avoiding in this way typical memory bottlenecks associated with storing the whole set of *T*_{3} amplitudes.

The discussed DUCC formalism also offers a possibility of integrating classical and quantum computations, where the CC/UCC calculations for *σ*_{ext} and forming $\chi QP$, $\chi RSPQ$, … amplitudes are performed on classical computers, while the diagonalization step takes advantage of quantum computing resources. For this reason, it is instructive to discuss the quantum resources as a function of the number of active orbitals (*N*_{act}), total number of spin orbitals (*N*_{S}), and rank of many-body effects included in the Γ operator expression (66).

The specific choice of the active-space or equivalently subsystem embedding algebra $h$ defines how efficient is the process of integrating out remaining degrees of freedom (i.e., the parameters/amplitudes defining the *T*_{ext}/*σ*_{ext} operator). In particular, the proper choice of $h$ will impact the accuracy of low-cost (perturbative) estimates of *T*_{ext}. Here, we will consider two cases: (1) choice of $h$ based on the energy threshold and (2) choice of $h$ based on locality criteria (or equivalently subsystem separability discussed in Sec. III).

In the first case, the active space is chosen in analogy to typical applications of multireference methods [complete active space self-consistent field (CASSCF), CASPT2,^{122,123} NEVPTn,^{124,125} multi-reference many-body perturbation theory (MRMBPT),^{126–133} and DMRG^{134,135} methods] where active spaces usually contain high- and low-lying occupied and virtual orbitals. In this situation, the first order contribution to *T*_{ext} can be, for example, written as

where in this specific example, *s*-amplitude contains one active occupied (J) and one active virtual (B) spin-orbital indices; we will also assume that in the above example, *i* and *a* represent inactive occupied and virtual spin-orbital indices. If active-space orbitals are well separated (energetically) from the remaining orbitals, one can expect that the perturbative denominators used to define *T*_{ext} are much larger than those corresponding to excitations within active space and which are determined in the diagonalization procedure. In this case, the use of perturbative techniques should provide reliable *T*_{ext} estimates.

In contrast to the energy separation criteria, when “spatial” arguments are invoked to define the active space, the “smallness” of *σ*_{ext}-amplitudes is determined by the decay law of the orbitals implicated in a specific excitation. Here, we will consider two specific situations shown in Fig. 6: (a) the “active” or strongly correlated subsystem is weakly interacting with subsystem B (which can be described by low-order contributions of the many-body perturbation theory) and (b) two “active” strongly correlated centers (*A*_{1} ad *A*_{2}) are embedded in a weakly correlated medium (for example, solution). In case (a), the matrix element in Eq. (70),

should be “small” since spin orbitals *i* and *a* vs *J* and *B* are spatially well separated [see Fig. 6(a)]. In the second case, one should utilize the joint active spaces (SES $h$) defined by active orbitals defining subsystems *A*_{1} [with corresponding SES $h1$] and *A*_{2} [with corresponding SES $h2$], i.e.,

where the downfolded Hamiltonian $H\xafextDUCC(h)$ is given by the expression

Once subsystems *A*_{1} and *A*_{2} are spatially separated and localized basis set is employed and $\sigma int(h1)$ and $\sigma int(h2)$ commute

then further downfolding of $H\xafextDUCC(h)$ is possible, i.e.,

The above considerations indicate that any system “separability” parameter can be used to define appropriate model space.

## VII. RESOURCE ESTIMATES FOR SIMULATION ON A QUANTUM COMPUTER

The discussed DUCC formalism also offers a possibility of integrating classical and quantum computations, where the CC/UCC calculations for *σ*_{ext} and forming $\chi QP$, $\chi RSPQ$, … amplitudes are performed on classical computers, while the diagonalization step takes advantage of quantum computing resources. For this reason, it is instructive to discuss the quantum resources as a function of the number of active orbitals (*N*_{act}) and total number of spin-orbitals (*N*_{S}) and rank of many-body effects included in the Γ operator. As an example of the improvement afforded by our downfolding approach, let us evaluate the worst-case resources required by the quantum phase estimation algorithm for obtaining an eigenvalue of the original Hamiltonian (involving only one- and two-body interactions) in full space of all *N*_{S} orbitals vs DUCC models involving *N*_{act} active spin orbitals for two situations: Γ operator contains (1) up to two-body interactions (Γ_{1} and Γ_{2}) and (2) up to three-body interactions (Γ_{1}, Γ_{2}, and Γ_{3}).

Given an arbitrary unitary *U* with eigenvalues $U|\psi j\u2009=ei\theta j|\psi j\u2009$ and an input state |*ψ*⟩ = *∑*_{j}*α*_{j}|*ψ*_{j}⟩, quantum phase estimation returns an eigenphase estimate of *θ*_{j} with error *δ* for a randomly chosen eigenvector |*ψ*_{j}⟩ that is sampled with probability |*α*_{j}|^{2}. More precisely, the estimate $\theta ^$ is drawn from the distribution

where *f*(*x*) is a function that depends on the choice of phase estimation algorithm and is sharply peaked about *x* = 0 with width *δ*. In common variants of phase estimation algorithms,^{136–138} controlled-*U* must be applied $O(1/\delta )$ times to obtain a single $\theta ^$ estimate and is the dominant cost.

One common example of *U* is the real-time evolution operator *e*^{−iHt}, in which case, *θ*_{j} is an eigenvalue of $H=\u2211jO(Nterms)hjPj$, for $O(Nterms)$ Pauli operators *P*_{j} acting on at most $O(NS)$ qubits, and positive coefficients *h*_{j}, scaled by a constant *t*. Although real-time evolution may be approximated using Trotter-Suzuki product formulas,^{93} it is difficult to obtain tight error bounds on scaling of the approximation error with *N*_{act}. Thus, we consider the case where *U* is a quantum walk^{94,139,140} with eigenvalues commensurate with those of $e\xb1i\u2061sin\u22121(H/\lambda )$, for *λ* = *∑*_{j}|*h*_{j}| ≥ ∥*H*∥. Unlike Trotter-Suzuki formulas, there is no approximation error in the eigenvalues apart from errors in the Hamiltonian representation. This walk-based approach has become popular of late under the name qubitization.^{39,139–142} Note that an eigenphase estimate $\theta ^$ of the quantum walk and its error *δ* may be related to the eigenvalue estimate of *H* by computing $\lambda \u2061sin(\theta ^)$. Thus, if we wish to learn the eigenvalue within error *ϵ*, we need to apply the walk operator $O(\lambda /\u03f5)$ times.

In qubitization, the walk operator is of the form *U* = (PREPARE^{†}⊗ I) · SELECT · (PREPARE ⊗ I) · ((1 − 2|0⋯0⟩⟨0⋯0|) ⊗ *I*) for unitary subroutines PREPARE and SELECT. The PREPARE subroutine prepares from the all zero state |0⋯0⟩, a state of the form $\u2211jhj/\lambda |\u2009j\u2009$ within error *ϵ*, using $O(Nterms)$ quantum gates. Here and in the following, we count the number of arbitrary single- and two-qubit quantum gates. The SELECT subroutine applies each *P*_{j} in the Hamiltonian selected by a register that stores the value of *j*. In a naive implementation, this requires $O(NSNterms)$ gates, but optimizations specific to the Jordan-Wigner representation of fermionic Hamiltonians reduce this to $O(NS)$,^{141} which is subdominant to the cost of the PREPARE circuit.

The overall cost of the algorithm is found by multiplying the cost of the walk operator by the number of iterations needed within phase estimation. Thus, the gate complexity of obtaining an estimate of an eigenvalue of *H* to error *δ* is in $O\u0303(\lambda Nterms/\delta )$. In the worst-case, we may assume that $Nterms=O(NS4)$ and $\lambda =O(|H|maxNS4)$, where |*H*|_{max} is the maximum absolute value of any entry in the Hamiltonian.

While |*H*|_{max} is difficult to estimate in general, there are cases when its scaling can be asymptotically estimated. If we assume that a local orbital basis is used consisting of atomic orbitals centered at each of the nuclei with charge *Z*_{i}, then it is straightforward to show that $|hpq|\u2208O(maxiZi2)$ and $|hpqrs|\u2208O(maxiZi)$.^{143} Therefore, we anticipate that $\lambda \u2208O(NS2maxiZi2+NS4maxiZi)$ for such problems. If we consider the nuclear charges to be fixed, we then expect the worst case scaling of a simulation of the total Hamiltonian within error *δ* to be in $O(NS8/\delta )$ in such cases.

If we consider simulating the downfolded Hamiltonian, one simply replaces $Nterms=O(Nact4)$ with the number of terms in the downfolded effective Hamiltonian, and similarly for the normalization constant $\lambda =O(|\Gamma |maxNact4)$. Importantly, the number of active space orbitals is much smaller than that of the full Hamiltonian, that is, *N*_{act} ≪ *N*_{S}. In common practice, one normally chooses active spaces such that amplitude corrections of downfolding are small. In other words, |Γ|_{max} is expected to be similar to that of the maximum absolute value of *H* restricted to the active space. If core electrons are also moved out of the active space, one may expect |Γ|_{max} ≪ |*H*|_{max}. The case of simulating three-body interactions is less-studied, but the essential idea is identical. With three-body terms, the worst-case $Nterms=\lambda =O(Nact6)$ although *λ* could be much smaller if |Γ_{3}|_{max} ≪ |Γ|_{max}, which would be expected of a small correction.

In practice, it will be essential to understand costs using realistic examples rather than building intuition with the worst-case analysis. In particular, realistic cases might be more accurately captured with low-rank approximations that reduce *N*_{terms} to as small as $O\u0303(Nact2)$^{43} for the two-body case—low-rank approximations for the three-body case are still not well understood. Moreover, various quantum circuit optimization techniques^{144} enable a further $O(Nact)$ reduction in non-Clifford gate complexity, which are the dominant expensive gates in fault-tolerant quantum computation. More advanced quantum simulation algorithms for structured Hamiltonians could also be applicable, such as those that negate the cost of simulating diagonal terms or exploit large separations in energy scales.^{145,146}

## VIII. NUMERICAL TESTS

We performed a series of numerical tests with the QDK simulator for two systems: H_{2} and Be. While the main goal of the H_{2} system was to test the feasibility of the DUCC techniques in studying ground-state potential energy surfaces (PESs), the Be system was employed to study the impact of the basis set size employed on the efficiency of the DUCC approach. For H_{2}, we used cc-pVTZ basis set,^{147} whereas for Be, we used a series of cc-pVDZ, cc-pVTZ, and cc-pVQZ basis sets.^{147} All tensor contractions from Tables I–V were implemented and integrated in NWChem^{148} with the Tensor Contraction Engine module. In all calculations, all electrons were correlated. To generate the format amenable for the QDK processing, we used recently published NWChem-QDK interface.^{51} Also, the discussion of the quality of the DUCC energies is mainly based on the averaged energy values reported (along with the standard deviations) in Tables VI and VII.

Method . | 0.8 a.u. . | 1.4008 a.u. . | 4.00 a.u. . | 10.00 a.u. . |
---|---|---|---|---|

FCI (active space) | −0.9830 | −1.1467 | −1.0070 | −0.9971 |

DUCC (active space)^{a} | −1.0086 | −1.1679 | −1.0152 | −1.0016 |

±0.0172 | ±0.0167 | ±0.0221 | ±0.0177 | |

FCI (full space) | −1.0157 | −1.1725 | −1.0149 | −0.9996 |

Method . | 0.8 a.u. . | 1.4008 a.u. . | 4.00 a.u. . | 10.00 a.u. . |
---|---|---|---|---|

FCI (active space) | −0.9830 | −1.1467 | −1.0070 | −0.9971 |

DUCC (active space)^{a} | −1.0086 | −1.1679 | −1.0152 | −1.0016 |

±0.0172 | ±0.0167 | ±0.0221 | ±0.0177 | |

FCI (full space) | −1.0157 | −1.1725 | −1.0149 | −0.9996 |

^{a}

The DUCC energies were obtained with the quantum phase estimation algorithm in combination with the Trotter-Suzuki simulation method. The trotter step size of 0.1 and bit precision of 9 were used in the calculations. Results exclude values for excited states.

Method . | 5 orbitals . | 6 orbitals . | 9 orbitals . | All orbitals . |
---|---|---|---|---|

cc-pVDZ | ||||

Full CI | −14.5952 | −14.5968 | −14.6169 | −14.6174 |

DUCC^{a} | −14.5984 ± 0.0089 | −14.5993 ± 0.0081 | −14.6156 ± 0.0083 | |

cc-pVTZ | ||||

Full CI | −14.5889 | −14.5902 | −14.6168 | −14.6238 |

DUCC^{a} | −14.5971 ± 0.0082 | −14.5994 ± 0.0079 | −14.6188 ± 0.0078 | |

cc-pVQZ | ||||

Full CI | −14.5852 | −14.5858 | −14.6133 | −14.6401 |

DUCC^{a} | −14.6096 ± 0.0079 | −14.6126 ± 0.0091 | −14.6328 ± 0.0063 |

Method . | 5 orbitals . | 6 orbitals . | 9 orbitals . | All orbitals . |
---|---|---|---|---|

cc-pVDZ | ||||

Full CI | −14.5952 | −14.5968 | −14.6169 | −14.6174 |

DUCC^{a} | −14.5984 ± 0.0089 | −14.5993 ± 0.0081 | −14.6156 ± 0.0083 | |

cc-pVTZ | ||||

Full CI | −14.5889 | −14.5902 | −14.6168 | −14.6238 |

DUCC^{a} | −14.5971 ± 0.0082 | −14.5994 ± 0.0079 | −14.6188 ± 0.0078 | |

cc-pVQZ | ||||

Full CI | −14.5852 | −14.5858 | −14.6133 | −14.6401 |

DUCC^{a} | −14.6096 ± 0.0079 | −14.6126 ± 0.0091 | −14.6328 ± 0.0063 |

^{a}

The DUCC energies were obtained with the quantum phase estimation algorithm in combination with the Trotter-Suzuki simulation method. The trotter step size of 0.1 and bit precision of 10 were used in the calculations. Results exclude values for excited states and outliers based on a recursive Grubb’s test (*α* = 0.05).

In Fig. 7 and Table VI, we show energies obtained by diagonalization of (1) bare Hamiltonian in the active space spanned by 1 occupied and 3 lowest-lying virtual orbitals, (2) DUCC effective Hamiltonian in the same active space as in point (1), and (3) bare Hamiltonian in the full space. One can observe that DUCC energies are significantly closer to the exact FCI ones than the ones obtained in FCI calculations for small active space. For example, for 1.4008 bohr, the DUCC procedure reduces 26 mhartree of error of the FCI procedure in the active space to 4.6 mhartree. Although for larger H–H separations (e.g., for 10 bohr), the energy discrepancy between diagonalization in the active space and full space is smaller and amounts to 2.5 mhartree, it is worth noticing that the DUCC approach further reduces this error to 2.0 mhartree. It is important to notice that the DUCC formalism provides good approximations to the exact FCI energies for various internuclear distances (see Fig. 7 and Table VI) irrespective of balance between static and dynamic correlation effects (the former class of correlation effects can be evaluated by discrepancies between FCI results in active space and exact ones).

The beryllium atom is often employed as a benchmark system for various coupled cluster formulations including single- and multireference coupled clusters. In this paper, we study the correlated effects in using DUCC formalism for basis sets and active spaces of various sizes. To evaluate the accuracy of the DUCC formalism as a function of active space size, we considered three cases of active spaces defined by 5, 6, and 9 lowest-lying restricted Hartree-Fock orbitals. Also, cc-pVDZ, cc-pVTZ, and cc-pVQZ basis sets contributing to 14, 30, and 55 orbitals were considered in the present studies. We collated results of our simulations in Table VII. For the cc-pVDZ basis set, one can observe that the DUCC energies are of the same quality as corresponding active-space FCI results. For 9 active orbitals, both active-space FCI and DUCC methods reach the quality of the exact FCI energy. For smaller active spaces, DUCC provides only minor improvements of the corresponding active-space FCI energies. For the cc-pVTZ basis set for small 5- and 6-orbital active spaces, the DUCC approach provides significant 8 and 9 mhartree corrections to the active-space FCI results toward the exact FCI energy for the cc-pVTZ basis set. The effectiveness and capacity of the DUCC formalism are particularly displayed with the cc-pVQZ basis set. The standard active-space FCI calculations give errors of about 55, 54, and 27 mhartree for 5-, 6-, and 9-orbital active spaces. On the other hand, the DUCC approach only gives errors of about 31, 28, and 7 mhartree, a reduction of 24, 27, and 20 mhartree for 5-, 6-, and 9-orbital active spaces, respectively, when compared to the corresponding active-space FCI results. These results suggest that dynamical (or out-of-active-space) correlation effects can be “downfolded” more effectively for larger basis sets, where the correlation effects are distributed more evenly over the virtual orbitals. We consider the DUCC for the Be system as satisfactory given the source of the *T* amplitudes, the inclusion of single commutators, and limiting many-body effects in DUCC Hamiltonians to single and pairwise interactions.

## IX. CONCLUSIONS

We have shown that the SES-CC methodology can be extended to unitary CC formalisms, which provides a procedure for downfolding many-body effects into a Hermitian effective Hamiltonian, in contrast to the earlier canonical SES-CC work in which the effective Hamiltonian was not Hermitian. We introduced the DUCC model which decouples the classes of excitations used to define the effective Hamiltonian from those obtained in the corresponding eigenvalue problem. These techniques may provide a convenient way of decoupling two types of degrees of freedom corresponding to parameters defining low- and high-energy components defining the electronic wave function of interest. Using the computational chemistry nomenclature, these two subsets can be identified with static and dynamical correlation effects. However, one can also envision slightly different scenarios for the DUCC formalism application where different types of effects (scales)—for example, short- vs long-range correlations effects—are decoupled using an appropriate form of the local DUCC ansatz or equivalently adequate definition of the corresponding active space. This development, which provides a rigorous scheme for obtaining Hermitian effective Hamiltonians, opens doors for obtaining the ground-state energy with previously unobtainable diagonalization techniques, such as those used in quantum algorithms or the DMRG method. This effort also provides a rigorous generalization of the genuine active-space CC and SES-CC ideas to a Hermitian formulation of the effective Hamiltonian approach, which is amenable to quantum computing. The pilot DUCC calculations with QDK quantum simulator for H_{2} and Be systems using larger cc-pVTZ (H_{2} and Be) and cc-pVQZ (Be) basis sets, respectively, can capture a significant portion of the out-of-active-space correlation effects even when a basic form of approximation [Eq. (62)] of the DUCC Hamiltonian and non-DUCC external amplitudes are used. For the H_{2} systems, we showed that the DUCC procedure significantly improves the quality of the total energies obtained by the diagonalization of bare Hamiltonians in the active space for several internuclear geometries corresponding to equilibrium and stretched bond lengths. Moreover, for both H_{2} and Be systems, we demonstrated that acceptable estimates of the FCI energies (by correlating all orbitals) can be obtained with 6–8 fold reduction in the size of orbital space, which only increases as one goes to larger basis sets. The last result holds for large-size basis sets. More tests are needed to evaluate: (1) the performance of the DUCC formalism for systems characterized by various types of correlation effects (especially concerning the balance between static and dynamical correlation effects), (2) the effects of higher-order commutators in expansion (57), (3) the many-body structure and rank of excitations in the effective Hamiltonian (58), and (4) the choice of approximate *T*_{ext} amplitudes.

Integrating-out the corresponding fermionic degrees of freedom, which leads to low-dimensionality second-quantized effective Hamiltonians, will open up the possibility of performing quantum simulations on existing quantum simulators, as well as on larger molecular systems. These problems will be tested in the forthcoming papers. An interesting development area is also associated with recent advances in compressing the second quantized form of the electronic Hamiltonian based on the composite Cholesky-Singular Value Decompositions (SVDs) of one- and two-electron integrals.^{43,52} One can envisage a further extension of the applicability of Cholesky-SVD decomposition to compress *χ*-amplitudes defining DUCC downfolded Hamiltonians.

From the classical computing viewpoint, DUCC downfolding techniques have broader implications and can be used in the context of the density renormalized group approach (DMRG)^{26–28,149} by providing a dressed form of the effective Hamiltonian. This approach may complement existing perturbative techniques used in the context of the DMRG theory to account for the dynamical correlation effects.^{150–152}

## ACKNOWLEDGMENTS

This work was supported by the “Embedding Quantum Computing into Many-body Frameworks for Strongly Correlated Molecular and Materials Systems” project, which is funded by the U.S. Department of Energy (DOE), Office of Science, Office of Basic Energy Sciences, the Division of Chemical Sciences, Geosciences, and Biosciences. A portion of this research was funded by the Quantum Algorithms, Software, and Architectures (QUASAR) Initiative at the Pacific Northwest National Laboratory (PNNL). It was conducted under the Laboratory Directed Research and Development Program at PNNL. Pacific Northwest National Laboratory is operated by Battelle for DOE under Contract DE-AC05-76RL01830.