Quantum self-consistent equation-of-motion method for computing molecular excitation energies, ionization potentials, and electron affinities on a quantum computer

Near-term quantum computers are expected to facilitate material and chemical research through accurate molecular simulations. Several developments have already shown that accurate ground-state energies for small molecules can be evaluated on present-day quantum devices. Although electronically excited states play a vital role in chemical processes and applications, the search for a reliable and practical approach for routine excited-state calculations on near-term quantum devices is ongoing. Inspired by excited-state methods developed for the unitary coupled-cluster theory in quantum chemistry, we present an equation-of-motion-based method to compute excitation energies following the variational quantum eigensolver algorithm for ground-state calculations on a quantum computer. We perform numerical simulations on H2, H4, H2O, and LiH molecules to test our quantum self-consistent equation-of-motion (q-sc-EOM) method and compare it to other current state-of-the-art methods. q-sc-EOM makes use of self-consistent operators to satisfy the vacuum annihilation condition, a critical property for accurate calculations. It provides real and size-intensive energy differences corresponding to vertical excitation energies, ionization potentials and electron affinities. We also find that q-sc-EOM is more suitable for implementation on NISQ devices as it is expected to be more resilient to noise compared with the currently available methods.


I. INTRODUCTION
Quantum chemistry is expected to be one of the first areas which can have demonstrable quantum advantage in the near term [1][2][3][4][5][6][7][8].This is owed to the fact that the computational effort required for exact evaluation of electron correlation on a classical computer-whose accurate calculation is essential for a reliable comparison with experimental values-scales factorially with the number of molecular orbitals.This unfavourable scaling is expected to reduce drastically when wavefunctions are instead prepared on quantum devices.
For estimation of molecular ground-state properties on noisy intermediate-scale quantum (NISQ) era devices, variational quantum eigensolver (VQE) based algorithms have gained popularity due to their relatively low circuit depth and resilience to noise [9,10].This has led to a series of successful demonstrations involving the computation of molecular ground-state energies of small molecules on present-day quantum devices and simulators [4,6,[11][12][13][14][15][16][17][18][19][20][21][22].However, estimation of just the molecular groundstate energy is not sufficient for describing many interesting chemical processes that involve electronic excitations in some form [23].For example, accurate modelling of chemical phenomena such as photochemical reactions, catalytic processes involving transition metal complexes, photosynthesis, solar cell operation, etc. require an accurate simulation of both molecular ground and excited states.The electronically excited states of such systems are generally strongly correlated and hence, require the use of sophisticated quantum chemical theories for their accurate description.A number of methods have been developed in this regard in the last few decades [24][25][26][27][28][29][30][31][32].and widely used through several software packages [33,34].The equation-of-motion coupled-cluster (EOM-CC) [26] approach, originally developed by Stanton and Bartlett, is a popular example that is routinely used to calculate molecular excited-state properties such as excitation energies and transition dipole moments [35][36][37][38][39]. EOM-CC has also been extended to calculate energies required to add or remove electrons from the ground-state electronic configuration [40][41][42][43][44].For example, IP-EOM-CC [40,42] and EA-EOM-CC [41] approaches have been developed which can compute accurate vertical ionization potentials (IPs) and vertical electron affinities (EAs), respectively.IPs/EAs are defined as the difference in energy between the ground state and the states obtained by a single electron detachment/attachment process.Some of the other advantages associated with the EOM-CC formalism are its theoretical rigour, the accuracy and correct scaling behavior of energy differences computed, and the ability to systematically improve the results.However, standard quantum chemistry methods like EOM-CC sometimes face challenges in a quantitative determi-nation of excited states and their properties, notably for same-symmetry conical intersections [45][46][47][48] and when the ground state has a prominent multi-reference character [49][50][51][52].Since VQE algorithms are expected to provide accurate ground-state wavefunctions, even in the case of strongly correlated systems, NISQ era devices can help address these challenging problems with practical computational expenses.
We would like to note that methods for the estimation of molecular excited states on a quantum computer based on other popular quantum algorithms have also been proposed.A number of approaches are based on quantum phase estimation algorithm [53][54][55][56][57][58][59] with new developments for efficient implementation on quantum computers [60][61][62].Methods based on Krylov subspace diagonalization [63,64] and quantum annealing [65,66] have also been proposed.While these methods are theoretically exact (in absence of any noise) and expected to provide a significant computational advantage over exact treatment on a classical computer, they will mostly be useful in fault-tolerant quantum computing and not suitable for NISQ era quantum computers due to their high quantum resource requirements and low tolerance to noise.
Significant effort has been made in developing methods for the calculation of molecular excitation energies within the framework of VQE in the last few years.These techniques can be broadly classified into circuit optimization and diagonalization-based approaches.In the former approach, optimal parametrized circuits are obtained for every excited state, usually by minimizing a cost function involving energies of one or multiple excited states.Subspace-search VQE (SS-VQE) [67], orthogonal state reduction variational eigensolver (OSRVE) [68], variational quantum deflation (VQD) [69,70] and the folded spectrum method [4] are some examples.These approaches, however, generally require increased quantum resources, specifically the gate depth.This makes them challenging for near term applications.Moreover, there is no guarantee for them to find the entire spectrum when the states are close in energy to one another.On the other hand, the diagonalization-based approaches use a classical computer to diagonalize the Hamiltonian in a subspace and can provide several excited states simultaneously.In this regard, methods like the Quantum Krylov subspace expansion [63,71,72], the Quantum subspace expansion (QSE) [73][74][75][76], and the quantum equation-of-motion (qEOM) [77] have been developed recently.QSE has had significant success in the last few years and has also been extended to capture the missing correlation from large virtual orbital spaces [75,78].However, it requires an estimate of higher than 2-body reduced density matrices (RDMs), prompting the use of cumulant approximations [75] inspired by developments in quantum chemistry [31,32,79].Furthermore, a significant drawback of the QSE approach is the lack of size-intensivity of the computed excitation energies.The property of size-intensivity ensures correct scaling of ex-citation energies computed by a method with increasing size of the system.The violation of this property can lead to errors and even non-physical predictions, for instance, the QSE excitation energies of a "super molecule" consisting of two non-interacting systems is not guaranteed to be the same as the excitation energies of the two systems calculated separately (see Fig. 7).This may become a severe limitation when QSE will be applied to larger systems in the future and the underlying groundstate wavefunction is imprecise.
In search of a size-intensive alternative, the EOM formalism based qEOM method was proposed by Ollitrault et al. [77] for electronic excitation energies (EEs).The qEOM method provides good agreement for EEs with the exact results obtained by the full configuration interaction (FCI) method.However, the qEOM formalism (in Ref. [77]) does not necessarily satisfy the vacuum annihilation condition (VAC), also known as the killer condition [80,81], which ensures that the ground-state wavefunction cannot be de-excited.This may result in the appearance of large errors when the formalism is extended to calculate properties like IPs and EAs.Moreover, the qEOM method, just like QSE, requires higherbody RDMs which significantly increases the measurement challenges.
In this work, we propose a generally applicable EOMbased formalism for the calculation of molecular properties like EEs, IPs, and EAs following a VQE ground-state calculation on a quantum computer.Our formalism, which we refer to as q-sc-EOM, satisfies the VAC; produces size-intensive and real energy differences between the ground state and the excited states/charged states; does not involve measurements of higher than 2-body RDM-type quantities; is expected to be more resilient to noise than the current diagonalization-based state-ofthe-art methods.
This paper is organized as follows: Section II A discusses the theoretical formalism of q-sc-EOM using selfconsistent operators, while the implementation details and circuit design are explained in Sec.II B. Section III provides the computational details for the simulation data in this paper.In Sec.IV, we discuss the results obtained in this work.Specifically, Sec.IV A analyses the performance of the q-sc-EOM method in calculating EEs, IPs, and EAs of H 2 , LiH, and H 2 O molecules, while Section IV B compares the performance of the q-sc-EOM method with the QSE and qEOM formalisms.Finally, the key conclusions from the paper are summarized in Sec.V.

II. THEORY
The excitation energy of a given excited state can be obtained by the action of the commutator of the Hamiltonian and the corresponding excitation operator acting on the exact ground-state wavefunction.For an arbitrary k th excited state, this can be expressed as where |Ψ gr is the ground-state wavefunction, and E gr and E k refer to the energies of the ground and the k th excited state, respectively.Ĥ is the molecular Hamiltonian, which in the second quantization formalism can be written as where h pq and pq||rs are the one-and two-electron elements of the Hamiltonian, respectively.â † p and âp refer to the fermionic creation and annihilation operators (with respect to physical vacuum), respectively.Following common notations, here, we use indices {p, q, r, s . . .} for arbitrary molecular orbitals while {a, b, . . .} and {i, j, . . .} refer to unoccupied and occupied orbitals, respectively in the Hartree-Fock (HF) wavefunction.The state-transfer operator Ôk is defined as where |Ψ k is the wavefunction of the k th excited state.These operators should ensure an important property, referred to as vacuum annihilation or Killer condition, which expresses that the ground state cannot be deexcited [80][81][82][83][84][85]: In the case of exact operators O † k and an exact groundstate (Ψ gr ), the above equation can be reformulated as where Ô † k = |Ψ gr Ψ k |.It can be seen that this condition is automatically satisfied for exact state-transfer operators acting on an exact ground-state due to the orthogonality of eigenstates.However, one needs to ensure that the VAC is satisfied when approximate statetransfer operators are used [80,85].The state-transfer operators that satisfy the VAC are referred to as the "selfconsistent" operators [80].
This work utilizes the framework of unitary coupledcluster (UCC) theory, where the ground-state wavefunction is given by where |Ψ 0 is the HF wavefunction, and e σ is a unitary operator.An approximate form of σ, which is a cluster operator, can be written using single and double excitations as where σ a i and σ ab ij are amplitudes for the associated excitations.In order to capture the dominant electroncorrelation effects contributing to electronic excitations, an excitation manifold can be constructed by including all possible single and double excitations and deexcitations, represented by { Ĝ † I } ∪ { ĜI }.Here, ĜI can refer to any single (â † a âi ) or double (â † a â † b âj âi ) excitation operator.However, to satisfy the VAC, the excitation operator manifold is rotated, forming the self-consistent manifold { Ŝ † I } ∪ { ŜI }, where ŜI = e σ ĜI e −σ .
Similar excitation manifolds can also be constructed using particle number non-conserving excitation operators that are needed for computation of IPs and EAs.This technique was developed by Mukherjee and co-workers in Refs.[80,86].This self-consistent operator manifold can now be used to develop excited-state methods that satisfy the VAC.Following Eq. ( 8), the state-transfer operator for electronic excitations, ( Ôk ) EE , can now be written as a linear combination of all possible operators from the selfconsistent excitation manifold, given by where ( Rk 1 ) EE and ( Rk 2 ) EE are single and double excitation operators defined as Here, (A k ) I and (B k ) I † are the amplitudes corresponding to the excitation (I) and de-excitation operators (I † ), respectively, for the k th excited state.Here I refers to all possible single and double excitations.State-transfer operators for singly charged states can also be defined in a similar manner, where ( Rk 1 ) IP/EA and ( Rk 2 ) IP/EA refer to the particle number non-conserving single and double excitation operators, defined as and A. q-sc-EOM method In the VQE algorithm, the unitary evolution operator U (θ) (e σ in UCC theory) is implemented on a quantum computer using a parameterized circuit.The parameters (θ) of the circuit are optimized to variationally minimize the molecular energy and obtain the molecular ground state |Ψ VQE , such that Projecting Eq. ( 1) onto the k th excited state wavefunction and using the state-transfer operator defined in Eq. ( 9), the q-sc-EOM excitation energy from the ground state to the k th excited state (E 0k ) is given by Expressions for IPs and EAs can also be derived using the associated state-transfer operators defined in Eq. (11).As discussed in Ref. [77], Eq. ( 15) provides sizeintensive energy differences.By inserting the expression for ( Ôk ) EE from Eq. ( 9) in Eq. ( 15), it can be seen that the final equation for the excitation energy for the k th excited state has a parametric dependence on the amplitudes (A k ) I and (B k ) I † , where I refers to all possible single and double excitations.A variational minimization of the resulting equation (δE 0k = 0) with respect to these amplitudes leads to the following secular equation where the matrix elements of matrices M, Q, V, W are defined as Upon careful inspection, one can see that the matrices W and Q are zero due to the use of self-consistent operators.This is a simplification (compared to the qEOM formalism of Ref. [77]) as it reduces the secular equation to The matrix elements of M and V can be further simplified as and Thus, in q-sc-EOM the overlap matrix (V) is guaranteed to be the identity matrix and does not need to be computed on a quantum computer.It has an two key benefits, namely, it leads to a Hermitian formulation for excitation energy and it converts the generalized eigenvalue equation into an eigenvalue equation.The formalism developed so far can be written in the form of a concise eigenvalue equation as where S and D refer to single and double excitations, respectively.Thus, M SS refers to the block of matrix M with two single excitation operators in the double commutator (see Eq. ( 17)), while I XX is an identity matrix of dimension X.It should be emphasized that, since the above formulation is Hermitian, the eigenvalues obtained using Eq. ( 21) are guaranteed to be real (unlike in EOM-CCSD or qEOM) [45,77,87].This also ensures that this formulation is free from problems related to different left and right eigenfunctions encountered in classical EOM-CC methods [26].Finally, each element of matrix M can be computed on a quantum computer using Eq.(19).The resulting eigenvalue equation can then be solved classically to obtain q-sc-EOM EE values.Here, Krylov subspace based formalisms, such as the Davidson algorithm [88][89][90], can be used to obtain the excitation energies of a few lowlying excited states while avoiding the explicit evaluation and diagonalization of the M matrix.It should be noted that this method is closely related to UCC-based excitedstate methods in quantum chemistry, and thus Eq. ( 21) resembles the equation for EEs for UCC-based methods as derived, for example, in Refs.[84,87].

B. Circuit design and implementation details
Here, we discuss our proposed implementation of the q-sc-EOM formalism on a quantum computer.The discussion can be divided into two parts: circuit details to evaluate the diagonal and the off-diagonal elements of the matrix M in Eq. (19).The state preparation for the diagonal elements involves the same circuit as the one optimized for the ground state, but it is now applied to a classical state that represents an excited Slater |0i X   determinant.We refer to this classical state as the reference state.Thus, the circuit can be prepared in two steps: the first step is the creation of a reference state, whereas the second step involves the action of the previously optimized ground-state VQE circuit on the newly formed reference state.The molecular Hamiltonian is then measured using this prepared state, as done for the ground-state energy evaluation.To give an example of a state preparation circuit, we consider the H 2 molecule in the STO-3G basis and use a singly excited determinant in Eq. (19).We choose the single-excitation operator for the excitation from 1s β to 2s β (represented by |0011 → |1001 notation in the qubit representation using the Jordon-Wigner mapping).The classical state that corresponds to such an excited Slater determinant can be created by the action of two NOT gates, as shown in the circuit in Fig. 1a.U(θ) in Fig. 1a refers to the optimized circuit prepared for the VQE ground-state evaluation, and the shaded region represents the circuit needed to create the reference state.Note that the orbitals 1s α , 1s β , 2s α , and 2s β are mapped onto the qubits in the bottom-to-top order, such that the lowest energy orbitals are at the bottom.The off-diagonal elements of the matrix M can be evaluated using an entangled state of two excited Slater determinants, as shown below.A representative offdiagonal element M I,J can be written in terms of diagonal elements as where the term M I+J,I+J is given by A similar expression holds for the imaginary part of M I,J .
The matrix elements M I,I and M J,J in Eq. ( 22) are diagonal elements that are evaluated using the method described previously (without the ground-state energy term).The element M I+J,I+J can be evaluated using Eq. ( 23), which involves the creation of the I + J state, application of the unitary U (θ), and finally, the measurement of the Hamiltonian using this prepared state.The state I + J is a superposition state of two classical states I ( ĜI |Ψ 0 ) and J ( ĜJ |Ψ 0 ).Notice that both of these states are excited Slater determinants, which can be represented through qubit states in the computational basis.
An entangled state, such as I + J, is commonly created by using an ancilla qubit (for example, see Ref. [13]).
In the case of q-sc-EOM, we can use a simpler method to create this entanglement without adding any ancilla qubit.Being classical states, I and J are trivial to create on a quantum computer using NOT gates.An entangled state can then be created using a Hadamard gate and a series of CNOT gates.For example, the I state can be created using NOT gates and then transformed to the I + J state using Hadamard and CNOT gates.We can take the example of two single excitations to demonstrate this, specifically 1s α → 2s α and 1s β → 2s β (denoted as |0011 → |1001 and |0011 → |0110 in the qubit representation, respectively).Here, we need to create the entangled state |I + J = √ (|1001 + |0110 ).The circuit for creating this state and then evaluating M I+J,I+J is shown in Fig. 1b.The shaded region represents the circuit used to create the entangled |I + J state.The gate sequences can be stored for each excitation operator and can be applied to the HF state at the start of the EOM procedure.A maximum of CNOT gates will be needed to create any entangled state for the off-diagonal terms, when using up to double excitations in Eq. ( 19).This is a very small number compared with the total number of CNOT gates required to prepare the ground state using VQE-based algorithms.It should be noted that this proposed circuit design loses on the phase factor of excitation operators in I and J classical states, but all our numerical simulations have consistently shown that these phase factors have no effect on the computed eigenvalues of the eigenvalue equation.To preserve the phase factor, an ancilla qubit based creation of entangled I + J states could be used instead.
Regarding the resource requirements for running the qsc-EOM method on a quantum computer, the number of qubits required for the evaluation of each matrix element is the same as that in the computation of the ground state.The circuit that is implemented in the groundstate VQE calculation remains unchanged in the generation of the EOM matrix elements as well.The only difference is that the reference state is changed from the HF determinant to an excited Slater determinant or an entangled state involving two excited Slater determinants, as discussed above.Thus, unlike qEOM and QSE, where the excitation operators are measured together with the Hamiltonian, the excitation operators act directly on the HF state in q-sc-EOM.This provides a notable advantage that q-sc-EOM methods do not need the estimation of higher than 2-body RDMs.Such a framework also makes it easier to include higher excitations when required, such as triples, whose inclusion becomes important in higher-accuracy charged-state calculations [91].The calculation of elements of the M matrix can also be run in parallel on a quantum computer.The evaluation of the M matrix on a quantum computer requires the measurement of O(o 4 v 4 ) matrix elements, where o and v are the number of occupied and unoccupied orbitals, respectively.Despite this scaling, generally, the matrix M is very sparse.The number of elements in the M matrix can be drastically reduced using this sparsity through the use of spatial symmetry, spin symmetry, etc., which will be a topic of future study.The major advantages afforded by the use of quantum computers come from the efficient evaluation of the q-sc-EOM matrix elements and from the accurate ground state wavefunction provided by the VQE-based algorithm.

III. COMPUTATIONAL DETAILS
All the computations in this work utilize the STO-3G basis set.One-and two-electron integrals are calculated using the PySCF [34] program with the HF reference state.The Jordon-Wigner mapping and the transformation of the second-quantized operators to the Pauli form are carried out using the OpenFermion [92] software package.A classical noise-free simulator, where exact unitary operations are applied to the state vector representing the wavefunction, is used to assess the accuracy of the formalism developed in this work.The ground-state wavefunction is calculated using the fermionic ADAPT-VQE method using the generalized singles and doubles (GSD) operator pool [11].We use the gradient convergence criterion with a threshold of 10 −3 for all the ground-state energy calculations.All the formalisms discussed in this work, namely q-sc-EOM, qEOM, and QSE, utilize the ground-state energy and wavefunction obtained from the ADAPT-VQE simulation.It should be noted that we have extended the qEOM formulation of Ref. [77], originally developed for the calculation of EEs, to calculate IPs and EAs for our theoretical investigation.The EE results obtained using the qEOM approach in this work are verified against Qiskit's qEOM implementation [93].The code used for generating the data in this work can be found in Ref. [94].

IV. RESULTS AND DISCUSSION
We test the accuracy of the q-sc-EOM approach for three small molecules: H 2 , LiH, and H 2 O, and compare the results obtained with the exact FCI values.Computations are carried out for EEs, IPs, and EAs for these molecules.The total energies of the electronically excited, single electron-detached and single electronattached states are computed by adding the EEs, IPs and EAs, calculated using q-sc-EOM, to the ground-state energy obtained using the ADAPT-VQE simulations.Since the ground-state energy computed using ADAPT-VQE is near-exact for the molecules considered in this study, the errors in the energies of the electronically excited states and the single electron attached/detached states, with respect to the FCI, are almost entirely due to the error in the post-VQE procedure.For LiH and H 2 O, we invoke the frozen-core approximation.Thus, the number of qubits required for the q-sc-EOM computation for H 2 , LiH, and H 2 O are 4, 10 and 12, respectively.We plot energy errors relative to the FCI values and use shading to indicate errors below 0.1 eV, as this value is generally the desired accuracy for these properties, so that they can be quantitatively compared with experiments.We also compare the performance of the q-sc-EOM formalism with QSE and qEOM in subsection IV B.
A. EE, IP, and EA calculations with q-sc-EOM In Fig. 2, we show the energies of the ground state along with first few electronically excited, single electrondetached and single electron-attached states of the H 2 molecule calculated using q-sc-EOM, as a function of the inter-hydrogen distance.The corresponding FCI results are shown as gray lines.The errors in the energies with respect to the exact FCI values are shown in the subgraph on top of each panel.It can be seen that the errors in q-sc-EOM computed energies, or in other words, q-sc-EOM evaluated EEs, IPs, and EAs are essentially identical to FCI, with errors of less than 10 −8 Hartree.This is because the q-sc-EOM formalism for the H 2 molecule using the STO-3G basis set is exact, as the singles and doubles excitations used in Eq. ( 19) span the full excitation space and the VAC is satisfied.Figure 3 shows the energies of ground state along with the first few electronically excited, single electron-detached and single electron-attached states for the LiH molecule as a function of the Li-H bond length in a similar layout as the previous figure.For both EEs and IPs, the q-sc-EOM method gives nearly exact results, and as expected, errors less than 10 −8 Hartree were obtained with respect to the FCI values.However, we start to see the appearance of non-negligible errors for the EA results.This is because the computation of EAs for the LiH molecule with q-sc-EOM is no longer exact due to the addition of an electron.Thus, triple excitation operators should be added to the excitation manifold for an exact treatment for EAs.However, q-sc-EOM is still able to produce EAs within the desired accuracy, as seen from the shaded region in the error plot at the top of Fig. 3c.For the H 2 O molecule, we investigate the performance of q-sc-EOM as a function of the O-H bond length where both O-H bonds are stretched symmetrically.From Fig. 4, one can see that the errors in EEs are within the desired accuracy up to the O-H bond length of 1.4 Å.The errors build up as the two O-H bonds are broken further, due to the appearance of strong correlation effects in the wavefunction.Classical EOM-based methods, such as EOM-CCSD, show similar trends in errors as well.The errors in IPs and EAs are larger compared to EEs and are above our desirable error limit of 0.1 eV.This is well-known in classical quantum chemistry, where EOM-based methods require at least an approximate treatment of triple excitations in the EOM framework to reach an accuracy of 0.1 eV relative to FCI values for IPs and EAs [91,95,96].For higher accuracy in charged excitations, carefully selected triple excitations can be added to the excitation manifold in q-sc-EOM.This will be a subject of future study.

B. Comparison with QSE and qEOM
In this subsection, we discuss the connections between the q-sc-EOM, QSE, and qEOM formalisms and compare them in the contexts of a) their accuracy in computing energy differences, and b) quantum resource requirements and sensitivity to noise.The EOM formulation in Eq. ( 15) with excitation operators taken from the excitation manifold represented as { Ĝ † I } ∪ { ĜI } (not the self-consistent excitation manifold) may lead to a violation of the VAC.There are two methods discussed in the quantum chemistry literature to impose the VAC in the EOM formalism: the projected operator approach (see Ref. [81]) and the self-consistent operator formalism (see Ref. [80]).If we use projected operators based on Ref. [81], we arrive at a QSE-type formalism (as done in Ref. [97]) to calculate excited-state properties.The q-sc-EOM formalism developed in this work is based on the use of the self-consistent excitation operators to impose the VAC.It should be noted that this concept has been utilized in the development of different excited-state methods in classical quantum chemistry [87,[98][99][100].(grey lines) and H2-H4 system with H2 not interacting with H4 (H2 and H4 separation taken as 100 Å) (triangles) using STO-3G basis set.An inexact ground state is used in the plot computed using ADAPT-VQE stopped after adding 3 operators.The difference in the dashed line with triangles represent the size-intensivity error for each excitation energy using QSE.

Accuracy of energy differences
The qEOM method, which uses the EOM formulation in Eq. ( 16) with the conventional excitation manifold, { Ĝ † I } ∪ { ĜI } (details can be found in Ref. [77]), may lead to large errors in calculated IPs and EAs due to the violation of the VAC.This can be seen through the qEOM-evaluated IPs and EAs added to almost exact ADAPT-VQE ground state of the H 2 molecule in the STO-3G basis in Fig. 5. Large deviation from FCI results can be observed in this image.It should be noted that although single and double excitations span all the possible excitations in the case of H 2 molecule in STO-3G basis, the VAC is still not satisfied in qEOM.This is because the excitation manifold is not complete with respect to the exact ground-state when we use the { Ĝ † I } ∪ { ĜI } operator manifold corresponding to electron detached/atached states.One way to solve this issue is by increasing the size of the operator manifold.However, this would significantly increase the computational cost.Figure 6 shows the energies of the three lowest excited states evaluated using qEOM and q-sc-EOM methods for a rectangle geometry H 4 molecular system as a function of H 2 • • • H 2 separation distance.The H-H bond distance is fixed at 1.5 Å in the two H 2 molecules.We observe from Fig. 6 that q-sc-EOM is more robust in strongly correlated situations, compared with the qEOM method.The first and third excited states of H 4 computed using the qEOM method show qualitatively wrong behaviour in the region shown in the figure .The EEs obtained using QSE are not necessarily sizeintensive for an inexact ground state.This is due to the inclusion of the identity operator in the operator manifold used in QSE.We illustrate this point using an H 2 • • • H 4 molecular system as an example.Fig 7 shows the difference in EEs computed for an isolated H 2 molecule and an H 2 • • • H 4 molecular supersystem with no interaction between the H 2 and H 4 subsystems (the distance between H 2 and H 4 is taken as 100 Å).An inexact ground state is taken using an ADAPT-VQE simulation that is stopped after adding 3 operators.An identity operator is added to the operator manifold of QSE which uses the operator manifold represented by { ĜI } in Section II A. QSE computations give an error of ∼81 mH in this test, while the q-sc-EOM method shows the correct behavior.In this scenario, the two EEs should be identical for a method that provides size-intensive EEs.The magnitude of this size-intensivity error in QSE will depend on the accuracy of ground state.Since it is expected that near-term quantum computers may not provide exact ground-state energies for all molecular cases, these size-intensivity errors may cause problems.We note here that, just like q-sc-EOM, qEOM method provides sizeintensive EEs, IPs and EAs as well.

Noise-sensitivity
Along with the above mentioned thoeretical benefits, q-sc-EOM also has important computational advantages.q-sc-EOM is expected to be more noise-resilient compared with qEOM and QSE.A prime reason for this is because q-sc-EOM requires upto 2-RDM type measurements while QSE and qEOM require measurement of higher-body RDMs.Estimation of higher-body RDMs Error in excitation energy (E h ) b) QSE (with exact overlap (V)) vs q-sc-EOM QSE with exact V (State 1) QSE with exact V (State 2) QSE with exact V (State 3) q-sc-EOM (State 1) q-sc-EOM (State 2) q-sc-EOM (State 3) FIG. 8. (a) Error in the first few excitation energies computed for H4 using the QSE versus q-sc-EOM frameworks plotted against the maximum magnitude of errors in the matrix elements before diagonalization.(b) Error in first few excitation energies computed for H4 using the QSE (with no error in overlap matrix (V ) versus q-sc-EOM frameworks plotted against the maximum magnitude of errors in the matrix elements before diagonalization.A rectangle geometry H4 molecular system is used for the calculations with 1.5 Å and 2.0 Å as H-H distances.This shows that most of the noise sensitivity in QSE arises from noise in the computed overlap matrix, while the overlap matrix in q-sc-EOM is exactly the identity matrix and thus noise-free.
can significantly increase the noise in measured matrix elements.This has been noted in Ref. [101] and shown using simple noise models that error in expectation value will scale by a factor exponential in the number of qubits measured together.q-sc-EOM strictly requires measurement of upto 2-RDMs or 4 qubit measurements, while QSE and qEOM require measurements of upto 12 qubits at a time.
Additionally, we have also carried out a theoretical study based on matrix perturbation theory to compare the noise-resilience of QSE with q-sc-EOM on a noisy quantum computer.The analysis is carried out for an interacting pair of H 2 molecules (H 4 structure with a rectangular geometry with 1.5 Å and 2.0 Å bond distances).We model the effect of noise by adding random errors directly to each matrix element utilized by the two methods.It should be noted that in real implementations where these matrix elements are measured on quantum computers, there will be multiple sources of errors (gate errors, measurement errors, etc.), which will finally create a net error in each matrix element.In our noise study, it is reasonable to add errors directly as we assume that all matrix elements, both in QSE and q-sc-EOM, will have errors of same magnitude when the same quantum resources are used.This is an appropriate (and rather conservative) assumption for two reasons: First, all measurements in both q-sc-EOM and QSE utilize the same ground-state quantum circuit (with the difference being that in q-sc-EOM the circuit is applied on different reference states).Second, q-sc-EOM makes use of at most 2-body RDMs, while QSE requires the estimation of higher-body RDMs as well, which is expected to generate errors of same or higher magnitude in the case of QSE.
A more detailed noise analysis performed on an actual quantum device will be the subject of future studies.
In Fig. 8a, we show the performance of QSE compared with q-sc-EOM, where to each matrix element we add random offsets sampled from a uniform distribution inline with matrix perturbation theory.The horizontal axes in the figure correspond to the upper bound of this distribution (i.e., the maximum allowed error).The EEs are calculated after solving the eigenvalue equation of qsc-EOM and the generalized eigenvalue equation of QSE (the latter resembles Eq. ( 18)).Errors are defined with respect to the EE values obtained from the associated method in the absence of noise.Figure 8a shows the error averaged over 100,000 calculations of the EEs with different random offsets in each case.We can observe in the figure that QSE produces much larger errors at the same level of noise compared to q-sc-EOM, thus showing that we can expect q-sc-EOM to be more resilient to noise than QSE.
A key difference to consider between QSE and q-sc-EOM is that the overlap matrix V in q-sc-EOM is exactly known to be the identity matrix and thus does not need to be measured, while in the case of QSE, the overlap matrix must be measured on the quantum computer.Note that the latter is also true for qEOM.Figure 8b shows that if the overlap matrix is exactly known in QSE (i.e., we compute it exactly without noise), the noise-resilience of QSE is similar to that of q-sc-EOM.This suggests that the knowledge of the exact overlap matrix in q-sc-EOM is critical in providing noise-resilience, whereas methods that measure the overlap matrix and thus solve the generalized eigenvalue problem (such as QSE and qEOM) are expected to be more sensitive to noise.This noise-sensivity, we believe, is a direct result of a noise-sensitive matrix inversion step in solving generalized eigenvalue problem.A detailed analysis of this problem for general quantum algorithms for ground and excited state estimation will be presented in a future work.

V. CONCLUSION
In this work, we propose a new method, named q-sc-EOM, for calculating molecular excitation energies using a quantum computer.The method can be implemented on top of any quantum variational algorithm used to obtain the ground state of the target molecule.Our approach is inspired by excited-state methods developed in quantum chemistry, specifically the ones based on the UCC theory.q-sc-EOM has several important benefits compared to current state-of-the-art excited-state methods for NISQ era, with theoretical benefits including (a) q-sc-EOM uses the self-consistent operators that satisfy the vacuum annihilation condition, thus it can be generalized to evaluate accurate vertical excitation energies, ionization potentials, and electron affinities; (b) Energy differences obtained using q-sc-EOM are strictly sizeintensive, an important property that ensures their correct scaling with the size of the molecular system; (c) q-sc-EOM is a Hermitian theory, providing guaranteed real energy differences.The major computational benefit of q-sc-EOM is that it can be expected to more resilient to noise because (a) q-sc-EOM does not require higher than 2-body RDMs; and (b) it requires classical solution of eigenvalue equation, bypassing the noise-sensitive step associated inversion of overlap matrix in solving generalized eigenvalue equation in QSE and qEOM.These benefits provide important theoretical and practical advantages for the computation of excitation energies on near-term quantum devices.
NISQ era devices are expected to be noisy with limited resources.Thus, to achieve an advantage through these devices in quantum chemistry problems over classi-cal computation, one should use methods that are meaningfully accurate, while at the same time resistant to errors and resource-efficient.The q-sc-EOM method proposed in this work is promising in this regard because it exhibits many of the crucial properties of the highly successful EOM-based quantum chemistry methods.At the same time, it remains resource-efficient and is expected to be resilient to noise compared to the currently available diagonalization-based methods.Our future studies will combine q-sc-EOM with the recently developed transcorrelated Hamiltonian formalism [102,103] to obtain quantitatively accurate excited-state properties with minimal utilization of quantum resources, which otherwise generally requires the use of large basis sets.
e e m 3 7 J 0 p X s q F 1 H D o g 5 v H m j W b 2 7 U a 5 4 B o 9 7 5 u 7 0 7 h z 9 9 7 9 3 Q f N h 4 8 e P 9 l r 7 T + 9 0 u J s 8 p 6 7 e 6 u a c P K f 4 5 / w T / g C t w R a 8 e k 6 Q c w k j X P b 2 f 2 z Y x m o 5 w z p T 3 v 6 4 q 7 2 r p x 8 9 b a 7 f a d u / f u P

FIG. 1 .
FIG.1.Proposed quantum circuit for the estimation of a representative element of the M matrix for the H2 molecule using (a) an excited Slater determinant as the reference state needed to compute diagonal elements and (b) an entangled state involving two excited Slater determinants as the reference required for the evaluation of off-diagonal elements.

FIG. 2
FIG. 2. q-sc-EOM energies of (a) electronically excited, (b) single electron-detached, and (c) single electron-attached states along with the ground state (black circles) of the H2 molecule as a function of bond length.The gray lines correspond to the FCI results.The deviations from the FCI results are shown in the upper panel, where the shaded region indicates errors below 0.1 eV.

FIG. 4
FIG. 4. q-sc-EOM energies of first few (a) electronically excited, (b) single electron-detached, and (c) single electron-attached states along with the ground state (black circles) of the H2O molecule as a function of the O-H bond lengths where both O-H bonds are stretched symmetrically.The gray lines correspond to the FCI results.The deviations from the FCI results are shown in the upper panel, where the shaded region indicates errors below 0.1 eV.

FIG. 5 .
FIG. 5. Energies of (a) single electron-detached, and (b) single electron-attached states along with the ground state (black circles) of the H2 molecule plotted as a function of the H-H bond length using the STO-3G basis.The FCI results are plotted in gray.

FIG. 6 .FIG. 7 .
FIG.6.Energies of the three lowest excited states along with the ground state (black circles) computed using a) q-sc-EOM and b) qEOM formalisms for the dissociation of a rectangular geometry of H4 into H2 • • • H2 as a function of the H2 • • • H2 separation distance.Both of the H2 molecules have a bond length of 1.5 Å.The FCI results are plotted in gray.