Downfolding of many-body Hamiltonians using active-space models: extension of the sub-system embedding sub-algebras approach to unitary coupled cluster formalisms

TL;DR

This paper extends the SES-CC formalism to unitary CC methods, enabling efficient downfolding of many-electron Hamiltonians into low-energy models with Hermitian effective Hamiltonians suitable for quantum algorithms.

Contribution

It introduces a unitary CC extension of the SES-CC formalism, providing a Hermitian effective Hamiltonian for active space models and separating dynamical from internal correlations.

Findings

Hermitian effective Hamiltonians suitable for FCI solvers

Separation of external and internal cluster amplitudes

Potential for hybrid classical-quantum algorithms

Abstract

In this paper we outline the extension of recently introduced the sub-system embedding sub-algebras coupled cluster (SES-CC) formalism to the unitary CC formalism. 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 formalis (DUCC) 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 complete active space 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 (FCI) 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.

Figures

Tables

Table 1. Table 1: The algebraic form of the γ Q P subscript superscript 𝛾 𝑃 𝑄 \gamma^{P}_{Q} and γ R ​ S P ​ Q subscript superscript 𝛾 𝑃 𝑄 𝑅 𝑆 \gamma^{PQ}_{RS} amplitudes in Eq.( 65 ) stemming from the H N subscript 𝐻 𝑁 H_{N} term.

Amplitude	Expression	Amplitude	Expression	Amplitude	Expression	Amplitude	Expression
$H_N$ term
${\gamma }_A^B=$	$f_A^B$	${\gamma }_I^J=$	$f_I^J$	${\gamma }_A^I=$	$f_A^I$	${\gamma }_I^A=$	$f_I^A$
${\gamma }_{IA}^{BC}=$	$v_{IA}^{BC}$	${\gamma }_{IJ}^{KA}=$	$v_{IJ}^{KA}$	${\gamma }_{AB}^{CI}=$	$v_{AB}^{CI}$	${\gamma }_{KA}^{IJ}=$	$v_{KA}^{IJ}$
${\gamma }_{IA}^{JB}=$	$v_{IA}^{JB}$	${\gamma }_{AB}^{CD}=$	$v_{AB}^{CD}$	${\gamma }_{IJ}^{KL}=$	$v_{IJ}^{KL}$	${\gamma }_{AB}^{IJ}=$	$v_{AB}^{IJ}$
${\gamma }_{IJ}^{AB}=$	$v_{IJ}^{AB}$

Table 2. Table 2: Algebraic form of γ Q P subscript superscript 𝛾 𝑃 𝑄 \gamma^{P}_{Q} and γ R ​ S P ​ Q subscript superscript 𝛾 𝑃 𝑄 𝑅 𝑆 \gamma^{PQ}_{RS} amplitudes in Eq.( 65 ) (continued). Amplitudes defining T ext , 1 subscript 𝑇 ext 1 T_{{\rm ext},1} and T ext , 2 subscript 𝑇 ext 2 T_{{\rm ext},2} operators are denoted by s i a subscript superscript 𝑠 𝑎 𝑖 s^{a}_{i} and s i ​ j a ​ b subscript superscript 𝑠 𝑎 𝑏 𝑖 𝑗 s^{ab}_{ij} ( i , j , … 𝑖 𝑗 … i,j,\ldots and a , b , … 𝑎 𝑏 … a,b,\ldots are generic occupied and virtual spin-orbitals, respectively) while amplitudes defining T ext , 1 † superscript subscript 𝑇 ext 1 † T_{{\rm ext},1}^{\dagger} and T ext , 2 † superscript subscript 𝑇 ext 2 † T_{{\rm ext},2}^{\dagger} are denoted as s a i superscript subscript 𝑠 𝑎 𝑖 s_{a}^{i} and s a ​ b i ​ j superscript subscript 𝑠 𝑎 𝑏 𝑖 𝑗 s_{ab}^{ij} , respectively). By definition of the external parts of T 𝑇 T and T † superscript 𝑇 † T^{\dagger} , all s 𝑠 s -amplitudes that carry active spin-orbital indices only disappear. These terms pertain to active spaces that contain all correlated occupied orbitals and a subset of the virtual ones. For simplicity we assume a restricted/unrestricted Hartree-Fock (RHF/UHF) reference | Φ ⟩ ket Φ |\Phi\rangle , where all non-diagonal Fock matrix elements disappear. The Einstein summation convention is invoked.

Amplitude	Expression
${\left(H_N{T}_{\mathrm{ext}}\right)}_{C,\mathrm{open}}$ term
${\gamma }_A^B+=$	$-f_A^M{s}_M^B+v_{eA}^{MB}{s}_M^e-\frac{1}{2}{v}_{eA}^{MN}{s}_{MN}^{eB}$
${\gamma }_I^J+=$	$f_e^J{s}_I^e+v_{eI}^{MJ}{s}_M^e+\frac{1}{2}{v}_{ef}^{MJ}{s}_{MI}^{ef}$
${\gamma }_A^I+=$	$v_{eA}^{MI}{s}_M^e$
${\gamma }_I^A+=$	$f_e^A{s}_I^e-f_I^M{s}_M^A+v_{eI}^{MA}{s}_M^e+f_e^M{t}_{MI}^{eA}-\frac{1}{2}{v}_{eI}^{MN}{s}_{MN}^{eA}+\frac{1}{2}{v}_{ef}^{MA}{s}_{MI}^{ef}$
${\gamma }_{IA}^{BC}+=$	$v_{eA}^{BC}{s}_I^e-v_{AI}^{MB}{s}_M^C+v_{AI}^{MC}{s}_M^B+f_A^M{s}_{MI}^{BC}-v_{eA}^{MB}{s}_{MI}^{eC}+v_{eA}^{MC}{s}_{MI}^{eB}+\frac{1}{2}{v}_{IA}^{MN}{s}_{MN}^{BC}$
${\gamma }_{IJ}^{KA}+=$	$v_{eJ}^{KA}{s}_I^e-v_{eI}^{KA}{s}_J^e+v_{IJ}^{MK}{s}_M^A+f_e^K{s}_{IJ}^{eA}+\frac{1}{2}{v}_{ef}^{KA}{s}_{IJ}^{ef}-v_{eJ}^{MK}{s}_{MI}^{eA}+v_{eI}^{MK}{s}_{MJ}^{eA}$
${\gamma }_{AB}^{CI}+=$	$-v_{AB}^{MI}{s}_M^C$
${\gamma }_{KA}^{IJ}+=$	$v_{eA}^{IJ}{s}_K^e$
${\gamma }_{IA}^{JB}+=$	$v_{eA}^{JB}{s}_I^e+v_{IA}^{MJ}{s}_M^B-v_{eA}^{MJ}{s}_{MI}^{eB}$
${\gamma }_{AB}^{CD}+=$	$v_{AB}^{MC}{s}_M^D-v_{AB}^{MD}{s}_M^C+\frac{1}{2}{v}_{AB}^{MN}{s}_{MN}^{CD}$
${\gamma }_{IJ}^{KL}+=$	$v_{eJ}^{KL}{s}_I^e-v_{eI}^{KL}{s}_J^e+\frac{1}{2}{v}_{ef}^{KL}{s}_{IJ}^{ef}$
${\gamma }_{IJ}^{AB}+=$	$v_{eJ}^{AB}{s}_I^e-v_{eI}^{AB}{s}_J^e+v_{IJ}^{MA}{s}_M^B-v_{IJ}^{MB}{s}_M^A+f_e^A{s}_{IJ}^{eB}-f_e^B{s}_{IJ}^{eA}+f_J^M{s}_{MI}^{AB}-f_I^M{s}_{MJ}^{AB}$
	$+\frac{1}{2}{v}_{ef}^{AB}{s}_{IJ}^{ef}+\frac{1}{2}{v}_{IJ}^{MN}{s}_{MN}^{AB}-v_{eJ}^{MA}{s}_{MI}^{eB}+v_{eI}^{MA}{s}_{MJ}^{eB}+v_{eJ}^{MB}{s}_{MI}^{eA}-v_{eI}^{MB}{s}_{MJ}^{eA}$
