Guoming Wang

Daniel Stilck França

Gumaro Rendon

Peter D Johnson

Faster ground state energy estimation on early fault-tolerant quantum computers via rejection sampling

A major thrust in quantum algorithm development over the past decade has been the search for the quantum algorithms that will deliver practical quantum advantage first. Today's quantum computers and even early fault-tolerant quantum computers will be limited in the number of operations they can implement per circuit. We introduce quantum algorithms for ground state energy estimation (GSEE) that accommodate this design constraint. The first estimates ground state energies and has a quadratic improvement on the ground state overlap parameter compared to other methods in this regime. The second certifies that the estimated ground state energy is within a specified error tolerance of the true ground state energy, addressing the issue of gap estimation that beleaguers several ground state preparation and energy estimation algorithms. We note, however, that the scaling of this certification technique is, unfortunately, worse than that of the GSEE algorithm. These algorithms are based on a novel use of the quantum computer to facilitate rejection sampling. After a classical computer is used to draw samples, the quantum computer is used to accept or reject the samples. The set of accepted samples correspond to draws from a target distribution. While we use this technique for ground state energy estimation, it may find broader application. Our work pushes the boundaries of what operation-limited quantum computers are capable of and thus brings the target of quantum advantage closer to the present.

Introduction

Estimating the ground state energy of quantum many-body systems is a fundamental problem in condensed matter physics, quantum chemistry, material science and quantum information. It serves as a basic subroutine in drug discovery and materials design [START_REF] Aspuru-Guzik | Simulated quantum computation of molecular energies[END_REF][START_REF] Cao | Quantum chemistry in the age of quantum computing[END_REF][START_REF] Vincent E Elfving | How will quantum computers provide an industrially relevant computational advantage in quantum chemistry[END_REF][START_REF] Goings | Reliably assessing the electronic structure of cytochrome p450 on today's classical computers and tomorrow's quantum computers[END_REF]. In this problem, the system of interest is characterized as an n-qubit Hamiltonian H with unknown spectral decomposition H = j E j |E j E j |, where E 0 < E 1 ≤ E 2 ≤ ... are the eigenvalues of H and the |E j 's are the orthonormal eigenstates of H, and we need to estimate the lowest eigenvalue E 0 of H within additive error ǫ with high probability. Due to its significance, this problem has been extensively studied in the past decades, and various methods have been proposed to tackle it. Among these methods, two main approaches to ground state energy estimation (GSEE) have been explored.

GSEE assuming a fault-tolerance cost model

The first approach to GSEE is based on the quantum phase estimation (QPE) algorithm [START_REF] Yu | Quantum measurements and the abelian stabilizer problem[END_REF][START_REF] Michael | Quantum computation and quantum information[END_REF] and often has rigorous performance guarantees. In this approach, we assume that there exists an efficient procedure to prepare a state |ψ that has non-trivial overlap with the ground state of H (i.e. | E 0 |ψ | = Ω(1/ poly (n))). This assumption is reasonable in many practical settings. For example, Refs. [START_REF] Babbush | Chemical basis of trotter-suzuki errors in quantum chemistry simulation[END_REF][START_REF] Norm M Tubman | Postponing the orthogonality catastrophe: efficient state preparation for electronic structure simulations on quantum devices[END_REF][START_REF] Sugisaki | Quantum chemistry on quantum computers: A method for preparation of multiconfigurational wave functions on quantum computers without performing post-hartree-fock calculations[END_REF][START_REF] Mcardle | Quantum computational chemistry[END_REF] provide several methods for creating such states for quantum chemistry Hamiltonians. Given the "ansatz state" |ψ , we can perform standard QPE on the time evolution of H using |ψ as the initial state. Assuming p 0 := | E 0 |ψ | 2 , we can obtain an ǫ-accurate estimate of E 0 by using O(p -1 0 ) QPE circuits, each of which uses O log ǫ -1 p -1 0 ancilla qubits and implements a time evolution of H for a duration at most O(ǫ -1 p -1 0 ). Subsequently, the polynomial scaling of the algorithm runtime is O(ǫ -1 p -2 0 ). Several variants of this basic algorithm have been proposed in which the overall evolution time is reduced to Õ(ǫ -1 p -1 0 ) [START_REF] Knill | Optimal quantum measurements of expectation values of observables[END_REF][START_REF] Nagaj | Fast amplification of qma[END_REF][START_REF] Poulin | Sampling from the thermal quantum gibbs state and evaluating partition functions with a quantum computer[END_REF][START_REF] Dong | Ground-state preparation and energy estimation on early fault-tolerant quantum computers via quantum eigenvalue transformation of unitary matrices[END_REF], or the maximal evolution time is reduced to O ǫ -1 log p -1 0 [START_REF] Knill | Optimal quantum measurements of expectation values of observables[END_REF][START_REF] Nagaj | Fast amplification of qma[END_REF][START_REF] Poulin | Sampling from the thermal quantum gibbs state and evaluating partition functions with a quantum computer[END_REF][START_REF] Lin | Heisenberg-limited ground-state energy estimation for early faulttolerant quantum computers[END_REF][START_REF] Wan | Randomized quantum algorithm for statistical phase estimation[END_REF][START_REF] Dong | Ground-state preparation and energy estimation on early fault-tolerant quantum computers via quantum eigenvalue transformation of unitary matrices[END_REF], or the number of ancilla qubits is reduced to O(1) [START_REF] Berry | How to perform the most accurate possible phase measurements[END_REF][START_REF] Higgins | Entanglement-free heisenberg-limited phase estimation[END_REF][START_REF] Lin | Heisenberg-limited ground-state energy estimation for early faulttolerant quantum computers[END_REF][START_REF] Dong | Ground-state preparation and energy estimation on early fault-tolerant quantum computers via quantum eigenvalue transformation of unitary matrices[END_REF]. Despite these improvements, all of these algorithms require implementing a time evolution of H up to O(ǫ -1 ) time. This implies that the number of operations per circuit on the quantum computer is quite large when the tolerable error ǫ is small. As an example, [START_REF] Isaac H Kim | Faulttolerant resource estimate for quantum chemical simulations: Case study on li-ion battery electrolyte molecules[END_REF] estimated that the number of T gates required to implement ground state energy estimation to sufficient accuracy for a class of small molecules is greater than 10 10 ; the error rate of each gate must be significantly less than 10 -10 . Today's quantum hardware and that of the near future are far more limited in the number of gates that can be implemented per circuit. As an example, IBM recently announced their 100 × 100 challenge [START_REF] Gambetta | Quantum-centric supercomputing: The next wave of computing[END_REF], in which, by 2024, they aim to realize accurate expectation value calculations using circuits of 100×100 = 10 4 quantum gates, each of error rate 10 -3 . There is large consensus in the community that the only feasible approach to increasing the number of operations from 10 4 to greater than 10 10 is by building large-scale quantum architectures that can implement low-error-rate fault-tolerant operations.

The GSEE algorithms above were developed with a fault-tolerant quantum computer in mind. It is by now standard to assume a cost model where the overhead of reducing error rates (through fault-tolerant protocols) does not change the polynomial scaling of a quantum algorithm; that is, the cost is logarithmic in the number of operations. Under this model, it is usually favorable to reduce algorithm runtime at the cost of an increased number of operations per circuit. For this reason, methods for ground state energy estimation have been developed [START_REF] Lin | Near-optimal ground state preparation[END_REF] that increase the number of ancilla qubits and operations to reduce the overall runtime.

GSEE assuming a pre-fault-tolerance cost model The second approach to GSEE is the variational quantum eigensolver (VQE) [START_REF] Peruzzo | A variational eigenvalue solver on a photonic quantum processor[END_REF][START_REF] Aspuru-Guzik | Simulated quantum computation of molecular energies[END_REF][START_REF] Jarrod R Mcclean | The theory of variational hybrid quantum-classical algorithms[END_REF]. This approach was initially conceived assuming a cost model where fault-tolerant protocols are unavailable and thus gate error rates are fixed. For a given circuit error rate, this limits the maximum number of quantum operations per circuit. The method uses an ansatz circuit of limited number of operations to generate an approximation of the ground state of H. While this algorithm is amenable to near-term implementation, it lacks performance guarantees for two reasons. First, the accuracy of its output depends on the representation power of the chosen variational ansatz which is often hard to understand. Second, in order to find the optimal ansatz parameters, one needs to solve the associated non-convex optimization problem which could be challenging [START_REF] Bittel | Training variational quantum algorithms is NP-hard[END_REF]. Furthermore, since VQE needs to repeatedly evaluate the energy of the ansatz state (via direct sampling or more sophisticated methods [START_REF] Wang | Minimizing estimation runtime on noisy quantum computers[END_REF][START_REF] Dax Enshan Koh | Foundations for bayesian inference with engineered likelihood functions for robust amplitude estimation[END_REF]), this algorithm becomes too time-consuming to be practical for chemical system sizes of industrial relevance, as shown in a recent work [START_REF] Jérôme F Gonthier | Identifying challenges towards practical quantum advantage through resource estimation: the measurement roadblock in the variational quantum eigensolver[END_REF] and are more sensitive to noise than previously expected [START_REF] Stilck | Limitations of optimization algorithms on noisy quantum devices[END_REF][START_REF] De | Limitations of variational quantum algorithms: A quantum optimal transport approach[END_REF][START_REF] Quek | Exponentially tighter bounds on limitations of quantum error mitigation[END_REF][START_REF] Wang | Noise-induced barren plateaus in variational quantum algorithms[END_REF].

GSEE assuming a limited-number-of-operations cost model The two modes of ground state energy estimation algorithms introduced so far were developed assuming two very different cost models: logarithmic-overhead error reduction through fault-tolerance and fixed error rate.

These models correspond to future fault-tolerant architectures and today's near-term intermediate scale devices, respectively. These cost models may not be sufficient to describe the capabilities and limitations of early fault-tolerant quantum computers [START_REF] Fellous-Asiani | Limitations in quantum computing from resource constraints[END_REF]. This motivates the development of cost models that interpolate between these two regimes. While rigorously developing such cost models is beyond the scope of our work, we will assume a proxy for such a model: any early fault-tolerant quantum computer will be limited in the number of logical operations it can implement per circuit while keeping the total error rate below some amount. This motivates the development of quantum algorithms in which the number of operations can be reduced below that allotted by the device, while paying a cost of an increase in runtime.

Several recent works have been developed under similar premises. The first was the so-called α-VQE method [START_REF] Wang | Accelerated variational quantum eigensolver[END_REF] in which a variable-depth amplitude estimation is introduced. Later, Wang et al. introduced a variable-depth amplitude estimation algorithm that is robust to substantial amounts of error [START_REF] Wang | Minimizing estimation runtime on noisy quantum computers[END_REF]. Similar methods were explored in the context of quantum algorithms for finance [START_REF] Alcazar | Quantum algorithm for credit valuation adjustments[END_REF][START_REF] Giurgica-Tiron | Low depth algorithms for quantum amplitude estimation[END_REF]. Some of these methods have been implemented on quantum hardware [START_REF] Katabarwa | Reducing runtime and error in vqe using deeper and noisier quantum circuits[END_REF][START_REF] Giurgica-Tiron | Low-depth amplitude estimation on a trapped-ion quantum computer[END_REF]. Other algorithms have been introduced for early fault-tolerant quantum computers with the aim of reducing the number of ancilla qubits compared to their traditional counterparts [START_REF] Suzuki | Amplitude estimation without phase estimation[END_REF][START_REF] Lin | Heisenberg-limited ground-state energy estimation for early faulttolerant quantum computers[END_REF][START_REF] Zhang | Computing Ground State Properties with Early Fault-Tolerant Quantum Computers[END_REF].

Recently, Wang et al. [START_REF] Wang | Quantum algorithm for ground state energy estimation using circuit depth with exponentially improved dependence on precision[END_REF] applied this methodology to the task of ground state energy estimation. They asked whether there exists a method for ground state energy estimation that uses fewer operations per circuit than previous approaches (i.e. the merit of VQE), while admitting a rigorous performance guarantee (i.e. the merit of traditional GSEE methods). Inspired by the recent work of Lin and Tong [START_REF] Lin | Heisenberg-limited ground-state energy estimation for early faulttolerant quantum computers[END_REF], they showed that, given a lower bound ∆ on the spectral gap of H around E 0 , one can estimate E 0 within additive error ǫ with high probability by using Õ ∆ 2 ǫ -2 p -2 0 quantum circuits, where each circuit evolves H up to O ∆ -1 polylog ∆ǫ -1 p -1 0 time and requires a single ancilla qubit. Note that the maximal evolution time is only poly-logarithmic in 1/ǫ, which is exponentially better than previous methods. Furthermore, the parameter ∆ can be chosen anywhere between Θ(ǫ) and Θ(∆ true ), where ∆ true := E 1 -E 0 , and tuning it leads to a class of GSEE algorithms that smoothly interpolate between an algorithm with circuit depth O ∆ -1 true polylog ∆ true ǫ -1 p -1 0 and runtime Õ ∆ true ǫ -2 p -2 0 and an algorithm with circuit depth O ǫ -1 polylog p -1 0 and runtime Õ(ǫ -1 p -2 0 ). In other words, these algorithms smoothly transition from having low circuit depth to achieving Heisenberg-limit scaling in the runtime with respect to ǫ. Later, Ding and Lin [START_REF] Ding | Even shorter quantum circuit for phase estimation on early faulttolerant quantum computers with applications to ground-state energy estimation[END_REF] recovered this result using a different approach.

Although Refs. [START_REF] Wang | Quantum algorithm for ground state energy estimation using circuit depth with exponentially improved dependence on precision[END_REF] and [START_REF] Ding | Even shorter quantum circuit for phase estimation on early faulttolerant quantum computers with applications to ground-state energy estimation[END_REF] significantly reduce the circuit depth for GSEE, one aspect of their algorithms is somewhat unsatisfactory -their runtimes are quadratic in 1/p 0 . This means that if the overlap between the input state |ψ and the ground state |E 0 is small, these algorithm could take a long time to finish. Is it possible to improve the dependence on 1/p 0 to linear while preserving the nice features of these algorithms? In this work, we answer this question in affirmative. Precisely, we prove that, given a lower bound ∆ on the spectral gap of H around E 0 , one can estimate E 0 within additive error ǫ with high probability by using Õ(∆ 2 ǫ -2 p -1 0 ) quantum circuits, where each circuit has O ∆ -1 polylog ∆ǫ -1 p -1 0 depth and requires O(1) ancilla qubits. This algorithm works with two common input models of H: 1. The block-encoding model, in which H is embedded as a sub-matrix of a larger unitary operator; 2. The Hamiltonian-evolution model, in which a controlled time evolution of H can be efficiently implemented. Moreover, as in Ref. [START_REF] Wang | Quantum algorithm for ground state energy estimation using circuit depth with exponentially improved dependence on precision[END_REF], we can tune the parameter ∆ anywhere between Θ(ǫ) and Θ(∆ true ), and doing so yields a class of GSEE algorithms that smoothly transition from having circuit depth O ∆ -1 true polylog ∆ true ǫ -1 p -1 0 and runtime Õ ∆ true ǫ -2 p -1 0 to having circuit depth O(ǫ -1 polylog p -1

0

) and runtime Õ ǫ -1 p -1 0 . Therefore, our algorithm is also able to trade circuit depth for runtime like the one in Ref. [START_REF] Wang | Quantum algorithm for ground state energy estimation using circuit depth with exponentially improved dependence on precision[END_REF] but has a shorter runtime than the latter.

The core component of our algorithm is a novel procedure for performing rejection sampling on a quantum computer. We apply it to draw random values from a neighborhood of E 0 . Specifically, suppose we want to draw a sample from a target Gaussian-like distribution ν(x) such that ν(x) peaks around E 0 . It might be difficult to do so directly. Instead, we first draw a random value x from another distribution µ(x), and construct a unitary operator U x (H) such that

| 0 k U x (H) 0 k |ψ | 2 ∝ ν(x)
µ(x) for some integer k. Then we accept x if a measurement on the first k qubits of U x (H) |ψ yields outcome 0 k , and discard x otherwise. This procedure is akin to the standard rejection sampling algorithm.

We show that one can utilize quantum singular value transformation (QSVT) [START_REF] Gilyén | Quantum singular value transformation and beyond: Exponential improvements for quantum matrix arithmetics[END_REF] or similar techniques [START_REF] Hao Low | Hamiltonian simulation by uniform spectral amplification[END_REF][START_REF] Dong | Ground-state preparation and energy estimation on early fault-tolerant quantum computers via quantum eigenvalue transformation of unitary matrices[END_REF] to implement U x (H) with a low-depth circuit. Furthermore, the circuit only involves a few ancilla qubits and the rejection sampling procedure is efficient (i.e. not many raw samples are discarded). Once we obtain sufficiently many samples from the desired distribution ν(x), we compute their mean which is close to E 0 with high probability.

