This Perspective focuses on the several overlaps between quantum algorithms and Monte Carlo methods in the domains of physics and chemistry. We will analyze the challenges and possibilities of integrating established quantum Monte Carlo solutions into quantum algorithms. These include refined energy estimators, parameter optimization, real and imaginary-time dynamics, and variational circuits. Conversely, we will review new ideas for utilizing quantum hardware to accelerate the sampling in statistical classical models, with applications in physics, chemistry, optimization, and machine learning. This review aims to be accessible to both communities and intends to foster further algorithmic developments at the intersection of quantum computing and Monte Carlo methods. Most of the works discussed in this Perspective have emerged within the last two years, indicating a rapidly growing interest in this promising area of research.

## I. INTRODUCTION

The solution to quantum many-body problems in chemistry and physics is one of the most anticipated applications of a quantum computer, as first proposed by Feynman.^{1} Over time, it has been proposed that many other classes of problems can benefit from quantum speedup, including cryptography, data science, machine learning, finance, linear algebra, and optimization.^{2} However, physics and chemistry remain among the main candidates for demonstrating practical quantum advantage over conventional methods because they contain classes of problems with the following characteristics: (i) they are very challenging for classical computation, and exponential quantum speed-ups are possible; and (ii) they are defined by a small number of variables, the Hamiltonian parameters, thus featuring a limited cost of data loading and reading.^{3}

Among all possible problems in physics, here we will focus on electronic structure and spin models (including classical spin models), as their implementation requires a relatively lower cost compared to models like high-energy physics. Excellent review articles on quantum algorithms for quantum chemistry^{4,5} and materials science^{6} have already been published a few years ago, the latest one in 2020. The purpose of this manuscript is not to duplicate such presentations but rather to concentrate on a more frontier topic that is becoming relevant due to several works appearing in recent months. However, notice that about 40% of references cited in this Perspective are pretty recent, i.e., from 2021 onwards. This indicates how fast the whole field of quantum computing is growing.

We will analyze points of contact between quantum computing and Monte Carlo (MC), quantum Monte Carlo (QMC) methods.^{7} There are many common themes between the two worlds. Shot noise arising from the measurements of the quantum register finds a parallel in the statistical root of Monte Carlo. Both methods require extracting and utilizing expectation values computed in the presence of statistical noise.

The existence of shot noise is one of the major issues for near-term simulations; in a variational setting, this implies a problematically large number of circuit repetitions.^{8} On the other hand, such uncorrelated wave function collapses can have computational value if used as important sampling in Monte Carlo. We will also report attempts at cross-fertilization between the two fields in designing variational *Ansätze* and optimization methods for ground-state and dynamical problems. Additionally, we discuss the requirements that a classical-quantum hybrid QMC algorithm, relying on a quantum computing subroutine, must meet. Regarding classical applications, we review several proposals for accelerating the Metropolis algorithm using quantum hardware and examine their practicality under realistic hardware constraints.

Therefore, the purpose of this manuscript is to review various Monte Carlo techniques that can be useful for creating new quantum algorithms or designing new applications of already-known quantum primitives. Conversely, this Perspective also aims to be an accessible presentation of the potential and limitations of quantum computing for Monte Carlo experts and, more broadly, computational physicists.

This Perspective is timely as (1) on the experimental side, the first steps toward fault-tolerant hardware have been made.^{9–11} Moreover, experiments at the threshold of quantum advantage for quantum dynamics are now possible given the existence of $\u223c100$ qubits of hardware, albeit noisy.^{12} (2) On the quantum variational algorithm side, we have increasing evidence that the quantum measurement noise—the focus of this Perspective—is the major, unavoidable, bottleneck of near-term quantum algorithms.^{8,13} (3) In the last year, thorough resource assessment papers for quantum chemistry have appeared,^{15–16} which clearly reaffirm the threshold for quantum advantage for ground state problems, at least with today’s algorithms, to be deep into the future fault-tolerant quantum computing regime and question the previously claimed exponential advantage for ground-state electronic structure application.^{17} (4) Finally, we observe the emergence of a new class of hybrid quantum algorithms revisiting classical and quantum Monte Carlo, opening completely new possibilities for quantum advantage in these areas. In this Perspective, we report about 20 recent (i.e., those that appeared in the last two years) works that aim to complement quantum computing and Monte Carlo in several sub-fields: hybrid quantum computing (QMC), variational circuit development, parameter optimization, time-dependent simulations, and classical sampling.

The manuscript is organized as follows. In Sec. II, we briefly mention the types of quantum hardware and their fundamental limitations, namely the hardware errors in the Noisy Intermediate-Scale Quantum (NISQ) regime and the fairly long [compared to conventional central processing unit (CPUs)] gate time of future, error-corrected machines. In Sec. III, we introduce general quantum algorithms for physics and the requirements for quantum advantage. In Sec. IV, we review the basics of the encoding of a fermionic Hamiltonian into a quantum computer. After these introductions, in Sec. V, we discuss the variational method, the simplest kind of algorithm for physics and chemistry, and its limitations due to the shot noise. In Sec. VI, we instead describe how the same calculation could be performed in the fault-tolerant era. The more technical Sec. VII introduces the local energy estimator that is central in QMC and explains why these algorithms, while still being stochastic, do not suffer from the severe variance problem of variational quantum algorithms. In Sec. VIII, we review attempts to create hybrid quantum–classical QMC algorithms, as well as other points of contact between QMC and quantum algorithms. Finally, in Sec. IX, we reverse our perspective and discuss quantum computing methods to speed-up classical sampling using digital machines, quantum simulators, and annealers.

The concept map of the Perspective is shown in Fig. 1.

## II. QUANTUM HARDWARES

Unlike conventional reviews on algorithms for quantum chemistry, it is necessary to briefly introduce the hardware on which these must be executed. Understanding the possibilities and limitations of the hardware is crucial to get an idea of the feasibility of current and future algorithms.

There are many types of quantum computers and quantum simulators. The difference between the two classes is that a quantum computer is built with the idea of being universal and, therefore, able to support any type of program. A quantum simulator is designed to perform a narrower range of tasks, such as optimizing classical cost-functions^{18,19} or simulating particular Hamiltonians.^{20}

To tend toward universality, a quantum computer must support the execution of basic local operations, called *quantum gates*, just like a regular computer. For example, an architecture capable of synthesizing the following set of gates, ${CNOT,H,S,T}$ or ${Rx(\theta ),Ry(\theta ),Rz(\theta ),S,CNOT}$, that act on a maximum of two qubits (see textbooks on quantum computing for the definition of these gates^{21}) is capable of approximating any possible unitary operation on *n* qubits.^{21} On the other hand, a special-purpose quantum simulator can implement the global operation of interest directly, such as *e*^{iHt}, where *H* is the quantum many-body Hamiltonian operator, without the problem of having to compile it using a gate set.^{22} Currently, the greatest engineering effort is focusing on gate-based “digital quantum computers,” although it is not excluded that algorithms of interest for chemistry and materials science can be executed on quantum simulators.

There is then a second distinction that is important to keep in mind. At present, the size of digital quantum computers is on the order of *n* ≈ 100 qubits.^{12,23} In principle, these computers have access to a $2O(100)$ dimensional Hilbert space. However, practical quantum advantage has not been achieved yet. The reason is that these machines are not properly digital but are subject to hardware noise. For example, each gate has a finite failure probability. As we will see, the circuits necessary to write a quantum algorithm for chemistry require a considerable number of gates and, therefore, even a small infidelity propagates devastatingly, and the total error accumulates until it completely compromises the success of the algorithm. Most current hardware is built with the idea of executing a universal gate set but is still affected by hardware noise. They generally go by the name of NISQ (Noisy Intermediate-Scale Quantum) machines.^{24}

The final step to achieving a true digital quantum computer is to realize hardware capable of executing gates without errors, just like their classical counterpart. Many detractors of quantum computers base their skepticism on the impossibility of maintaining a macroscopic coherent wavefunction for an arbitrary number of operations.^{25} Fortunately, there is a theorem that does not exclude the possibility of a digital universal computer: below a certain noise threshold, it is possible to correct this hardware noise faster than it can accumulate during runtime (the same would not hold for analog computers).^{21} Practical proposals to realize this idea include using multiple physical qubits to realize a logical qubit and more operations to realize a single “error-corrected” logical gate.^{26,27}

It is important to understand which types of algorithms have the hope of being executed on a NISQ machine and which will require a fault-tolerant machine to properly contextualize the ever-growing literature on quantum computing for quantum chemistry. At present, sub-communities have formed that are dedicated to developing NISQ algorithms, and others that are increasingly growing in number are developing algorithms for the fault-tolerant era. This situation is unprecedented. In classical computing, to attempt a comparison, it would be as if in the 1960s there was a community developing algorithms for chemistry on punched cards and another preparing for exascale computing without a clear idea of whether and how a high-performance computing (HPC) facility would be built.

However, the technological progress toward a fault-tolerant machine is steady. Several experiments published in 2022–2023 have already demonstrated some building blocks necessary for quantum error correction.^{9–11,28–30} Clearly, we are still in the infancy of fault-tolerant hardware, and it is not yet clear when a large-scale error-corrected machine, able to accommodate electronic structure calculations, will appear. This also depends on progress on the algorithmic side of compressing memory and runtime resources.