${\left(T_{\mathrm{ext}}^{†}{H}_N\right)}_{C,\mathrm{open}}$ term
${\gamma }_B^A+=$	$-f_M^A{s}_B^M+v_{MB}^{eA}{s}_e^M-\frac{1}{2}{v}_{MN}^{eA}{s}_{eB}^{MN}$
${\gamma }_J^I+=$	$f_J^e{s}_e^I+v_{MJ}^{eI}{s}_M^e+\frac{1}{2}{v}_{MJ}^{ef}{s}_{ef}^{MI}$
${\gamma }_I^A+=$	$v_{MI}^{eA}{s}_M^e$
${\gamma }_A^I+=$	$f_A^e{s}_e^I-f_M^I{s}_A^M+v_{MA}^{eI}{s}_M^e+f_M^e{s}_{eA}^{MI}-\frac{1}{2}{v}_{MN}^{eI}{s}_{eA}^{MN}+\frac{1}{2}{v}_{MA}^{ef}{s}_{ef}^{MI}$
${\gamma }_{BC}^{IA}+=$	$v_{BC}^{eA}{s}_e^I-v_{MB}^{AI}{s}_C^M+v_{MC}^{AI}{s}_B^M+f_M^A{s}_{BC}^{MI}-v_{MB}^{eA}{s}_{eC}^{MI}+v_{MC}^{eA}{s}_{eB}^{MI}+\frac{1}{2}{v}_{MN}^{IA}{s}_{BC}^{MN}$
${\gamma }_{KA}^{IJ}+=$	$v_{KA}^{eJ}{s}_e^I-v_{KA}^{eI}{s}_e^J+v_{MK}^{IJ}{s}_A^M+f_K^e{s}_{eA}^{IJ}+\frac{1}{2}{v}_{KA}^{ef}{s}_{ef}^{IJ}-v_{MK}^{eJ}{s}_{eA}^{MI}+v_{MK}^{eI}{s}_{eA}^{MJ}$
${\gamma }_{CI}^{AB}+=$	$-v_{MI}^{AB}{s}_C^M$
${\gamma }_{IJ}^{KA}+=$	$v_{IJ}^{eA}{s}_e^K$
${\gamma }_{JB}^{IA}+=$	$v_{JB}^{eA}{s}_e^I+v_{MJ}^{IA}{s}_B^M-v_{MJ}^{eA}{s}_{eB}^{MI}$
${\gamma }_{CD}^{AB}+=$	$v_{MC}^{AB}{s}_D^M-v_{MD}^{AB}{s}_C^M+\frac{1}{2}{v}_{MN}^{AB}{s}_{CD}^{MN}$
${\gamma }_{KL}^{IJ}+=$	$v_{KL}^{eJ}{s}_e^I-v_{KL}^{eI}{s}_e^J+\frac{1}{2}{v}_{KL}^{ef}{s}_{ef}^{IJ}$
${\gamma }_{AB}^{IJ}+=$	$v_{AB}^{eJ}{s}_e^I-v_{AB}^{eI}{s}_e^J+v_{MA}^{IJ}{s}_B^M-v_{MB}^{IJ}{s}_A^M+f_A^e{s}_{eB}^{IJ}-f_B^e{s}_{eA}^{IJ}+f_M^J{s}_{AB}^{MI}-f_M^I{s}_{AB}^{MJ}$
	$+\frac{1}{2}{v}_{AB}^{ef}{s}_{ef}^{IJ}+\frac{1}{2}{v}_{MN}^{IJ}{s}_{AB}^{MN}-v_{MA}^{eJ}{s}_{eB}^{MI}+v_{MA}^{eI}{s}_{eB}^{MJ}+v_{MB}^{eJ}{s}_{eA}^{MI}-v_{MB}^{eI}{s}_{eA}^{MJ}$

Table 3. Table 3: Algebraic form of γ Q P subscript superscript 𝛾 𝑃 𝑄 \gamma^{P}_{Q} and γ R ​ S P ​ Q subscript superscript 𝛾 𝑃 𝑄 𝑅 𝑆 \gamma^{PQ}_{RS} amplitudes in Eq.( 65 ) (continued). Amplitudes defining T ext , 1 subscript 𝑇 ext 1 T_{{\rm ext},1} and T ext , 2 subscript 𝑇 ext 2 T_{{\rm ext},2} operators are denoted by s i a subscript superscript 𝑠 𝑎 𝑖 s^{a}_{i} and s i ​ j a ​ b subscript superscript 𝑠 𝑎 𝑏 𝑖 𝑗 s^{ab}_{ij} ( i , j , … 𝑖 𝑗 … i,j,\ldots and a , b , … 𝑎 𝑏 … a,b,\ldots are generic occupied and virtual spin-orbitals, respectively) while amplitudes defining T ext , 1 † superscript subscript 𝑇 ext 1 † T_{{\rm ext},1}^{\dagger} and T ext , 2 † superscript subscript 𝑇 ext 2 † T_{{\rm ext},2}^{\dagger} are denoted as s a i superscript subscript 𝑠 𝑎 𝑖 s_{a}^{i} and s a ​ b i ​ j superscript subscript 𝑠 𝑎 𝑏 𝑖 𝑗 s_{ab}^{ij} , respectively). By definition of the external parts of T 𝑇 T and T † superscript 𝑇 † T^{\dagger} , all s 𝑠 s -amplitudes that carry active spin-orbital indices only disappear. These terms pertain to active spaces that contain all correlated occupied orbitals and a subset of the virtual ones. For simplicity we assume a restricted/unrestricted Hartree-Fock (RHF/UHF) reference | Φ ⟩ ket Φ |\Phi\rangle , where all non-diagonal Fock matrix elements disappear. The Einstein summation convention is invoked.