Although good estimates of the spectral gap are available for many physically relevant models [START_REF] Lee | Is there evidence for exponential quantum advantage in quantum chemistry[END_REF], this is not always the case. Furthermore, our algorithm will yield good estimates of the ground state energy for evolution times t ≪ ∆ -1 if the initial state has more structure. For instance, if it has large overlap with the ground state and little overlap with other low-energy states. Thus, it is highly desirable to ensure that our algorithm can be proven to give reliable estimates when we do not know ∆ or have reasons to believe we have a structured initial state. To address this, we devise a test that does not rely on knowledge of the spectral gap and only accepts if the ground state energy estimate produced by the algorithm is within some tolerated error of the true ground state energy 1 . Furthermore, the test will always accept if we evolve for times proportional to the inverse of the spectral gap. Although our algorithm, together with the certification result, is the first to deliver certifiably correct energy estimates in the regime of depths ≪ ǫ -1 , its sample complexity of Õ(ǫ -4 p -3 0 ) might make its implementation prohibitive in some contexts. The remainder of this paper is organized as follows. In Section 2, we formally define our problems and summarize our main results. Then in Section 3, we describe our method for performing rejection sampling on a quantum computer and utilize it to design a quantum algorithm for GSEE. The implementation details of certain parts of the algorithm depend on the input model of H and are deferred to Sections 4 and 5, where we develop the full-fledged versions of this algorithm in the block-encoding model and the Hamiltonian-evolution model, respectively, and analyze the costs of the resulting algorithms. After that, in Section 6, we present a method for certifying the correctness of the outputs of our algorithms. Finally, Section 7 concludes this paper.

2 Overview of the main results

Problem formulation

Let us start with a formal definition of the GSEE problem:

Definition 2.1 (GSEE). Suppose H is an n-qubit Hamiltonian with unknown spectral decomposi- tion H = N -1 j=0 E j |E j E j |, where N = 2 n , -1 ≤ E 0 < E 1 ≤ E 2 ≤ ... ≤ E N -1 ≤ 1, E 1 -E 0 ≥ ∆ for known ∆ > 0,
and the |E j 's are orthonormal eigenstates of H. Moreover, suppose we can prepare an n-qubit state |ψ = N -1 j=0 γ j |E j such that |γ 0 | 2 ≥ η for known η ∈ (0, 1). Our goal is to estimate E 0 within accuracy (i.e. additive error) ǫ with probability at least 1δ, for given ǫ > 0 and δ ∈ (0, 1). Furthermore, we want to achieve this by using quantum circuits that depend on H and |ψ and classical post-processing of their outcomes.

The above definition of the GSEE problem has omitted one important detail -how H is accessed. In this paper, we will consider two common input models of H: the block-encoding model [START_REF] Gilyén | Quantum singular value transformation and beyond: Exponential improvements for quantum matrix arithmetics[END_REF] and the Hamiltonian-evolution model. Definition 2.2 (block-encoding). Let n, m ∈ N and α, γ ∈ R + be arbitrary. We say that an (n + m)-qubit unitary operator U is an (α, m, γ)-block-encoding of an n-qubit linear operator A if

α ( 0 m | ⊗ I) U (|0 m ⊗ I) -A ≤ γ. (1)
For simplicity, we will only consider (1, m, 0)-block-encodings of H in this work. In fact, our results in Sections 4 can be readily extended to general block-encodings of H, but doing so would complicate the notation and obscure our main message. Definition 2.3 (GSEE in the block-encoding model). For GSEE in the block-encoding model, we assume that there exists a method to implement a (1, m, 0)-block-encoding U of H for some integer m ≥ 0. Our quantum circuits will consist of multiple queries to U and U † and other primitive quantum gates, and act on 0 k |ψ for some k ≥ m. The other assumptions are the same as in Definition 2.1.

Definition 2.4 (GSEE in the Hamiltonian-evolution model).

For GSEE in the Hamiltonianevolution model, we assume that there exists a method to implement the controlled time evolution of H, i.e. the unitary operation controlled-e iHt := |0 0| ⊗ I + |1 1| ⊗ e iHt . Each of our quantum circuits will consist of multiple queries to controlled-e iHt and controlled-e -iHt for some t = Θ(1) and other primitive quantum gates, and act on 0 k |ψ for some k ≥ 1. The other assumptions are the same as in Definition 2.1.

Main results

Our main results are threefold: a new way of performing rejection sampling on a quantum device, how to use it to obtain state-of-the-art early fault-tolerant algorithms for GSEE, and techniques for certifying ground state energy estimates. We will now informally summarize these technical contributions.

First, we start with a novel way to perform rejection sampling on a quantum computer:

Theorem 2.5 (Rejection sampling from the spectral measure). Let H be an n-qubit Hamiltonian with spectral decomposition

H = i E i |E i E i |,
and let |ψ be an n-qubit state. Define the probability measure p on R as

p(x) = i p i δ(x -E i ), p i = | E i |ψ | 2 , (2)
Furthermore, let ν, µ be two probability density functions on R such that supp(p * ν) ⊂ supp(µ), where * denotes the convolution. Finally, assume that there exists a constant M ∈ C such that for every x ∈ supp(µ) we have access to a (1, m, 0)-block-encoding U x of an operator h x (H) which satisfies

|h x (H)| 2 = ν(H -xI) |M | 2 µ(x) , h x (H) ≤ 1, (3)
and that we can generate samples distributed according to µ. Then we can generate a sample distributed according to p * ν from an expected number |M | 2 uses of such U x 's.

Using this rejection sampling routine, we are then able to obtain the following GSEE algorithm:

Theorem 2.6 (GSEE algorithm, informal). Suppose H is an n-qubit Hamiltonian with unknown spectral decomposition

H = N -1 j=0 E j |E j E j |, where N = 2 n , -1 ≤ E 0 < E 1 ≤ • • • ≤ E N -1 ≤ 1, E 1 -E 0 ≥ ∆ for known ∆ > 0.
Moreover, suppose we can prepare an n-qubit state |ψ such that | ψ|E 0 | 2 ≥ η for known η ∈ (0, 1). Then in both the block-encoding model and the Hamiltonianevolution model, there exists an algorithm that, for any given ǫ > 0 and δ ∈ (0, 1), estimates E 0 within additive error ǫ with probability at least 1δ by making use of

Õ η -1 (∆ 2 ǫ -2 + log ∆ -1 ) log δ -1
quantum circuits, where each circuit has depth proportional to Õ ∆ -1 log(η -1 ∆ǫ -1 ) .

We emphasize that in Theorem 2.6, ∆ is merely a lower bound on the true spectral gap ∆ true := E 1 -E 0 of H, not necessarily ∆ true itself. In fact, we can tune it between Θ(ǫ) and Θ(∆ true ), and doing so leads to a class of GSEE algorithms: Corollary 2.7 (GSEE algorithms, informal). Suppose H is an n-qubit Hamiltonian with unknown spectral decomposition

H = N -1 j=0 E j |E j E j |, where N = 2 n , -1 ≤ E 0 < E 1 ≤ • • • ≤ E N -1 ≤ 1, ∆ true = E 1 -E 0 is
the spectral gap of H. Moreover, suppose we can prepare an n-qubit state |ψ such that | ψ|E 0 | 2 ≥ η for known η ∈ (0, 1). Then in both the block-encoding model and the Hamiltonian-evolution model, for every β ∈ [0, 1], there exists an algorithm that, for any given ǫ > 0 and δ ∈ (0, 1), estimates E 0 within additive error ǫ with probability at least 1δ by making use of

Õ η -1 ǫ -2+2β ∆ 2-2β true + log ǫ -β ∆ β-1 true log δ -1
quantum circuits, where each circuit has depth proportional to

Õ(ǫ -β ∆ β-1 true log(η -1 ∆ 1-β true ǫ β-1 )).
Note that this result interpolates between a standard limit (β = 0) and a Heisenberg limit (β = 1) by exchanging the circuit depth for the number of samples. Importantly, rather than scaling as η -2 like previous methods [START_REF] Wang | Quantum algorithm for ground state energy estimation using circuit depth with exponentially improved dependence on precision[END_REF][START_REF] Ding | Even shorter quantum circuit for phase estimation on early faulttolerant quantum computers with applications to ground-state energy estimation[END_REF], the number of repetitions here is quadratically improved to an η -1 scaling. This is achieved without increasing (i.e. worsening the scaling of) the number of quantum operations per circuit, which is what standard approaches to quadratic speedups (e.g. amplitude amplification) require.

Figure 1 compares the performance of our GSEE algorithm with others in the literature. In particular, one can see that our algorithm has the best performance in terms of η in the regime of circuit depths between Õ(∆ -1 true ) and Õ(ǫ -1 ). So far, we have only considered algorithms for GSEE with circuit depths proportional to ∆ -1 , where ∆ > ǫ is a lower bound on the spectral gap ∆ true = E 1 -E 0 of H. However, it is desirable to develop algorithms that do not require this knowledge to accurately estimate the ground state energy. Thus, we develop a method for certifying the ground state energy estimate, albeit with a higher sample complexity:

Circuit depth (T max ) Runtime (T tot ) η -2 ǫ -1 η -1 ǫ -1 η -1/2 ǫ -1 ∆ -1 true η -1 ǫ -2 ∆ true η -2 ǫ -2 ∆ true ǫ -1 η -1/2 ǫ -1 [40, 41]
Cor. 2.7 [START_REF] Lin | Heisenberg-limited ground-state energy estimation for early faulttolerant quantum computers[END_REF][START_REF] Wan | Randomized quantum algorithm for statistical phase estimation[END_REF] [14] [START_REF] Dong | Ground-state preparation and energy estimation on early fault-tolerant quantum computers via quantum eigenvalue transformation of unitary matrices[END_REF] Figure 1: This figure shows the landscape of early fault-tolerant GSEE algorithms plotted according to their runtimes and circuit depths. The loosely dashed green curve represents our previous interpolation result [START_REF] Wang | Quantum algorithm for ground state energy estimation using circuit depth with exponentially improved dependence on precision[END_REF] (which is also implied by Ref. [START_REF] Ding | Even shorter quantum circuit for phase estimation on early faulttolerant quantum computers with applications to ground-state energy estimation[END_REF]), while the lower densely dashed violet curve depicts the result in Corollary 2.7 (which encompasses the results in Corollary 4.8 and Corollary 5.8). For the purposes of this plot, we have assumed that ǫ < η∆ true , leading us to place η -1 ǫ -2 ∆ true above η -2 ǫ -1 , though this assumption is not necessary for the main results of our paper.

Theorem 2.8 (Ground state energy certification, informal). Under the same conditions as Corollary 2.7, there exists an algorithm that is given, as an input, arbitrary σ > 0 and Ê0 ∈ [E 0 -σ, E 0 +σ] (but no lower bound on ∆ true ). It makes

Õ(η -3 ǫ -4 log δ -1 ) (4)
uses of quantum circuits of depth Õ(σ -1 ) and produces a new estimate Ê′ 0 of E 0 that it either accepts or rejects. If | Ê′ 0 -E 0 | > 4ǫ, it rejects Ê′ 0 with probability at least 1-δ. If it accepts, then | Ê′ 0 -E 0 | ≤ ǫ with probability of incorrect acceptance at most δ. Furthermore, the maximal σ for which the algorithm outputs an estimate Ê′ 0 it accepts satisfies

| Ê′ 0 -E 0 | ≤ ǫ is σ = Ω ∆ true /log ǫ -1 η -1 .
Thus, with a higher sample complexity, we can also produce a certificate that the estimate was correct, even if we do not know the spectral gap ∆ true of H. Furthermore, the algorithm described above will always accept with high confidence if we go to depths proportional to the true spectral gap, the minimal depth required in Cor. 2.7. Note, however, that in some situations the algorithm can accept and, thus, give a correct estimate, even if σ ≫ ∆ true . As we discuss in more detail in Section 6, this can happen if, e.g. the initial state has little support on low-energy excitations. We see that our certification algorithm does not certify the spectral gap itself, only our ground state energy estimate. To the best of our knowledge, this is the first method to achieve certifiable GSEE in the regime where the circuit depth is O(∆ -1 true ). Finally, note that the assumption that we are given an estimate that satisfies |E 0 -Ê0 | ≤ σ does not increase the overall depth of circuits required for the certification method. This is because we will show how to obtain such an estimate with circuits of depth O(σ -1 ) later in Section 3.3. Finally, let us briefly comment on the role of δ. Our algorithm will be based on empirically estimating the variance of a random variable up to a certain precision. As usual with such protocols, the empirical estimate will be within a tolerated error up to a given failure probability. It is the case that our protocol always rejects correctly if it is given estimates of the variance within the tolerated error. Thus, the δ quantifies the probability of the estimates deviating more than required for the verification and, then, we do not posess a guarantee on its performance.

GSEE via rejection sampling

In the section, we will discuss our GSEE algorithm which will be based on a novel way of performing rejection sampling on a quantum computer. But since we believe that this rejection sampling routine will find applications elsewhere, we will first describe it in a higher level of generality.

Rejection sampling for the spectral measure

Given the spectral measure p associated with an initial state |ψ ψ| and a Hamiltonian H =

i |E i E i |, i.e. p(x) = i p i δ(x -E i ), p i = | E i |ψ | 2 , (5)
we will show how to sample from a random variable with distribution p * ν on a quantum computer, where ν is another probability measure and * denotes the convolution. It is then not difficult to see that the density of p * ν is given by:

(p * ν)(x) = i p i ν(x -E i ). (6)
Sampling from such distributions for various choices of ν gives us access to information about the spectrum of the Hamiltonian. For instance, if ν is itself close to a delta distribution at 0, then sampling from this distribution is close to sampling from p itself. This will, roughly speaking, correspond to our algorithm later, where we will pick ν to be a Gaussian with small variance. Let us now discuss the sampling routine. We will assume that we can classically sample from a random variable X distributed according to a density µ such that supp(p * ν) ⊂ supp(µ). Furthermore, for every x ∈ supp(µ), we can implement a (1, m, 0)-block-encoding U x of the operator ν(H-xI) M μ(x) , where ν and μ satisfy:

ν(H -xI) M μ(x) 2 = ν(H -xI) |M | 2 µ(x) , (7)
and M ∈ C is a constant such that for all x ∈ supp(µ):

i p i ν(E i -x) |M | 2 µ(x) ≤ 1. ( 8
)
There are many cases where explicit block encoding circuits have been developed [START_REF] Lee | Even more efficient quantum computations of chemistry through tensor hypercontraction[END_REF][START_REF] Loaiza | Reducing molecular electronic hamiltonian simulation cost for linear combination of unitaries approaches[END_REF]. Let us now explain the quantum algorithm that allows us to sample from p * ν in this setting with an expected number O(|M | 2 ) runs of such U x 's. First, we generate a sample x from the density µ. After that, we implement the block-encoding U x of ν(H-xI) M μ(x) and apply it to |0 m |ψ . Then we measure the m ancilla qubits of the resulting state. If the outcome is 0 m , we accept the sample. Let us see why, conditioned on the acceptance, x will be distributed according to p * ν.

Note that:

U x |0 m |ψ = |0 m ν(H -xI) M μ(x) |ψ + j∈{0,1} m , j =0 m |j ψ ′ j , (9)
where the ψ ′ j 's are unnormalized states. By expanding ν(H-xI) M μ(x) |ψ in the eigenbasis of H, we see that the probability of measuring 0 m on the first register is expressed by:

i ν(E i -x) M μ(x) 2 | E i |ψ | 2 = i p i ν(E i -x) |M | 2 µ(x) = (p * ν)(x) |M | 2 µ(x) . ( 10
)
Thus the probability we accept is exactly the same as if we ran the rejection sampling algorithm for the densities p * ν and µ. Furthermore, by standard properties of rejection sampling, we will, in expectation, need O(|M | 2 ) trials to generate a single sample from the target distribution. From the discussion above, we obtain our new rejection sampling routine:

Theorem 3.1. Let H = i E i |E i E i | be a
Hamiltonian on n qubits and |ψ an n-qubit state.

Definite the probability measure p on R as

p(x) = i p i δ(x -E i ), p i = | E i |ψ | 2 , (11)
Furthermore, let ν, µ be two probability density functions on R such that supp(p * ν) ⊂ supp(µ), where * denotes the convolution. Finally, assume that there exists a constant M ∈ C such that for every x ∈ supp(µ) we have access to a (1, m, 0)-block-encoding of a unitary U x of an operator h x which satisfies

|h x (H)| 2 = ν(H -xI) |M | 2 µ(x) , h x (H) ≤ 1, (12)
and that we can generate samples distributed according to µ. Then we can generate a sample distributed according to p * ν from an expected number |M | 2 uses of such U x 's.