The last concept that needs to be presented in this brief account of hardware is related to the clock frequency of a quantum computer, which will always be slower than classical gates. This is because the execution of a quantum gate requires the manipulation of a wave function by an external control; the quantum gate can never be faster than the classical electronic apparatus that controls it. At present, the execution time of a noisy CNOT i.e., controlled not, or noisy single qubit rotation is of the order of 100 ns in the NISQ era for superconducting qubits, corresponding to 100 MHz.^{23} The expected logical clock rate in the fault-tolerant regime is much slower, of the order of 10–100 kHz, because every error-corrected gate requires a large number of elementary operations that involve a large number of physical qubits.^{31} Often, in literature, one hears about T-depth as a proxy for the complexity of an algorithm.^{16,32,33} The reason is that truly digital hardware can only operate using a set of discrete gates that can be error-corrected. A rotation of an arbitrary angle does not belong to this category, and therefore every continuous rotation must be compiled into a series of discrete operations, such as S, H, and T gates. The T gates are the most expensive to synthesize (i.e., 100–10000 more costly than a CNOT)^{31} and, therefore, their number and how many can be executed in parallel determine the runtime. Since in chemistry, rotations of an arbitrary angle may be necessary for orbital basis rotations or to realize infinitesimal Trotter step operations, which are ubiquitous, fault-tolerant quantum algorithms generally include a large number of T gates. It is interesting to note that in the NISQ era, the opposite is true: rotations, being one-qubit operations, are comparatively simple gates to perform, while efforts are made to reduce the number of CNOT gates, which are currently the most noisy, as they require the control of two physical qubits. Recent proposals include the possibility of retaining analog rotation gates and error-corrected Clifford gates, which are easier to synthesize.^{34,35} This hybrid approach is interesting but has yet to be demonstrated in practical algorithms. Moreover, the gate times depend on the specific hardware architecture. Other platforms such as spin qubits, trapped ions, or photonic hardware will imply different hardware constraints.

## III. ALGORITHMS AND QUANTUM ADVANTAGE

After this essential overview of hardware, we are now in a position to introduce more concretely the most popular algorithms for the quantum many-body problem. The application where quantum advantage appears most clear and easy to justify is quantum dynamics. In this case, an exponential advantage can be obtained by virtue of the exponential compression of the Hilbert space into a linear memory with the number of particles.^{36} If we consider, for simplicity, spin-1/2 lattice models composed of *n* spins, it is easy to see that an *exact* simulation becomes unfeasible as soon as *n* ∼ 50. Storing a quantum state of 50 qubits with double-precision coefficients for each of the 2^{50} possible components requires 16 PB of memory. To perform arbitrary discrete-time evolution, we would need to manipulate such an array thousands of times. For instance, direct matrix exponentiation *e*^{iHt} for a typical many-body quantum Hamiltonian, *H*, and 50 spin-1/2 particles would require an array made of 10^{15} entries undergoing a matrix-vector multiplication of size 2^{50} × 2^{50} = 10^{30}.

The memory requirement for evolving a system of *N* electrons in second quantization using *n* orbitals is the same. Second quantization offers optimal memory usage if a Fock encoding is used, where each binary string represents a possible configuration of orbital occupation (see Sec. IV).^{37} As we will see, the choice of second quantization introduces non-locality in the Hamiltonian, which is translated as a sum of tensor products of Pauli matrices (qubits are not fermionic particles) and, therefore, requires very long circuits. From an asymptotic point of view, first quantization can be advantageous for Hamiltonian evolution tasks, as it preserves the locality of interactions at the cost of introducing expensive arithmetic operations to calculate the Coulomb term. In addition, an antisymmetric function in real space must be provided.^{38}

Simulating lattice spin models, therefore, appears to be the most obvious choice in the search for a quantum advantage in Hamiltonian simulations. Beverland *et al.*^{16} show that a fault-tolerant simulation of a 10 × 10 plaquette of the two-dimensional quantum Ising model requires on the order of 10^{5}–10^{6} physical superconducting qubits using the digital Trotter algorithm. This example is instructive, as it clearly shows that, although the system requires *a priori* a memory of *n* = 100 logical qubits (1 logical qubit = 1 spin), error correction and extra resources for distilling T gates (often called T-factories) push the total computation of qubits up to a million. This system sets a lower bound on the resources needed to simulate fermionic systems that are more complex than spin-1/2 models for the same number of particles/spin-orbitals. Simulating lattice models is also possible with quantum simulators, although with calibration errors, and it is, therefore, likely that there will be competition between the two quantum computing paradigms toward the first simulation beyond what is possible classically.^{22,39}

Indeed, it is also possible that the first demonstration of practical quantum advantage in real-time simulations will be achieved before the fault-tolerant regime, for example, in simulations of dynamical phase transition.^{40} At the time of writing this review, researchers at IBM showcased the record-sized real-time dynamics of a 2D heavy-hexagon lattice Ising model using trotterization, 127 qubits, and error mitigation techniques.^{12} Just ten days later, approximate tensor network simulations achieved the same result,^{41,42} raising the bar once again to declare a quantum advantage, similar to what happened with the first claim of quantum supremacy on random circuit sampling.^{23,43,44} Demonstrating a definitive quantum advantage in quantum dynamics tasks is a less well-defined goal, as classical limits are still not thoroughly explored. One can expect that increasingly sophisticated classical methods will be adopted to counter the new claims of quantum advantage. The case of chemistry is different since classical methods are fairly established, and it is much clearer which electronic structure Hamiltonians are beyond the reach of classical computing.

Coming back to our chemistry problems, if the required answer requires precision that can only be achieved with a long circuit, then we must prepare ourselves for the fact that our algorithm will only be available in the digital era, which could take a decade or more when a technology capable of controlling millions of qubits will be available. As mentioned, because fault-tolerant hardware has a lower clock frequency than a conventional computer, a non-exponential asymptotic speed-up may not be sufficient to guarantee actually shorter runtimes for reasonable system sizes. Babbush *et al.*^{45} recently discussed how a quadratic speed-up is insufficient for practical advantage in many applications.

To conclude, quantum advantage in chemistry problems can be obtained either in several years using a fault-tolerant algorithm with a superquadratic speed-up or in a heuristic way using NISQ or analog hardware. In this case, only short-depth algorithms can be used, which often require a classical feedback loop, such as optimization of variational parameters that define the circuit and repeated executions of the circuit. Variational methods fall exactly into this class of solutions.^{46,47} We will see in the next chapter how these can also be used and their limitations, especially their relationship with a type of noise that cannot be eliminated even in a fault-tolerant regime.

## IV. FERMIONIC HAMILTONIANS

^{4–6}for a general introduction and more details. The most general fermionic Hamiltonian reads

*a*

_{p}are the fermionic creation and destruction operators, which create (annihilate) a particle in the spin-orbital

*p*.

To proceed, the next step is to define an encoding from fermionic Hilbert space to qubit space so that a base vector of the latter, for example, 1100, has a unique correspondence in the fermionic space.

In second quantization, this is trivial using Fock states. A bit-string denotes the occupation numbers of spin-orbitals in a chosen ordering. For instance, the string 1100 can represent the Hartree–Fock state of an H_{2} molecule described with two spatial molecular orbitals. Here, there are two electrons of opposite spin, occupying the first two spin-orbitals, {|*ψ*_{g,↑}⟩, |*ψ*_{g,↓}⟩}, while leaving empty the two higher energy ones, {|*ψ*_{u,↑}⟩, |*ψ*_{u,↓}⟩}.

*H*

_{2}case, the exact ground state (in this small atomic basis) can be expressed as a linear combination of these two strings alone, but in the general case, the ground state of a molecule with

*N*electrons and

*n*spin-orbitals is written as a linear combination,

*c*

_{i}are complex coefficients and |

*i*⟩ is a basis vector whose binary format denotes the occupation number of the spin-orbitals in a chosen ordering. Although many

*c*

_{i}’s are zero due to particle and spin conservation, as is known, an exponentially large number of them remains finite when

*N*and

*n*increase.

The main advantage of a quantum computer is, therefore, clear. With only *n* qubits, it is possible to store in memory a wavefunction with 2^{n} complex coefficients. The fundamental question now is whether it is possible to devise a quantum algorithm with polynomial complexity capable of manipulating these 2^{n} coefficients to find the ground state of a fermionic problem, or at least with a better approximation than the best classical method.

*minus*and

*plus*operators change the occupation number of the target mode

*p*, while the string of

*Z*operators is needed to enforce antisymmetrization. Therefore, each fermionic operator in Eq. (1) translates into a combination of the tensor product of Pauli operators.

*h*

_{j}is a real scalar coefficient, and

*N*

_{P}is of order $O(n4)$ since there are

*n*

^{4}terms to transform in Eq. (1). Each term

*P*

_{j}in the Hamiltonian is typically referred to as a

*Pauli string*and is a tensor product $Pj=\u2297i=1n\sigma i\alpha $ of

*n*Pauli matrices

*σ*

^{α}∈ {

*I*,

*Z*,

*X*,

*Y*}. What is important for our discussion is that (1) the coefficient

*h*

_{p}can take very large values (in modulus). Just to give an idea,

*∑*

_{j}|

*h*

_{j}| is of the order of 40 Ha, already a moderate example of an H

_{2}O molecule described in the STO-6G basis.

^{8}Then (2) the operators

*P*

_{j}can have a number of non-identity gates, which is $O(n)$ due to the non-locality introduced by the Jordan–Wigner transformation. This implies that the circuit for real-time dynamics is longer compared to local quantum spin Hamiltonians,

^{37}and when measured, the expectation value ⟨

*P*

_{j}⟩ is exponentially susceptible to bit-flip measurement errors.

^{48}Other fermion-to-qubit mappings exist, such as the parity mapping and the Bravyi–Kitaev mapping, which are the most popular ones after the Jordan–Wigner transformation. In particular, the Bravyi–Kitaev mapping improves asymptotically the locality of the resulting qubit Hamiltonian. However, its effectiveness on finite-size problems should be evaluated on a case-by-case basis.

^{4}The form of Eq. (4) is, however, more general than the Jordan–Wigner workflow, so we assume it as the starting point of our discussion.

## V. VARIATIONAL QUANTUM ALGORITHMS

In short, variational methods aim to use shallow parametrized circuits that can be optimized to minimize the calculated energy.^{47,49,50} The *N*_{par} variational parameters can be the angles *θ*’s of the rotation gates defined earlier. This strategy takes the name of Variational Quantum Eigensolver (VQE), but it is nothing more than a common variational calculation on a quantum computing platform. First of all, it is important to notice that even shallow circuits, i.e., those featuring constant depth vs *n*, can display quantum advantage, although on quite artificial tasks,^{51} or cannot be efficiently simulated classically.^{52} Therefore, variational algorithms are reasonable candidates for quantum advantage in the near term. As of the present, the literature features many small-scale hardware demonstrations, still away from the quantum advantage threshold. The most notable either use heuristic circuits^{53} or more structured physically inspired circuit *Ansätze*.^{54,55} The current largest variational simulation of a chemical system reaches a system size of about 20 qubits.^{56} Performing variational calculations of many-body quantum systems has advantages in principle but also many limitations in practice.

