A Full Quantum Eigensolver for Quantum Chemistry Simulations

Quantum simulation of quantum chemistry is one of the most compelling applications of quantum computing. It is of particular importance in areas ranging from materials science, biochemistry, and condensed matter physics. Here, we propose a full quantum eigensolver (FQE) algorithm to calculate the molecular ground energies and electronic structures using quantum gradient descent. Compared to existing classical-quantum hybrid methods such as variational quantum eigensolver (VQE), our method removes the classical optimizer and performs all the calculations on a quantum computer with faster convergence. The gradient descent iteration depth has a favorable complexity that is logarithmically dependent on the system size and inverse of the precision. Moreover, the FQE can be further simplified by exploiting a perturbation theory for the calculations of intermediate matrix elements and obtaining results with a precision that satisfies the requirement of chemistry application. The full quantum eigensolver can be implemented on a near-term quantum computer. With the rapid development of quantum computing hardware, the FQE provides an efficient and powerful tool to solve quantum chemistry problems.


Introduction
Quantum chemistry studies chemical systems using quantum mechanics. One primary focus of quantum chemistry is the calculation of molecular energies and electronic structures of a chemical system which determine its chemical properties. Molecular energies and electronic structures are calculated by solving the Schrodinger equation within chemical precision. However, the computational resources needed scale exponentially with the system size on a classical computer, making the calculations in quantum chemistry intractable in high dimension.
Quantum computers, originally envisioned by Benioff, Manin, and Feynman [1][2][3], have emerged as promising tools for tackling this challenge with polynomial overhead of computational resources. Efficient quantum simulations of chemistry systems promise breakthroughs in our knowledge for basic chemistry and revolutionize research in new materials, pharmaceuticals, and industrial catalysts.
The universal quantum simulation method [4] and the first quantum algorithm for simulating fermions [5] have laid down the fundamental block of quantum chemistry simulation. Based on these techniques and quantum phase estimation algorithm [6], Aspuru-Guzik et al. presented a quantum algorithm for preparing ground states undergoing an adiabatic evolution [7], and many theoretical and experimental works [8][9][10][11][12][13][14][15][16][17][18][19][20][21][22][23][24] have been developed since then. In 2002, Somma et al. proposed a scalable quantum algorithm for the simulation of molecular electron dynamics via Jordan-Wigner transformation [25]. The Jordan-Wigner transformation directly maps the fermionic occupation state of a particular atomic orbital to a state of qubits, which enables the quantum simulation of chemical systems on a quantum computer. Then, the Bravyi-Kitaev transformation [26][27][28][29][30] encodes both locality of occupation and parity information onto the qubits, which is more efficient in operation complexity. In 2014, Peruzzo et al. developed the variational quantum eigensolver (VQE) [18,31], which finds a good variational approximation to the ground state of a given Hamiltonian for a particular choice of ansatz. Compared to quantum phase estimation and Trotterization of the molecular Hamiltonian, the VQE requires a lower number of controlled operations and shorter coherence time. However, the VQE is a classical and quantum hybrid algorithm; the optimizer is performed on a classical machine.
Meanwhile, implementations of quantum chemistry simulation have been developing steadily. Studies in present-day quantum computing hardware have been carried out, such as nuclear magnetic resonance system [32,33], photonic system [34][35][36], nitrogen-vacancy center system [37], trapped ion [38,39], and superconducting system [40][41][42]. Rapid development in quantum computer hardware with even the claims of quantum supremacy greatly stimulates the expectation of its real applications. Quantum chemistry simulation is considered a real application in Noisy Intermediate-Scale Quantum (NISQ) computers [22,43,44]. The FQE is an effort on this background. In the FQE, not only calculation of Hamiltonian matrix part is done on quantum computer, but also the optimization by gradient descent is performed on quantum computer. FQE can be used in near-term NISQ computers and in future fault-tolerant large quantum computers.