Let us illustrate the rejection sampling routine with one of the leading examples of this work: conditioned normal random variables. Suppose that we wish to sample from p * n σ, [a,b] , where n σ, [a,b] is the density of a normal random variable with mean 0 and variance σ 2 conditioned on an interval

[a, b]. That is n σ,[a,b] corresponds to the density: n σ,[a,b] (x) = e -x 2 2σ b a e -y 2 2σ dy χ [a,b] (x). ( 13
)
Here χ [a,b] (x) is the indicator function. One possible strategy we will follow to sample from this distribution is to perform rejection sampling with the density µ given by the uniform distribution on [a, b]. In this case, if we implement a block-encoding of the operator

h σ (H -xI) = e -(H-xI) 2 4σ 2 , (14)
we see that it satisfies the requirement in Eq. ( 7) with |M | 2 given by

|M | 2 = σ √ 2π (b -a) b a (p * n σ,[-∞,∞] )(x)dx , (15)
as with this choice of |M | 2 , the ratio of the density of p * n σ,[a,b] and the uniform distribution on [a, b] coincide with the probability of accepting. Furthermore, the convolution p * h σ is clearly upper bounded by 1. Note that we do not need to know the constant in Eq. ( 15) to run the algorithm, only h σ . But |M | 2 will determine the number of repetitions. Building upon recent work on GSEE, we will show how to estimate the ground state energy by first sampling from the convolution of the spectral measure with a uniform distribution on [0, 1] (the approach of Ref. [START_REF] Lin | Heisenberg-limited ground-state energy estimation for early faulttolerant quantum computers[END_REF]) and then refining the estimate by sampling with the convolution of the same measure with a normal random variable (the approach of Ref. [START_REF] Wang | Quantum algorithm for ground state energy estimation using circuit depth with exponentially improved dependence on precision[END_REF]). However, the previous algorithms worked by estimating the density ν * p through generalized Hadamard tests. As here we generate samples instead of estimating densities, we will obtain a better scaling in the ground-state overlap.

We believe that the rejection sampling subroutine developed in this work is of independent interest and will find application elsewhere. Indeed, to the best of our knowledge, it is the first quantum algorithm for sampling from distributions associated with the spectral measure with a constant number of ancilla qubits.

Basic idea for GSEE

We will now discuss our algorithm for GSEE based on rejection sampling. This algorithm consists of two stages. The first stage generates a coarse estimate Ẽ0 of E 0 by using the strategy of Ref. [START_REF] Dong | Ground-state preparation and energy estimation on early fault-tolerant quantum computers via quantum eigenvalue transformation of unitary matrices[END_REF]. The second stage refines this estimate and improves its accuracy to arbitrarily small ǫ. To achieve this goal, we make the following observation about Gaussian functions. Let f σ (x) := e -x 2 /(2σ 2 ) be (an unnormalized) Gaussian. If we truncate f 2 σ (x) = e -x 2 /σ 2 to a sufficiently large interval around 0, then the mean of this truncated function is extremely close to 0, as long as the middle point of the interval is not far from 0:

-x 2 /σ 2 dx s+w s-w e -x 2 /σ 2 dx ≤ ǫ 2.24e . ( 16
)
If we perturb f 2 σ (x) slightly, then the resulting function still possesses a similar property:

Lemma 3.3. Let σ > 0 and ǫ ∈ (0, σ] be arbitrary, and let w = 2σ ln (eσ/ǫ). Then for any

s ∈ [-w/2, w/2], if g : [s -w, s + w] → R + satisfies that |g(x) -e -x 2 /σ 2 | ≤ 0.03ǫ σ ln(eσ/ǫ) for all x ∈ [s -w, s + w],
then we have

The proofs of Lemma 3.2 and Lemma 3.3 are given in Appendix A. As a corollary of Lemma 3.3, we obtain the following strategy for estimating an unknown quantity E with high accuracy: Corollary 3.4. Under the same assumptions as in Lemma 3.3, if X is a random variable such that it takes values from [E + sw, E + s + w] for some unknown E ∈ R and its probability density function is proportional to g(x -E), then its expected value E [X] is O(ǫ)-close to E. Now we apply this strategy to GSEE in which the target quantity is E 0 . We first choose appropriate σ = O(∆/ log (η -1 ∆ǫ -1 )) and w = O(σ log (σ/ǫ)). Then we find a coarse estimate Ẽ0 of E 0 such that s := Ẽ0 -E 0 ⊆ [-w/2, w/2]. Next, we pick x uniformly at random from [ Ẽ0w, Ẽ0 + w], and construct a quantum circuit that implements f σ (H -xI) approximately and probabilistically, and run this circuit on the initial state 0 k |ψ for appropriate k ≥ 1. Then the probability of this circuit approximately performing

f σ (H -xI) on |ψ is roughly proportional to f 2 σ (x -E 0 )
. Informally, this is because this probability is close to

f σ (H -xI) |ψ 2 = N -1 j=0 p j f 2 σ (E j -x) ≈ p 0 f 2 σ (E 0 -x), (18)
since Gaussian functions decrease exponentially from their peaks and x ∈ [ Ẽ0w, Ẽ0 + w] is quite far from the E j 's for all j ≥ 1. We can determine whether this event happens by measuring the k ancilla qubits and checking if particular outcomes are obtained in the end. If this event happens, then x is accepted and appended to an array S; otherwise, x is rejected and discarded. Clearly, the entries of S are i.i.d. random variables. Let Y be such a random variable. By Corollary 3.4, its expected value E [Y ] is O(ǫ)-close to E 0 . Furthermore, we can prove that the probability of a random x ∈ [ Ẽ0w, Ẽ0 + w] getting accepted by this procedure is Ω(η/ log (∆ǫ -1 )), i.e. this rejection sampling method is efficient. Finally, we repeat the above procedure Õ(η -1 ∆ 2 ǫ -2 ) times and compute the mean of the accepted samples. Then by Chernoff's bound, this quantity will be O(ǫ)-close to E [Y ] with high probability, which implies that it is also an O(ǫ)-accurate estimate of E 0 with high probability, as desired.

Formal description of the GSEE algorithm

We will now describe our GSEE algorithm in a way that is agnostic to the input model of H. This algorithm can be executed in both the block-encoding and Hamiltonian-evolution models. It requires the block-encodings of certain functions of H whose implementation details depend on the specific input model and will be explained in Sections 4 and 5. The other components of the algorithm remains identical in these two models. We first propose a basic GSEE algorithm (i.e. Algorithm 1) which generalizes the one from Ref. [START_REF] Dong | Ground-state preparation and energy estimation on early fault-tolerant quantum computers via quantum eigenvalue transformation of unitary matrices[END_REF]. This algorithm will be used to generate a coarse estimate of E 0 .

Before we give a formal description of the algorithm, let us consider a simplified and idealized version to build up some intuition. Suppose that for λ ∈ (-1, 1), we can implement a (c, m, 0)block-encoding V λ of the operator f λ (H), where

f λ : [-1, 1] → [0, 1] satisfies f (x) = 1 for x ∈ [-1, λ] and f (x) = 0 for x ∈ (λ, 1]. Let us now suppose that |ψ is an initial state such that | ψ|E 0 | 2 ≥ η
for some known η ∈ (0, 1). If we pick λ > E 0 , then f λ (H) |ψ 2 ≥ η, as the operator f λ (H) will act as the identity on the ground space of H. On the other hand, if we pick λ < E 0 we have f λ (H) |ψ 2 = 0, as it will send all the eigenspaces of H to 0. Furthermore, by the definition of block-encoding, if we measure the m ancilla qubits of the state V λ |0 m |ψ , the probability of getting outcome 0 m is given by c 2 f λ (H) |ψ 2 . Thus, if we could implement such a block-enconding, we could also infer what the ground state energy is by performing a binary search. Indeed, we only need to find the smallest λ such that f λ (H) |ψ 2 ≥ η. As for a given λ, whether f λ (H) |ψ 2 ≥ η or f λ (H) |ψ 2 = 0 can be determined with success probability of 2/3 by repeating the aforementioned experiment O(η -1 ) times and checking how often we obtain measurement outcome 0 m . So the block-encoding V λ would give rise to an algorithm with O(η -1 log ǫ -1 )) runs to estimate the ground state energy with accuracy ǫ. However, implementing block-encodings of non-smooth functions like f λ cannot be accomplished directly and it is necessary to perform suitable approximations of f λ by smooth functions. This motivates the introduction of the function f a,b;ǫ below. The fact that we only approximate f λ also introduces subtleties on how we post-process the samples and perform the binary search, but the strategy is comparable to the one described above.

Algorithm 1 Basic GSEE algorithm

1: procedure Basic GSEE(ǫ, δ, η, c, m) 2: l ← -1; 3: r ← 1; 4: L ← ⌈log 3/2 ǫ -1 ⌉; 5: M ← ⌈12cη -1 ln Lδ -1 ⌉; 6: ǫ ′ ← min( √ 0.1η, 0.05); 7:
while rl > 2ǫ do 8: return (l + r)/2. 25: end procedure Theorem 3.5. Suppose H is an n-qubit Hamiltonian with unknown spectral decomposition

K ← 0; 9: a ← (2l + r)/3; 10: b ← (l + 2r)/3; 11: Let V a,b;ǫ ′ be a (c, m, 0)-block-encoding of f a,b;ǫ ′ (H), where f a,b;ǫ ′ : [-1, 1] → R satisfies that 1 -ǫ ′ ≤ f a,b;ǫ ′ (x) ≤ 1 for all x ∈ [-1, a], |f a,b;ǫ ′ (x)| ≤ 1 for all x ∈ [a, b], and |f a,b;ǫ ′ (x)| ≤ ǫ ′ for all x ∈ [b, 1];
H = N -1 j=0 E j |E j E j |, where N = 2 n , -1 ≤ E 0 < E 1 ≤ E 2 ≤ ... ≤ E N -1 ≤ 1,
and the |E j 's are orthonormal eigenstates of H. Moreover, suppose we can prepare an n-qubit state |ψ = N -1 j=0 γ j |E j such that |γ 0 | 2 ≥ η for known η ∈ (0, 1). Finally, suppose for arbitrary a, b ∈ (-1, 1) and ǫ ′ ≥ 0, we can implement a (c, m, 0)-block-encoding of f a,b;ǫ ′ (H), where f a,b;ǫ

′ : [-1, 1] → R satisfies that 1 -ǫ ′ ≤ f a,b;ǫ ′ (x) ≤ 1 for all x ∈ [-1, a], |f a,b;ǫ ′ (x)| ≤ 1 for all x ∈ [a, b], and |f a,b;ǫ ′ (x)| ≤ ǫ ′ for all x ∈ [b, 1].
Then for arbitrary ǫ > 0 and δ ∈ (0, 1), the output of Algorithm 1 (i.e. BASIC GSEE(ǫ, δ, η, c, m)) is ǫ-close to E 0 with probability at least 1δ.

Proof. We claim that in each iteration of the while-loop, if E 0 ≤ a = (2l+r)/3, then K > 0.5c -1 ηM with probability at least 1δ/L; if E 0 ≥ b = (l + 2r)/3, then K < 0.5c -1 ηM with probability at least 1δ/L. If this is true, then after each iteration, E 0 remains between l and r with probability at least 1δ/L. Since there are at most L iterations in the while-loop, E 0 is between the final l and r with probability at least 1δ. Then since lr ≤ 2ǫ in the end, we know that (l + r)/2 is ǫ-close to E 0 with probability at least 1δ, as desired.

To prove the above claims, we need to bound the probability of obtaining measurement outcome 0 m in each iteration of the for-loop in the cases E 0 ≤ a and E 0 ≥ b respectively. If E 0 ≤ a, then this probability is

c -1 N -1 j=0 p j f 2 a,b;ǫ ′ (E j ) ≥ c -1 p 0 f 2 a,b;ǫ ′ (E 0 ) ≥ c -1 η(1 -ǫ ′ ) 2 ≥ 0.9c -1 η. ( 19
)
On the other hand, if E 0 ≥ b, then this probability is

c -1 N -1 j=0 p j f 2 a,b;ǫ ′ (E j ) ≤ c -1 max 0≤j≤N -1 f 2 a,b;ǫ ′ (E j ) ≤ c -1 (ǫ ′ ) 2 ≤ 0.1c -1 η. (20)
Next, we invoke the Chernoff-Hoeffding Theorem which states that if X 1 , . . . , X n are i.i.d. random variables which take values in {0, 1}. Let p = E [X 1 ]. Then for any q ∈ (p, 1] and r ∈ [0, p), we have

P 1 n n i=1 X i ≥ q ≤ e -nD(q p) , (21)
P 1 n n i=1 X i ≤ r ≤ e -nD(r p) , (22)
where D(x y) := x ln (x/y)+(1-x) ln ((1x)/(1y)) is the Kullback-Leibler divergence between Bernoulli random variables with parameters x and y respectively. We also use the fact that

D(x y) ≥ (x -y) 2 2 max(x, y)
, ∀x, y ∈ (0, 1). ( 23)

If E 0 ≤ a, then we have p ≥ 0.9c -1 η, and hence

P K/M ≤ 0.5c -1 η ≤ e -M D(0.5c -1 η p) (24)
in which

D(0.5c -1 η p) ≥ (p -0.5c -1 η) 2 2p ≥ 0.088c -1 η, (25)
as (p-0.5c -1 η) 2 2p increases monotonically for p ≥ 0.5c -1 η. Then by our choice of M , we know that P K ≤ 0.5c -1 ηM ≤ δ/L.

On the other hand, if E 0 ≥ b, we have p ≤ 0.1c -1 η, and hence

P K/M ≥ 0.5c -1 η ≤ e -M D(0.5c -1 η p) (26)
in which

D(0.5c -1 η p) ≥ (0.5c -1 η -p) 2 2 • 0.5c -1 η ≥ 0.16c -1 η. (27)
Then by our choice of M , we obtain that P K ≥ 0.5c -1 ηM ≤ δ/L. The theorem is thus proved.

Next, we present our main GSEE algorithm (i.e. Algorithm 2) which utilizes Algorithm 1 to obtain a coarse estimate of E 0 and then refines this estimate by the strategy mentioned in Section 3.2:

Algorithm 2 Advanced GSEE algorithm 1: procedure ADV GSEE(ǫ, δ, ∆, η, c 1 , m 1 , c 2 , m 2 ) 2: if ǫ ≥ ∆/8 then 3: return BASIC GSEE (ǫ, δ, η, c 1 , m 1 ); 4: end if 5: σ ← ∆ 5 √ ln(ec 2 ) √ ln(eη -1 ∆ǫ -1 ) ; 6: ǫ 1 ← min( ǫ 1.1 , σ); 7:
w ← 2σ ln (eσ/ǫ 1 );

8: Ẽ0 ← BASIC GSEE (w/2, δ/3, η, c 1 , m 1 ); 9:
ǫ ′′ ← 0.0075ǫ 1 η σ ln(eσ/ǫ 1 ) ; 10: for i = 1, 2, . . . , M do 13:

M ← ⌈32c 2 2 η -1 ln (eσ/ǫ 1 ) • max(2w 2 ǫ -2 ln 6δ -1 , ln 3δ -1 )⌉;
Choose ξ i uniformly at random from [ Ẽ0w, Ẽ0 + w];

14:

Let W σ,ξ i ;ǫ ′′ be a (c 2 , m 2 , 0)-block-encoding of g σ,ξ i ;ǫ ′′ (H), where g σ,ξ i ;ǫ

′′ : [-1, 1] → R satisfies that |g σ,ξ i ;ǫ ′′ (x) -e -(x-ξ i ) 2 /(2σ 2 ) | ≤ ǫ ′′ and |g σ,ξ i ;ǫ ′′ (x)| ≤ 1 for all x ∈ [-1, 1]. 15:
Measure the first m 2 qubits of W σ,ξ i ;ǫ ′′ |0 m 2 |ψ in the standard basis;

= N -1 j=0 E j |E j E j |, where N = 2 n , -1 ≤ E 0 < E 1 ≤ E 2 ≤ ... ≤ E N -1 ≤ 1, E 1 -E 0 ≥ ∆ for known ∆ > 0,
and the |E j 's are orthonormal eigenstates of H. Moreover, suppose we can prepare an n-qubit state |ψ = N -1 j=0 γ j |E j such that |γ 0 | 2 ≥ η for known η ∈ (0, 1). In addition, suppose for arbitrary a, b ∈ (-1, 1) and ǫ ′ ≥ 0, we can implement a (c