Amplitude	Expression
$\frac{1}{2}{\left(T_{\mathrm{ext}}^{†}{\left(F_N{T}_{\mathrm{ext}}\right)}_{C,\mathrm{open}}\right)}_{C,\mathrm{open}}$ term
${\gamma }_I^J+=$	$\frac{1}{2}{s}_e^J{f}_f^e{s}_I^f-\frac{1}{2}{s}_e^M{f}_M^J{s}_I^e-\frac{1}{4}{s}_{ef}^{MJ}{f}_I^N{s}_{MN}^{ef}+\frac{1}{2}{s}_{eg}^{MJ}{f}_f^e{s}_{MI}^{fg}+\frac{1}{4}{s}_{ef}^{IM}{f}_M^N{s}_{NJ}^{ef}$
${\gamma }_B^A+=$	$\frac{1}{2}{s}_B^M{f}_M^N{s}_N^A-\frac{1}{2}{s}_e^M{f}_B^e{s}_M^A+\frac{1}{4}{s}_{fB}^{MN}{f}_e^A{s}_{MN}^{ef}+\frac{1}{4}{s}_{eB}^{MN}{f}_f^e{s}_{MN}^{Af}-\frac{1}{2}{s}_{eB}^{MK}{f}_K^N{s}_{NM}^{eA}$
${\gamma }_I^A+=$	$\frac{1}{2}{s}_e^M{f}_f^A{s}_{MI}^{ef}-\frac{1}{2}{s}_e^M{f}_I^N{s}_{MN}^{eA}-\frac{1}{2}{s}_e^M{f}_M^N{s}_{IN}^{Ae}+\frac{1}{2}{s}_e^M{f}_f^e{s}_{IM}^{Af}$
${\gamma }_{KL}^{IJ}+=$	$\frac{1}{4}{s}_{ef}^{IJ}{f}_L^M{s}_{MK}^{ef}-\frac{1}{4}{s}_{ef}^{IJ}{f}_K^M{s}_{ML}^{ef}+\frac{1}{2}{s}_{eg}^{IJ}{f}_f^e{s}_{KL}^{fg}$
${\gamma }_{JB}^{IA}+=$	$\frac{1}{2}{s}_{eB}^{MI}{f}_J^N{s}_{MN}^{eA}-\frac{1}{2}{s}_{eB}^{MI}{f}_f^A{s}_{MJ}^{ef}-\frac{1}{2}{s}_{eB}^{MI}{f}_f^e{s}_{MJ}^{fA}+\frac{1}{2}{s}_{eB}^{MI}{f}_M^N{s}_{NJ}^{eA}$
${\gamma }_{CD}^{AB}+=$	$\frac{1}{4}{s}_{CD}^{MN}{f}_e^A{s}_{MN}^{eB}-\frac{1}{4}{s}_{CD}^{MN}{f}_e^B{s}_{MN}^{eA}-\frac{1}{2}{s}_{CD}^{MK}{f}_M^N{s}_{NK}^{AB}$
${\gamma }_{KA}^{IJ}+=$	$-\frac{1}{2}{s}_{eA}^{IJ}{f}_K^M{s}_M^e$
${\gamma }_{AB}^{CI}+=$	$-\frac{1}{2}{s}_{AB}^{MI}{f}_e^C{s}_M^e$
${\gamma }_{IJ}^{KA}+=$	$\frac{1}{2}{s}_e^K{f}_f^e{s}_{IJ}^{fA}-\frac{1}{2}{s}_e^K{f}_J^M{s}_{IM}^{eA}+\frac{1}{2}{s}_e^K{f}_I^M{s}_{JM}^{eA}+\frac{1}{2}{s}_e^K{f}_f^A{s}_{IJ}^{ef}$
${\gamma }_{CI}^{AB}+=$	$\frac{1}{2}{s}_C^M{f}_M^N{s}_{NI}^{AB}-\frac{1}{2}{s}_C^M{f}_e^B{s}_{IM}^{eA}+\frac{1}{2}{s}_C^M{f}_e^A{s}_{IM}^{eB}+\frac{1}{2}{s}_C^M{f}_I^N{s}_{MN}^{AB}$
$\frac{1}{2}{\left({\left(T_{\mathrm{ext}}^{†}{F}_N\right)}_{C,\mathrm{open}}{T}_{\mathrm{ext}}\right)}_{C,\mathrm{open}}$ term
${\gamma }_J^I+=$	$\frac{1}{2}{s}_J^e{f}_e^f{s}_f^I-\frac{1}{2}{s}_M^e{f}_J^M{s}_e^I-\frac{1}{4}{s}_{MJ}^{ef}{f}_N^I{s}_{ef}^{MN}+\frac{1}{2}{s}_{MJ}^{eg}{f}_e^f{s}_{fg}^{MI}+\frac{1}{4}{s}_{IM}^{ef}{f}_N^M{s}_{ef}^{NJ}$
${\gamma }_A^B+=$	$\frac{1}{2}{s}_M^B{f}_N^M{s}_A^N-\frac{1}{2}{s}_M^e{f}_e^B{s}_A^M+\frac{1}{4}{s}_{MN}^{fB}{f}_A^e{s}_{ef}^{MN}+\frac{1}{4}{s}_{MN}^{eB}{f}_e^f{s}_{Af}^{MN}-\frac{1}{2}{s}_{MK}^{eB}{f}_N^K{s}_{eA}^{NM}$
${\gamma }_A^I+=$	$\frac{1}{2}{s}_M^e{f}_A^f{s}_{ef}^{MI}-\frac{1}{2}{s}_M^e{f}_N^I{s}_{eA}^{MN}-\frac{1}{2}{s}_M^e{f}_N^M{s}_{Ae}^{IN}+\frac{1}{2}{s}_M^e{f}_e^f{s}_{Af}^{IM}$
${\gamma }_{IJ}^{KL}+=$	$\frac{1}{4}{s}_{IJ}^{ef}{f}_M^L{s}_{ef}^{MK}-\frac{1}{4}{s}_{IJ}^{ef}{f}_M^K{s}_{ef}^{ML}+\frac{1}{2}{s}_{IJ}^{eg}{f}_e^f{s}_{fg}^{KL}$
${\gamma }_{IA}^{JB}+=$	$\frac{1}{2}{s}_{MI}^{eB}{f}_N^J{s}_{eA}^{MN}-\frac{1}{2}{s}_{MI}^{eB}{f}_A^f{s}_{ef}^{MJ}-\frac{1}{2}{s}_{MI}^{eB}{f}_e^f{s}_{fA}^{MJ}+\frac{1}{2}{s}_{MI}^{eB}{f}_N^M{s}_{eA}^{NJ}$
${\gamma }_{AB}^{CD}+=$	$\frac{1}{4}{s}_{MN}^{CD}{f}_A^e{s}_{eB}^{MN}-\frac{1}{4}{s}_{MN}^{CD}{f}_B^e{s}_{eA}^{MN}-\frac{1}{2}{s}_{MK}^{CD}{f}_N^M{s}_{AB}^{NK}$
${\gamma }_{IJ}^{KA}+=$	$-\frac{1}{2}{s}_{IJ}^{eA}{f}_M^K{s}_e^M$
${\gamma }_{CI}^{AB}+=$	$-\frac{1}{2}{s}_{MI}^{AB}{f}_C^e{s}_e^M$
${\gamma }_{KA}^{IJ}+=$	$\frac{1}{2}{s}_K^e{f}_e^f{s}_{fA}^{IJ}-\frac{1}{2}{s}_K^e{f}_M^J{s}_{eA}^{IM}+\frac{1}{2}{s}_K^e{f}_M^I{s}_{eA}^{JM}+\frac{1}{2}{s}_K^e{f}_A^f{s}_{ef}^{IJ}$
${\gamma }_{AB}^{CI}+=$	$\frac{1}{2}{s}_M^C{f}_N^M{s}_{AB}^{NI}-\frac{1}{2}{s}_M^C{f}_B^e{s}_{eA}^{IM}+\frac{1}{2}{s}_M^C{f}_A^e{s}_{eB}^{IM}+\frac{1}{2}{s}_M^C{f}_N^I{s}_{AB}^{MN}$

Table 4. Table 4: Translation of particle-hole normal product forms for typical strings of creation/annihilation operators into the physical-vacuum normal product form (all creation operators are to the left with respect to the annihilation operators).

normal product form:		normal product form:
particle-hole formalism		physical vacuum formalism
$N{\left[a_B^{†}{a}_A\right]}^{\mathrm{ph}}$	=	$a_B^{†}{a}_A$
$N{\left[a_J^{†}{a}_I\right]}^{\mathrm{ph}}$	=	$a_J^{†}{a}_I-{\delta }_{IJ}$
$N{\left[a_I^{†}{a}_A\right]}^{\mathrm{ph}}$	=	$a_I^{†}{a}_A$
$N{\left[a_A^{†}{a}_I\right]}^{\mathrm{ph}}$	=	$a_A^{†}{a}_I$
$N{\left[a_B^{†}{a}_C^{†}{a}_A{a}_I\right]}^{\mathrm{ph}}$	=	$-a_B^{†}{a}_C^{†}{a}_I{a}_A$
$N{\left[a_K^{†}{a}_A^{†}{a}_J{a}_I\right]}^{\mathrm{ph}}$	=	$-a_A^{†}{a}_K^{†}{a}_J{a}_I-{\delta }_{IK}{a}_A^{†}{a}_J+{\delta }_{JK}{a}_A^{†}{a}_I$
$N{\left[a_C^{†}{a}_I^{†}{a}_B{a}_A\right]}^{\mathrm{ph}}$	=	$a_C^{†}{a}_I^{†}{a}_B{a}_A$
$N{\left[a_I^{†}{a}_J^{†}{a}_A{a}_K\right]}^{\mathrm{ph}}$	=	$-a_I^{†}{a}_J^{†}{a}_K{a}_A-{\delta }_{IK}{a}_J^{†}{a}_A+{\delta }_{JK}{a}_I^{†}{a}_A$
$N{\left[a_J^{†}{a}_B^{†}{a}_A{a}_I\right]}^{\mathrm{ph}}$	=	$a_B^{†}{a}_J^{†}{a}_I{a}_A-{\delta }_{IJ}{a}_B^{†}{a}_A$
$N{\left[a_C^{†}{a}_D^{†}{a}_B{a}_A\right]}^{\mathrm{ph}}$	=	$a_C^{†}{a}_D^{†}{a}_B{a}_A$
$N{\left[a_K^{†}{a}_L^{†}{a}_J{a}_I\right]}^{\mathrm{ph}}$	=	$a_K^{†}{a}_L^{†}{a}_J{a}_I-{\delta }_{JL}{a}_K^{†}{a}_I+{\delta }_{IL}{a}_K^{†}{a}_J+{\delta }_{JK}{a}_L^{†}{a}_I-{\delta }_{IK}{a}_L^{†}{a}_J+{\delta }_{JL}{\delta }_{IK}-{\delta }_{JK}{\delta }_{IL}$
$N{\left[a_I^{†}{a}_J^{†}{a}_B{a}_A\right]}^{\mathrm{ph}}$	=	$a_I^{†}{a}_J^{†}{a}_B{a}_A$
$N{\left[a_A^{†}{a}_B^{†}{a}_J{a}_I\right]}^{\mathrm{ph}}$	=	$a_A^{†}{a}_B^{†}{a}_J{a}_I$

Table 5. Table 5: The algebraic form of χ Q P subscript superscript 𝜒 𝑃 𝑄 \chi^{P}_{Q} and χ R ​ S P ​ Q subscript superscript 𝜒 𝑃 𝑄 𝑅 𝑆 \chi^{PQ}_{RS} amplitudes as functions of γ Q P subscript superscript 𝛾 𝑃 𝑄 \gamma^{P}_{Q} and γ R ​ S P ​ Q subscript superscript 𝛾 𝑃 𝑄 𝑅 𝑆 \gamma^{PQ}_{RS} ones.