Preparing the Hamiltonian for Quantum Chemistry
Simulation. A molecular system contains a collection of nuclear charges Z i and electrons. The fundamental task of quantum chemistry is to solve the eigenvalue problem of the molecular Hamiltonian. The eigenstates of the manybody Hamiltonian determine the dynamics of the electrons as well as the properties of the molecule. The corresponding Hamiltonian of the system includes kinetic energies of nuclei and electrons and the Coulomb potentials of nuclei-electron, nuclei-nuclei, and electron-electron, and it can be expressed in first quantization as in atomic units ðℏ = 1Þ, where R i , Z i , M i and r i are the positions, charges, and masses of the nuclei and the positions of the electrons, respectively. Under the Born-Oppenheimer approximation which assumes the nuclei as a fixed classical point, this Hamiltonian is usually rewritten in the particle number representation in a chosen basis: where ⋯ denotes higher-order interactions and a † i and a j are the creation and annihilation operators of particle in orbital i and j, respectively. The parameters h ij and h ijkl are the onebody and two-body integrations in the chosen basis functions fψ i g. In Galerkin formulation, the scalar coefficients in Equation (2) can be calculated by In order to perform calculations on a quantum computer, we need to map fermionic operators to qubit operators. We choose the Jordan-Wigner transformation to achieve this task due to its straightforward expression.
The Jordan-Wigner transformation maps Equation (2) into a qubit Hamiltonian form: where Roman indices i and j denote the qubit on which the operator acts and Greek indices α and β refer to the type of Pauli operators; i.e., σ i x means Pauli matrix σ x acting on a qubit at site i. Apparently, H in Equation (2) is a linear combination of unitary Pauli matrices. The methods used in this paper finding the molecular ground state and its energy are all based on it.
In this work, we present the FQE to find the molecular ground-state energy by gradient descent iterations. Gradient descent is one of the most fundamental ways for optimization that looks for the target energy value along the direction of the steepest descent. Here, it is performed in a quantum computer with the help of linear combination of unitary operators. We analyze the relationships between the gradient descent iteration depth and the precision of the groundstate energy. The explicit quantum circuit to implement the algorithm is constructed. As illustrative examples, the ground-state energies and electronic structures of four molecules, H 2 , LiH, H 2 O, and NH 3 are presented. Taking H 2 O and NH 3 as examples, a comparison between the FQE and VQE, a representative hybrid method, is given. FQE can be accelerated further by harnessing a perturbation theory in chemical precision. Finally, we analyze the computation complexity of FQE and summarize the results.