1 , m 1 , 0)-block-encoding of f a,b;ǫ ′ (H), where f a,b;ǫ ′ : [-1, 1] → R satisfies that 1 -ǫ ′ ≤ f a,b;ǫ ′ (x) ≤ 1 for all x ∈ [-1, a], |f a,b;ǫ ′ (x)| ≤ 1 for all x ∈ [a, b], and |f a,b;ǫ ′ (x)| ≤ ǫ ′ for all x ∈ [b, 1]. Finally, suppose for arbitrary σ, ǫ ′′ ∈ (0, 1) and ξ ∈ [-2, 2], we can implement a (c 2 , m 2 , 0)-block-encoding of g σ,ξ;ǫ ′′ (H), where c 2 ≥ 1 and g σ,ξ;ǫ ′′ : [-1, 1] → R satisfies that |g σ,ξ;ǫ ′′ (x) -e -(x-ξ) 2 /(2σ 2 ) | ≤ ǫ ′′ and |g σ,ξ;ǫ ′′ (x)| ≤ 1 for all x ∈ [-1, 1].
Then for arbitrary ǫ > 0 and δ ∈ (0, 1), the output of Algorithm 2 (i.e.

ADV GSEE(ǫ, δ, ∆, η, c 1 , m 1 , c 2 , m 2 )) is ǫ-close to E 0 with probability at least 1 -δ.
Proof. Let Ê0 denote the output of Algorithm 2. If ǫ ≥ ∆/8, then Ê0 = BASIC GSEE (ǫ, δ, η, c 1 , m 1 ) is ǫ-close to E 0 with probability at least 1δ by Theorem 3.5. So we only need to focus on the case ǫ < ∆/8.

For convenience, we introduce the following notation. Let A i = 1 if ξ i is added to S, and A i = 0 otherwise, for i = 1, 2, . . . , M . Then since ξ 1 , ξ 2 , . . . , ξ M are i.i.d. random variables, A 1 , A 2 , . . . , A M are i.i.d. Bernoulli variables. Next, let K = M i=1 A i be the size of S at the end of the algorithm, and let X k be the k-th number added to S, for k = 1, 2, . . . , K. Then X 1 , X 2 , . . . , X K are i.i.d. random variables, and Ê0 = 1 K K k=1 X k has the same expected value as X 1 . We claim that:

P [|E [X 1 ] -E 0 | ≥ ǫ/2] ≤ δ/3, (28)
P K ≥ 8w 2 ǫ -2 ln 6δ -1 ≥ 1 -δ/3, (29)
and

P | Ê0 -E [X 1 ] | ≥ ǫ/2 | K ≥ 8w 2 ǫ -2 ln 6δ -1 ≤ δ/3. ( 30
)
Let G denote the event K ≥ 8w2 ǫ -2 ln 6δ -1 . If the above claims are true, then by Eqs. ( 29), [START_REF] Quek | Exponentially tighter bounds on limitations of quantum error mitigation[END_REF] and union bound, we get

P | Ê0 -E [X 1 ] | ≥ ǫ/2 ≤ P | Ê0 -E [X 1 ] | ≥ ǫ/2 | G + P [¬G] ≤ δ/3 + δ/3 = 2δ/3. (31)
Then by Eqs. ( 28), [START_REF] Wang | Noise-induced barren plateaus in variational quantum algorithms[END_REF], triangle inequality and union bound, we know that

P | Ê0 -E 0 | ≥ ǫ ≤ P | Ê0 -E [X 1 ] | ≥ ǫ/2 + P [|E [X 1 ] -E 0 | ≥ ǫ/2] ≤ 2δ/3 + δ/3 = δ, (32)
as desired.

It remains to prove Eqs. ( 28), ( 29) and [START_REF] Quek | Exponentially tighter bounds on limitations of quantum error mitigation[END_REF]. Consider each iteration in the for-loop. Suppose ξ ∈ [ Ẽ0w, Ẽ0 + w] is chosen 2 . Let W σ,ξ;ǫ ′′ be the (c 2 , m 2 , 0)-block-encoding of g σ,ξ;ǫ ′ (H). That is,

W σ,ξ;ǫ ′′ |0 m 2 |ψ = c -1 2 |0 m 2 g σ,ξ;ǫ ′′ (H) |ψ + s∈{0,1} m 2 , s =0 m 2 |s φ s,ξ,ǫ ′′ , (33)
where the φ s,ξ,ǫ ′′ 's are unnormalized states. Consequently, the probability of obtaining measurement outcome 0 m 2 when measuring the first m 2 qubits of this state (in the standard basis) is

q(ξ) := P [A i = 1|ξ i = ξ] = c -2 2 ψ| g 2 σ,ξ;ǫ ′′ (H) |ψ = c -2 2 N -1 j=0 p j g 2 σ,ξ;ǫ ′′ (E j ). ( 34
)
Meanwhile, let f σ (x) := e -x 2 /(2σ 2 ) be Gaussian, and let

r(ξ) := c -2 2 ψ| f 2 σ (H -ξI) |ψ = c -2 2 N -1 j=0 p j f 2 σ (E j -ξ). (35)
Then by the assumptions on g σ,ξ;ǫ ′′ , we have

|g 2 σ,ξ;ǫ ′′ (x) -f 2 σ (x -ξ)| = |g σ,ξ;ǫ ′′ (x) -f σ (x -ξ)| • |g σ,ξ;ǫ ′′ (x) + f σ (x -ξ)| (36) ≤ 2|g σ,ξ;ǫ ′′ (x) -f σ (x -ξ)| (37) ≤ 2ǫ ′′ , ∀x ∈ [-1, 1], (38)
and hence |q(ξ)r(ξ)| ≤ 2c -2 2 ǫ ′′ . Now since ξ i is uniformly distributed in [ Ẽ0w, Ẽ0 + w], we get that

P [A i = 1] = Ẽ0 +w Ẽ0 -w p(ξ i = ξ) • P [A i = 1|ξ i = ξ] dξ = 1 2w Ẽ0 +w Ẽ0 -w q(ξ)dξ. (39)
This implies each X k has probability density function

p(X k = ξ) = p(ξ i = ξ|A i = 1) (40) = p(ξ i = ξ) • P [A i = 1|ξ i = ξ] P [A i = 1] (41)
= q(ξ)

Ẽ0 +w Ẽ0 -w q(ξ)dξ , ∀ξ ∈ [ Ẽ0 -w, Ẽ0 + w]. (42)
Furthermore, the expected value of X k is

E [X k ] = Ẽ0 +w Ẽ0 -w ξ • p(X k = ξ)dξ = Ẽ0 +w Ẽ0 -w ξq(ξ)dξ Ẽ0 +w Ẽ0 -w q(ξ)dξ . (43)
Next, it is convenient to make a change of variables. Let s := Ẽ0 -E 0 . Then by Theorem 3.5, we know that s ∈ [-w/2, w/2] with probability at least 1δ/3. From now on, we will focus on the case when this event happens. Moreover, let x := ξ -E 0 and h(x) := q(E 0 + x) = q(ξ). Then we have x ∈ [sw, s + w] ⊆ [-3w/2, 3w/2]. Moreover, we have

P [A i = 1] = 1 2w s+w s-w h(x)dx, (44)
and

E [X k ] = E 0 + s+w s-w xh(x)dx s+w s-w h(x)dx . (45)
In Appendix B, we prove that

E j -ξ ≥ σ ln 0.5c 2 2 /ǫ ′′ (46)
for all j ≥ 1. Consequently,

f 2 σ (E j -ξ) ≤ 2ǫ ′′ c -2
2 for all j ≥ 1. Then by Eq. ( 35), we obtain

c 2 2 p -1 0 r(ξ) -f 2 σ (E 0 -ξ) = c 2 2 p -1 0 N -1 j=1 p j f 2 σ (E j -ξ) (47)
≤ c 2 2 p -1 0 max 1≤j≤N -1 f 2 σ (E j -ξ) (48)
≤ c 2 2 p -1 0 • 2ǫ ′′ c -2 2 (49) = 2p -1 0 ǫ ′′ , (50)
Meanwhile, recall that |q(ξ)r(ξ)| ≤ 2c -2 2 ǫ ′′ . Therefore, we have

c 2 2 p -1 0 q(ξ) -f 2 σ (E 0 -ξ) ≤ c 2 2 p -1 0 (q(ξ) -r(ξ)) + c 2 2 p -1 0 r(ξ) -f 2 σ (E 0 -ξ) (51) ≤ 2p -1 0 ǫ ′′ + 2p -1 0 ǫ ′′ (52) = 4p -1 0 ǫ ′′ (53) ≤ 4η -1 ǫ ′′ (54) = 0.03ǫ 1 σ ln (eσ/ǫ 1 ) (55)
where the fourth step follows from p 0 ≥ η and the last step follows from the definition of ǫ ′ . Now since h(x) = q(ξ) and

f 2 σ (E 0 -ξ) = f 2 σ (x), we get c 2 2 p -1 0 h(x) -f 2 σ (x) ≤ 0.03ǫ 1 σ ln (eσ/ǫ 1 ) , ∀x ∈ [s -w, s + w]. (56)
Thus, by Lemma 3.3 we know that

s+w s-w xh(x)dx s+w s-w h(x)dx ≤ 0.55ǫ 1 ≤ 0.5ǫ. ( 57
)
Combining Eqs. ( 45) and (57) yields that |E [X k ] -E 0 | ≤ ǫ/2 whenever |s| ≤ w/2, which happens with probability at least 1δ/3. So Eq. ( 28) holds. Furthermore, Lemma 3.3 also implies that s+w s-w c 2 2 p -1 0 h(x)dx ≥ σ. Combining this fact and Eq. (44) yields

P [A i = 1] ≥ p 0 σ 2c 2 2 w ≥ ησ 2c 2 2 w = η 4c 2 2 ln (eσ/ǫ 1 ) , ( 58
)
where the second step follows from p 0 ≥ η and the last steps follows from the definition of w. Then we have

E [K] = M P [A 1 = 1] ≥ M η 4c 2 2 ln (eσ/ǫ 1 ) = max(16w 2 ǫ -2 ln 6δ -1 , 8 ln 3δ -1 ), (59)
where the last step follows from the definition of M . Now we invoke the multiplicative Chernoff bound:

P K ≥ E [K] 2 ≥ 1 -e -E[K]/8 , (60)
and get

P K ≥ 8w 2 ǫ -2 ln 6δ -1 ≥ 1 -δ/3. (61)
So Eq. ( 29) is proved. Finally, since Ê0 = 1 K K k=1 X k where X 1 , X 2 , . . . , X k are i.i.d. and each X k takes values in [ Ẽ0w, Ẽ0 + w], by Hoeffding's inequality, we get

P | Ê0 -E [X 1 ] | ≥ ǫ/2 ≤ 2e -ǫ 2 K/(8w 2 ) . ( 62
)
If K ≥ 8w 2 ǫ -2 ln 6δ -1 , then the RHS of Eq. ( 62) is at most δ/3. So Eq. ( 30) is proved as well.

The theorem is thus proved.

Note that Algorithm 2 relies on block-encodings of f a,b;ǫ ′ (H) and g σ,ξ;ǫ ′′ (H) (which are approximate threshold and Gaussian functions of H respectively), and their implementation depends on how H is accessed. In the next two sections, we will give concrete implementation of these unitary operations in the block-encoding and Hamiltonian-evolution models respectively, and analyze the costs of the resulting GSEE algorithms in these models.

GSEE in the block-encoding model

In this section, we focus on GSEE in the block-encoding model. That is, we assume that there exists a method to implement an (1, m, 0)-block-encoding U of H for some integer m ≥ 0. Our quantum circuits will consist of multiple queries to U and U † as well as other primitive quantum gates, and act on 0 k |ψ for some integer k ≥ m.

We will develop a full-fledged version of Algorithm 2 in the block-encoding model. This means that we need to give concrete implementation of the block-encodings of f a,b;ǫ ′ (H) and g σ,ξ;ǫ ′′ (H) in Theorem 3.6 in this input model. To accomplish this, we utilize the technique of quantum singular value transformation (QSVT) [START_REF] Gilyén | Quantum singular value transformation and beyond: Exponential improvements for quantum matrix arithmetics[END_REF]. This technique allows us to implement any bounded low-degree polynomial of a linear operator A by using a block-encoding of A and its inverse a small number of times. Specifically, we will invoke Theorem 31 of Ref. [START_REF] Gilyén | Quantum singular value transformation and beyond: Exponential improvements for quantum matrix arithmetics[END_REF] which states that: We will use QSVT to approximately implement the threshold and Gaussian functions of the Hamiltonian H. To this end, we need low-degree polynomial approximations of such functions in the relevant regions. Precisely, for the threshold functions, we have the following result which follows directly from Corollary 6 of Ref. [START_REF] Hao Low | Hamiltonian simulation by uniform spectral amplification[END_REF] or Lemma 14 of Ref. [START_REF] Gilyén | Quantum singular value transformation and beyond: Exponential improvements for quantum matrix arithmetics[END_REF]: Lemma 4.2. For every a, b ∈ [-1, 1], a < b and ǫ ′ ∈ (0, 1), there exists an efficiently-computable