Amplitude	Expression		Amplitude	Expression
${\chi }_A^B=$	${\gamma }_A^B-{\sum }_M{\gamma }_{MA}^{MB}$		${\chi }_I^J=$	${\gamma }_I^J-{\sum }_M{\gamma }_{MI}^{MJ}$
${\chi }_A^I=$	${\gamma }_A^I-{\sum }_M{\gamma }_{MA}^{MI}$		${\chi }_I^A=$	${\gamma }_I^A-{\sum }_M{\gamma }_{MI}^{MA}$
${\chi }_{IA}^{BC}=$	${\gamma }_{IA}^{BC}$		${\chi }_{IJ}^{AK}=$	${\gamma }_{IJ}^{AK}$
${\chi }_{AB}^{CI}=$	${\gamma }_{AB}^{CI}$		${\chi }_{AK}^{IJ}=$	${\gamma }_{AK}^{IJ}$
${\chi }_{AI}^{BJ}=$	${\gamma }_{AI}^{BJ}$		${\chi }_{AB}^{CD}=$	${\gamma }_{AB}^{CD}$
${\chi }_{IJ}^{KL}=$	${\gamma }_{IJ}^{KL}$		${\chi }_{AB}^{IJ}=$	${\gamma }_{AB}^{IJ}$
${\chi }_{IJ}^{AB}=$	${\gamma }_{IJ}^{AB}$

Equations

$$|\Psi ⟩=e^T|\Phi ⟩,$$

$$T=\sum _{k=1}^m{T}_k.$$

$$T_k=\frac{1}{(k!{)}^2}\sum _{\begin{array}{c}{i}_1,\dots ,i_k;\\ a_1,\dots ,a_k\end{array}}{t}_{i_1\dots i_k}^{a_1\dots a_k}{E}_{i_1\dots i_k}^{a_1\dots a_k},$$

$$E_{i_1\dots i_k}^{a_1\dots a_k}=a_{a_1}^{†}\dots a_{a_k}^{†}{a}_{i_k}\dots a_{i_1}.$$

$$|{\Phi }_{i_1\dots i_k}^{a_1\dots a_k}⟩=E_{i_1\dots i_k}^{a_1\dots a_k}|\Phi ⟩.$$

$$(P+Q)H{e}^T|\Phi ⟩=E(P+Q)e^T|\Phi ⟩,$$

$$Q=\sum _{k=1}^m\sum _{\begin{array}{c}{i}_1<\dots <i_k;\\ a_1<\dots <a_k\end{array}}|{\Phi }_{i_1\dots i_k}^{a_1\dots a_k}⟩⟨{\Phi }_{i_1\dots i_k}^{a_1\dots a_k}|.$$

$$Q{e}^{-T}H{e}^T|\Phi ⟩=Q\bar{H}|\Phi ⟩=0,$$

$$E=⟨\Phi |\bar{H}|\Phi ⟩.$$

$$H=\sum _{p,q}{h}_q^p{a}_p^{†}{a}_q+\frac{1}{4}\sum _{p,q,r,s}{v}_{rs}^{pq}{a}_p^{†}{a}_q^{†}{a}_s{a}_r,$$

$$T=T_{\mathrm{int}}(\mathfrak{h})+T_{\mathrm{ext}}(\mathfrak{h}).$$

$$|\Psi ⟩$$

$$|\Psi (\mathfrak{h})⟩=e^{T_{\mathrm{int}}(\mathfrak{h})}|\Phi ⟩.$$