Quantum Gradient Descent
Iteration. The classical gradient descent algorithm is usually employed to obtain the minimum of an target function f ðXÞ. One starts from an initial point X ð0Þ = ðx 0 1 , x 0 2 , ⋯, x 0 N Þ ∈ ℝ N , then moves to the next point along the direction of the gradient of the target function, namely, where γ 0 is a positive learning rate that determines the step size of the iteration. In searching the minimum energy of a Hamiltonian, the target function can be expressed as a quadratic optimization problem in the form f ðXÞ = X T HX. At where H g is a linear combination of unitary operators (LCU) which was proposed in [45] in designing quantum algorithms and studied extensively [46][47][48][49][50][51][52]. This nonunitary evolution can be implemented in a unitary quantum circuit by adding ancillary qubits that transform it into unitary evolution in a larger space [53]. The realization of LCU can be viewed as a quantum computer wave function passing through M-slits and operated by a unitary operation in each slit, and then, the wave functions are combined and the result of the calculation is read out by a measurement [49]. As shown in Figure 1, we perform the evolution described by Equation (9) with the following four steps. Wave division: the register is a composite system which contains a work system and an ancillary register. Firstly, the initial point X = ðx 1 , ⋯, x N Þ T is efficiently mapped as an initial state jX ðtÞ i of the work system. In quantum chemistry, Hartree-Fock (HF) product state is usually used as an initial state. And the ancillary register is initialized from j0i m , where m = log 2 M, to a specific superposition state jψ s i, where is a normalization constant and jii is the computational basis. This is equivalent to let the state jx ðtÞ i pass through M-slits. β i is a factor describing the properties of the slit, which is determined by the forms of the Hamiltonian in Equation (8). This can be done by the initialization algorithm in [54]. Moreover, the quantum random access memory (qRAM) approach can be used to prepare jx ðtÞ i and jψ s i, which consume Oðlog NÞ and Oðlog MÞ basic steps or gates, respectively, after qRAM cell is established. We denote the whole state of the composite system as jΦi = jψ s ijx ðtÞ i.
Entanglement: then, a series of ancillary systemcontrolled operations ∑ M−1 i=0 jiihij ⊗ H g i are implemented on the work qubits. The work qubits and the ancilla register are now entangled, and the state is transformed into The corresponding physical picture is that different unitary operations are implemented simultaneously in different subspaces, corresponding to different slits.
Wave combination: we perform m Hadamad gates on ancillary register to combine all the wave functions from the M different subspaces. We merely focus on the component in a subspace where the ancillary system is in state j0i. The state of the whole system in this subspace is Measurement: then, we measure the ancillary register. If we obtain j0i, our algorithm succeeds and we obtain the state iÞ, where the work system is in jx ðt+1Þ i = H g jx ðtÞ i. And then, this will be used as input for the next iteration in the quantum gradient descent process. The probability of obtaining j0i for the state is The successful probability after n measurements is 1 − ð1 − kH g jx ðtÞ ik 2 /ℂ 2 MÞ n , which is an exponential function of n. The number of measurements is ℂ 2 M/kH g jx ðtÞ ik 2 .
The measurement complexity will grow exponentially with respect to the number of iteration steps [55]. Alternatively, one can use the oblivious amplitude amplification [51] to amplify the amplitude of the desired term (ancillary qubits in state j0i) up to a deterministic order with Oð ffiffiffiffi ffi M p Þ repetitions before the measurement. Then, the measurement complexity will be the product of iteration depth k and Oð ffiffiffiffi ffi M p Þ, linearly dependent on the number of iteration steps. After obtaining j0i, we can continue the gradient descent process by repeating the above four steps, with jx ðtÞ i replaced by j x ðt+1Þ i in a wave-division step. We can preset a threshold defined as as criterion for stopping the iteration. Thus, we judge if the iterated state satisfies the criterion by measuring the expectation value of Hamiltonian around the expected number of iteration, which is easier than constructing the tomography. If the next iterative state jx ðt+1Þ i does not hit our preset threshold, this output jx ðt+1Þ i will be regarded as the new input state jx ðtÞ i and run the next iteration. Otherwise, the iteration can be terminated and the state jx ðt+1Þ i is the final result jx f i, as one good approximation of the ground state. The ground-state energy can be calculated by Measuring the expectation values during the iteration procedure will destroy the state of the work system, stopping the quantum gradient descent process. So, determining the iteration depth k in advance is essential. After k times iterations, the approximation error is limited to (ignoring constants) where λ 1 and λ 2 are the two largest absolute values of the eigenvalues of Hamiltonian H (see Supplemental Material (available here) for proof). The iteration depth is logarithmically dependent on the system size and the inverse of precision. The algorithm may be terminated at a point with a preset precision ε. It can be seen that the choise of γ has little impact on converge rate when γ is large. This makes this algorithm very robust to this parameter. The rate of convergence primarily depends upon the ratio of λ 1 and λ 2 . The gap between the iterative result and the ground state depends on the choice of initial point. If we choose an ansatz state with a large overlap with the exact ground state, the iterative process will converge to the ground state in fewer iterations. Usually, the mean-field state which represents a good classical approximation to the ground state of Hamiltonian H, such as a Hartree-Fock (HF) product state, is chosen as an initial state. Compared to the VQE, the FQE does not need to make measurements of the expectation values of Hamiltonian during each iteration procedure and this substantially reduces the computation resources.

Perturbation Theory.
The FQE involves multitime iterations to obtain an accurate result, which is difficult to implement in the present-day quantum computer hardware. Here, we present an approximate method to find the ground state and its energy by using the gradient descent algorithm and perturbation theory. The perturbation theory is widely used and plays an important role in describing real quantum systems, because it is impossible to find exact solutions to the Schrodinger equation for Hamiltonians even with moderate complexity. The Hamiltonian described by Equation (4) can be divided into two classes, H 0 and H′. H 0 consists of a set of Pauli terms containing only σ i α=z and the identity matrices, and Pauli terms σ i α=x,y belong to H′. H 0 is a diagonal matrix with exact solutions that can be regarded as a simple system. H ′ usually is smaller compared to H 0 and is treated as a "perturbing" Hamiltonian. The energy levels and eigenstates associated with the perturbed system can be expressed as "corrections" to those of the unperturbed system. We begin with the time-independent Schrodinger equation: where E n and jψ n i are the nth energy and eigenstate, respectively. Unperturbed Hamiltonian H 0 satisfies the timeindependent Schrodinger equation: H 0 jni = E ð0Þ n jni. Our goal is to express E n and jψ n i in terms of E 0 n and jni. Denote the expectation value of H ′ as H ′ nn = hn | H ′ | ni, and it is easy to see that hn | H ′ | ni is zero because H ′ only contains Pauli terms σ i α=x,y . In the first-order approximation, the energies and eigenstates are expressed as In the second-order approximation, they are The matrix elements in the firstand second-order approximations can be obtained by one iteration of the quantum circuit in Figure 1. Here, we let H ′ be equal to H g . Explicitly, the first-order approximation only involves H ′ mn , a series of transition probabilities of the state after H′ implemented on state jni, and they can be obtained by performing the quantum circuit of Figure 1 directly. For the second-order approximation, matrix elements such as value jH ′ mn j 2 and H ′ mn H ′ kn can be calculated by H ′ mn . Then, the approximate ground energy and ground state up to the second order are obtained. We will show the performance of the FQE and perturbation theory in the next section.  Figure 3), and the details are given in Supplemental Material. The converge rates of the four molecules depend on the system size and the ratio of the two largest absolute eigenvalues of the Hamiltonian H, which are consistent with the theoretical analysis above. We also studied the influence of noises which is also shown in Figure 2. The noise term is chosen in the form of  (1:6 × 10 − 3 a:u). This indicates that our method is robust to a certain type of noise, which is important in the implementation of quantum simulation on near-term quantum devices. For more noisy situations, see Supplemental Material for details, where the parameters of noise are 10 times the above values (see Figure 4). The convergence deteriorates and some oscillations occur as the number of iterations increases. In Figure 5, a comparison with the VQE is shown for H 2 O and NH 3 . In VQE calculation, the initial state jx 0 i is mapped to an ansatz state by a parameterized unitary opera-