real polynomial f a,b;ǫ ′ (x) of degree O log(1/ǫ ′ ) b-a such that 1 -ǫ ′ ≤ f a,b;ǫ ′ (x) ≤ 1 for all x ∈ [-1, a], |f a,b;ǫ ′ (x)| ≤ 1 for all x ∈ [a
(2, m + 2, 0)-block-encoding V of f a,b;ǫ ′ (A) with O log(1/ǫ ′ ) b-a queries of U , U † , and O( (m+1) log(1/ǫ ′ ) b-a
) other primitive quantum gates.

Meanwhile, for the Gaussian functions, we have: Lemma 4.4. For every σ, ǫ ′′ ∈ (0, 1), there exists an efficiently-computable real polynomial

P σ;ǫ ′′ (x) of degree O(max(σ -1 , log (1/ǫ ′′ ))• log (1/ǫ ′′ )) such that |P σ;ǫ ′′ (x)| ≤ 1 and |P σ;ǫ ′′ (x)-e -x 2 /(2σ 2 ) | ≤ ǫ ′′ for all x ∈ [-π, π].
The proof of Lemma 4.4 is given in Appendix C. This lemma implies: Corollary 4.5. For every σ, ǫ ′′ ∈ (0, 1) and ξ ∈ [-2, 2], there exists an efficiently-computable real polynomial g σ,ξ;ǫ ′′ (x) of degree O(max(σ

-1 , log (1/ǫ ′′ )) • log (1/ǫ ′′ )) such that |g σ,ξ;ǫ ′′ (x)| ≤ 1 and |g σ,ξ;ǫ ′′ (x) -e -(x-ξ) 2 /(2σ 2 ) | ≤ ǫ ′′ for all x ∈ [-1, 1].
Combining Lemma 4.1 and Corollary 4.5 yields: Corollary 4.6. Suppose U is a (1, m, 0)-block-encoding of a Hermitian matrix A. For every σ, ǫ ′′ ∈ (0, 1) and ξ ∈ [-2, 2], let g σ,ξ;ǫ ′′ (x) be the polynomial in Corollary 4.5. Then we can implement a (2, m + 2, 0)-block-encoding W of g σ,ξ;ǫ ′′ (A) with O(max(σ -1 , log (1/ǫ ′′ )) • log (1/ǫ ′′ )) queries of U and U † , and O((m + 1) • max(σ -1 , log (1/ǫ ′′ )) • log (1/ǫ ′′ )) other primitive quantum gates.

Equipped with concrete implementation of f a,b;ǫ ′ (H) and g σ,ξ;ǫ ′′ (H), we are now ready to state the main result about GSEE in the block-encoding model: Theorem 4.7. Suppose H is an n-qubit Hamiltonian with unknown spectral decomposition

H = N -1 j=0 E j |E j E j |, where N = 2 n , -1 ≤ E 0 < E 1 ≤ E 2 ≤ ... ≤ E N -1 ≤ 1, E 1 -E 0 ≥ ∆ for known ∆ > 0,
and the |E j 's are orthonormal eigenstates of H. Moreover, suppose we can prepare an n-qubit state |ψ = N -1 j=0 γ j |E j such that |γ 0 | 2 ≥ η for known η ∈ (0, 1). Finally, suppose U is a (1, m, 0)-block-encoding of H. Then there exists a quantum algorithm that, for any given ǫ > 0 and δ ∈ (0, 1), estimates E 0 within additive error ǫ with probability at least 1δ by making use of

Õ(η -1 (∆ 2 ǫ -2 + log ∆ -1 ) log δ -1 )
quantum circuits, where each circuit uses

Õ(∆ -1 log(η -1 ∆ǫ -1 ))
queries of U and U † , and Õ((m + 1)∆ -1 log(η -1 ∆ǫ -1 ))

other primitive quantum gates.

Proof. We estimate the ground state energy of H by running Algorithm 2 with appropriate parameters. Specifically, Corollary 4.

′ : [-1, 1] → R satisfies that 1-ǫ ′ ≤ f a,b;ǫ ′ (x) ≤ 1 for all x ∈ [-1, a], |f a,b;ǫ ′ (x)| ≤ 1 for all x ∈ [a, b], and |f a,b;ǫ ′ (x)| ≤ ǫ ′ for all x ∈ [b, 1].
Moreover, Corollary 4.6 implies that for arbitrary σ, ǫ ′′ ∈ (0, 1) and ξ ∈ [-2, 2], a (2, m + 2, 0)-blockencoding of g σ,ξ;ǫ ′′ (H) can be implemented with O(max(σ -1 log (1/ǫ ′′ ), log (1/ǫ ′′ ))) queries of U and U † , and O((m + 1) max(σ -1 log (1/ǫ ′′ ), log (1/ǫ ′′ ))) other primitive quantum gates, where

g σ,ξ;ǫ ′′ : [-1, 1] → R satisfies that |g σ,ξ;ǫ ′′ (x)| ≤ 1 and |g σ,ξ;ǫ ′′ (x) -e -(x-ξ) 2 /(2σ 2 ) | ≤ ǫ ′′ for all x ∈ [-1, 1]
. Consequently, we can run Algorithm 2 with parameters (ǫ, δ, ∆, η, 2, m + 2, 2, m + 2), and by Theorem 3.6, its output is ǫ-close to E 0 with probability at least 1δ, as desired.

It remains to analyze the resource cost of ADV GSEE(ǫ, δ, ∆, η, 2, m + 2, 2, m + 2). If ǫ ≥ ∆/8, then we only need to run BASIC GSEE (ǫ, δ, η, 2, m + 2). This procedure makes use of

Õ(η -1 log ǫ -1 log δ -1 ) = Õ(η -1 log ∆ -1 log δ -1 )
quantum circuits, where each circuit implements a (2, m+2, 0)-block-encoding of f a,b;ǫ ′ (H) in which ba = Ω(ǫ) and ǫ ′ = Ω( √ η). We can construct this circuit by using

O log (1/ǫ ′ ) b -a = O(ǫ -1 log η -1 ) = O(∆ -1 log η -1 )
queries of U and U † , and

O (m + 1) log (1/ǫ ′ ) b -a = O((m + 1)ǫ -1 log η -1 ) = O((m + 1)∆ -1 log η -1 )
other primitive quantum gates. Otherwise, ǫ < ∆/8. In this case, we need to first run BASIC GSEE (w/2, δ/3, η, 2, m + 2), where

w ≥ 2σ = Ω ∆ log (η -1 ∆ǫ -1 ) . ( 63
)
As a result, BASIC GSEE (w/2, δ/3, η, 2, m + 2) makes use of

Õ(η -1 log w -1 log δ -1 ) = Õ(η -1 log ∆ -1 log η -1 ∆ǫ -1 log δ -1 )
quantum circuits, where each circuit implements a (2, m+2, 0)-block-encoding of f a,b;ǫ ′ (H) in which ba = Ω(w) and ǫ ′ = Ω( √ η). We can construct this circuit by using

O log (1/ǫ ′ ) b -a = O(w -1 log η -1 ) = O(∆ -1 log (η -1 ∆ǫ -1 ) log η -1 )
queries of U and U † , and

O (m + 1) log (1/ǫ ′ ) b -a = O((m + 1)w -1 log η -1 ) = O((m + 1)∆ -1 log (η -1 ∆ǫ -1 ) log η -1 )
other primitive quantum gates. Moreover, we also need to run

M = Õ(η -1 ∆ 2 ǫ -2 log δ -1 ) (64)
quantum circuits in the for-loop, where each circuit implements a (2, m + 2, 0)-block-encoding of g σ,ξ;ǫ ′′ (H) in which σ -1 = O(∆ -1 log (η -1 ∆ǫ -1 )) and (ǫ ′′ ) -1 = Õ(η -1 ∆ǫ -1 ). We can build each of these circuits by using

O(max(σ -1 log (1/ǫ ′′ ), log 1/ǫ ′′ )) = Õ(∆ -1 log(η -1 ∆ǫ -1 ))
queries of U and U † , and

O((m + 1) max(σ -1 log (1/ǫ ′′ ), log 1/ǫ ′′ )) = Õ((m + 1)∆ -1 log(η -1 ∆ǫ -1 ))
other primitive quantum gates. Overall, regardless of the value of ǫ, Algorithm 2 always makes use of

Õ(η -1 (∆ 2 ǫ -2 + log ∆ -1 ) log δ -1 )
quantum circuits, where each circuit uses Õ(∆ -1 log(η -1 ∆ǫ -1 ))

queries of U and U † , and Õ((m + 1)∆ -1 log(η -1 ∆ǫ -1 ))

other primitive quantum gates, as claimed.

Note that in Theorem 4.7, ∆ is merely a lower bound on the true spectral gap ∆ true := E 1 -E 0 of the Hamiltonian H. In fact, we can tune it between Θ(ǫ) and Θ(∆ true ), and doing so leads to a smooth interpolation between an algorithm with low circuit depth and an algorithm with Heisenberg-limit scaling of runtime with respect to ǫ: Corollary 4.8. Suppose H is an n-qubit Hamiltonian with unknown spectral decomposition

H = N -1 j=0 E j |E j E j |, where N = 2 n , -1 ≤ E 0 < E 1 ≤ E 2 ≤ ... ≤ E N -1 ≤ 1, ∆ true = E 1 -E 0 is
the spectral gap of H, and the |E j 's are orthonormal eigenstates of H. Moreover, suppose we can prepare an n-qubit state |ψ = N -1 j=0 γ j |E j such that |γ 0 | 2 ≥ η for known η ∈ (0, 1). Finally, suppose U is a (1, m, 0)-block-encoding of H. Then for every β ∈ [0, 1], there exists an algorithm that, for any given ǫ > 0 and δ ∈ (0, 1), estimates E 0 within additive error ǫ with probability at least 1δ by making use of

Õ η -1 ǫ -2+2β ∆ 2-2β true + log ǫ -β ∆ β-1 true log δ -1
quantum circuits, where each circuit uses

Õ(ǫ -β ∆ β-1 true log(η -1 ∆ 1-β true ǫ β-1 ))
queries of U and U † , and

Õ((m + 1)ǫ -β ∆ β-1 true log(η -1 ∆ 1-β true ǫ β-1 ))
other primitive quantum gates.

Proof. This corollary is the direct consequence of setting ∆ = Θ(ǫ β ∆ 1-β true ) for β ∈ [0, 1] in Theorem 4.7.

In particular, setting β = 0 or 1 in Corollary 4.8 yields:

• ∆ = Θ(∆ true ), for which we get an algorithm that makes use of Õ η -1 ǫ -2 ∆ 2 true + log ∆ -1 true log δ -1 quantum circuits, where each circuit uses

Õ(∆ -1 true log(η -1 ∆ true ǫ -1 ))
queries of U and U † , and

Õ((m + 1)∆ -1 true log(η -1 ∆ true ǫ -1 ))
other primitive quantum gates; or

• ∆ = Θ(ǫ)
, for which we get an algorithm that makes use of Õ η -1 log ǫ -1 log δ -1 quantum circuits, where each circuit uses Õ(ǫ -1 log(η -1 ))

queries of U and U † , and Õ((m + 1)ǫ -1 log(η -1 ))

other primitive quantum gates.

GSEE in the Hamiltonian-evolution model

In this section, we focus on GSEE in the Hamiltonian-evolution model. That is, we assume that there exists a method to implement the controlled time evolution of H, i.e. the unitary operation controlled-e iHt := |0 0| ⊗ I + |1 1| ⊗ e iHt . Each of our quantum circuits will consist of multiple queries to controlled-e iHt and controlled-e -iHt for some t = Θ(1) and other primitive quantum gates, and act on 0 k |ψ for some k ≥ 1.

We will develop a full-fledged version of Algorithm 2 in the Hamiltonian-evolution model. This means that we need to give concrete implementation of the block-encodings of f a,b;ǫ ′ (H) and g σ,ξ;ǫ ′′ (H) in Theorem 3.6 in this input model. To this end, we will utilize the technique of quantum eigenvalue transformation of unitary matrices with real polynomials (QET-U) [START_REF] Dong | Ground-state preparation and energy estimation on early fault-tolerant quantum computers via quantum eigenvalue transformation of unitary matrices[END_REF]. This technique allows us to implement certain low-degree trigonometric polynomials of H/2 by using controlled-e iH and controlled-e -iH a small number of times. Specifically, we will invoke Theorem 1 of Ref. [START_REF] Dong | Ground-state preparation and energy estimation on early fault-tolerant quantum computers via quantum eigenvalue transformation of unitary matrices[END_REF] which states that: It turns out that both the threshold and Gaussian functions can be approximated by such trigonometric polynomials. Specifically, for the threshold functions, Lemma 6 of Ref. [START_REF] Lin | Heisenberg-limited ground-state energy estimation for early faulttolerant quantum computers[END_REF] implies: Lemma 5.2. For every a, b ∈ [0, 1], a < b and ǫ ′ ∈ (0, 1), there exists an efficiently-computable even real polynomial Meanwhile, for the Gaussian functions, we have the following result: Lemma 5.4. For every σ, ǫ ′′ ∈ (0, 1), there exist an efficiently-computable even real polynomial

F a,b;ǫ ′ (x) of degree Õ log(1/ǫ ′ ) b-a such that 1 -ǫ ′ ≤ F a,b;ǫ ′ (cos(x/2)) ≤ 1 for all x ∈ [0, a], |F a,b;ǫ ′ (cos(x/2))| ≤ 1 for all x ∈ [a, b],
G σ;ǫ ′′ (x) of degree O(σ -1 T log (1/ǫ ′′ )), where T = max(2π, Θ(σ log (1/ǫ ′′ ))), such that |G σ;ǫ ′′ (cos(πx/T ))| ≤ 1 and |G σ;ǫ ′′ (cos(πx/T )) -e -x 2 /(2σ 2 ) | ≤ ǫ ′′ for all x ∈ [-π, π].
The proof of Lemma 5.4 is given in Appendix C. Combining Lemma 5.1 and Lemma 5.4 yields: Corollary 5.5. For every σ, ǫ ′′ ∈ (0, 1), let T and G σ,ǫ ′′ (x) be defined as in Lemma 5.4. Then we can implement a (1, 1, 0)-block-encoding of G σ;ǫ ′′ (cos(πH/T )) with O(σ -1 T log (1/ǫ ′′ )) queries of controlled-e i2πH/T and controlled-e -i2πH/T , and O(σ -1 T log (1/ǫ ′′ )) other primitive quantum gates.

Note that for any ξ, t ∈ R, we can implement controlled-e -i(H-ξI)t by using one controlled-e iHt and one Z-rotation gate on the control qubit. Combining this fact and Corollary 5.5, we conclude that:

Corollary 5.6. For every σ, ǫ ′′ ∈ (0, 1), let T and G σ,ǫ ′′ (x) be defined as in Lemma 5.4. Then for every ξ ∈ [-2, 2], we can implement a (1, 1, 0)-block-encoding of G σ;ǫ ′′ (cos(π(H -ξI)/T )) with O(σ -1 T log (1/ǫ ′′ )) queries of controlled-e i2πH/T and controlled-e -i2πH/T , and O(σ -1 T log (1/ǫ ′′ )) other primitive quantum gates.

Corollaries 5.3 and 5.6 provide the methods to approximately implement the threshold functions and Gaussian functions of H. This enables us to obtain the following result about GSEE in the Hamiltonian-evolution model: Theorem 5.7. Suppose H is an n-qubit Hamiltonian with unknown spectral decomposition

H = N -1 j=0 E j |E j E j |, where N = 2 n , 0 ≤ E 0 < E 1 ≤ E 2 ≤ ... ≤ E N -1 ≤ 1, E 1 -E 0 ≥ ∆ for known ∆ > 0,
and the |E j 's are orthonormal eigenstates of H. Moreover, suppose we can prepare an n-qubit state |ψ = N -1 j=0 γ j |E j such that |γ 0 | 2 ≥ η for known η ∈ (0, 1). Then there exists a quantum algorithm that, for any given ǫ > 0 and δ ∈ (0, 1), estimates E 0 within additive error ǫ with probability at least 1δ by making use of

Õ(η -1 (∆ 2 ǫ -2 + log ∆ -1 ) log δ -1 )
quantum circuits, where each circuit uses

Õ(∆ -1 log(η -1 ∆ǫ -1 ))
queries of controlled-e iHt and controlled-e -iHt for some t = Θ(1), and

Õ(∆ -1 log(η -1 ∆ǫ -1 ))
other primitive quantum gates.

Proof. We estimate the ground state energy of H by running Algorithm 2 with appropriate parameters. Specifically, Corollary 5.3 implies that for arbitrary a, b ∈ [0, 1], a < b, and ǫ ′ ≥ 0, a (1, 1, 0)-block-encoding of f a,b;ǫ ′ (H) := F a,b;ǫ ′ (cos(H/2)) can be implemented with Õ log(1/ǫ ′ ) b-a queries of controlled-e iH and controlled-e -iH , and Õ( log(1/ǫ ′ ) b-a ) other primitive quantum gates, where

f a,b;ǫ ′ satisfies that 1 -ǫ ′ ≤ f a,b;ǫ ′ (x) ≤ 1 for all x ∈ [0, a], |f a,b;ǫ ′ (x)| ≤ 1 for all x ∈ [a, b], and |f a,b;ǫ ′ (x)| ≤ ǫ ′ for all x ∈ [b, 1]
. Moreover, Corollary 5.6 implies that for arbitrary σ, ǫ ′′ ∈ (0, 1) and ξ ∈ [-2, 2], there exists T = max(2π, Θ(σ log (1/ǫ ′′ ))) such that a (1, 1, 0)-block-encoding of g σ,ξ;ǫ ′′ (H) := G σ,ǫ ′′ (cos(π(H -ξI)/T )) can be implemented with O(σ -1 T log (1/ǫ ′′ )) queries of controlled-e i2πH/T , controlled-e -i2πH/T , and O(σ -1 T log (1/ǫ ′′ )) other primitive quantum gates, where g σ,ξ;ǫ ′′ satisfies that |g σ,ξ;ǫ

′′ (x)| ≤ 1 and |g σ,ξ;ǫ ′′ (x) -e -(x-ξ) 2 /(2σ 2 ) | ≤ ǫ ′′ for all x ∈ [-1, 1].
Consequently, we can run Algorithm 2 with parameters (ǫ, δ, ∆, η, 1, 1, 1, 1), and by Theorem 3.6, its output is ǫ-close to E 0 with probability at least 1δ, as desired.

It remains to analyze the resource cost of ADV GSEE(ǫ, δ, ∆, η, 1, 1, 1, 1). If ǫ ≥ ∆/8, then we only need to run BASIC GSEE (ǫ, δ, η, 1, 1). This procedure makes use of

Õ(η -1 log ǫ -1 log δ -1 ) = Õ(η -1 log ∆ -1 log δ -1 )
quantum circuits, where each circuit implements a (1, 1, 0)-block-encoding of f a,b;ǫ ′ (H) in which ba = Ω(ǫ) and ǫ ′ = Ω( √ η). We can build this circuit by using

Õ log (1/ǫ ′ ) b -a = Õ(ǫ -1 log η -1 ) = Õ(∆ -1 log η -1 )
queries of controlled-e iH and controlled-e -iH , and

Õ log (1/ǫ ′ ) b -a = Õ(ǫ -1 log η -1 ) = Õ(∆ -1 log η -1 )
other primitive quantum gates. Otherwise, ǫ < ∆/8. In this case, we need to first run BASIC GSEE (w/2, δ/3, η, 1, 1), where

w ≥ 2σ = Ω ∆ log (η -1 ∆ǫ -1 ) . ( 65
)
As a result, BASIC GSEE (w/2, δ/3, η, 1, 1) makes use of

Õ(η -1 log w -1 log δ -1 ) = Õ(η -1 log ∆ -1 log η -1 ∆ǫ -1 log δ -1 )
quantum circuits, where each circuit implements a (1, 1, 0)-block-encoding of f a,b;ǫ ′ (H) in which ba = Ω(w) and ǫ ′ = Ω( √ η). We can build this circuit by using

Õ log (1/ǫ ′ ) b -a = Õ(w -1 log η -1 ) = Õ(∆ -1 log (η -1 ∆ǫ -1 ) log η -1 )
queries of controlled-e iH and controlled-e -iH , and

Õ log (1/ǫ ′ ) b -a = Õ(w -1 log η -1 ) = Õ(∆ -1 log (η -1 ∆ǫ -1 ) log η -1 )
other primitive quantum gates. Moreover, we also need to run

M = Õ(η -1 ∆ 2 ǫ -2 log δ -1 ) (66)
quantum circuits in the for-loop, where each circuit implements a (1, 1, 0)-block-encoding of g σ,ξ;ǫ ′′ (H) in which σ -1 = O(∆ -1 log (η -1 ∆ǫ -1 )) and (ǫ ′′ ) -1 = Õ(η -1 ∆ǫ -1 ). Furthermore, one can verify that σ log (1/ǫ ′′ ) = O(∆) = O(1) by our choice of parameters. This implies that we can construct each circuit in the for-loop by using

O(σ -1 T log (1/ǫ ′′ )) = Õ(∆ -1 log(η -1 ∆ǫ -1 ))
queries of controlled-e i2πH/T and controlled-e -i2πH/T for T = max(2π, Θ(σ log (1/ǫ ′′ ))) = Θ(1), and

O(σ -1 T log (1/ǫ ′′ )) = Õ(∆ -1 log(η -1 ∆ǫ -1 ))
other primitive quantum gates. Overall, regardless of the value of ǫ, Algorithm 2 always makes use of

Õ(η -1 (∆ 2 ǫ -2 + log ∆ -1 ) log δ -1 )
quantum circuits, where each circuit uses

Õ(∆ -1 log(η -1 ∆ǫ -1 ))
queries of controlled-e iHt and controlled-e -iHt for some t = Θ(1), and

Õ(∆ -1 log(η -1 ∆ǫ -1 ))
other primitive quantum gates, as claimed.

Note again that ∆ in Theorem 5.7 is just a lower bound on the true spectral gap ∆ true = E 1 -E 0 of the Hamiltonian H. By tuning it between Θ(ǫ) and Θ(∆ true ), we obtain a class of GSEE algorithms that smoothly transition from having low circuit depth to having Heisenberglimit scaling of runtime with respect to ǫ: Corollary 5.8. Suppose H is an n-qubit Hamiltonian with unknown spectral decomposition

H = N -1 j=0 E j |E j E j |, where N = 2 n , 0 ≤ E 0 < E 1 ≤ E 2 ≤ ... ≤ E N -1 ≤ 1, ∆ true = E 1 -E 0 is
the spectral gap of H, and the |E j 's are orthonormal eigenstates of H. Moreover, suppose we can prepare an n-qubit state |ψ = N -1 j=0 γ j |E j such that |γ 0 | 2 ≥ η for known η ∈ (0, 1). Then for every β ∈ [0, 1], there exists an algorithm that, for any given ǫ > 0 and δ ∈ (0, 1), estimates E 0 within additive error ǫ with probability at least 1δ by making use of

Õ η -1 ǫ -2+2β ∆ 2-2β true + log ǫ -β ∆ β-1 true log δ -1
quantum circuits, where each circuit uses

Õ(ǫ -β ∆ -1+β true log(η -1 ∆ 1-β true ǫ β-1 ))
queries of controlled-e iHt and controlled-e -iHt for some t = Θ(1), and

Õ(ǫ -β ∆ -1+β true log(η -1 ∆ 1-β true ǫ β-1 ))
other primitive quantum gates.

Proof. This corollary is the direct consequence of setting ∆ = Θ(ǫ β ∆ 1-β true ) for β ∈ [0, 1] in Theorem 5.7.

In particular, setting β = 0 or 1 yields:

• ∆ = Θ(∆ true ), for which we get an algorithm that makes use of Õ(η -1 ǫ -2 ∆ 2 true log δ -1 ) quantum circuits, where each circuit uses Õ(∆ -1 true log(η -1 ∆ true ǫ -1 )) queries of controllede iHt and controlled-e -iHt for some t = Θ(1), and Õ(∆ -1 true log(η -1 ∆ true ǫ -1 )) other primitive quantum gates; or • ∆ = Θ(ǫ), for which we get an algorithm that makes use of Õ(η -1 log ǫ -1 log δ -1 ) quantum circuits, where each circuit uses Õ(ǫ -1 log(η -1 )) queries of controlled-e iHt and controllede -iHt for some t = Θ(1), and Õ(ǫ -1 log(η -1 )) other primitive quantum gates.

Certifying the correctness of the energy estimate

In this section we introduce a method for certifying that the energy estimate output by Algorithm 2 is correct. This may be favorable in cases where we are not confident in the correctness of ∆, the lower bound on the energy gap. Furthermore, in the case that the input state has additional structure (i.e. has sufficiently small weight on low energy excited states) our algorithm can provide a correct estimate for evolution times that are arbitrarily smaller than the spectral gap. To illustrate this point, consider that the input state is given by the ground state. In this case, any constant evolution time can be used to obtain a reliable estimate. In both of these cases, it may be desirable to assume that the algorithm settings are correct, run Algorithm 2, and then certify the result. This way we do not need to necessarily trust our estimate of the spectral gap and can profit from better input states.

We describe a procedure that achieves exactly that. That is, given a lower-bound on the overlap of the input state with the ground state, access to samples from the distribution p * n σ , an error tolerance ǫ > 0 and an estimate of the ground state Ê0 , we will be able to certify that |E 0 -Ê0 | ≤ ǫ. It is also important to stress once more what our test will not achieve, namely certifying a spectral gap estimate. It only certifies the ground state energy estimate. However, we will show that if the evolution time is proportional to the inverse of the spectral gap, our test will accept the estimate with high probability. It is also worth mentioning that our test comes at the cost of a higher sample complexity than the estimation procedure itself, namely Õ(ǫ -4 η -3 ).

We give an intuitive explanation of how the verification works in Fig. 6 and explain the protocol in detail in Appendix D. But it essentially relies on two observations. The first, formalized in Lemma D.3, is that if the density of p * n σ has a well-defined peak around our ground state energy estimate, then p * n σ conditioned on an interval of size ∼ σ around Ê0 will again approximately have the moments of a convex combination of Gaussians. Furthermore, we then show that the variance of the random variable conditioned on this interval will exceed a given threshold if the estimate is incorrect. Finally, we show that if we evolve by a time proportional to the inverse of the spectral gap, the test will accept with high probability. We summarize the certification algorithm below: Theorem 6.1. Under the same conditions as Theorem 3.6, there exists an algorithm that, given arbitrary σ > 0 and Ê0 ∈ [E 0σ, E 0 + σ] (but no bound on the spectral gap ∆ true ), makes

Õ(η -3 ǫ -4 log δ -1 ) (67)
uses of the block-encoding of g σ,ξ;ǫ ′′ (H) in expectation and generates a new estimate Ê′

0 of E 0 . If | Ê′ 0 -E 0 | ≥ 4ǫ
, it and rejects Ê′ 0 with probability at least 1δ. If it accepts, then | Ê′ 0 -E 0 | ≤ ǫ, with probability of incorrect acceptance at most δ. Furthermore, there exists some

σ 0 = Ω ∆ true log (ǫ -1 η -1 ) (68)
such that if σ ≤ σ 0 , then this algorithm outputs an estimate Ê′ 0 that it accepts with probability 1δ and is such that

| Ê′ 0 -E 0 | ≤ ǫ.
We defer all the proofs and descriptions of the algorithm to Appendix D. Note that the minimal time for which the test accepts, in Eq. (68), is of the same order as when our algorithm is guaranteed to succeed. Unfortunately, the sample complexity of the test is higher than that of the algorithm at this depth. Thus, the test is only advantageous when we do not wish to trust our estimate of the spectral gap or have reasons to believe that the input is such that we can obtain a good estimate with smaller circuit depths. Our algorithm will be based on empirically estimating the variance of a random variable up to a certain precision. As is common with such protocols, depending on the number of samples the empirical estimate will have a tolerated error up to a certain failure probability. It is noteworthy that our protocol consistently provides accurate rejections when given variance estimates within the tolerated error. Therefore, the parameter δ quantifies the probability of the estimates deviating beyond the acceptable threshold for verification, and as a result, we cannot give guarantees on its performance when given wrong estimates.

Furthermore, the test can be combined with a binary search procedure to determine σ 0 , the largest σ for which the test accepts. Although this would yield reliable ground state energy estimates at depths that are potentially ≪ ∆ true , the sample complexity of such an algorithm would be Õ(ǫ -4 η -3 ). However, this might still be advantageous if we wish to optimize the maximal circuit depth. In this case, the density quickly goes to 0 around our estimate Ê0 . We then continue to check the variance conditioned around Ê0 to see if we accept the estimate.

Conclusion and outlook

To summarize, we have introduced a quantum algorithm for estimating the ground state energy of a Hamiltonian using fewer operations per circuit than previous methods. This makes the algorithm more suitable for early fault-tolerant quantum computation. This algorithm relies on a new way to perform rejection sampling on quantum computers. We have shown how it works with both the block-encoding and Hamiltonian-evolution input models of the Hamiltonian. The circuit depth required by our algorithm is governed by the spectral gap of the Hamiltonian rather than the target accuracy, and the runtime of our algorithm is smaller than that of previous algorithms with similar circuit depths. Furthermore, this algorithm can be tuned to trade circuit depth for runtime, and its runtime can be also reduced by parallelizing the sample collection across multiple quantum devices (i.e. trading time resources for space resources). All of these features make our algorithm a great candidate for realizing quantum advantage for industrially-relevant problems on early fault-tolerant quantum computers [START_REF] Wang | Accelerated variational quantum eigensolver[END_REF][START_REF] Wang | Minimizing estimation runtime on noisy quantum computers[END_REF][START_REF] Lin | Heisenberg-limited ground-state energy estimation for early faulttolerant quantum computers[END_REF][START_REF] Dong | Ground-state preparation and energy estimation on early fault-tolerant quantum computers via quantum eigenvalue transformation of unitary matrices[END_REF][START_REF] Zhang | Computing Ground State Properties with Early Fault-Tolerant Quantum Computers[END_REF][START_REF] Wang | State Preparation Boosters for Early Fault-Tolerant Quantum Computation[END_REF][START_REF] Wang | Quantum algorithm for ground state energy estimation using circuit depth with exponentially improved dependence on precision[END_REF][START_REF] Wang | Qubit-efficient randomized quantum algorithms for linear algebra[END_REF].

In addition, we also devised the first method for certifying ground state energy estimates in the regime of circuit depth Õ(∆ -1 true ). Although the sample complexity of the current method is high, it can provide certifiably correct estimates of the ground state energy from arbitrarily small circuit depth depending on the properties of the initial state.

There are several related research directions that are worth further exploration:

• In this work, we have focused on GSEE in the regime of circuit depth ≥ 1/∆ true (recall that ∆ true is the spectral gap of H around E 0 ). It remains unclear how hard this problem is to solve when the circuit depth (accounting for both the state preparation and estimation subroutines) is ≪ 1/∆ true ? Is it possible to devise an algorithm with reasonable runtime (e.g. poly (1/ǫ, 1/p 0 )) in this setting?

• In many applications, one needs to estimate ground state properties beyond their energies. So far, this problem was not extensively studied in the low-depth regime. One exception is Ref. [START_REF] Zhang | Computing Ground State Properties with Early Fault-Tolerant Quantum Computers[END_REF] which develops a low-depth algorithm for estimating the expectation value of an observable O with respect to the ground state of a Hamiltonian H. This problem is referred to as ground state property estimation (GSPE). In order to estimate E 0 | O |E 0 within additive error ǫ with high probability, the algorithm of Ref. [START_REF] Zhang | Computing Ground State Properties with Early Fault-Tolerant Quantum Computers[END_REF] requires circuit depth Õ(∆ -1 polylog ǫ -1 p -1 0

) and runtime Õ(∆ -1 ǫ -2 p -2 0 ) (where ∆ and p 0 are defined as before). Note that this runtime is worse than our runtime for ground state energy estimation. Since GSEE can be viewed as a special case of GSPE in which O = H, can we apply the techniques in this paper to improve the efficiency of the latter?

• Can we design more efficient methods for certifying the correctness of ground state energy estimates? Our current certifcation method has a quadratically worse sample complexity than our GSEE algorithm. If we could devise a certification method with the same sample complexity as the GSEE algorithm, we could run binary search to find the minimal evolution time required to obtain reliable estimates on the ground state energy without extra assumptions.

• Finally, we expect that our rejection sampling method could find more applications beyond GSEE. In many applications, one needs to generate samples from a continuous distribution X with probability density ν(x). The traditional approach is to prepare a quantum state that approximates Ω ν(x) |x dx and then to measure it. But this approach has two shortcomings. First, the quantum circuit for creating this state might be difficult to implement. Second, one needs to discretize the target distribution, and the number of grid points (= the dimension of the resulting system) depends on the desired accuracy. On the other hand, our approach does not rely on this discretization, as the sample value is encoded into the circuit parameter. Furthermore, our circuit might be easier to implement than the one in the direct approach and might require O(1) qubits regardless of the target accuracy. Provided that rejection sampling is efficient (i.e. not too many raw samples are discarded), our method could accelerate this central component of many algorithms or make it more amenable to near-term implementation.

We have contributed to the growing body of work developing quantum algorithms for early fault-tolerant quantum computers. Our hope is that these results will inspire further exploration to discover the quantum algorithms that will deliver advantage first.

A Proofs of Lemma 3.2 and Lemma 3.3

Proof of Lemma 3.2. Let us focus on the case s ≥ 0 here. The other case s < 0 can be handled similarly.

We will use the following bound on the complementary error function (erfc). For any z > 0,

erfc z = 2 √ π ∞ z e -t 2 dt = 2 π π/2 0 e -z 2 / sin 2 θ dθ ≤ 2 π π/2 0 e -z 2 dθ = e -z 2 . ( 69
)
It implies that

z 0 e -t 2 dt = ∞ 0 e -t 2 dt - ∞ z e -t 2 dt ≥ √ π 2 1 -e -z 2 . ( 70
)
Now by the definition of w and the assumption 0 ≤ s ≤ w/2, we have w + s ≥ ws ≥ w/2 > 0. Consequently, we get

≥ √ πσ 1 -e -w 2 /(4σ 2 ) (74) ≥ √ πσ 1 - ǫ eσ (75) ≥ √ πσ 1 -e -1 (76) ≥ 1.12σ, (77)
where the sixth step follows from 0 < ǫ ≤ σ. Meanwhile, a direct calculation indicates

s+w s-w xe -x 2 /σ 2 dx = σ 2 2 e -(w-s) 2 /σ 2 -e -(w+s) 2 /σ 2 (78) ≤ σ 2 2 e -(w-s) 2 /σ 2 (79) = σ 2 2 e -w 2 /(4σ 2 ) (80) ≤ σǫ 2e . ( 81
)
The lemma is thus proved.

Proof of Lemma 3.3. As in the proof of Lemma 3.2, we will focus on the case s ≥ 0 here, and the other case s < 0 can be handled similarly. On one hand, we have

s+w s-w g(x)dx ≥ s+w s-w e -x 2 /σ 2 dx - s+w s-w |e -x 2 /σ 2 -g(x)|dx (82) ≥ 1.12σ -2w • 0.03ǫ σ ln (eσ/ǫ) (83) = 1.12σ - 0.12ǫ ln (eσ/ǫ) (84) ≥ 1.12σ -0.12ǫ (85) ≥ σ, (86)
where the first step follows from triangle inequality, the second step follows from Lemma 3.2 and the assumption on g(x), the third step follows from the definition of w, and the fourth and fifth steps follow from 0 < ǫ ≤ σ.

On the other hand, we also have x(e -x 2 /σ 2g(x))dx (87)

≤ σǫ 2e + 2w • 3w 2 • 0.03ǫ σ ln (eσ/ǫ) (88)
≤ σǫ 2e + 0.36σǫ (89)
≤ 0.55σǫ, (90)
where the first step follows from triangle inequality, the second step follows from Lemma 3.2, the assumption on g(x) and the fact that |x| ≤ w + s ≤ 3w/2 for all x ∈ [sw, s + w], the fourth step follows from the definition of w.

The lemma is thus proved.

B Proof of Eq. [START_REF] Loaiza | Reducing molecular electronic hamiltonian simulation cost for linear combination of unitaries approaches[END_REF] We continue to use the notation in the proof of Theorem 3.6. Recall that E 0ξ ∈ [-3w/2, 3w/2] and E j -E 0 ≥ ∆ for all j ≥ 1. So we have

E j -ξ = (E j -E 0 ) + (E 0 -ξ) (91) ≥ ∆ - 3w 2 (92) = ∆ -3σ ln (eσ/ǫ 1 ) (93)
= ∆ 1 -3 ln (eσ/ǫ 1 ) 5 ln (ec 2 ) ln (eη -1 ∆ǫ -1 ) .

Meanwhile, by the definition of σ and ǫ ′′ , we have σ ln 0.5c 2 2 /ǫ ′′ = ∆ ln 200c 2 2 σ ln (eσ/ǫ 1 ) /(3ǫ 1 η) 5 ln (ec 2 ) ln (eη -1 ∆ǫ -1 ) .

Therefore, in order to prove Eq. ( 46), it is sufficient to show that 5 ln (ec 2 ) ln (eη -1 ∆ǫ -1 ) ≥ 3 ln (eσ/ǫ 1 ) + ln 200c 2 2 σ ln (eσ/ǫ 1 ) /(3ǫ 1 η) .

To achieve this, we consider two cases separately:

• Case 1: ǫ 1 = σ. In this case, we have 3 ln (eσ/ǫ 1 ) + ln 200c 2 2 σ ln (eσ/ǫ 1 ) /(3ǫ 1 η) = 3 + ln 200c 2 2 /(3η) (97)

≤ 5 ln 3ec 2 2 η -1 (98)
≤ 5 ln (ec 2 ) ln (3eη -1 ) (99)

≤ 5 ln (ec 2 ) ln (eη -1 ∆ǫ -1 ), ( 100
)
where the second step follows from the facts that c 2 , η -1 ≥ 1 and 3 + ln (200t/3) ≤ 5 ln (3et) for all t ≥ 1, the third step follows from the facts that c 2 , η -1 ≥ 1 and ln a 2 t ≤ ln (ea) • ln (t) for all a ≥ 1 and t ≥ 3e, and the last step follow from the assumption ∆ǫ -1 > 8. So Eq. ( 96) holds.

• Case 2:

ǫ 1 = ǫ 1.1 < σ. In this case, let A := η -1 ∆ǫ -1 ≥ ∆ǫ -1 ≥ 8.
where the second step follows from ∆ǫ -1 ≤ A and c 2 ≥ 1, and ln 200c

C Polynomial and trigonometric approximations of Gaussian functions

In this appendix, we prove Lemma 4.4 and Lemma 5.4 which state that there exist low-degree polynomial and trigonometric approximations of Gaussian functions, respectively.

Proof of Lemma 4.4. By Theorem 4.1 of Ref. [START_REF] Sachdeva | Faster algorithms via approximation theory[END_REF], for every b > 0 and γ ∈ (0, 1), there exists a real polynomial R b,γ (x) of degree O( max(b, log (1/γ)) • log (1/γ)) such that |e -x -R b,γ (x)| ≤ γ for all x ∈ [0, b]. Furthermore, Ref. [START_REF] Sachdeva | Faster algorithms via approximation theory[END_REF] also gives an explicit construction of R b,γ (x). Now we replace x with x 2 /(2σ 2 ) and set b = π 2 /(2σ 2 ), γ = ǫ ′′ /3 in the above result, and obtain that: for every σ > 0 and ǫ ′′ ∈ (0, 1), there exists an efficiently-computable real polynomial Pσ;ǫ ′′ (x) of degree O(max(σ -1 , log (1/ǫ ′′ )) • log (1/ǫ ′′ )) such that |e -x 2 /(2σ 2 ) -Pσ;ǫ ′′ (x)| ≤ ǫ ′′ /3 for all x ∈ [-π, π]. Note that this implies | Pσ;ǫ ′′ (x)| ≤ 1 + ǫ ′′ /3 for all x ∈ [-π, π]. Now let P σ;ǫ ′′ (x) := (1ǫ ′′ /3) Pσ;ǫ ′′ (x). Then we have that for every x ∈ [-π, π],

|P σ;ǫ ′′ (x)| = (1 -ǫ ′′ /3)| Pσ;ǫ ′′ (x)| ≤ (1 -ǫ ′′ /3)(1 + ǫ ′′ /3) ≤ 1, (108)
and P (x)e -x 2 /(2σ 2 ) ≤ P (x) -P (x) + P (x)e -x 2 /(2σ 2 ) (109)

≤ ǫ ′′ 3 • | P (x)| + ǫ ′′ 3 (110) ≤ ǫ ′′ 3 1 + ǫ ′′ 3 + ǫ ′′ 3 (111) ≤ ǫ ′′ . ( 112
)
Thus P σ;ǫ ′′ (x) satisfies all the requirements of Lemma 4.4.

Proof of Lemma 5.4. We claim that for every σ, ǫ ′′ ∈ (0, 1), there exist T = max(2π, Θ(σ log (1/ǫ ′′ ))) and N = O(σ -1 T log (1/ǫ ′′ )) and a 0 , a 1 , . . . , a N ∈ R + such that N j=0 a j ≤ 1 and

N j=0 a j cos(2πjx/T ) -e -x 2 /(2σ 2 ) ≤ ǫ ′′ , ∀x ∈ [-π, π]. (113)
If this claim holds, then since cos(2πjx/T ) = T 2j (cos(πx/T )), where T 2j (x) is the 2j-th Chebyshev polynomial of the first kind, N j=0 a j cos(2πjx/T ) can be written as G σ;ǫ ′′ (cos(πx/T )) for some even real polynomial G σ;ǫ (x) of degree 2N = O(σ -1 T log (1/ǫ ′′ )). So we get

G σ;ǫ ′ (cos(πx/T )) -e -x 2 /(2σ 2 ) ≤ ǫ ′′ , ∀x ∈ [-π, π]. (114)
Furthermore, since a 0 , a 1 , . . . , a N ≥ 0, we have

G σ;ǫ ′′ (cos(πx/T )) = N j=0 a j cos(2πjx/T ) ≤ N j=0 a j ≤ 1, ∀x ∈ R. ( 115
)
Therefore G σ;ǫ ′′ (x) satisfies all the requirements of Lemma 5.4. The proof of the above claim is inspired by the Nyquist-Shannon sampling theorem. Specifically, let f (x) := ∞ k=-∞ δ(x -kT ) be the Dirac comb with period T . It has Fourier transform f (ξ) =

1 T ∞ n=-∞ δ(ξ -n/T ).
Meanwhile, let g(x) := e -x 2 /(2σ 2 ) be Gaussian with mean 0 and variance σ 2 . It has Fourier transform ĝ(ξ) = √ 2πσe -2π 2 σ 2 ξ 2 . Convolving f and g yields

(f * g)(x) = ∞ k=-∞ g(x -kT ). ( 116
)
On the other hand, the Fourier transform of f * g is

f * g(ξ) = f (ξ)ĝ(ξ) = 1 T ∞ n=-∞ ĝ(n/T )δ(ξ -n/T ). ( 117
)
Applying inverse Fourier transform to both sides of this equation leads to

(f * g)(x) = 1 T ∞ n=-∞ ĝ(n/T )e i2πxn/T . (118)
Comparison of Eqs. ( 116) and ( 118) indicates that

∞ k=-∞ g(x -kT ) = 1 T ∞ n=-∞ ĝ(n/T )e i2πxn/T . ( 119
)
Now we claim that for some T = max(2π, Θ(σ log (1/ǫ ′′ ))) and

N = O(σ -1 T log (1/ǫ ′′ )), -1 k=-∞ g(x -kT ) + ∞ k=1 g(x -kT ) ≤ ǫ ′′ 6 , ∀x ∈ [-π, π], (120) and 1 T
-N -1 n=-∞ ĝ(n/T )e i2πxn/T + 1 T ∞ n=N +1 ĝ(n/T )e i2πxn/T ≤ ǫ ′′ 6 , ∀x ∈ R. (121)
If these claims are true, then combining them and Eq. (119) yields

g(x) - 1 T N n=-N ĝ(n/T )e i2πxn/T = g(x) -ĝ(0) - 2 T N n=1 ĝ(n/T ) cos(2πxn/T ) ≤ ǫ ′′ 3 , ∀x ∈ [-π, π]. (122)
This also implies that

ĝ(0) + 2 T N n=1 ĝ(n/T ) ≤ 1 + ǫ ′′ 3 . (123)
Now let a 0 = (1ǫ ′′ /3)ĝ(0) and a j = 2T -1 (1ǫ ′′ /3)ĝ(j/T ) for j = 1, 2, . . . , N . Then we have a 0 , a 1 , . . . , a N > 0 and

N j=0 a j = (1 -ǫ ′′ /3) ĝ(0) + 2 T N n=1 ĝ(n/T ) ≤ (1 -ǫ ′′ /3)(1 + ǫ ′′ /3) < 1, (124)
and for any

x ∈ [-π, π], N j=0 a j cos(2πjx/T ) -g(x) ≤ N j=0 a j cos(2πjx/T ) -ĝ(0) - 2 T N n=1 ĝ(n/T ) cos(2πxn/T ) + ĝ(0) + 2 T N n=1 ĝ(n/T ) cos(2πxn/T ) -g(x) (125)
≤ ǫ ′′ 3 • ĝ(0) + 2 T N n=1 ĝ(n/T ) cos(2πxn/T ) + ǫ ′′ 3 (126) ≤ ǫ ′′ 3 ĝ(0) + 2 T N n=1 ĝ(n/T ) + ǫ ′′ 3 (127) ≤ ǫ ′′ 3 1 + ǫ ′′ 3 + ǫ ′′ 3 (128) ≤ ǫ ′′ , (129)
as desired.

It remains to prove Eqs. ( 120) and (121). To prove the former, we impose the constraint T ≥ 2π. This means that |x| ≤ T /2 for all

x ∈ [-π, π]. Then it follows that -1 k=-∞ g(x -kT ) + ∞ k=1 g(x -kT ) ≤ 2 ∞ l=0 g((l + 1/2)T ) (130) = 2 ∞ l=0 exp -(l + 1/2) 2 T 2 /(2σ 2 ) (131) = 2 exp -T 2 /(8σ 2 ) ∞ l=0 exp -(l 2 + l)T 2 /(2σ 2 ) (132) ≤ 2 exp -T 2 /(8σ 2 ) ∞ l=0 exp -lT 2 /σ 2 (133) = 2 exp -T 2 /(8σ 2 ) 1 -exp (-T 2 /σ 2 ) . (134)
By picking some T = Θ(σ log (1/ǫ ′′ )), we can ensure that the RHS of Eq. ( 134) is at most ǫ ′′ /6. Overall, Eqs. (120) holds as long as we pick some T = max(2π, Θ(σ log (1/ǫ ′′ ))).

To prove Eq. ( 121), we note that

1 T -N -1 n=-∞ ĝ(n/T )e i2πxn/T + 1 T ∞ n=N +1 ĝ(n/T )e i2πxn/T ≤ 1 T -N -1 n=-∞ ĝ(n/T ) + 1 T ∞ n=N +1 ĝ(n/T ) (135) = 2 T ∞ n=N +1 ĝ(n/T ) (136) = 2 √ 2πσ T ∞ n=N +1 exp -2π 2 σ 2 n 2 /T 2 (137) ≤ 2 √ 2πσ T ∞ N exp -2π 2 σ 2 y 2 /T 2 dy (138) = 2 √ π ∞ √ 2πσN/T exp -z 2 dz (139) = erfc( √ 2πσN/T ) (140) ≤ exp -2π 2 σ 2 N 2 /T 2 . ( 141
)
By picking some N = O(σ -1 T log (1/ǫ ′′ )), we can guarantee that exp -2π 2 σ 2 N 2 /T 2 ≤ ǫ ′′ /6 and hence Eq. ( 121) holds for all x ∈ R.

The lemma is thus proved.

D Certification of the ground state energy estimate

Our algorithm has the desirable feature that it only requires us to consider circuit depths that are proportional to the spectral gap of the Hamiltonian. However, in many applications of practical interest we may not know the spectral gap in advance and would like to check if the output of the algorithm was correct. Furthermore, as we will explain soon, depending on the initial state, it might be the case that the algorithm still outputs the right estimate even if we evolve for times that are smaller than the spectral gap. If we suspect we have such a structured initial state, it is desirable to be able to verify that the output was correct.

We will now outline a test such that, given σ > 0, an estimate Ê0 of the ground state energy E 0 up to σ, a precision tolerance ǫ, a failure probability δ > 0 and a lower bound on the overlap p 0 ≥ η, accepts only if the estimate satisfies | Ê0 -E 0 | ≤ ǫ, with a failure probability of wrong certification at most δ. The quantum circuits required to run the algorithm are of depth O(σ -1 ). Furthermore, we will also show that if σ is proportional to the spectral gap ∆ true , it is guaranteed to accept with probability at least 1δ (although it might also accept at a larger value of σ/smaller circuit depth for more structured states).

Note that we will not certify that the spectral gap of the Hamiltonian is at least some ∆, but only that the ground state energy estimate is correct up to ǫ under the promise that the initial state had an overlap with the ground space of at least η > 0. Certifying the spectral gap would require additional assumptions on the initial state and its support on the first excited state. To make this point clear, consider the case in which the initial state is ρ = |E 0 E 0 |, i.e. the ground state. In this case, it is not difficult to see that our algorithm already outputs the correct estimate with access to controlled evolution times of order O(1), as E 0 will be the maximum of p * n σ for all σ > 0, where n σ is again a Gaussian with standard deviation σ. Thus, given such an initial state and enough samples, our algorithm will always output the correct estimate regardless of the circuit depth. On the other hand, it is clear that the state ρ cannot give us any information about the spectral gap of the Hamiltonian. In short, our test will only certify that the estimate of the energy was within a tolerated error.

Intuitively speaking, our algorithm outputs a correct estimate whenever σ is small enough to make sure that the contribution of the other energies in n σ at E 0 is sufficiently suppressed to make sure that we have a peak close to the ground state. If we have pollution from other energies, then they will make the peak around E 0 wider than if not. Thus, our test will consist on first making sure that we have a peak close to E 0 and then measuring the width of the peak.

We will first collect an auxilliary result about mixtures of Gaussians, as p * n σ is a mixture of Gaussians, i.e. a convex combination of Gaussian random variables. The Lemma below formalizes the idea that if we have pollution from other energies, then the peak around E 0 will be wider than it should be, where we measure the width by the variance.

Lemma D.1. Let X = i p i N i , where N i ∼ N (E i , σ 2 ) are normal random variables with the same variance σ 2 but different means E i . Furthermore, assume E 0 < E 1 < • • • < E m and that p 0 ≥ η. If |E(X) -E 0 | ≥ cǫ for some constants c > 1, ǫ > 0, then: E(|X -E(X)| 2 ) ≥ σ 2 + ηc 2 ǫ 2 (142) Proof. Recall the standard identity E(|X -E(X)| 2 ) = E(X 2 ) -E(X) 2 . Let us compute these expectation values. Clearly, E(X) = i p i E i . Furthermore, E(X 2 ) = i p i E(N 2 i ) = i p i (σ 2 + E 2 i ). (143)
Thus,

E(X 2 ) -E(X) 2 = σ 2 + i p i E 2 i -( i p i E i ) 2 . Now note that i p i E 2 i -( i p i E i ) 2
is nothing but the variance of a discrete random variable Y that takes value E i with probability p i and E(Y ) = E(X). Then

E(|Y -E(Y )| 2 ) ≥ ηc 2 ǫ 2 . ( 144
)
To see this, note that by our assumption that |E(X) -E 0 | ≥ cǫ, with probability p 0 the difference |Y -E(Y )| 2 will be at least c 2 ǫ 2 , which yields the claim.

The result above shows that if the mean of a Gaussian mixture deviates significantly from the minimum mean of its components, then this will result in a variance that is larger than their individual variances σ 2 .

Our test will be based on "zooming into" the peak around the estimate we found with our sampling algorithm and testing whether the variance is larger than expected. That is, we will condition the distribution around Ê0 and estimate its variance.

But it is not a given that conditioning around Ê0 will also yield a distribution that is a mixture of Gaussians. Thus, we will show that this is approximately the case if we have a well-defined peak in the distribution.

Before we prove our theorem showing how conditioning on an interval yields another mixture of Gaussians, let us again give some intuition. Given some peak Ê0 that we identified, we will truncate to some interval [ Ê0 -L, Ê0 + L], with

L = 2c log (σ -1 τ -1 )σ (145)
for some constants c > 1, 1 > τ > 0 to be specified later. The mass of energies that are located to the right of E 0 + 3L is negligible on the interval [ Ê0 -L, Ê0 + L], as it is at most O(Lτ 2c ). We see that we can safely discard these random variables when truncating to the interval. Furthermore, we will also approximate well all Gaussians whose means are located in [ Ê0 -L/2, Ê0 + L/2], as their mass outside this interval is O(τ c 4 ). Thus, the only energies that "pose a threat" of pollution are those that have substantial weight on the edge of our interval [ Ê0 -L, Ê0 + L]. This is where we use the assumption that we found a peak: we will impose that the mass going around the edges of the intervals is small. Note that only the right edge of the interval has to be considered, as we assume in our algorithm that we know the ground state energy up to σ and, thus, there are no energies to the left of the interval. More formally, we want:

Ê0 +2L Ê0 + L 2 (p * n σ )(x)dx ≤ c ′ ǫ 2 η 2 . ( 146
)
for some constant c ′ we will choose later. We will call distributions that satisfy this property peaked around Ê0 .

Definition D.2 (Distribution peaked around value). Given a point distribution p on R, ǫ > 0, η > 0 and n σ the probability density function of a Gaussian random variable with standard deviation σ, we say that p * n σ is peaked around a point Ê0 up to ǫ with parameters c, c ′ , τ > 0 if for L = σ 2c log (σ -1 τ -1 ) we have that:

Ê0 +2L Ê0 + L 2 (p * n σ )(x)dx ≤ c ′ ǫ 2 η 2 . ( 147
)
To further simplify notation we will define n σ,i to be the density of the random variable N (E i , σ 2 ). We will start by proving that whenever Eq. ( 146) holds, the mass on [ Ê0 -L, Ê0 + L] is mostly given by energies in the interval [ Ê0 -L/2, Ê0 + L/2].

Lemma D.3. Let X = i p i N i , where N i ∼ N (E i , σ 2 ) are normal random variables with the same variance σ 2 but different means E i and L be defined as in Eq. (145). For some Ê0 that is whithin σ of E 0 , group the means of the N i into three subsets: (152)

G = {i : |E i -Ê0 | ≤ L 2 }, (148)
B = {i : L 2 ≤ |E i -Ê0 | ≤ 2L}, (149)
We will estimate each term in the sum separately. Let us start with energies in F . As, by definition, all energies in that set are at least 2L away from Ẽ0 , they are at least L away from the upper end of the integration limits. Thus, by our choice of L we have for i ∈ F : To see this, note that by definition of B the mean of these Gaussians are greater than Ê0 + L 2 and the integral on the r.h.s. is over a larger interval that is closer to their mean. Thus, as the mass of the Gaussian is maximal around the mean, the inequality follows. By our assumption in Eq. ( 146) we conclude that:

i∈B p i Ê0 +L Ê0 -L n σ,i (x)dx ≤ 2c ′ ǫ 2 η 2 . ( 158
)
Putting together inequalities in Eqs. ( 154), ( 156) and (158), Eq. (151) follows.

We can show in a similar fashion that the following holds: Proof. We can apply the same estimates as we used to show Lemma D.3 to show that the contributions from energies in B and F are small. This is because we assume that [ Ê0 -L, Ê0 + L] ⊂ [-1, 1] and for x ∈ [-1, 1]:

|xn σ,i | ≤ n σ,i , x 2 n σ,i ≤ n σ,i (161)
Furthermore, it is also not difficult to show that for i ∈ G the value of the first two moments are approximated if we integrate over [ Ê0 -L, Ê0 + L].

This shows that under Eq. ( 146) the moments of the random variable conditioned on the interval [ Ê0 -L, Ê0 + L] are indeed well-approximated by those of the energies in the interval [ Ê0 -L/2, Ê0 + L/2]. Thus, if we check that Eq. ( 146) holds and we estimate the moments conditioned on the interval, then we can approximate the quantity in Eq. ( 142) and check that the variance is larger than expected if we got the energy right up to cǫ.

D.1 Certification algorithm and its correctness

With these technical results at hand, let us give a description of our energy certification in Algorithm 3, which we refer to as GSEE CERT . We deliberately formulated it in a way independent Algorithm 3 Energy certification algorithm algorithm GSEE CERT . return Accept estimate and return M . return Reject estimate. return Reject estimate.

16:

end if 17: end procedure of how we obtained the estimates for the relevant quantities. This is because different algorithms for GSEE have different access models to the function p * n σ , be it samples, as in this work, or to estimates of the density, as in Ref. [START_REF] Wang | Quantum algorithm for ground state energy estimation using circuit depth with exponentially improved dependence on precision[END_REF]. Furthermore, we assume that we are certain about the correctness of the estimates used in the algorithm for simplicity. But it is straightforward to extend the analysis to the case were the estimates are only correct up to some given failure probability. We then have: Theorem D.5. If Ê0 is such that | Ê0 -E 0 | ≤ σ, Algorithm 3 accepts, and all the estimates used were correct, then we have:

|M -E 0 | ≤ 4ǫ. ( 162
)
Proof. If the algorithm accepted, then we already know that the density p * n σ is peaked around Ê0 . Otherwise we would have already rejected in the first at Line 15. Let us assume that

| Ê0 -E 0 | > 4ǫ. ( 163
)
We will now show that the test will reject. As before, let X = i p i N (E , σ 2 ) and X be the same random variable conditioned on the interval [ Ê0 -L, Ê0 -L]. Furthermore, also let as before G be the set of energies E i ∈ [ Ê0 -L/2, Ê0 -L/2]. By Lemma D.3 and Cor. D.4, we have that the variance of X approximates up to ǫ 2 η 160 that of the random variable Y = i∈G p i N (E i , σ 2 ) / i∈G p i . By Lemma D.1, if Eq. (163) holds, then the variance of X conditioned on [ Ê0 -L, Ê0 + L] must exceed

E(|X -Ê0 | 2 |X ∈ [ Ê0 -L, Ê0 -L]) ≥ σ 2 + 16ηǫ 2 - ηǫ 2 160 > σ 2 + 2ǫ 2 η, (164)
which means that the test will have rejected. Thus, if we accept, then Eq. (162) must hold.

Let us now show that the test will always accept if we evolve for a time proportional to the gap over the precision: Lemma D.6. Let ∆ true be the spectral gap of H, p 0 ≥ η, a desired precision ǫ > 0 and | Ê0 -E 0 | ≤ σ. Let σ 0 be the minimal σ > 0 such that Algorithm 3 accepts for a given p. Then:

σ 0 ≤ ∆ true 10 log (2η -1 ǫ -1 ) . ( 165
)
Proof. Let us first show that for our value of σ 0 in Eq. ( 165), the distribution will be peaked around E 0 . By our choice of σ 0 , it is easy to see that the contribution of the other energies on [E 0 -2L, E 0 + 2L] will be of order O((ǫη) 50 ), which immediately implies that the distribution is peaked. Furthermore, the probability distribution on [E 0 -2L, E 0 + 2L] is well-approximated by a single Gaussian up to corrections of order O((ǫη) 50 ). Thus, the variance conditioned on [E 0 -L, E 0 + L] will deviate from σ by less than ǫ 2 η. As our estimate S of the variance deviates from the true one by at most ǫ 2 η 160 , we have |Sσ 2 | ≤ 2ǫ 2 η, which means we will accept.

We conclude that the test is sound: in the worst case, if we pick σ of similar magnitude as the spectral gap, we are guaranteed to accept. But it can be the case that σ 0 ≪ ∆ true depending on the initial state. To illustrate this, take once again the case where p 0 = 1. Then the test will always accept, irrespective of the value of the gap. This is because we will always only have one Gaussian component in n σ and, thus, the variance will always be σ 2 .

Lemma 3 . 2 .

 32 Let σ > 0 and ǫ ∈ (0, σ] be arbitrary, and let w = 2σ ln (eσ/ǫ). Then for any s ∈ [-w/2, w/2], we have s+w s-w e -x 2 /σ 2 dx ≥ 1.12σ and s+w s-w xe -x 2 /σ 2 dx ≤ σǫ 2e , and hence s+w s-w xe

  g(x)dx≤ 0.55ǫ.

12 :

 12 for j = 1, 2, . . . , M do 13: Measure the first m qubits of V a,b;ǫ ′ |0 m |ψ in the standard basis; 14:if the measurement outcome is 0 m then 15:K ← K + 1; if K < 0.5c -1 ηM then 19:

11 :S

 11 ← empty array; ⊲ S may contain duplicate values 12:

Lemma 4 . 1 .

 41 Suppose U is a (1, m, 0)-block-encoding of a Hermitian matrix A. Let P ∈ R[x] be a degree-d real polynomial such that |P (x)| ≤ 1 for all x ∈ [-1, 1]. Then we can implement a (2, m + 2, 0)-block-encoding Ũ of P (A) with 2d queries of U and U † , and O((m + 1)d) other primitive quantum gates.

  , b], and |f a,b;ǫ ′ (x)| ≤ ǫ ′ for all x ∈ [b, 1]. Combining Lemma 4.1 and Lemma 4.2 yields: Corollary 4.3. Suppose U is a (1, m, 0)-block-encoding of a Hermitian matrix A. For every a, b ∈ [-1, 1], a < b and ǫ ′ ∈ (0, 1), let f a,b;ǫ ′ (x) be the polynomial in Lemma 4.2. Then we can implement a

  3 implies that for arbitrary a, b ∈ (-1, 1) and ǫ ′ ≥ 0, a (2, m+2, 0)-block-encoding of f a,b;ǫ ′ (H) can be implemented with O log(1/ǫ ′ ) b-a queries of U and U † , and O( (m+1) log(1/ǫ ′ ) b-a ) other primitive quantum gates, where f a,b;ǫ