Current technology allows the execution of circuits with more than 100 qubits and a depth of about 60 two-qubit gates.^{12} While error correction is not yet available, there are error mitigation methods that enable an unbiased estimation of the expectation value of operators.^{57,58} It seems, therefore, that all the ingredients for enabling NISQ variational methods are present, as a circuit can define a variational *Ansatz* that could outperform the best classical *Ansatz* for a given problem. A central point of this Perspective is thus what is missing to translate this potential into practical variational computation and, hopefully, achieve quantum advantage.

The first point to establish is what kind of circuit can be used to create the ground state of our target many-body quantum Hamiltonian. As we have seen, the advantage of a quantum computer is the mere possibility of storing an exponentially large wavefunction with a linear number of qubits in memory [cf. Eq. (2)]. However, this gives us no guarantee that (1) a quantum circuit with a finite, and possibly small, depth can give us a better approximation than the best classical method, and (2), even more importantly from a conceptual point of view, that it is possible to optimize the parameters even assuming that the ground state is contained in the variational *Ansatz*.

An obvious drawback of variational calculations in the NISQ regime is the presence of noise. A gate error of 0.1% propagating through a circuit with a depth of 100, composed of 100 qubits, produces a state with fidelity on the order of 10^{−3}–10^{−4} compared to the noiseless case.^{12} However, there are error mitigation methods that can be applied to obtain unbiased estimates of expectation values. Assuming, as a working hypothesis, that the noiseless version of a quantum circuit can generate a variational *Ansatz* that is better than the best classically available, or even the exact one, everything becomes a game of exponents, as error mitigation incurs an exponential cost in circuit repetitions.^{57–58} It remains to be determined whether the exponent of this post-processing step is mild enough to guarantee reasonable runtimes for classically intractable molecules.

In this Perspective, we will not focus on hardware error mitigation, which is introduced in the recent Ref. 56, but on another complementary issue. As we will see, a major problem with quantum variational methods is even more surprising: even assuming that we have prepared exactly the target state, computing its energy is, so far, an inefficient procedure. Although this inefficiency is not the same as in the definition of complexity theory, where an algorithm is said to be inefficient if it is exponentially scaling, it is from a practical point of view. This is a completely new condition that does not happen, for example, in variational calculations using QMC.

### A. The variance problem

The fundamental (and obvious) concept at the root of this is that a quantum computer obeys the postulates of quantum mechanics: we cannot access the state that is created by the circuit exactly, but only through measurements. We can measure the expectation value of the Hamiltonian by post-processing the read-out of each measurement, and then prepare exactly the same state and repeat the measurement, and so on.

*ψ*⟩. Following Eq. (4), we see that the mean value of the Hamiltonian is the linear combination of the expectation values of the

*N*

_{p}Pauli strings,

*P*

_{j}’s are obtained from

*N*

_{P}independent sets of measurements; the error on the estimate is then given by

*M*

_{j}repeated measurements or

*shots*. Depending on the magnitude of Var[

*P*

_{j}], several shots, thus the circuit’s repetitions, need to be allocated to measure the expectation value of the

*j*th Pauli string to reach the target additive precision.

Since we need to evaluate each ⟨*P*_{j}⟩ independently, their statistical fluctuations are not correlated, so one needs to reach chemical accuracy by resolving each expectation value with very high accuracy as it gets multiplied by a possibly very large prefactor *h*_{j}. In Wecker *et al.*,^{8} it is estimated that *M* = *∑*_{j}*M*_{j} = 10^{9} measurements would be needed to reach *ϵ* ∼ 10^{−3} Ha for H_{2}O and up to 10^{13} for Fe_{2}S_{2} described with 112 spin-orbitals and STO-3G basis.

This happens even in the case where we prepare the exact ground state, thus violating the *zero-variance* property of the ground-state if the energy is calculated this way. This issue is often called *the variance problem* and is one of the most overlooked issues in the VQE community, which seems to be more active in new circuit *Ansätze* development. However, several works aiming to mitigate this problem have been put forward. The simplest one consists of grouping all *P*_{j}’s that commute qubit-wise and that can be measured simultaneously.^{53} Other methods aim to find better grouping schemes, introducing general commutativity rules at the expense of longer measurement circuits.^{48,60–62}

For the evaluation of two-body reduced density matrices for fermionic systems, it is possible to devise an asymptotically optimal measurement scheme^{62,63} with a number of unique measurement circuits scaling as $Ngroups\u223cO(n2)$. However, there is still the problem that all these correlators, estimated stochastically, need to be summed up to ⟨*H*⟩. Therefore, *N*_{groups} is not a faithful representation of the number of total measurements to be performed to achieve chemical accuracy since each partition needs to be measured *M*_{per groups to ϵ} times, as for the sum to achieve a total error of *ϵ*.

There are also methods based on a different philosophy, namely to systematically *approximate* the electronic Hamiltonian [Eq. (1)] to reduce the number of terms in the sum using a low-rank representation of the Hamiltonian.^{64,65} This technique finds application in real-time dynamics simulations and may considerably reduce the runtime of error-corrected algorithms.^{15} While it can certainly mitigate the variance problem in the NISQ era, it does not qualitatively solve the problem for the very same reason outlined earlier. We observe that several works aim to reduce the number of basis states required to represent the Hamiltonian. However, this may not necessarily improve the number of measurements, as these fewer terms may have a larger variance. An interesting example of this is seen in variational quantum algorithms applied to classical cost functions. This scenario is important for optimization, and in this case, the most popular variational method is called the Quantum Approximate Optimization Algorithm (QAOA). Despite the fact that the cost function, by definition, only needs to be measured on one basis—the computational basis—the impact of shot noise is still quite detrimental to the overall performance, even in this case.^{66}

Finally, we also report the celebrated *shadow tomography*^{67} method, which may be useful for estimating local qubit operators but has an exponential scaling for non-local ones, such as our *P*_{j}’s.

*N*

_{iter}optimization steps, is

Let us consider a concrete example of an H_{2}O molecule on a very minimal basis consisting of 12 spin-orbitals. Following a state-of-the-art variance reduction method,^{48} we quote a number of 10^{8} circuit repetitions to compute a single point energy within 10^{−3} Ha accuracy (notice that this is the error from the exact value obtained with this minimal basis set, not the exact value using a converged basis set). Assuming a circuit execution time with measurements of 1 *µ*s, we are limited to about hundreds of optimization steps per day (totally neglecting classical communications and reset times).

Notice that NISQ hardware means that hardware noise will always be present, and this noise usually varies with time. That is why an optimization run that lasts for more than one day will likely never converge: the optimal parameters ** θ**, which yield the lowest energy with the hardware configuration at present, may not be optimal the day after. The fact that hardware errors are device and time-dependent realistically excludes the possibility to parallelize the shots using several machines, as is customary in conventional QMC. In this case, one would claim the possibility of always preparing the same trial state on different noisy machines.

### B. The noisy optimization problem

The second step of any numerical variational calculation is optimization. In this Perspective, we aim to maintain a high-level tone and will not be concerned with details such as which optimization method is better or worse depending on the type of circuit. In general, optimizing the parameters is a complicated task that is delegated to the classical part of the algorithm. Heuristic circuits are very short but have a much larger number of variational parameters than those inspired by chemistry, such as the unitary coupled cluster, or physics, such as the Hamiltonian variational *Ansatz*.^{4} The latter have longer circuits but allow for a more stable optimization.

A wealth of literature focuses on theoretical roadblocks, such as the existence of barren plateaus^{68} and the fact that optimization itself is an NP-hard (non-deterministic polynomial-time) problem.^{69} However, we observe that even in the conventional case, the optimization of parameters occurs in a corrugated landscape. Nevertheless, it is almost routine to optimize thousands of variational parameters in variational Monte Carlo (VMC).^{70} Moreover, barren plateaus are a concept borrowed from the quantum machine learning community and are likely not relevant, or at least not the real bottleneck, in the case of using structured *Ansätze*^{13,71,72} (i.e., non-random or heuristic), which should be the norm for studying physical systems. The reason why this concept has never arisen in conventional VMC is that no one has ever tried to optimize molecular systems from quasi-random trial states.

^{7}which is not possible in this case. The other possibility is to calculate the expectation value of the generalized forces

*f*

_{i}, defined as the derivative of the energy with respect to the variational parameters. For simple circuits, calculating the force is possible thanks to a technique called the parameter shift rule,

^{5}and extensions are possible for more structured circuits.

^{13}Due to the no-free-lunch theorem, the statistical error that we had using energy translates into a statistical error in the forces. The optimization effectively becomes a stochastic gradient descent, which resembles a discrete Langevin equation at a finite effective temperature,

*f*

_{i}= −

*∂*

_{i}⟨

*H*⟩, $\eta ishot$ is a Gaussian distributed random variable, and

*δ*is a finite integration step. Astrakhantsev

*et al.*

^{13}have shown that the statistical error defines an effective temperature,

*T*

^{shot}, proportional to the variance of the random variable $\eta ishot$. Below a certain number of samples,

*M**, therefore, above a certain effective noise temperature, the search is unsuccessful. Above, optimization becomes possible. Moreover, in the

*M*≫

*M** regime, the infidelity of the state preparation seems to scale as 1/Δ

^{2}, where Δ is the energy gap between the ground and the first excited state. These numerical results have been obtained on challenging

*j*

_{1}−

*j*

_{2}Heisenberg models (cf. Ref. 73), and it would be interesting to check how they generalize to chemistry problems. On the optimistic side, the critical number of samples

*M** seems to not scale exponentially with the system’s size, though allowing in principle efficient VQE optimization in the presence of quantum measurement noise. Moreover, in this state-of-the-art VQE study, barren plateaus are not observed.

*∂*