$$b_{\ell }=\{\begin{array}{ll}{a}_{\ell }& \ell \in V\\ a_{\ell }^{†}& \ell \in O\end{array},b_{\ell }^{†}=\{\begin{array}{ll}{a}_{\ell }& \ell \in O\\ a_{\ell }^{†}& \ell \in V\end{array},$$

$$Q_{\mathrm{int}}{\bar{H}}_N|\Phi ⟩$$

$$Q_{\mathrm{ext}}{\bar{H}}_N|\Phi ⟩$$

$$⟨\Phi |\bar{H}|\Phi ⟩$$

$$(P+Q_{\mathrm{int}}){\bar{H}}_{\mathrm{ext}}{e}^{T_{\mathrm{int}}}|\Phi ⟩$$

$$Q_{\mathrm{ext}}{\bar{H}}_N|\Phi ⟩$$

$$T_{\mathrm{int}}$$

$$T_{\mathrm{ext}}$$

$$Q_{\mathrm{int}}$$

$$Q_{\mathrm{ext}}$$

$${\bar{H}}_{\mathrm{ext}}^{\mathrm{eff}}(\mathfrak{h})=(P+Q_{\mathrm{int}}(\mathfrak{h})){\bar{H}}_{\mathrm{ext}}(\mathfrak{h})(P+Q_{\mathrm{int}}(\mathfrak{h})),$$

$${\mathfrak{g}}^{(N)}\to \dots {\mathfrak{g}}^{(N)}(i_{\mathrm{R}})\to \dots {\mathfrak{g}}^{(N)}(2_{\mathrm{R}})\to {\mathfrak{g}}^{(N)}(1_{\mathrm{R}}),$$

$$|\Psi ⟩=e^{\sigma }|\Phi ⟩$$

$${\sigma }^{†}=-\sigma .$$

$$\sigma =T-T^{†},$$

$$e^X{e}^Y=e^{X+Y+{\sum }_{k=1}^{\infty }{C}_{\mathrm{BCH}}^{(k)}(X,Y)},$$

$$e^{X+Y}$$

$$R_{\mathrm{Z}}(X,Y)=\prod _{\infty }^{k=2}{e}^{C_{\mathrm{Z}}^{(k)}(X,Y)},$$

$$R_{\mathrm{Z}}^{-1}(X,Y)=\prod _{k=2}^{\infty }{e}^{-C_{\mathrm{Z}}^{(k)}(X,Y)}.$$

$$H{e}^{\sigma }|\Phi ⟩=E{e}^{\sigma }|\Phi ⟩.$$

Full text

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 Paldus and Li (1999); Bartlett and Musiał (2007), multi-reference concepts Mukherjee et al. (1975); Lindgren and Mukherjee (1987); Jeziorski and Monkhorst (1981); Meissner and Bartlett (1990); Kaldor (1991); Mahapatra et al. (1998); Meissner (1998); Li and Paldus (2004), and symmetry-breaking mechanisms have resulted in a myriad of approximations trying to describe static and dynamic correlations effects. For example, coupled cluster (CC) Coester (1958); Coester and Kummel (1960); Čížek (1966); Paldus et al. (1972); Purvis and Bartlett (1982), density matrix renormalization group (DMRG) White (1992); Schollwöck (2005); Chan and Sharma (2011), and density matrix methods Mazziotti and Erdahl (2001); Mazziotti (2012) have already demonstrated their efficiency in coping with complicated electron correlation effects in modest size molecular and materials systems. 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 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. Bravyi and Kitaev (2000); Seeley et al. (2012); Poulin et al. (2017); Setia and Whitfield (2018); Wecker et al. (2015); Babbush et al. (2018a); Motta et al. (2018); McClean et al. (2016); Shen et al. (2017a)

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. Kowalski (2018); Kowalski et al. (2018) 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 $A$ set is obtained from predefined classes of excitations (or sub-algebra(s)), and the $B$ set contains all the remaining parameters needed to describe the whole system. This technique, known as sub-system embedding sub-algebras 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. Inglesfield (1981)) and Löwdin partitioning techniques Löwdin (1965, 1982) 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 also introduce 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.

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 $|\Phi\rangle$ that can provide an adequate starting point for the exponential (CC) expansion of the ground-state electronic wave function $|\Psi\rangle$

[TABLE]

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

[TABLE]

In the second quantized form, each $T_{k}$ component can be expressed as

[TABLE]

where indices $i_{1},i_{2},\ldots$ ($a_{1},a_{2},\ldots$) refer to occupied (unoccupied) spin orbitals in the reference function $|\Phi\rangle$, $t^{i_{1}\ldots i_{k}}_{a_{1}\ldots a_{k}}$ are cluster amplitudes, and $E^{a_{1}\ldots a_{k}}_{i_{1}\ldots i_{k}}$ are excitation operators. The excitation operators are defined through strings of standard creation ($a_{p}^{\dagger}$) and annihilation operators ($a_{p}$)

[TABLE]

In the particle-hole formalism, the excitation operators are expressed in terms of particle/hole creation operators only. When acting on the reference function $|\Phi\rangle$, the $E^{a_{1}\ldots a_{k}}_{i_{1}\ldots i_{k}}$ operators produce the so-called excited configurations $|\Phi_{i_{1}\ldots i_{k}}^{a_{1}\ldots a_{k}}\rangle$ defined as

[TABLE]

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

[TABLE]

where $P$ and $Q$ are projection operators onto the reference function ($P=|\Phi\rangle\langle\Phi|$) and onto the excited configurations generated by the $T$ operator when acting on the reference function,

[TABLE]

A careful diagrammatic analysis Paldus and Li (1999) leads to an equivalent (at the solution), energy-independent form of the CC equations for the cluster amplitudes

[TABLE]

where the $\bar{H}=e^{-T}He^{T}$ is referred to as the similarity transformed Hamiltonian. It can also be shown that the $\bar{H}$ is expressible in terms of connected diagrams, i.e., $\bar{H}=(He^{T})_{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

[TABLE]

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

[TABLE]

where $p$, $q$, $r$, $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 $h^{p}_{q}$ and $v^{pq}_{rs}$ represent one- and two-electron integrals (in the above representation, the $v^{pq}_{rs}$ 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}$ Babbush et al. (2018a). 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. Peng and Kowalski (2017); Motta et al. (2018) for details)

From the point of quantum computing applications, the net effect of the number of basis set functions and the number of non-vanishing 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 Foster and Boys (1960); Edmiston and Ruedenberg (1963); Pipek and Mezey (1989). 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 Brueckner and Wada (1956); Nesbet (1958); Stolarczyk and Monkhorst (1984); Kobayashi et al. (1994); Nooijen and Lotrich (2000); Crawford et al. (1997) that maximize the overlap $|\langle\Phi_{B}|\Psi\rangle|$ between normalized lowest energy Slater determinant $|\Phi_{B}\rangle$ and the correlated wave function $|\Psi\rangle$. 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\ll N_{e}$. In this way, one can define standard approximations such as CCSD ($m=2$) Purvis and Bartlett (1982), CCSDT ($m=3$) Noga and Bartlett (1987, 1988); Scuseria and Schaefer (1988), CCSDTQ ($m=4$) Kucharski and Bartlett (1991); Oliphant and Adamowicz (1991), etc. 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 Dean and Hjorth-Jensen (2004); Kowalski et al. (2004); Hagen et al. (2008, 2016); Hirata et al. (2001, 2004); Booth et al. (2013); McClain et al. (2017). 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 sub-system embedding sub-algebras

Certain properties of CC equations are inextricably associated with the possibility of partitioning of cluster operators in CC wave function into components corresponding to various sub-algebras of excitation Lie algebra denoted here as $\mathfrak{g}^{(N)}$, which is generated by all excitation operators $E^{a_{l}}_{i_{l}}=a_{a_{l}}^{\dagger}a_{i_{l}}$. In a recent paper Kowalski (2018), we have analyzed properties of the CC equations stemming from the presence of CC-approximation-specific sub-algebras 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 sub-algebras provide a means to re-cast 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 non-linear eigenvalue problem (either in the context of dressed configuration interaction Hamiltonian Sánchez-Marín et al. (1996), inclusion of high-order excitations Kállay and Surján (2000), or the analysis of multiple solutions of CC equations Živković and Monkhorst (1978)), 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 actives-space indices. In typical applications the number of active spin orbitals ($N_{\rm act}$) is significantly smaller compared to the total number of spin orbitals $N_{\rm S}$, i.e., $N_{\rm act}\ll N_{\rm S}$.

Let us start the discussion by introducing basic notions defining sub-algebras. An important class of sub-algebras of the $\mathfrak{g}^{(N)}$ excitation algebra is closely related to ideas underlying the active-space concepts in quantum chemistry, where one can define sub-algebras $\mathfrak{h}$ through all possible excitations $E^{a_{1}\ldots a_{m}}_{i_{1}\ldots i_{m}}$ that excite electrons from a subset of active occupied orbitals (denoted as R) to a subset of active virtual orbitals (denoted as S). These sub-algebras will be denoted as $\mathfrak{g}^{(N)}({\rm R,S})$. The $\mathfrak{g}^{(N)}({\rm R,S})$ sub-algebras can also be viewed as generators of various complete active spaces (CAS(R,S)) spanned by the reference function $|\Phi\rangle$ and all excited configurations obtained by acting with elements of $\mathfrak{g}^{(N)}({\rm R,S})$ onto $|\Phi\rangle$ (see Fig. 2). Specific examples of these sub-algebras contain sub-algebras 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 $\mathfrak{g}^{(N)}({\rm S})$ and $\mathfrak{g}^{(N)}({\rm R})$, respectively. In an alternative notation, we will denote $\mathfrak{g}^{(N)}({\rm R},{\rm S})$ as $\mathfrak{g}^{(N)}(x_{\rm R},y_{\rm S})$ where numbers $x$ and $y$ refer to the number of orbitals included in sets R and S, respectively. Special sub-algebras $\mathfrak{g}^{(N)}({\rm R})$ and $\mathfrak{g}^{(N)}({\rm S})$ will be denoted as $\mathfrak{g}^{(N)}(x_{\rm R})$ and $\mathfrak{g}^{(N)}(y_{\rm S})$. In this paper, we will entirely focus on active-space excitations sub-algebras for the closed-shell single-reference CC formulations.

An important property of the excitation sub-algebras is the fact, that in general, an arbitrary cluster operator $T$ can be decomposed in a part that belongs to sub-algebra of interest $\mathfrak{h}$ (denoted as an internal part, $T_{\rm int}(\mathfrak{h})$) and a part that belongs to $\mathfrak{g}^{(N)}-\mathfrak{h}$ (denoted as an external part, $T_{\rm ext}(\mathfrak{h})$), i.e.,

[TABLE]

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

[TABLE]

where the CAS-type wave function $|\Psi(\mathfrak{h})\rangle$ is defined as

[TABLE]

In (12), we use the fact that $[T_{\rm ext},T_{\rm 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 $E_{i_{1},\ldots,i_{k}}^{a_{1},\ldots,a_{k}}=b_{a_{1}}^{\dagger}\cdots b_{a_{k}}^{\dagger}b^{\dagger}_{i_{k}}\cdots b^{\dagger}_{i_{1}}$, for anti-commuting $b_{\ell}^{\dagger}$, where $b_{\ell}$/$b_{\ell}^{\dagger}$ operators are defined as

[TABLE]

and $O$ and $V$ represent sets of occupied and unoccupied spin orbitals in $|\Phi\rangle$, respectively.

The sub-algebras $\mathfrak{h}$ that satisfy two important requirements:

The $|\Psi(\mathfrak{h})\rangle=e^{T_{\rm int}(\mathfrak{h})}|\Phi\rangle$ is characterized by the same symmetry properties as $|\Psi\rangle$ and $|\Phi\rangle$ vectors (for example, spin and spatial symmetries). 2. 2.

The $e^{T_{\rm int}(\mathfrak{h})}|\Phi\rangle$ Ansatz generates the FCI expansion for the subsystem defined by the CAS corresponding to the $\mathfrak{h}$ sub-algebra,

play an important role in further analysis of the structure of the CC equations and we will refer all sub-algebras satisfying requirements (1) and (2) as sub-system embedding sub-algebras (SESs).

In Ref. Kowalski (2018) it has been shown that $\mathfrak{g}^{(N)}(1_{\rm R},y_{\rm S})$ along with $\mathfrak{g}^{(N)}(x_{\rm R},1_{\rm S})$ form SESs for CCSD and $\mathfrak{g}^{(N)}(2_{\rm R},y_{\rm S})$ along with $\mathfrak{g}^{(N)}(x_{\rm R},2_{\rm S})$ are SESs for CCSDTQ. The largest sub-algebras in these classes are $\mathfrak{g}^{(N)}(1_{\rm R})$ and $\mathfrak{g}^{(N)}(2_{\rm R})$ sub-algebras. In Ref. Kowalski (2018), we have also demonstrated that for any SES $\mathfrak{h}$ (corresponding to some CC approximation) the standard

[TABLE]

and hybrid

[TABLE]

representations of CC equations are equivalent at the solution. In Eqs. (15)-(19) we used simplified notations: (1) the projection operators $Q_{\rm int}(\mathfrak{h})$ and $Q_{\rm ext}(\mathfrak{h})$ ($Q=Q_{\rm int}(\mathfrak{h})+Q_{\rm ext}(\mathfrak{h})$) project onto subspaces spanned by all excited configurations generated by acting $T_{\rm int}(\mathfrak{h})$ and $T_{\rm ext}(\mathfrak{h})$ onto the reference function, respectively, (2) $\bar{H}_{N}=\bar{H}-\langle\Phi|\bar{H}|\Phi\rangle$ is the normal product form of the similarity transformed Hamiltonian $\bar{H}$, (3) the $\bar{H}_{\rm ext}$ operator is defined as $\bar{H}_{\rm ext}\equiv\bar{H}_{\rm ext}(\mathfrak{h})=e^{-T_{\rm ext}(\mathfrak{h})}He^{T_{\rm ext}(\mathfrak{h})}$, and (4) for notational simplicity we used the following notational convention,

[TABLE]

The above mentioned equivalence means that cluster amplitudes corresponding to excitations included in $\mathfrak{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 $\mathfrak{h}=\mathfrak{g}^{(N)}(0)$ (where $\mathfrak{g}^{(N)}(0)$ contains no excitations) - in this case $T_{\rm int}(\mathfrak{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 $\bar{H}_{\rm ext}^{\rm eff}(\mathfrak{h})$,

[TABLE]

is the corresponding complete active space. By construction, the cluster amplitudes $T_{\rm ext}$, used to define $\bar{H}_{\rm ext}(\mathfrak{h})$, are decoupled from cluster amplitudes $T_{\rm int}$ that define the components of the corresponding eigenvector. One should also notice that the many-body expansion of $\bar{H}_{\rm ext}(\mathfrak{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 re-cast 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 non-interacting 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

[TABLE]

which results in separate eigenvalue problems corresponding to sub-systems of various sizes. This observation can be used to define a new CC approximations where instead of referring to adding higher and higher ranks of excitations as a design principle (used to defined standard approximations such as CCD, CCSD, CCSDT, CCSDTQ, etc.) 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. Kowalski (2018) we discussed an approximation (sub-algebra flow CCSD(2) approximation (SAF-CCSD(2)) that employs all amplitudes contained in all $\mathfrak{g}^{(N)}(2_{\rm R})$ SESs, which leads to a CC model that contains all singly and double 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. Kowalski (2018)). 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.

To summarize this section, techniques based on the utilization of the sub-system embedding sub-algebras can be used to downfold the full electronic Hamiltonian to an arbitrary active space corresponding to some subsystem embedding sub-algebra. 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) Luis and Peřina (1996); Cleve et al. (1998); Berry et al. (2007); Childs (2010); Wecker et al. (2015); Häner et al. (2016); Poulin et al. (2017) or variational quantum eigensolver (VQE) Peruzzo et al. (2014); McClean et al. (2016, 2017); Fontalvo et al. (2017); Shen et al. (2017a); Kandala et al. (2018); Colless et al. (2018). Unfortunately, the main caveat is related to the fact that effective Hamiltonians $\bar{H}_{\rm ext}(\mathfrak{h)}$ are non-Hermitian, which is a major obstacle in using these algorithms. The main goal of the following section 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 $\bar{H}_{\rm ext}(\mathfrak{h)}$ discussed in this section.

IV Unitary CC method

The unitary CC (UCC) formulations have been introduced in Refs. Pal (1984); Hoffmann and Simons (1988); Bartlett et al. (1989); Kutzelnigg (1991); Taube and Bartlett (2006); Sur et al. (2008); Evangelista (2011a); Cooper and Knowles (2010); Harsha et al. (2018) and have been intensively studied in the recent years in the context of various quantum algorithms Peruzzo et al. (2014); McClean et al. (2016); Romero et al. (2018); Shen et al. (2017b); Motta et al. (2018).

The standard UCC Ansatz (the generalization of standard UCC Ansatz has recently been discussed in the literature, for details see Ref. Lee et al. (2018))

[TABLE]

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

[TABLE]

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

[TABLE]

where $T$ has exactly the same structure as discussed in Eq. (2). In this and following Sections 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:

1.) The Baker-Campbell-Hausdorff (BCH) formula:

[TABLE]

where commutators $C_{\rm BCH}^{(k)}$ are defined as $C_{\rm BCH}^{(1)}(X,Y)=\frac{1}{2}[X,Y]$, $C_{\rm BCH}^{(2)}(X,Y)=\frac{1}{12}([X,[X,Y]]+[Y,[Y,X]])$, etc.

2.) The transposed variant of the Zassenhaus formula Wilcox (1967); Scholz and Weyrauch (2006):

[TABLE]

where $C_{\rm Z}^{(2)}=\frac{1}{2}[X,Y]$, $C_{\rm Z}^{(3)}=\frac{1}{3}[Y,[[X,Y]]+\frac{1}{6}[X,[X,Y]]$, etc. The $R_{\rm Z}(X,Y)$ function is defined as

[TABLE]

where its inverse is given by the formula

[TABLE]

The equations for $\sigma$ operator in the exact limit can be obtained by substituting ansatz (26) into the Schrödinger equation

[TABLE]

The equivalent representation, which provides the explicitly connected form of the equations, can be obtained by multiplying both sides by $e^{-\sigma}$ operator and decoupling equations for amplitudes from the equation for energy (see Ref. Taube and Bartlett (2006); Evangelista (2011b) for details), i.e.,

[TABLE]

In analogy to the standard single-reference approach let us explore if partitioning of $\sigma$ into part belonging to some SES $\mathfrak{h}$ ($\sigma_{\rm int}\equiv\sigma_{\rm int}(\mathfrak{h})$) and its complement ($\sigma_{\rm ext}\equiv\sigma_{\rm ext}(\mathfrak{h})$) leads to a downfolding of Hamiltonian in the form discussed in Section III. For this purpose we will (1) apply Zassenhaus formula (30) to factorize $e^{\sigma_{\rm int}+\sigma_{\rm ext}}$, (2) premultiply (33) from the left by $e^{-\sigma_{\rm ext}}R_{Z}^{-1}(\sigma_{\rm int},\sigma_{\rm ext})$, and (3) project resulting equations onto $P+Q_{\rm int}$ space. This procedure leads to the following equations

[TABLE]

For simplicity, the above equation can be re-written as

[TABLE]

where the $\bar{H}_{\rm ext}^{\rm UCC}$ operator is defined as

[TABLE]

or equivalently

[TABLE]

where the effective Hamiltonian $\bar{H}_{\rm ext}^{\rm eff(UCC)}$ in Eq. (39) is defined as

[TABLE]

In the above equation we use the fact that

[TABLE]

It can be shown that the form of the active-space eigenvalue-type problem (39) is equivalent to the $(P+Q_{\rm 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^{-\sigma_{\rm int}}e^{\sigma_{\rm int}}$ to the left of $\bar{H}_{\rm ext}^{\rm UCC}$ in Eq. (37) and noticing that $e^{-\sigma_{\rm int}}\bar{H}_{\rm ext}^{\rm UCC}e^{\sigma_{\rm int}}=e^{-\sigma}He^{\sigma}$. Then the resulting equation takes the form

[TABLE]

Using matrix representation of the $\sigma_{\rm int}$ operator in the CAS space denoted as $\bm{\sigma}_{\rm int}$ this equation can be re-written as

[TABLE]

where the first component of vector $\bm{x}$ corresponds to $\langle\Phi|e^{-\sigma}He^{\sigma}|\Phi\rangle-E$ while the remaining components correspond to projections of $e^{-\sigma}He^{\sigma}|\Phi\rangle$ onto excited (with respect to the reference determinant $|\Phi\rangle$) configurations $\langle\Phi_{\Delta}^{\rm(CAS)}|$ belonging to $\mathfrak{h}$-induced CAS. The $[e^{\bm{\sigma}_{\rm int}}]$ matrix is also non-singular, which is a consequence of the algebraic structure of $\sigma_{\rm int}$ rather than a particular values of cluster amplitudes defining $T_{\rm int}$ and $T_{\rm int}^{\dagger}$. To see it let us calculate its determinant

[TABLE]

where $\bm{T}_{\rm int}$ and $\bm{T}^{\dagger}_{\rm int}$ are matrix representations of $T_{\rm int}$ and $T^{\dagger}_{\rm int}$ in CAS. Since $T_{\rm int}$ and $T^{\dagger}_{\rm int}$ contains either excitations or de-excitations we have

[TABLE]

which means that ${\rm det}(e^{\bm{\sigma}_{\rm int}})=1$ and $(e^{\bm{\sigma}_{\rm int}})$ is a non-singular operator. Therefore, the only solution of Eq. (43) corresponds to $\bm{x}=0$, which proves the equivalence of Eq. (37) with $P$ and $Q_{\rm int}$ projections of Eqs. (34) and (35).

Although $\bar{H}_{\rm ext}^{\rm UCC}$ ($\bar{H}_{\rm ext}^{\rm eff(UCC)}$) is Hermitian, in contrast to $\bar{H}_{\rm ext}$ from Eq. (18) the $\bar{H}_{\rm ext}^{\rm UCC}$ ($\bar{H}_{\rm ext}^{\rm eff}$) operator does not decouple $\sigma_{\rm ext}$ amplitudes from the $\sigma_{\rm int}$ ones. If $\mathfrak{h}$ is chosen to contain highest/lowest lying occupied/virtual orbitals we can view this property of $\bar{H}_{\rm ext}^{\rm UCC}$ 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 the next section 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 $\mathfrak{h}$

[TABLE]

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 Stolarczyk and Monkhorst (1985); Meissner and Nooijen (1995); Nooijen (1996); Nooijen and Bartlett (1997); Meissner (1998). Moreover, ansatz (46) can also be viewed as a unitary generalization of the tailored CC formalism (TCC) Kinoshita et al. (2005); Hino et al. (2006); Melnichuk and Bartlett (2012, 2014), where equations for $\sigma_{\rm int}$ are represented in the form of eigenvalue problem discussed in previous sections.

In analogy to the previous section, we will focus on the exact formulation of the DUCC approach that includes all possible excitations in the $\sigma_{\rm ext}$ and $\sigma_{\rm int}$ operators. Since using BCH expansion (29) the DUCC expansion can be transformed to the alternative single-reference ansatz

[TABLE]

where $D$ is anti-Hermitian ($D^{\dagger}=-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

[TABLE]

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

[TABLE]

One can show that the equations (48) corresponding to projections onto $(P+Q_{\rm ext})$ sub-space can be written in equivalent form as an eigenvalue problem

[TABLE]

or equivalently, using effective Hamiltonian language,

[TABLE]

where

[TABLE]

and

[TABLE]

To show this fact it suffices to introduce the resolution of identity $e^{\sigma_{\rm int}}e^{-\sigma_{\rm int}}$ to the left of the $\bar{H}_{\rm ext}^{\rm DUCC}$ operator in Eq. (51) and notice that $e^{-\sigma_{\rm int}}\bar{H}_{\rm ext}^{\rm DUCC}e^{\sigma_{\rm int}}=e^{-\sigma_{\rm int}}e^{-\sigma_{\rm ext}}He^{\sigma_{\rm int}}e^{\sigma_{\rm ext}}$. Next, in analogy to Eqs. (42) and (43), Eq. (53) can be represented as

[TABLE]

where the first component of the $[\bm{y}]$ vector is equivalent to $\langle\Phi|e^{-\sigma_{\rm int}}e^{-\sigma_{\rm ext}}He^{\sigma_{\rm int}}e^{\sigma_{\rm ext}}|\Phi\rangle-E$ while the remaining components correspond to projections of $e^{-\sigma_{\rm int}}e^{-\sigma_{\rm ext}}He^{\sigma_{\rm int}}e^{\sigma_{\rm ext}}|\Phi\rangle$ onto excited configurations belonging to $Q_{\rm int}$. Given the non-singular character of the $[e^{\bm{\sigma}_{\rm 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 $\sigma_{\rm ext}$ amplitudes only, providing in this way in Eq. (52) a rigorous decoupling of degrees of freedom corresponding to ${\sigma_{\rm ext}}$ and ${\sigma_{\rm int}}$ in the sense of discussion of Section III. Since amplitudes defining $\sigma_{\rm ext}$ are characterized by larger perturbative denominator compared to the $\sigma_{\rm int}$, where small denominators may occur, it is much safer to determine $\sigma_{\rm ext}$ using perturbative techniques. The Eq. (52) also offers a possibility of downfolding the Hamiltonian to the $(P+Q_{\rm int})$ space where the correlation effects from $Q_{\rm ext}$ can be included through the $\sigma_{\rm 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 $\bar{H}_{\rm ext}^{\rm eff(DUCC)}$ operator (denoted as the $\Gamma$ 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,\ldots$ and some small subset of virtual spin orbitals $a,b,\ldots$ are deemed active). In this case amplitudes defining the $\sigma_{\rm 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 $\bar{H}_{\rm ext}^{\rm eff(DUCC)}$ in powers of $\sigma_{\rm ext}$ operator and (2) approximation of $\sigma_{\rm ext}$ operator. The $\Gamma$ operator is defined by strings of creation/annihilation operators that carry only active spin orbitals and can be written as

[TABLE]

where subscript ${\rm act}$ designates terms of a given operator expression that contains creation/annihilation operators carrying active-space spin-orbital labels. Using Baker-Campbell-Hausdorff formula one can re-cast the above equation in the form of infinite expansion:

[TABLE]

Using particle-hole (ph) formalism one can also expand the $\Gamma$ operator into the sum of its many-body components $\Gamma_{i}$

[TABLE]

where $N_{e,{\rm act}}$ designates the number of active electrons and

[TABLE]

where $N[...]^{\rm ph}$ designates the particle-hole normal product form of a given operator expression, $P,Q,R,S,\ldots$ designate general active spin-orbital indices, and $\Gamma_{\rm scalar}$ designates full contracted (scalar) part of the $\Gamma$ operator with respect to the particle-hole vacuum $|\Phi\rangle$. In the above expansion we also assume that all multidimensional tensors $\gamma^{P_{1}\ldots P_{i}}_{Q_{1}\ldots Q_{i}}$ 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 $\Gamma$ expansion given by Eq.(61) as well as excitation level included in the $\sigma_{\rm ext}$ operator. Below, as a specific example we describe a variant of $\Gamma$-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. Bartlett et al. (1989), 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-\langle\Phi|H|\Phi\rangle$) 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.,

[TABLE]

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

[TABLE]

For the simplest approximation of $T_{{\rm ext},1}$ and $T_{{\rm 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 $\Gamma$ 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:

V.1 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_{\rm ext}$ and $T_{\rm ext}^{\dagger}$ given by Eq.(63))

[TABLE]

and retaining terms only through two-body interactions one obtains

[TABLE]

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,\ldots$ and $A,B,C,\ldots$, respectively), and $(...)_{C,{\rm open}}$ designates connected and open (i.e., having external lines in diagrammatic representation) part of a given operator expression. In Eq.(65) $(\Gamma)_{\rm scalar}$ corresponds to the scalar part of $\Gamma$ in particle-hole representation, i.e. $\Gamma_{\rm scalar}=\langle\Phi|\Gamma|\Phi\rangle$. 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, $(H_{N}T_{\rm ext})_{C,{\rm open}}=(T_{\rm ext}^{\dagger}H_{N})_{C,{\rm open}}^{\dagger}$). 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^{P}_{Q}$ and $\gamma^{PQ}_{RS}$ amplitudes are given in Tables 1, 2, and 3. For these tables, we assume that a restricted or unrestricted Hartree-Fock (RHF/UHF) reference is employed (i.e. all non-diagonal elements of the Fock matrix are zero), and as a result the third and sixth terms in Eq. 64 vanish.

V.2 Determination of the physical-vacuum representation of the $\Gamma$ operator.

In order to find an equivalent characterization of the $\Gamma$ operator given by Eq.(65) using the physical vacuum parametrization we will employ the set of identities from Table 4 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 $\Gamma$, in one- and two-body interaction approximation takes the form:

[TABLE]

where all $\chi^{P}_{Q}$ and $\chi^{PQ}_{RS}$ coefficients are listed in Table 5. 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).

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 $\mathfrak{h}$,

• -

the accuracy of $\sigma_{\rm ext}$ estimates represented by single and double excitations,

• -

the role of missing higher-rank many-body effects in the $\sigma_{\rm ext}$ and $\Gamma$ operators.

As stated earlier, the size of the active space may effect the accuracy of $\sigma_{\rm ext}$ amplitudes, which is a consequence of the fact that utilization of larger spaces, whose choice is driven by the value of orbital energies, prevent 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 $\sigma_{\rm ext}$ is proportionally smaller. Otherwise one needs to include higher many-body components of $\sigma_{\rm ext}$ and $\Gamma$ operator, for example, three- and/or four-body components, i.e., $\sigma_{{\rm ext},3}$, $\sigma_{{\rm ext},4}$, $\ldots$, and $\Gamma_{3}$, $\Gamma_{4}$, $\ldots$. 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 $\sigma_{{\rm ext},3}$ and/or $\sigma_{{\rm ext},4}$ operators. Instead one should resort to using genuine UCC formulations. For example, the $T_{{\rm ext},3}$ amplitudes can be extracted from the UCC(4) model discussed in Ref. Bartlett et al. (1989), where the sufficiency conditions for cluster amplitudes

[TABLE]

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 $\sigma_{\rm ext}$ and forming $\chi^{P}_{Q}$, $\chi^{PQ}_{RS}$, $\ldots$ 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_{\rm act}$) and total number of spin orbitals ($N_{S}$) and rank of many-body effects included in the $\Gamma$ operator expression 66.

The specific choice of the active-spave or equivalently sub-system embedding algebra $\mathfrak{h}$ defines how efficient is the process of integrating out remaining degrees of freedom (i.e. the parameters/amplitudes defining the $T_{\rm ext}$/$\sigma_{\rm ext}$ operator). In particular, the proper choice of $\mathfrak{h}$ will impact the accuracy of low-cost (perturbative) estimates of $T_{\rm ext}$. Here we will consider two cases: (1) choice of the $\mathfrak{h}$ based on the energy threshold and (2) choice of the $\mathfrak{h}$ based on locality criteria (or equivalently sub-system separability discussed in previous Sections).

In the first case the active space is chosen in analogy to typical applications of multi-reference methods (CASSCF, CASPT2 Andersson et al. (1990, 1992), NEVPTn Angeli et al. (2001a, b), MRMBPT Hirao (1992); Finley et al. (1995, 1996); Chaudhuri et al. (2005); Nakano et al. (1998); Rintelman et al. (2005); Witek et al. (2002, 2003), and DMRG Sokolov et al. (2017); Guo et al. (2018) methods) where active spaces usually contain high- and low-lying occupied and virtual orbitals. In this situation the first order contribution to $T_{\rm ext}$ can be for example written as

[TABLE]

where in this specific example $s$-amplitude contains one active occupied (J) and one active virtual (B) spinorbital indices, we will also assume that in the above example $i$ and $a$ represent inactive occupied and virtual spinorbital indices. If active-space orbitals are well separated (energetically) from the remaining orbitals, one can expect that the perturbative denominators used to define $T_{\rm 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_{\rm ext}$ estimates.

In contrast to the energy separation criteria, when "spatial" arguments are invoked to define the active space, the "smallness" of $\sigma_{\rm 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) "active" or strongly correlated sub-system is weakly interacting with sub-system 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 the case (a), the matrix element in Eq.(70)

[TABLE]

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 $\mathfrak{h}$) defined by active orbitals defining sub-systems $A_{1}$ (with corresponding SES $\mathfrak{h}_{1}$) and $A_{2}$ (with corresponding SES $\mathfrak{h}_{2}$), i.e.,

[TABLE]

where the downfolded Hamiltonian $\bar{H}_{\rm ext}^{\rm DUCC}(\mathfrak{h})$ is given by expression

[TABLE]

Once sub-system $A_{1}$ and $A_{2}$ are spatially separated and localized basis set is employed and $\sigma_{\rm int}(\mathfrak{h}_{1})$ and $\sigma_{\rm int}(\mathfrak{h}_{2})$ commute

[TABLE]

then further downfolding of $\bar{H}_{\rm ext}^{\rm DUCC}(\mathfrak{h})$ is possible, i.e.,

[TABLE]

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 $\sigma_{\rm ext}$ and forming $\chi^{P}_{Q}$, $\chi^{PQ}_{RS}$, $\ldots$ 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 the number of active orbitals ($N_{\rm act}$) and total number of spinorbitals ($N_{\rm S}$) and rank of many-body effects included in the $\Gamma$ 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_{\rm S}$ orbitals versus DUCC models involving $N_{\rm act}$ active spin orbitals for two situations (2) $\Gamma$ operator contains (1) up to two-body interactions ($\Gamma_{1}$ and $\Gamma_{2}$) and (2) up to three-body interactions ($\Gamma_{1}$, $\Gamma_{2}$, and $\Gamma_{3}$).

Given an arbitrary unitary $U$ with eigenvalues $U|\psi_{j}\rangle=e^{i\theta_{j}}|\psi_{j}\rangle$ and an input state $|\psi\rangle=\sum_{j}{\alpha_{j}}|{\psi_{j}}\rangle$, quantum phase estimation returns an eigenphase estimate of $\theta_{j}$ with error $\delta$ for a randomly chosen eigenvector $|{\psi_{j}}\rangle$ that is sampled with probability $|\alpha_{j}|^{2}$. More precisely, the estimate $\hat{\theta}$ is drawn from the distribution

[TABLE]

where $f(x)$ is a function that depends on the choice of phase estimation algorithm and is sharply peaked about $x=0$ with width $\delta$. In common variants of phase estimation algorithms Nielsen and Chuang (2011); Kimmel et al. (2015); Wiebe and Granade (2016), controlled-$U$ must be applied $\mathcal{O}(1/\delta)$ times to obtain a single $\hat{\theta}$ estimate, and is the dominant cost.

One common example of $U$ is the real-time evolution operator $e^{-iHt}$, in which case $\theta_{j}$ is an eigenvalue of $H=\sum^{\mathcal{O}(N_{\rm terms})}_{j}h_{j}P_{j}$, for ${\mathcal{O}(N_{\rm terms})}$ Pauli operators $P_{j}$ acting on at most $\mathcal{O}(N_{\rm S})$ qubits, and positive coefficients $h_{j}$, scaled by a constant $t$. Though real-time evolution may be approximated using Trotter-Suzuki product formulas Berry et al. (2007), it is difficult to obtain tight error bounds on scaling of the approximation error with $N_{\rm act}$. Thus we consider the case where $U$ is a quantum walk Childs (2010); Low and Chuang (2016, 2017) with eigenvalues commensurate with those of $e^{\pm i\sin^{-1}(H/\lambda)}$, for $\lambda=\sum_{j}|h_{j}|\geq\|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 Low and Chuang (2016, 2017); Babbush et al. (2018b); Berry et al. (2018); Poulin et al. (2018). Note that an eigenphase estimate $\hat{\theta}$ of the quantum walk and its error $\delta$ may be related to the eigenvalue estimate of $H$ by computing $\lambda\sin(\hat{\theta})$. Thus if we wish to learn the eigenvalue within error $\epsilon$ we need to apply the walk operator $\mathcal{O}(\lambda/\epsilon)$ times.

In qubitization, the walk operator is of the form $U={(\rm PREPARE^{\dagger}\otimes I)}\cdot{\rm SELECT}\cdot{(\rm PREPARE\otimes I)}\cdot((1-2|0\cdots 0\rangle\langle 0\cdots 0|)\otimes I)$ for unitary subroutines PREPARE and SELECT. The PREPARE subroutine prepares from the all zero state $|0\cdots 0\rangle$ a state of the form $\sum_{j}\sqrt{h_{j}/\lambda}|j\rangle$ within error $\epsilon$, using $\mathcal{O}(N_{\rm terms})$ 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 $\mathcal{O}(N_{\rm S}N_{\rm terms})$ gates, but optimizations specific to the Jordan-Wigner representation of fermionic Hamiltonians reduce this to $\mathcal{O}(N_{\rm S})$ Babbush et al. (2018b), 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 $\delta$ is in $\widetilde{\mathcal{O}}(\lambda N_{\mathrm{terms}}/\delta)$. In the worst-case, we may assume that $N_{\mathrm{terms}}=\mathcal{O}(N_{\rm S}^{4})$, and $\lambda=\mathcal{O}(|H|_{\mathrm{max}}N_{\rm S}^{4})$, where $|H|_{\mathrm{max}}$ is the maximum absolute value of any entry in the Hamiltonian.

While $|H|_{\mathrm{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 straight forward to show that $|h_{pq}|\in\mathcal{O}(\max_{i}Z_{i}^{2})$ and $|h_{pqrs}|\in\mathcal{O}(\max_{i}Z_{i})$ Babbush et al. (2015). Therefore, we anticipate that $\lambda\in\mathcal{O}(N_{\rm S}^{2}\max_{i}Z_{i}^{2}+N_{\rm S}^{4}\max_{i}Z_{i})$ 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 $\delta$ to be in $\mathcal{O}(N_{\rm S}^{8}/\delta)$ in such cases.

If we consider simulating the downfolded Hamiltonian, one simply replaces $N_{\mathrm{terms}}=\mathcal{O}(N^{4}_{\rm act})$ with the number of terms in the downfolded effective Hamiltonian, and similarly for the normalization constant $\lambda=\mathcal{O}(|\Gamma|_{\mathrm{max}}N^{4}_{\rm act})$. Importantly, the number of active space orbitals is much smaller than that of the full Hamiltonian, that is $N_{\rm act}\ll N_{\rm S}$. In common practice one normally chooses active spaces such that that amplitude corrections of downfolding are small. In other words, $|\Gamma|_{\mathrm{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 $|\Gamma|_{\mathrm{max}}\ll|H|_{\mathrm{max}}$. The case of simulating three-body interactions is less-studied, but the essential idea is identical. With three-body terms, the worst-case $N_{\mathrm{terms}}=\lambda=\mathcal{O}(N_{\rm act}^{6})$, though $\lambda$ could be much smaller if $|\Gamma_{3}|_{\mathrm{max}}\ll|\Gamma|_{\mathrm{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_{\mathrm{terms}}$ to as small as $\widetilde{\mathcal{O}}(N_{\rm act}^{2})$ Motta et al. (2018) for the two-body case – low-rank approximations for the three-body case are still not well understood. Moreover, various quantum circuit optimization techniques Low et al. (2018) enable a further $\mathcal{O}(N_{\rm act})$ 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 Low and Wiebe (2018); Low (2018).

VIII 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 coninical 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 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.

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-SVD decompositions of one- and two-electron integrals Peng and Kowalski (2017); Motta et al. (2018). One can envisage a further extension of the applicability of Cholesky-SVD decomposition to compress $\chi$-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) White (1992); Schollwöck (2005); Legeza and Sólyom (2003); Chan and Sharma (2011) 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 Yanai et al. (2010); Kurashige and Yanai (2011); Guo et al. (2016).

IX acknowledgement

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 Pacific Northwest National Laboratory (PNNL). It was conducted under the Laboratory Directed Research and Development Program at PNNL.