Lemma 5 . 1 .

 51 For every even real polynomial F (x) of degree d satisfying |F (x)| ≤ 1 for all x ∈ [-1, 1], we can implement a (1, 1, 0)-block-encoding of F (cos(H/2)) with O(d) queries of controllede iH and controlled-e -iH , and O(d) other primitive quantum gates.

Corollary 5 . 3 .

 53 and |F a,b;ǫ ′ (cos(x/2))| ≤ ǫ ′ for all x ∈ [b, 1]. Combining Lemma 5.1 and Lemma 5.2 yields: For every a, b ∈ [0, 1], a < b and ǫ ′ ∈ (0, 1), let F a,b;ǫ ′ (x) be the polynomial in Lemma 5.2. Then we can implement a (1, 1, 0)-block-encoding of F a,b;ǫ ′ (cos(H/2)) with Õ log(1/ǫ ′ ) b-a queries of controlled-e iH and controlled-e -iH , and O log(1/ǫ ′ ) b-a other primitive quantum gates.

Figure 2 :

 2 Figure 2: Example of density of mixture of Gaussians.In this case, the density quickly goes to 0 around our estimate Ê0 . We then continue to check the variance conditioned around Ê0 to see if we accept the estimate.

Figure 3 :Figure 4 :

 34 Figure 3: Example of density of mixture of Gaussians. In this case, the density quickly does not go to 0 quickly around Ê0 . We reject the estimate.