_{j}

*ψ*⟩ is the derivative of |

*ψ*⟩ as a function of the

*j*th variational parameter. Although most of the community believes that Eq. (9) comes from machine learning, where it is used in the

*natural gradient*,

^{74}it has actually been used for more than twenty years in VMC to optimize trial wave functions for chemistry and condensed matter.

^{7}It was introduced by Sorella in the stochastic reconfiguration method

^{75,76}and later given the geometric meaning of the metric of the space of variational parameters.

^{77}A weak regularization of the diagonal is sufficient to obtain stable and effective optimizations, as has also been shown in the quantum case.

^{78}

However, the measurement problem also heavily affects the calculation of algorithm efficiency in this case. In VMC, the matrix *S* can be evaluated with negligible overhead using the same samples, *x*’s, distributed as |*ψ*(*x*)|^{2} (which can be evaluated classically there), that were already generated for the energy calculation. In the quantum case, each matrix element must be statistically evaluated using (uncorrelated) repetitions of a specific circuit for the pair *i*, *j*. Moreover, it is not trivial to obtain a circuit for each element *S*_{ij}, and only a block-diagonal approximation of *S*, where *i* and *j* belong to the same block, is the most feasible solution.^{78}

At the moment, an interesting development to overcome this problem is a heuristic combination of stochastic reconfiguration with the Simultaneous Perturbation Stochastic Approximation (SPSA) optimizer. In this case, the matrix S, which would require $Npar2$ circuits, is approximated by a Hessian calculated using only two random directions in the parameter space.^{79} However, the numerical benchmarks proposed to validate the method are too small to fully understand the real possibilities of this simplified optimizer.

Finally, another connection between VQE and QMC arises in the context of using noisy ionic forces in molecular dynamics (MD) simulations. In such cases, it becomes impractical to follow an energy-conserving trajectory. Sokolov *et al.*^{80} utilize a similar technique proposed in QMC-powered MD simulations from many years ago.^{81} They use shot noise to define an effective Langevin MD, enabling unbiased simulations at constant temperature.

## VI. FAULT-TOLERANT QUANTUM CHEMISTRY

At this point, we would face the conundrum that a quantum computer can in principle store an exact wavefunction, but we cannot practically evaluate its energy or other expectation values using the variational methods introduced so far. However, quantum computing admits an efficient method for calculating energy, even more efficient than Monte Carlo methods. These methods are based on the quantum phase estimation (QPE) algorithm or its successive variants. The significant property of this approach is the ability to project an initial state onto the ground state of the Hamiltonian, thus accessing the ground state energy. This method is, therefore, more precise and accurate than a variational approach and is often considered the standard quantum computing approach to quantum chemistry.^{32} However, as we will see, the probability of success depends on many factors, including good initial state preparation, which needs to have a good overlap with the ground state. However, let us proceed step-by-step.

The standard version of this algorithm,^{21} which finds applications far beyond chemistry, works as follows. Suppose we have a unitary operator *U* and one of its eigenstates, |*ψ*_{n}⟩. We can evaluate its eigenvalue $\lambda n=ei2\pi \varphi n$ (expressible as a function of its phase since *U* is unitary) using an extra register of *r* qubits and *r* controlled operations $cU$, where the first operation controls the application of *U*, the second of *U*^{2}, the third *U*^{4}, and so on, until the final $U2r\u22121$. Finally, it is sufficient to apply a quantum Fourier transform to read the phase value in binary form, truncated to *r* bits.

To arrive at an implementation that interests us, we simply identify this generic unitary operator *U* with *e*^{iHt}, i.e., a Hamiltonian evolution operator, and the starting state as an eigenstate of *H*, for example, the ground state. In this case, the phase we read is *E*_{n}*t*. The operation $U2r\u22121$ translates into an evolution of time *t* 2^{r−1}.

It is possible to show that the error (due to truncation in *r* bits) we make on the phase *E*_{n}*t* scales as 1/*r*, where *r* is the total number of applications of *U*. This is because the discretization error scales as 2^{−r}, but each application of the controlled unitary doubles the length of the circuit. In the literature, this is often quoted as a quadratic speed-up compared to a Monte Carlo evaluation, whose error scales as $1/M$, where *M* is the number of samples and thus the number of function calls of the function to be evaluated on the generated distribution. If we interpret *r* as the number of “function calls,” i.e., $cU$ operations applied to the state |*ψ*_{n}⟩, the asymptotic comparison with Monte Carlo can be made, keeping these specifications in mind.

There are now two complications to consider. The controlled operation *e*^{iHt} cannot be implemented exactly but requires approximations. The most established is the Trotter step decomposition, which has the advantage of not requiring additional qubits.^{36,82} Recently, other methods that have better asymptotic scaling but require additional qubits, such as the linear combination of unitaries^{83} and qubitization,^{84} have surpassed Trotter’s method in popularity. For example, the first runtime and resource estimation works for chemistry by Reiher and co-workers^{32} assumed that QPE with Trotter had to be used, while now more recent works use qubitization (plus many other tricks to shorten the circuit).^{15} Notice, however, that the empirical performance of Trotter methods can be better than predicted upper bounds, and this is still an active area of research.^{82,85,86}

The first crucial observation is that now even a ground-state calculation requires a black-box that implements real-time dynamics or closely matching objects. This brings us back to the initial discussions about the complexity of implementing Hamiltonian simulations of fermionic systems, much more complex than their spin-lattice counterparts. Given that we require quantitative precision on the time evolution, the Hamiltonian evolution algorithm requires a full fault-tolerant implementation.

The second observation is that these algorithms require the challenging assumption of having the ground state as their input. What happens if the state to which we apply QPE estimation is not the ground state? Here is the good news and bad news. Let us start with the good news: unlike the classical case, if we input a generic state Φ (classically, the equivalent would be preparing a generic *Ansatz* and sampling with Metropolis from it) when we read the phase register, the state must collapse onto an eigenstate of |*ψ*_{n}⟩ and the read-out phase is *ϕ*_{n}. Therefore, the energy readout in the auxiliary register determines the collapse onto an eigenstate of *H* of the previously initialized state Φ. The good news is, therefore, that we will not read a random number but one of the possible eigenvalues. The bad news is that we do not know which one. In general, we will read the eigenvalue *n* with a probability given by the overlap |⟨Φ|*ψ*_{n}⟩|^{2}. It then becomes crucial that, if we are interested in the ground state, the initial state is not completely random but has a sizable overlap with the ground state.

Generally, in chemistry and materials science, we are interested in the runtime scaling with size, i.e., number of electrons, basis set, etc. If the overlap vanishes exponentially with size, the entire procedure becomes exponentially long, nullifying the exponential advantage that could be initially envisioned, given the compression in memory and the possibility of efficiently reading the energy. A recent study focused on this aspect, showing that this issue, though known in principle but often forgotten in practice, could seriously undermine the claim of exponential advantage for electronic structure.^{17} It should be noted, however, that even a polynomial advantage could be sufficient to solve problems that are still intractable, as can be seen from how Density Functional Theory has revolutionized chemistry and materials science thanks to its improved *N*^{3} scaling compared to the *N*^{6}, *N*^{7} scalings of coupled-cluster.

To conclude, the absence of an exponential speed-up does not rule out the existence of a practical quantum advantage, which is more difficult to identify *a priori*, but on a case-by-case basis. In this context, resource estimation studies focusing on particular molecular systems are of great importance. Goings *et al.*^{15} perform resource estimates to simulate a challenging system for classical methods, the cytochrome P450 enzyme. The estimates depend on the hardware noise that needs to be corrected. To simulate the ground state using an accurate active space, $\u223c5\xd7106$ (5 × 10^{5}) physical qubits with error rates of 0.1% (0.001%) would be needed. Concerning materials science systems, state-of-the-art studies are represented by Rubin *et al.*^{87} and Ivanov *et al.*,^{88} which move away from the plane-wave basis set and combine Bloch or Wannier orbitals, respectively, with most recent techniques such as sparse qubitization or tensor hyper-contraction. Resource estimates applied to the lithium nickel oxide battery cathode^{87} and transition metal oxides^{88} indicate longer fault-tolerant runtimes compared to molecular systems such as P450. Clearly, such estimates are based on state-of-the-art algorithms, including the most efficient way to encode fermionic Hamiltonians for phase estimation and the current state-of-error correction algorithm. Furthermore, algorithmic developments will improve the cost of the simulations.^{14} Orders of magnitude in efficiency have been gained compared to just ten years ago^{37} and, therefore, the threshold for quantum advantage could shift in one direction or another, approaching when new quantum strategies are invented or perhaps moving away thanks to the constant progress of “conventional” methods such as density matrix renormalization group (DMRG) or QMC.

## VII. THE LOCAL ENERGY IN VARIATIONAL MONTE CARLO

*ψ*⟩

*∑*

_{x}|

*x*⟩⟨

*x*| in the denominator and numerator. Notice that here we use the notation for a discrete Hilbert space, but the formula can be generalized to continuous models by replacing the sum with the integral (

*∑*

_{x}→

*∫dx*). Some steps are necessary to transform Eq. (10) into the typical Monte Carlo format, where we integrate the product of a probability distribution from which we can sample form,

*p*(

*x*), and an objective function. This is achieved formally by dividing and multiplying by

*ψ*(

*x*) = ⟨

*x*|

*ψ*⟩. Equation (10) then becomes

*local energy*is defined as

*x*(see below). The probability distribution is defined as

*p*(

*x*) = |⟨

*ψ*|

*x*⟩|

^{2}/

*∑*

_{x}|⟨

*ψ*|

*x*⟩|

^{2}. Now, if we assume that we can sample configurations

*x*∼

*p*(

*x*), i.e., using a Markov-chain algorithm, the energy can be evaluated stochastically as

*M*

_{VMC}decorrelated samples taken from a Markov-chain algorithm, such as Metropolis.

^{89}This technique is called Variational Monte Carlo (VMC)

^{7}and is efficient as long as (i) computing

*E*

_{L}(

*x*) is efficient and (ii) it is possible to run the Metropolis algorithm also efficiently. Since classical trial functions

*ψ*(

*x*) can be evaluated with numerical precision for each

*x*, it is also the ratio |

*ψ*(

*x*′)|

^{2}/|

*ψ*(

*x*)|

^{2}for each pair

*x*,

*x*′, which is needed to perform a Metropolis update.

^{90}

One of the most important features of the local energy is that its variance is zero when |*ψ*⟩ is an eigenstate of *H*. Indeed, if *H*|*ψ*_{0}⟩ = *E*_{0}|*ψ*_{0}⟩, then *E*_{L}(*x*) = ⟨*x*|*E*_{0}|*ψ*⟩_{0}/⟨*x*|*ψ*_{0}⟩ = *E*_{0}. In practice, this means that the local energy function will be closer to a constant value as the trial state |*ψ*⟩ approaches the ground state. This results in reduced statistical fluctuation in Eq. (13).

At the same time, it is possible to use the variance of the local energy as a cost function for optimization. This allows, in principle, to certify the success of the minimization, as the ground state is signaled by zero statistical error.

### A. The local energy in practice

*Ansatz*(for fermionic systems), usually complemented with an explicit correlation operator like the Jastrow factor.

^{7,91}In this case, evaluating the local energy reduces to applying the Laplacian operator to the function in real space and dividing by the function itself. For the sake of clarity, let us consider a toy example. The local energy for a (unnormalized) Gaussian trial

*Ansatz*$\psi (x)=e\u2212\theta x2$, in continuous space, and a typical one-dimensional Hamiltonian

*H*= −1/2(

*∂*

^{2}/

*∂x*

^{2}) +

*V*(

*x*) reads

*x*and the variational parameter

*θ*, which can be optimized in an outer loop. In this case, it can be observed that if the external potential is a harmonic oscillator

*V*(

*x*) =

*ω*/2