Calculations of Four
We performed numerical simulations of the VQE for the two molecules. When the learning rate γ ≥ 10 −3 , the VQE does not converge to the ground state. So, in order to compare with each other, we choose the proper learning rates in two methods separately. In both cases, the initial ansatz state is prepared as the HF product state. In H 2 O and NH 3 , the VQE converges the fastest with the learning rate γ = 10 −3 . The FQE converges faster with larger and larger learning rate until a fixed speed is reached. As shown in Figure 5, the FQE generally converges faster than the VQE and the advantage will be more obvious in complex molecules.
The above examples are calculated in fixed interatomic distance of the molecules. If we want to calculate the interatomic distance corresponding to the most stable structure, the variation of interatomic distances is necessary. In Figure 6, four examples are given to illustrate the performance of the perturbation theory. To obtain the potential-energy surfaces for H 2 , LiH, H 2 O, and NH 3 molecules, we studied the dependence of ground-state energy of their molecules on the variating interatomic distances, between the two atoms in H 2 and LiH, the distance between the oxygen atom and one hydrogen atom (the two hydrogen atoms are symmetric with respect to the oxygen atom) in H 2 O, and the distance between the nitrogen atom and the plane formed by the three hydrogen atoms in NH 3 . The lowest energy in potential-energy surfaces corresponds to the most stable structure of the molecules. As shown in the picture, the ground-state energy of each molecule calculated under the second-order approximation is already quite close to its exact value, which is obtained from Hamiltonian diagonalizations. The energy values up to the secondorder correction are compared with their exact values at the most stable interatomic distance corresponding to the lowest        Figure 6: Theory results (blue lines) and zero-order (orange lines), first-order (green lines), and second-order (red lines) energy plots of outcomes from numerical simulations, for several interatomic distances for H 2 , LiH, H 2 O (between the oxygen atom and one hydrogen atom), and NH 3 (between the nitrogen atom and the plane formed by the three hydrogen atoms). respectively. In some applications, the perturbation theory only requires one iteration, and an approximate result in chemical precision can be obtained.

Summary
An efficient quantum algorithm, full quantum eigensolver (FQE), for calculating the ground-state wave function and the ground energy using gradient descent (FQE) was proposed, and numerical simulations are performed for four molecules. In the FQE, the complexity of basic gate operations is polylogarithmical to the number of single-electron atomic orbitals. It achieves an exponential speedup compared with its classical counterparts. It has been shown that the FQE is robust against noises of reasonable strengths. For very noisy situations that do not allow many iterations, the FQE can be combined with the perturbation theory that gives the ground state and its energy in chemical precision with onetime iteration. The FQE is exceptionally useful in quantum chemistry simulation, especially for the near-term NISQ applications. The FQE is a full quantum algorithm, not only applicable for NISQ computers but also directly applicable for future large-scale fault-tolerant quantum computers.

Data Availability
The data that support the findings of this study are available from the corresponding authors on reasonable request.

Conflicts of Interest
The authors declare no competing interests.