Figure 5 :

 5 Figure 5: We estimate the variance of the random variable conditioned around Ê0 . If it is close to σ 2 , we accept it. If it differs from σ 2 by ∼ ηǫ 2 , we reject it.

πσ 1 -

 1 e -(w-s) 2 /σ 2 (73)

  xe -x 2 /σ 2 dx + s+w s-w

  Then we have ln (eσ/ǫ 1 ) = ln e0.22∆ǫ -1 /( ln (ec 2 ) ln (eA)) ≤ ln e0.22A/ ln (eA) ,

F

  = {i : 2L ≤ |E i -Ê0 |}. (150) If p * n σ is peaked around Ê0 with c = 4 and τ = c ′ ǫ 2 η 2 n σ )(x)dx -i∈G p i ≤ 4c ′ ǫ 2 η 2 . (p * n σ )(x)dx = i∈G p i Ê0 +L Ê0 -L n σ,i (x)dx + i∈B p i Ê0 +L Ê0 -L n σ,i (x)dx + i∈F p i Ê0 +L Ê0 -L n σ,i (x)dx.

2 )By our choice of L we have τ c 2 = c ′ ǫ 2 η 2

 222 i (x)dx ≤ 2L e c log(τ σ) σ √ 2π = O(τ c log τ -1 ) i (x)dx = O(τ c log τ -1 ). (154)Let us now estimate the contribution of terms in G. We have for i ∈ G: the property that all energies in G are L 2 away from Ê0 , we see that both integrals above are at most O(τ c

4 2. 2 n

 42 Thus, both Eq. (156) and Eq. (154) are bounded byc ′ ǫ 2 η 2 .To finish the proof, we need to consider the contribution of the energies in B. Note that for i ∈ B: σ,i (x)dx ≤ 2c ′ ǫ 2 η 2 . (157)

Corollary D. 4 .x 2

 42 Under the same conditions as Lemma D.3, assume further that σ is small enough to ensure that | Ê0 ± L| ≤ 1. If p * n σ is peaked around Ê0 with c = 4 and τ = c ′ ǫ 2 η 2 4 then:Ê0 +L Ê0 -L x(p * n σ )(x)dx -i∈G p i E i ≤ 4c ′ ǫ 2 η 2 . (p * n σ )(x)dx -i∈G p i (σ 2 + E 2 i ) ≤ 4c ′ ǫ 2 η 2 .(160)

1 : 2 :τ ← ǫ 2 η 160 3 : 2 p 6 :τ < ǫ 2 η 80 then 7 : 8 :S 9 :

 12326789 procedure GSEE CERT(ǫ, Ê0 , σ, η,p)Let X ∼ i p i N (E i , σ 2 ) with pdf p * n σ * n σ (x)dx up to ǫ 2 η 160 . if M ← estimate E(X|X ∈ [ Ê0 -L, Ê0 -L]) up to ǫ 2 η 160 ← estimate E((X -M ) 2 |X ∈ [ Ê0 -L, Ê0 -L]) up to ǫ 2 η 160 . if |Sσ 2 | ≤ 2ǫ 2 η then 10:

  ) ≤ σ/ǫ 1 , the second step follows from η -1 ≥ 1, the fourth step follows from the definition of σ, ǫ 1 and A, and the last step follows from c 2 ≥ 1.

	2 2 σ ln (eσ/ǫ 1 ) /(3ǫ 1 η) ≤ ln 200c 2 2 σ 2 /(3ǫ 2 1 η) 2 σ 2 /(3ǫ 2 1 η 2 ) ≤ ln 200c 2 2 + 2 ln η -1 σǫ -1 1 ) ≤ ln 67c 2 = ln 67c 2 2 + 2 ln 0.22A/( ln (ec 2 ) ln (eA))	(102) (103) (104) (105)
	≤ ln 67c 2 2 + 2 ln 0.22A/ ln (eA)	(106)
	where the first step follows from ln (eσ/ǫ 1 Now we combine Eqs. (101), (106) and the fact that
	5 ln (ce) ln (et) ≥ 3 ln e0.22t/ ln (et) + ln (67c 2 ) + 2 ln 0.22t/ ln (et)	(107)
	for all c ≥ 1 and t ≥ 8, and obtain Eq. (96).

Note that our test certifies that a ground state energy estimate is ǫ-close to the correct one, not that the spectral gap of the Hamiltonian was at least some ∆0.

One can prove w ≤ 0.8 by the problem assumptions and our choice of parameters. Then since Ẽ0 ∈ [-1, 1], we have -1.8 ≤ Ẽ0 -w < Ẽ0 + w ≤ 1.8. This means that we might choose ξi > 1 or < -1 occasionally. This decision is made to ease the proof.

D.2 Sample complexity of GSEE CERT with rejection sampling

Let us analyze the sample complexity of the test above with access to samples and with a success probability of correct verification of 1δ.

Corollary D.7. Under the same conditions as Theorem 3.6, assume further that |E 0 -Ê0 | ≤ σ. Then we can perform the run Algorithm 3 with probability of correct rejection or acceptance 1δ with access to

uses of the block-encoding of g σ,ξ;ǫ ′′ (H).

Proof. There are two steps that need to be performed for running Algorithm 3. First, we need to determine if the distribution is peaked around Ê0 . If it is peaked, we need to determine the variance up to an error ηǫ 2 on the interval [ Ê0 -L, Ê0 + L]. Let us start by showing how to determine if the distribution is peaked. For that we will perform rejection sampling with the Gaussian condtioned on [ Ê0 + L 2 , Ê0 + 2L]. As explained in Section 3 in Eq. ( 15), the probability of accepting is given by

for an explicitly known constant c. Thus, by running O(ǫ -2 η -1 ) rounds of the rejection sampling, we can decide if the distribution is peaked around Ê0 or not with probability of failure at most δ 2 by checking how many times we accepted the sample. If it is peaked, then we can perform rejection sampling again, but now on the interval [ Ê0 -L, Ê0 + L]. With access to O(η -2 ǫ -4 log δ -1 ) samples from the conditioned distribution, we can estimate the variance of the conditioned random variable up to a precision of O(ǫ 2 η) up to a failure probability of at most δ 2 . The probability of accepting a sample from the interval is given by

As before, note that for our choice of L we have that this is lower-bounded by

for some constant C > 0. Now, as E 0 ∈ [ Ê0 -L, Ê0 + L] by our choice of L and the assumption that |E 0 -Ê0 | ≤ σ, we have that Ê0 +L Ê0 -L (p * n σ )(x)dx = Ω(η). Thus, we need to run the rejection sampling an expected Õ(η -1 ) number of times to accept a sample. This gives a total complexity for this step of Õ(η -3 ǫ -4 log δ -1 ), which dominates the overall sample complexity of the certification algorithm. Furthermore, we see that the probability of obtaining wrong estimates is at most δ. As the algorithm always outputs the correct answer if given correct estimates, this gives the claim.

Note that the assumption that we are given an estimate that satisfies |E 0 -Ê0 | ≤ σ does not increase the overall depth of circuits required for the certification method. This is because we can always obtain such an estimate with circuits of depth O(σ -1 ) using e.g. Algorithm 1.