*x*

^{2}, the local energy becomes

*x*when the variational parameter takes the value

*θ*=

*ω*/2, for which the variational

*Ansatz*becomes exact. Moreover, it also takes the value $ELh.o./opt.=\omega /2$ with zero statistical fluctuations in Eq. (13) as

*E*

_{L}does not depend on the sampled point

*x*

_{i}anymore. Modern codes for solving the many-electron Schrödinger equation in chemistry or materials science feature sophisticated trial

*Ansatz*, which are in turn functions of atomic orbitals.

^{70}While in the past, introducing a new

*Ansatz*required coding new functions for evaluating the derivatives, now the evaluation of the local energy can be delegated to algorithmic differentiation routines. This allows for the adoption of fairly sophisticated

*Ansätze*in VMC.

^{70,92}

*H*

_{x,x′}= ⟨

*x*′|

*H*|

*x*⟩, where

*x*,

*x*′ can represent a specific spin configuration or an occupation state of fermions or bosons on a lattice. In this discrete basis, the local energy is written as follows:

*x*′ such that the Hamiltonian matrix elements |

*H*

_{x,x′}| ≠ 0, at fixed

*x*, is only polynomially increasing.

### B. Pauli measurements vs local energy

*σ*

^{α}are Pauli matrices, and consider the critical transition point at

*J*= Γ = 1. We also consider a short chain of

*L*= 10. We denote a generic computational basis configuration as |

*x*⟩ = (

*s*

_{1}, …,

*s*

_{L}), where

*s*

_{k}are eigenvalues {1, −1} of the $\sigma jz$ operator.

In this case, the spin Hamiltonian is already expressed in Pauli terms (one just needs to re-define the eigenvalue of the spin-*z* operator from {1, −1} to {0, 1}). Regarding the VQE approach, the energy can be measured in only two bases: the computational basis and the “*XX*…*X*” basis, obtained by applying a Hadamard gate, H, to each qubit at the end of the circuit that prepares the variational state.

In this numerical experiment, we use the variational Hamiltonian form with a sufficiently long circuit of up to 24 layers, resulting in up to 48 variational parameters (see the Appendix). By optimizing *Ansätze* characterized by different circuit depths (without shot noise for simplicity), it is possible to obtain trial states systematically closer to the exact ground state of the model.^{73} In Fig. 2, we use depths ranging from 12 to 24, and we can reach a relative error on the energy of 10^{−5} compared to the exact ground state energy, *E*_{0}.

However, the statistical error on the energy, which is evaluated with Eq. (5), does not improve. In fact, if we had tried to optimize the circuit using the noisy energy estimator, we would not have been able to obtain such accurate optimized trial states. This clearly demonstrates that the estimator does not possess the zero variance property, as opposed to the VMC calculation.

To obtain the standard deviation in Fig. 2, we repeat the estimation of the variational energy, Eq. (5) [Eq. (13) for the VMC case described below], 100 times to obtain a population of variational energies that could be obtained with the given variational setting, *M*_{j} (*M*_{VMC} for the VMC case) setups. We use a number of shots *M*_{j}, *M*_{VMC}, which is smaller (10^{2}), equal (10^{3}), and larger (10^{5}) than the Hilbert space of the model, i.e., 2^{10} = 1024.

For the VMC comparison, we deliberately use a fairly simple classical *Ansatz*, a long-range Jastrow state, which features only five variational parameters for *L* = 10 (see the Appendix). Although this classical *Ansatz* only reaches a moderate relative accuracy of 10^{−3}, at best, the statistical error on the energy consistently improves, outperforming the statistical error obtained with the quantum circuit. Notice that this is an easy model for VQE: the number of measurement basis is the minimum possible for a genuine quantum many-body system. Electronic structure Hamiltonians unfold into thousands of Pauli operators, which in turn require similar numbers of basis. This numerical example demonstrates the power of the local energy-based estimator compared to the Pauli measurement one. From this example, we can also understand the following lesson: even finding the smallest possible set of basis to measure *H* will not solve all our problems, as this estimator still lacks the zero-variance property.

A hybrid solution has been proposed by Torlai *et al.*^{93} They use quantum state tomography, using neural networks^{94} and a tomographically incomplete basis set, to obtain a classical reconstruction of the quantum state. Classical VMC can then be applied to this classical approximation to calculate, precisely, the energy. This method solves the variance problem, but it introduces a bias stemming from a possibly, and likely, imperfect reconstruction of the quantum state. Moreover, it raises the question of finding the range of applicability of the method. If the quantum states can indeed be represented by a classical *Ansatz*, then one could directly reach the ground state by optimizing that without the need for a quantum computer.

### C. Beyond variational: Projective methods

Local energy is a central concept in quantum Monte Carlo beyond the simplest VMC method because, in practice, every projective QMC method requires a trial wave function *ψ* to alleviate the sign problem in fermionic simulations or, more generally, to reduce statistical fluctuations. These projective methods, such as diffusion Monte Carlo^{7} or Auxiliary Field Monte Carlo (AFQMC),^{95} improve upon the VMC energy but still rely on a variational state for importance sampling. Therefore, local energy resurfaces in these contexts as well. An accurate QMC simulation is rarely seen without a good variational starting point.^{96}

One can, therefore, see a similarity between the importance of VMC, which is foundational for a more accurate calculation with projective QMC, and the significance of the initial state preparation for a successful execution of QPE in quantum computing. It is highly likely that this duality between variational and projective methods (in imaginary time in the classical case and in real-time in the quantum case) will extend to quantum computing. In that case, algorithms like VQE or its variants, despite being considered old-fashioned by some, will remain central even in the fault-tolerant regime as state preparation subroutines.

## VIII. QUANTUM COMPUTING MEETS QUANTUM MONTE CARLO

### A. The local energy in quantum computing

The existence of a local energy estimator in quantum computing would eliminate any variance problem in VQE. However, it is not as straightforward to apply this trick in quantum computing simply because evaluating the ratio ⟨*x*′|*ψ*⟩/⟨*x*|*ψ*⟩ becomes extremely demanding in general.^{97}

*ψ*(

*x*) = ⟨

*x*|

*ψ*⟩ needs to be evaluated from quantum measurements and, hence, is affected by statistical noise. While evaluating

*ψ*(

*x*) with additive precision is possible, the local energy involves ratios of amplitudes. Maintaining a fixed precision on the ratio

These statistical fluctuations are different compared to those found in a standard Monte Carlo calculation. In VMC, the local energy can always be computed with numerical precision, and the fluctuations arise from a finite number of samples, *M*_{VMC} in Eq. (13) (in the presence of an approximate trial state). Here, uncontrolled statistical fluctuations arise solely from the estimation of the local energy at a fixed *x*_{i}.

We are witnessing an increase in work that aims to combine quantum computing and quantum Monte Carlo. Huggins and co-workers proposed an interesting combination of quantum computing and AFQMC.^{98} In this work, they use a circuit to generate the trial wave function, from which samples are drawn (in this representation, the configuration *x* is a Slater determinant). The AFQMC algorithm then proceeds unchanged, and the supposed advantage of the method lies in using a circuit to generate an *Ansatz* that could be inaccessible classically. Mazzola and Carleo^{99} showed that the procedure, when adapted to many-body lattice models at criticality, thus using Green’s function Monte Carlo instead of AFQMC, exhibits exponentially scaling behavior with a hard exponent. This is due to Eq. (18) and the fact that strongly correlated states have vanishing overlaps on the configuration’s basis, necessitating an exponentially increasing number of samples to compute the local energy. It is estimated that a reasonably accurate ground state calculation of a 40-site transverse-field Ising model [Eq. (17)] requires an order of 10^{13} measurements. Inserting a gate frequency of 10 kHz (i.e., assuming a fault-tolerant implementation; see Sec. II) and a circuit depth of $O(10)$ layers to generate an accurate trial state implies a runtime of a few thousand years for a system in reach of exact classical diagonalization.

Other works appeared almost simultaneously last year on this topic. Zhang *et al.*^{100} introduced a quantum computing adaptation of the Full Configuration Interaction Quantum Monte Carlo (FCIQMD).^{101} In this work, a quantum circuit *U* is used to create a “quantum” walker $|x\u0303\u3009=U|x\u3009$, i.e., a linear combination of Slater determinants |*x*⟩, undergoing the FCIQMC subroutine. The idea is interesting as it could counter the exponential explosion of the determinants/walkers during the imaginary time projection by compressing logarithmically the memory to store them. A possible major drawback of this method is that the Hamiltonian $Hx\u0303,x\u0303\u2032$ on this new basis is not sparse anymore. Xu and Li^{102} proposed to use Bayesian inference to reduce the number of shots required to compute the local energy. Kanno *et al.*^{103} further combine the ideas of Ref. 98 with tensor networks. Yang *et al.*^{104} propose instead a way to speed-up *real*-time path-integral MC already on NISQ hardware. Tan *et al.*^{105} instead propose integration with stochastic series expansion, another flavor of QMC used for spin models.

Finally, two recent works propose to use quantum data in a conventional VMC framework. In this case, the local energy is calculated using conventional hardware. Montanaro and Stanisic^{106} propose the use of a VQE circuit as an important sampler to speed-up the first iteration of a VMC simulation. Moss *et al.*^{107} use quantum data from Rydberg atom simulators to train a classical neural-network *Ansatz* (as in Ref. 93) and further optimize it in a VMC fashion.

Overall, it is likely that an efficient way to estimate the local energy is possible only for sparse states, i.e., for which the number of non-zero overlaps ⟨*x*|*ψ*⟩ ≠ 0 grows only polynomially with the system’s size. However, it remains to be understood if a quantum computer is really needed to tackle such systems at this point.^{108} Furthermore, if a suitable basis transformation *U* can be found to reduce the support of such states, then (1) this transformation should not spoil the sparsity of the Hamiltonian *H*_{x,x′} to keep the evaluation of Eq. (16) efficient. (2) Moreover, if this transformation exists, it can be used to efficiently diagonalize the systems in a reduced sub-space without the need for QMC.

On a more positive note, it is not excluded that, despite exhibiting exponential scaling, the aforementioned approaches could yield a better exponent than the best classical method for some specific fermionic systems. To achieve this, it will be crucial to start with a classically intractable trial state to justify the subsequent imaginary-time projection. Further research and methodological advancements are required to assess the true potential of the method in the presence of shot noise.

Overall, the pursuit of an efficient method for calculating energy, inspired by the local energy in VMC, is a field of research that we hope will yield numerous fruitful results. It is necessary for the quantum computing and QMC communities to clearly understand the limitations and potentialities of their respective techniques in order to invent new hybrid algorithms at the interface of these two worlds.

### B. Classical-inspired circuits for VQE, quantum-inspired *Ansätze* for VMC

The techniques and methods that have been used for decades in QMC are so numerous that many have been (and many are waiting to be) exported to quantum computing. Trial functions play a central role in VMC. The use of explicitly correlated non-separable *Ansätze* has brought great success to VMC and is basically a clever solution to compress the electronic wavefunction, which, when described in the space of determinants, requires otherwise an exponential number of coefficients. The latest iteration of this concept is the introduction of neural network quantum states by Carleo and Troyer in 2017,^{109} which can be seen as more general forms of Jastrow,^{91,94} back-flow,^{110,111} and tensor network states.^{112}

As mentioned earlier, compressing the Hilbert space within a polynomial scaling size quantum memory enables the manipulation of linear combinations of arbitrarily large Slater determinants. However, when considering the variational approach, we are constantly seeking shorter quantum circuits that can capture as much entanglement as possible, within the coherence time limitation of NISQ systems.

Several works have already proposed ways to implement a Gutzwiller operator, which is essentially the simplest form of a Jastrow operator, as a quantum circuit. Murta and Fernandes-Rossier^{113} propose a method based on post-selection. Typically, the way to create non-unitary operators in quantum computing is by embedding them in a larger system that undergoes unitary evolution, a method also known as “block encoding.” This involves introducing ancillary qubits, and it can be certified that the non-unitary operator has been successfully applied to the quantum state if and only if the ancillary register is measured and read in a given state. However, the problem with this approach is that the success probability decreases with the system size, requiring many repetitions.

Seki and co-workers^{114} also propose a similar approach based on the linear combinations of unitaries and, therefore, also affected by a finite success probability.

Using a different approach, Mazzola and co-workers^{115} implicitly defined a hybrid quantum–classical wavefunction with a Jastrow operator in post-processing. The approach has since improved in its scalability (Ref. 119). There, a quantum circuit is used as an important sampler, and the measured configurations undergo post-processing by a neural network. Benfenati and co-workers^{117} instead implemented a Jastrow operator, moving it from the wavefunction to the Hamiltonian. This approach also does not require additional circuits compared to a VQE calculation. However, the re-defined Hamiltonian operator features much more Pauli terms to measure. Motta and co-workers devised an imaginary time evolution (QITE) operator without ancilla and post-selections.^{120} The original formulation of the method formally incurs an exponential dependence on the correlation length in the general case because it requires quantum state tomography. However, if truncated, it can generate heuristic trial states for variational calculations and is still subject to improvement.^{119}

Finally, it is interesting to note that the flow of information is not always from the older discipline to the newer one. Some circuit *Ansatz* used in quantum computing can be adapted to VMC. Inspired by the Hamiltonian variational circuit *Ansatz*^{8} (see the Appendix), Sorella devised a method called variational AFQM capable of obtaining state-of-the-art ground state energies of the Hubbard model for various *U*/*t* parameters and dopings in the thermodynamic limit.^{120}

### C. Variational real-time dynamics and updates in parameters’ space

Variational states are not only used for ground state calculations but they can also be used to study dynamics. The price to pay is that the variational state must be flexible enough to accurately describe excited states, and this can be a demanding constraint, while the advantage is the ability to use much shorter circuits compared to those used, for example, for trotterization.

From a classical perspective, this area of research has been very active in recent months, as it allows for countering quantum advantage experiments in the quantum dynamics application space (cf. Sec. III). Obviously, the use of a variational state does not allow for exact evolution, but it is also true that the errors of a NISQ machine do not allow for it either. The balance between classical and quantum advantages for real-time dynamics will be shifted in favor of the latter when fidelity enables the simulation of sufficiently large systems for a sufficiently long time, rendering them inaccessible to classical approximation methods.^{22}

The subfield of variational real-time dynamics also offers interesting parallels between quantum computing and (time-dependent) VMC.^{121} The formalism based on the time-dependent variational principle is the same. In practice, even the fundamental ingredients that allow for the update of variational parameters are the same: the matrix *S* defined in Sec. V B (cf. Ref. 121 with Ref. 124). As we have seen in the case of optimization, the fact that the elements of the *S* matrix are subject to statistical noise is a common issue in both implementations. In this case, as well, it is reasonable to expect cross-fertilization between the two techniques regarding both variational forms and efficient ways to evaluate the *S* matrix. The field of variational algorithms for real-time simulations is very active. In this area, concepts borrowed from tensor-network simulations are also useful to shorten the circuit.^{123–125}

Generally speaking, variational parameters can be updated using different pseudo-dynamics, ** θ**′ =

**+**

*θ**δ*

**, to achieve various objectives. While pure energy minimization is the most popular goal and real-time evolution following the time-dependent variational principle is the second, there are other possibilities. Patti**

*θ**et al.*

^{126}devised an iteration scheme to perform Markov chain Monte Carlo in the quantum circuit’s parameter space, i.e., to sample from

*p*(

**) ∼ exp[−**

*θ**β*⟨

*ψ*(

**)|**

*θ**H*|

*ψ*(

**)⟩]. The resulting equation is a generalization of stochastic gradient descent that ensures detailed balance. This approach could assist in escaping local minima during VQE optimization.**

*θ*Similar ideas have been proposed in the VMC context earlier. Mazzola *et al.*^{77} showed that one can obtain an upper bound for the free energy by sampling from $p(\theta )\u223c|S(\theta )|exp[\u2212\beta \u27e8\psi (\theta )|H|\psi (\theta )\u27e9]$. In VMC, this can be achieved either using a modified Langevin equation for ** θ** or a modified Metropolis acceptance, also known as the “penalty method.”

^{127}

In conclusion, manipulating trial states in the presence of statistical noise is a common feature of VMC in all its formulations and scopes. Many ideas have been proposed to achieve stable parameter updates. The VQE community could profit from this established knowledge but also share its own developments and ideas to advance both fields.

## IX. CLASSICAL MONTE CARLO MEETS QUANTUM COMPUTING

In this section, we completely shift our perspective. Not all chemistry problems are quantum many-body ones; for example, understanding protein folding is already a daunting task in its classical force-field formulation. Likewise, not all problems that a quantum computer can solve are genuine quantum mechanical problems. In fact, in many cases, the opposite is true: the most famous quantum algorithms, which have made the field of quantum information renowned, are focused on solving “classical” problems. For instance, Shor’s algorithm provides exponential speed-up for factoring integers, and Grover’s algorithm enables quadratic speed-up for searching in databases.^{21} Other examples include algorithms for linear algebra, optimization, and machine learning. Philosophically speaking, solving a purely classical problem with a quantum machine can be even more intellectually rewarding than simulating a quantum system, where the distinction between computation and simulation becomes less clear.

Up to this point, we have been exploring whether and how well-known techniques in quantum Monte Carlo can be adapted to quantum computing to simulate many-body quantum systems. Now we are questioning the opposite: Can a quantum computer be useful in speeding up a classical Monte Carlo algorithm where the Hamiltonian is defined solely using classical variables, e.g., classical spins? In addition, more specifically, can we achieve this already on NISQ machines?

### A. Autocorrelation of a Markov chain

Markov chain Monte Carlo (MC) algorithms are of fundamental importance in both science and technology to understand models that lack a simple analytical solution.^{128,129}

MC methods aim to generate statistically independent, representative configurations *x*_{i}, belonging to the computational space, distributed as a target Boltzmann distribution, *ρ*(*x*) = exp(−*βV*(*x*)), at a finite inverse temperature *β* = 1/*T*, where *V*(*x*) is a classical potential energy. A Markov chain MC algorithm sequentially generates these representative configurations through a transition probability matrix *P*(*x*, *x*′), which defines which states *x*, *x*′ can be connected along the chain and the relative probability of the transition *x* → *x*′ (each row of the matrix *P* is normalized to one).^{130}

Among the family of Markov Chain MC algorithms, the Metropolis algorithm is certainly the most popular one.^{89} Here, the transition process takes the form *P*(*x*, *x*′) = *T*(*x*, *x*′)*A*(*x*, *x*′), where *T*(*x*, *x*′) and *A*(*x*, *x*′) are, respectively, the *proposal* and the *acceptance* probability matrices. The algorithm works as follows: when at state *x*, a candidate trial configuration *x*′ is generated from the distribution *T*(*x*, ·). The trial configuration is accepted with probability *A*(*x*, *x*′). If accepted, the next element of the chain becomes *x*′; otherwise, it remains *x*. If *T* is a symmetric matrix, then the Metropolis acceptance is given by *A*(*x*, *x*′) = min[*e*^{−β(V(x}′^{)−V(x))}, 1].^{89} It is clear that the efficiency of the algorithm strongly depends on the choice of *T*(*x*, *x*′).^{131–133} The efficiency is given by the relaxation or mixing time, which quantifies the speed of convergence toward the equilibrium distribution *ρ*(*x*) and is formally given by the inverse of the gap, *δ*, between the largest and the second-largest (in modulus) eigenvalues of *P*. Two limiting cases exist; the first involves a *local* update scheme based, for instance, on some physical intuition about the system (e.g., a single spin-flip). This usually produces a new configuration, *x*′, that is similar to the parent *x*. This choice increases the acceptance rate because *V*(*x*) ∼ *V*(*x*′) but also results in a long sequence of statistically correlated samples, such that long simulations are needed to thoroughly explore the configuration space. On the contrary, a *non-local* update scheme is more effective in producing uncorrelated samples, but usually at the expense of a vanishing acceptance rate.

Interestingly, it took about 30 years after the invention of the Metropolis algorithm in 1953 before the introduction of efficient non-local update schemes for lattice models, the Swendsen-Wang^{134} and the Wolff^{135} algorithms. These *cluster* updates solved the *critical slowing down* of MC simulations at phase transitions in ferromagnetic Ising models, but unfortunately, they are not as effective for frustrated models.^{136,137}

### B. Sampling from a classical Boltzmann distribution using wavefunction collapses

Recently, it has been proposed to use a quantum computer, digital or NISQ, to generate efficient non-local Metropolis updates *T*(*x*, *x*′) for spin systems. The theoretical framework was first introduced by Mazzola^{138} in 2021. Shortly after, Layden *et al.*^{139} demonstrated on real quantum hardware a quantum-enhanced Markov chain algorithm.

Following Ref. 138, the general idea is rooted in the Fokker–Plank formalism of non-equilibrium statistical mechanics in continuous systems.^{140} The Fokker–Plank operator *H*_{FP} is a parent quantum Hamiltonian of the physical potential *V*(*x*), whose spectrum is closely connected with the number of local-minima of *V*(*x*). For instance, in a double-well model, *H*_{FP} has two lowest-lying eigenvalues and eigenstates corresponding to the symmetric (antisymmetric) combination of the two Gaussian localized states in the two wells, |*ψ*_{L}(*x*)⟩, |*ψ*_{R}(*x*)⟩.

This idea can be ported to lattice models. Let us consider, for simplicity, the ferromagnetic Ising model, defined by *V* = *H*_{1} in Eq. (17), as our classical potential. The task is to sample from the classical Boltzmann distribution exp[−*βH*_{1}(*x*)]. Here, the autocorrelation time of a local spin-flip update scheme is dominated, at low temperatures, by the rate of the *rare-event* processes that drive the system from the *all-up* (left) (↑↑⋯↑)≔|*ψ*_{L}⟩ to the *all-down* (right) (↓↓⋯↓)≔|*ψ*_{R}⟩ classical states. These processes are necessarily characterized by a nucleation event exponentially suppressed with *β* and the subsequent diffusion of the domain wall separating the ↑ and the ↓ regions.^{141} If, however, one can construct a quantum Hamiltonian *H* such that its low-lying eigenstates are mostly localized on the states (↑↑⋯↑) and (↓↓⋯↓), one could then prepare and sample configurations from these states with optimal autocorrelation times through repeated collapses of these wavefunctions.

This quantum Hamiltonian could be the quantum transverse field Hamiltonian in Eq. (17), where we add a quantum driver *H*_{2} to the classical potential, *H*_{1}. Here, in the small Γ limit, the two (unnormalized) lowest-lying eigenstates of *H* are |*ψ*_{0}⟩ ≈ |*ψ*_{L}⟩ + |*ψ*_{R}⟩ and |*ψ*_{1}⟩ ≈ |*ψ*_{L}⟩ − |*ψ*_{R}⟩.^{142} While the gap *E*_{01} between these states is exponentially vanishing with the system size, the gap between the rest of the states remains $O(1)$. It is clear that by sampling configurations *x*^{143} from either |*ψ*_{0}⟩ or |*ψ*_{1}⟩, we can achieve optimal autocorrelation times in the large *β* limit, as the states |*ψ*_{L}⟩, |*ψ*_{R}⟩ are sampled with equal probability.

At intermediate temperatures, the probability distribution obtained via the eigenstates projection is generally different from the classical Boltzmann distribution $\rho (x;\beta )=e\u2212\beta H1(x)$ one aims to achieve. For this reason, a standard Metropolis acceptance step needs to be performed.

One can define a valid Markov chain out of this physical intuition, rooted in (1) the preparation of a localized state, (2) a quantum propagation to prepare a linear combination of the low-energy eigenstates of *H*, and (3) a measurement to collapse into a new string state. In Ref. 138, it is proposed to use a QPE subroutine to prepare such low-energy eigenstates. Layden *et al.*^{139} simplify the idea and use a Hamiltonian simulation subroutine, *e*^{−iHt}, with randomized *t* and Γ values at each step. Crucially, they observe that such a quantum proposal update is symmetric, thus enabling a fast and practical evaluation of the acceptance step. The algorithm is sketched in Fig. 3.

A superquadratic speedup for spin glass instances is observed, i.e., a polynomial speed-up of order 3.6, compared with the best possible classical update strategy. Interestingly, the procedure is demonstrated on hardware, where it is found that hardware errors only impact the efficiency of the chain while the sampling remains unbiased, exactly due to the existence of a classical acceptance step after the quantum, noisy proposal move.

This quantum-enhanced Markov chain MC needs a quantum dynamics subroutine. This, in turn, can be implemented both in the fault-tolerant NISQ regime and in analog simulators, depending on their architecture constraints. Clearly, to reach a possible quantum advantage, one needs to deal with the fact that, as explained in Sec. II, quantum gate times are much slower than classical CPUs, and this may cancel a scaling advantage for reasonable system sizes. In particular, lattice MC simulations can also be executed on special-purpose classical hardware, such as Field Programmable Gate Arrays (FPGA) as demonstrated in Refs. 144, which may enjoy an even faster logical clock speed. Therefore, more work will be needed to assess whether this idea can bring a real benefit in this application space.

### C. Quantum walks and quantum Metropolis algorithms

For the sake of completeness and clarity, it is important to mention here another family of algorithms, known as quantum walks, which share the same objective as the quantum-enhanced Markov chain algorithms described earlier: accelerating the convergence of classical Markov chains. While the justification for the quantum-enhanced Markov chain algorithm of Sec. IX B is based on physical intuition^{138} and the potential gains are assessed heuristically, quantum walks come with a more rigorous guarantee: if they can be implemented, they provide a quadratic speed-up in autocorrelation times. This quadratic speed-up is related to concepts such as amplitude amplification or Grover’s algorithm.^{145}

*x*, and reuse it if the trial move leading to

*x*′ is not accepted. However, in quantum computing, the

*no-cloning*theorem prohibits the direct copying of a quantum state.

^{146}Furthermore, the acceptance step also involves arithmetic operations that are computationally more demanding in the quantum context. It is easy to imagine that a unitary “walk” operator must include rotations by the following arbitrary angle:

^{147}Now, it is important to note that quantum arithmetic is much more expensive in quantum computing because it must be reversible. For instance, the most efficient way to perform the addition is still a subject of ongoing research.

^{148}

The most commonly used definition of a quantum walk originates from Szegedy.^{149} For the sake of brevity, we have to refer the reader to the comprehensive review^{150} or to Ref. 147, where Szegedy’s algorithm is revisited from a more practical perspective, for details. In short, for any classical Markov chain defined by the transition matrix *P*(*x*, *x*′) (see Sec. IX A), a quantum walk represents a quantized version of it that offers a quadratic speed-up in mixing time. Formally, it enhances the gap from *δ* to $\delta $, such that the mixing time decreases quadratically from $O(1/\delta )$ to $O(1/\delta )$.

*W*of the form

*W*would require digital rotation of angles such as Eq. (19), which are in turn evaluated by a sequence of costly arithmetic operators. Lemieux

*et al.*

^{147}analyze the cost of quantum walks, showing that the quadratic speed-up that they can offer is overshadowed by the cost of implementing the

*W*operator, assuming even optimistic estimates for the gate time of fault-tolerant hardware. The quantum-enhanced method of Sec. IX B completely avoids performing the acceptance step on the quantum hardware, which is the main reason for its hardware feasibility.

We note that it is increasingly easy to get confused with the names of the methods and the combinations of words such as “quantum,” “Monte Carlo,” and “Metropolis.” In the traditional literature as well as in this Perspective article, “quantum Monte Carlo” refers to the family of Monte Carlo algorithms that are executed on conventional computers but aim to solve many-body quantum problems. However, there is a community in quantum computing for which this combination of words indicates a Monte Carlo algorithm executed on quantum hardware, including Montanaro’s algorithm for computing expectation values of multidimensional integrals using quantum amplitude estimation.^{151} This type of application finds application in finance,^{33} and while interesting, it is not discussed in this manuscript.

Finally, we also note that a quantum computing-based method to speed up a Monte Carlo simulation could, in principle, be used to accelerate a QMC algorithm as well. In this case, the leap would be twofold: using quantum hardware to accelerate a classical algorithm, such as a path-integral MC, to simulate quantum Hamiltonians.

Two philosophically more interesting algorithms can be mentioned for this purpose. Temme *et al.*^{152} proposed a “quantum Metropolis algorithm” for studying quantum many-body Hamiltonians. This method suggests performing a walk in the eigenstates of the quantum Hamiltonian, thus overcoming the sign problem.^{153} The algorithm includes performing and undoing QPE, an ancilla register that stores the energy, and the ability to perform measurements that only reveal one bit of information, enabling the acceptance step to overcome the no-cloning principle.

Yung and Aspuru-Guzik^{154} chose a different strategy for their “quantum–quantum Metropolis algorithm,” finding a way to extend Szegedy’s walk to quantum Hamiltonians. The runtime performance and scaling of such approaches are still yet to be assessed.

### D. Sampling with quantum annealers or simulators

The algorithms presented in Sec. IX C require a fault-tolerant computer. However, it cannot be ruled out that an advantage in the sampling problem could come from hardware that falls on the opposite spectrum, namely noisy quantum simulators or quantum annealers.^{18} First of all, let us observe that the quantum-enhanced Markov chain method of Sec. IX B can be implemented not only using trotterization but also through real-time dynamics in a quantum simulator.^{139} Furthermore, optimization and sampling tasks are closely connected.^{155} Special-purpose quantum simulators, called quantum annealers, have been built with the aim of optimizing large-scale spin-glass problems, but they have also been reconsidered as thermal samplers^{156} with some specific applications in machine learning.^{157,158}

The possibility of using a quantum annealer as a sampler arises from its deviation from adiabaticity. The existence of vanishing gaps during annealing implies that at the end of the experiment, the wave function does not localize in the classical minimum of the cost function but remains delocalized, producing a distribution of read-outs. The presence of hardware noise amplifies this effect even further.

This residual distribution could resemble a thermal Boltzmann distribution of some classical Hamiltonian, close to the problem Hamiltonian originally meant to be optimized, and at some effective temperature, which is difficult to determine.^{156,159} However, given all the possible hardware and calibration errors, it is unlikely to be hoped that this approach can generate an unbiased sampling from a target distribution.

Recently, Ghamari *et al.*^{160} proposed the use of this annealing process as an important sampler. Similarly to the quantum-enhanced Markov chain method, detailed balance is restored using a classical acceptance step. In this case, as well, a control parameter is the annealing runtime, which generates a more-or-less-localized final distribution.

Finally, Wild *et al.*^{161} proposed an adiabatic state preparation of Gibbs states that can also bring quantum speedup over classical Markov chain algorithms and could also be implemented on NISQ Rydberg atom devices.

## X. CONCLUSIONS

In this Perspective, we investigate many intersections between quantum algorithms and Monte Carlo methods. We begin with a brief review of quantum computing applications in many-body quantum physics. We outline the consensus that is emerging after these years in which quantum computing has become mainstream. With the availability of quantum computers with $\u223c100$ qubits and the ability to implement gate sets, albeit noisy,^{12} the field is taking on a more practical connotation beyond the traditional boundaries of quantum information theory.

We observe that different hardware platforms imply different gate frequencies, which must be taken into consideration from the perspective of achieving quantum advantage. Quantum advantage for high-accuracy many-body ground state calculations is likely to be deferred to the fault-tolerant era due to the existence of hardware noise in today’s machines and the existence of highly developed classical competitors.^{89} Although the possibility of obtaining a quantum advantage through variational methods is possible, especially in systems where classical methods struggle,^{73} here we face the additional challenge posed by the presence of quantum measurements shot noise.

We then list several points of contact between quantum and classical variational methods. First, we explain the difference between the statistical noise present in conventional QMC algorithms and the noise arising from quantum measurements. Classical QMC methods feature an energy estimator—the *local energy*—that enjoys the zero-variance property. This, along with stable optimizers,^{75,162} enables the optimization of wave functions featuring thousands of parameters. This is not currently possible in variational quantum computing. Even with access to an exact ground state preparation circuit, obtaining the energy with sufficient precision requires a costly number of circuit repetitions.^{8} It is clear that this problem arises even before the circuit optimization stage. In the current literature, this aspect is often overlooked, as several new algorithms or circuits are tested without realistic shot noise conditions.^{66}

We suggest that finding the quantum equivalent of the *local energy* should be a priority in variational algorithm development. Along the same lines, we reviewed attempts and challenges in “quantizing” QMC methods.^{101–100,102} Other areas where we expect to see cross-fertilization between the quantum and classical worlds include the development of variational forms: classically inspired circuits for VQE and quantum-inspired *Ansätze* for VMC. Several essential ingredients for variational real-time evolution and parameter optimization under noisy conditions have already been put forward in the VMC community and will be instrumental for their quantum counterparts.

Finally, after discussing how knowledge of QMC methods can provide new momentum for the development of quantum algorithms, we take the opposite direction, showing that quantum hardware can bring advantages to Monte Carlo itself. In this space, quantum walks have been present in the literature for several years and achieve a quadratic speed-up in autocorrelation times through the quantization of a classical Markov chain.^{149} Their scaling is discussed, as is typical in quantum information, using an oracular form, which assumes the existence of key subroutines without delving into the details of their concrete gate-level implementation. Recently, it has been shown that these oracles require fairly long circuits. In their necessary fault-tolerant implementation, this implies absolute runtimes that are still slower than classical Monte Carlo, even for large-scale classical-spin simulations.^{147}

A more hardware-friendly possibility is represented by a family of methods that use a quantum computer as an important sampler or perform only the proposal part of a Metropolis update on the quantum hardware.^{138,139,160} Physically speaking, one simply leverages the fact that quantum measurements are uncorrelated, making them an efficient engine for sampling. In this case, shot noise is no longer a limitation but rather becomes the computational resource for quantum advantage.^{138}

Overall, the purpose of this Perspective is to further connect two communities: the quantum algorithms and the Monte Carlo ones. As mentioned, many methods developed in QMC can be repurposed in quantum computing. On the other hand, QMC can be a formidable competitor that can hinder or delay the quantum advantage. This is true not only for quantum chemistry applications but also for optimization and beyond. For instance, QMC can reproduce the scaling of quantum annealing machines for classical optimization purposes under certain conditions.^{142,163,164}

However, the two communities can be complementary, and we hope that new impactful algorithms, either quantum or classical, will emerge thanks to this interaction to solve important problems in chemistry and condensed matter.

## ACKNOWLEDGMENTS

I acknowledge discussions on these topics in the past years with S. Sorella, G. Carleo, P. Faccioli, A. Zen, K. Nakano, F. Tacchino, P. Ollitrault, S. Woerner, M. Motta, D. Layden, and M. Troyer. I dedicate this manuscript to the memory of Sandro Sorella, who invented many techniques mentioned earlier and inspired the whole QMC community with his creativity and enthusiasm. I also acknowledge financial support from the Swiss National Science Foundation (Grant No. PCEFP2_203455).

## AUTHOR DECLARATIONS

### Conflict of Interest

The author has no conflicts to disclose.

### Author Contributions

**Guglielmo Mazzola**: Conceptualization (equal); Writing – original draft (equal).

## DATA AVAILABILITY

The data that support the findings of this study are available from the corresponding author upon reasonable request.

### APPENDIX: HAMILTONIAN VARIATIONAL AND JASTROW *ANSÄTZE*

**Quantum circuit**. We use the Hamiltonian variational (HV)

*Ansatz*,

^{8}which is a very powerful and transferable

*Ansatz*that is inspired by the adiabatic principle. The unitary operator defining the HV

*Ansatz*is made of

*d*blocks, and each block is a product of

*ℓ*operators $U\u0302j=exp(i\theta jkH\u0302j)$, with

*j*= 1, …,

*ℓ*indexing the non-commuting terms of the Hamiltonian. For the transverse field Ising model, we only need

*ℓ*= 2, cf. Eq. (17). In this case, the full unitary operator is

^{⊗L}by placing one Hadamard gate on each qubit. The total number of parameters is

*ℓd*. In our numerical experiment, we use a state-vector emulation of the operator based on linear algebra operations, and we do not compile our operator into a real circuit. Parameters are optimized using COBYLA and BFGS optimizers. These results are compatible with Ref. 73.

Once we have obtained the optimized state, we emulate the noisy energy estimation using the Pauli measurement method. We sample *M*_{1} spin configuration from |*ψ*(*θ*)|^{2} in the computational basis and *M*_{2} = *M*_{1} in the rotated basis, obtained using the $H\u2297H\cdots H$ operator (cf. Ref. 49). The total cost of the energy evaluation is therefore 2*M*_{j}, where *M*_{j} = *M*_{1} = *M*_{2} values are reported in the text.

**Classical**. For the case of VMC, we use a long-range Jastrow state of the form

*Ansatz**E*

_{0}. To generate different

*Ansätze*of different quality and showcase the zero-variance property of the local energy, we simply act on the first parameter

*λ*

_{1}, pulling it away from its minimum value and toward smaller values (this is performed to create a more challenging, delocalized state), while keeping the other fixed. When

*λ*

_{1}= −0.15, the variational energy degraded up to a 10% systematic error.

We extract *M*_{VMC} configurations on a computational basis since we only need this basis to compute the local energy, Eq. (16), as the name suggests. To calculate the standard deviation of the estimator, we simply repeat the numerical experiment for each *Ansatz* and *M* setup 100 times.

## REFERENCES

*Quantum Monte Carlo Approaches for Correlated Systems*

*Quantum Computation and Quantum Information: 10th Anniversary Edition*

*Reinventing Physics from the Bottom Down*

*CCZ>*to 2|

*T>*transformation

*ab initio*electronic simulations by quantum Monte Carlo

*ab initio*molecular dynamics simulations on quantum computers

*ab initio*molecular dynamics calculation

*x*′)/p(

*x*) and perform the acceptance test.

^{4}

*Handbook of Markov Chain Monte Carlo*

*Markov Chains and Mixing Times*

*x*denote a lattice multi-spin configuration.