Quantum Computing for Quantum Tunnelling

We demonstrate how quantum field theory problems can be embedded on quantum annealers. The general method we use is a discretisation of the field theory problem into a general Ising model, with the continuous field values being encoded into Ising spin chains. To illustrate the method, and as a simple proof of principle, we use a (hybrid) quantum annealer to recover the correct profile of the thin-wall tunnelling solution. This method is applicable to many nonperturbative problems.


I. INTRODUCTION
There has been increasing interest in the possibility of simulating Quantum Field Theory (QFT) on quantum computers [1], with the development of efficient algorithms to compute scattering probabilities in simple theories of scalars and fermions [2][3][4][5][6][7][8][9][10][11][12][13][14][15][16][17]. In particular it is known that by latticizing field theories, quantum computers should be able to compute scattering probabilities in QFTs with a run time that is polynomial in the desired precision, and in principle to a precision that is not bounded by the limits of perturbation theory. However a particularly difficult aspect of this programme is the preparation of scattering states [4-6, 8, 9, 14-17], with several works having proposed a hybrid classical/quantum approach to solving this problem [11,[17][18][19]. A complementary approach is to map field theory equations to discrete quantum walks [20][21][22][23] which can be simulated on a universal quantum computer.
In this paper we point out that certain nonperturbative quantum processes do not suffer from this difficulty, and lend themselves much more readily to study on quantum computers in the short term. These are the tunnelling and related processes, which are of fundamental importance for the explanation of quantum mechanical and quantum field theoretical phenomena, for example transmission rates of electron microsopes, first-order phase transitions during baryogenesis, or the potential initiation of stochastic gravitational wave spectra in the early Universe and many more.
Typically in tunnelling, the system begins in a false vacuum state that is non-dynamical and virtually triv-ial. The initial state can be very long lived, with tunnelling to a lower "true" vacuum state taking place via non-perturbative instanton configurations. In principle in such a process, the confinement of the initial state to a false vacuum prepares the state for us, so that the analytically straighforward perturbative phenomena are paradoxically the quantum computationally more difficult ones.
As opposed to quantum computing realised by a series of discrete "gate" operations, quantum annealers [24,25] perform continuous time quantum computations, and therefore they are well-suited to the study of tunnelling problems by direct simulation (although our discussion could ultimately be adapted to gate-model quantum computers as well) [26][27][28][29][30][31][32][33][34][35][36]. In particular these devices, produced by D-Wave Systems [37], can be seeded with initial conditions using the "reverse annealing" feature, [38] allowing the simulation of dynamics. In contrast with the quantum-gate devices, they are already quite large, 2048 qubits in the current generation, with work ongoing to develop much more connected 5000 qubit machines. Moreover they operate in a dissipative rather than fully coherent regime, which is likely to be realistic for many real theories in which there are interactions with matter. In the present context this would be relevant for studies of so-called thermal tunnelling rather than (or in addition to) quantum tunnelling. D-Wave devices have been able successfully to simulate condensed matter systems, sometimes showing advantages over classical counterparts [39][40][41].
The main objective of this work is to demonstrate how a field theory problem can be successfully encoded on a quantum annealing device, and to do this we will focus on the classic problem of obtaining tunnelling rates for a system stuck in a metastable minimum (a.k.a. false vacuum).  The thin-wall solution computed using the hybrid quantumclassical techniques as discussed later is overlaid on the right panel.

II. SET-UP OF A SIMPLE PROBLEM
A useful potential to focus on is the following quartic one: The potential is shown in Fig.1. On the left we show the "thick-wall" regime where is large. This limit is when the barrier is close to disappearing (or has disappeared altogether) and the walls become comparable in size to the bubble itself. For numerics we choose a = λ = 1 and = 0.3. The opposite "thin-wall" regime (for which we choose = 0.01) is the limit in which is small and is approximately the difference in vacuum energy density between the false and true minima.
We are interested in the situation where the system starts in the false vacuum, and our objective is to study the rate per unit volume of tunnelling out of it. The analytic calculation of this rate is a classic problem, but it is worth briefly recapping it in order to recast the result in a form that can easily be compared with the results from a quantum simulation. It proceeds as follows.
First let us remove the extraneous constant term by working with U (φ) = V (φ) − V (φ + ), which has U (φ + ) = 0. Using the well-known technique of [42][43][44][45], the bubble profile is given by finding a "bounce solution" to the following differential equation: where in four dimensions, c takes the value 2 or 3 for a finite temperature O(3) symmetric bubble, or a purely quantum tunnelling O(4) symmetric instanton, respectively. The required "bounce" is subject to the boundary condition that dφ/dρ = 0 as ρ → 0, ∞, which determines the starting value φ(0), which is the field-value at the centre of the radially symmetric bubble or instanton (also called the escape-point). The resulting φ(ρ) profile for our particular choice of parameters is shown in Fig. 2.
Once such a solution is determined, the tunnelling rate per unit volume can be estimated from its classical action: respectively. The quantum determinant prefactors A 4 , A 3 are notoriously difficult to calculate, but for our purposes it will be sufficient to focus on the influence of the classical action.
The expressions for the action can be expressed in simple analytic terms in the two limits. In the thick wall limit the bounce action can be accurately approximated by expanding around the value = 0 , above which the barrier disappears (i.e. when the discriminant vanishes), which gives a cubic potential about the false vacuum. This critical value corresponds Then following the rescaling procedure of [45], the tunnelling actions for the O(4) and O(3) symmetric solutions can be written in terms of standard actions: The thin-wall regime is somewhat easier to study numerically, and semi-analytically the actions can be expressed in terms of the action S 1 for the one-dimensional c = 0 problem 1 : These limiting regimes give simple power-law behaviour for the tunnelling actions, against which the scaling of the (logarithm of) tunnelling rates could be tested, providing a useful laboratory for directly studying quantum annealing results.
As we stated in the introduction, the purpose if this study is not to recover these classical instanton solutions for the tunnelling per se, as they are well-known, but rather to demonstrate that the corresponding field-theory configuration can be suitably encoded into a quantum annealer. Once we have established this as a working principle, one could even envisage testing for the above behaviour directly. Therefore we will in what follows focus on using a quantum annealer to recover the simple c = 0 solution required for the thin-wall regime, as a proof of principle. We will therefore set ourselves the task of minimising the corresponding action integral, which should yield a solution of the form shown in Fig.2b.

III. ENCODING THE FIELD THEORY
Let us start with the central problem, which is how to formulate a continuous scalar field theory on quantum annealers. A quantum annealer is based on the adiabatic theorem of quantum mechanics, which implies that a physical system will remain in the ground state if a given perturbation acts slowly enough, and if there is a gap between the ground state and the rest of the system's energy spectrum [24]. For the annealer to provide a solution to a mathematical problem, e.g. the calculation of φ(ρ) for Eq. 7, we have to find a mapping such that the expectation value of its Hamiltonian can be identified with its solution, i.e. that it allows in this example to identify The form of the Hamiltonian available to a quantum annealer is that of a general Ising model, in addition to a time-dependent transverse field: where the Pauli Z operator, with the subscript indicating which spin it acts upon, and σ X is its friend pointing in the Xdirection. The gradual decrease of ∆(t) → 0 from a large value should drive the system into the ground state of the time-independent part of the Hamiltonian, and this is where we will put the field theory: It is worth noting that the couplings J ij and h i could also be adiabatically adjusted in the annealing process, and this could ultimately be used to adjust the potential U (φ) of a system in the quantum annealer so as to observe tunnelling, assuming it can be encoded. We will further split the Hamiltonian into three generic pieces, as Here, H (QFT) is the Hamiltonian corresponding to the minimisation of the action in Eq. 7 and H (BC) is a Hamiltonian that we add to enforce the boundary conditions 2 . However our first task is to encode continuous field values over a continuous domain, with only the discrete Ising model to hand: this is what H (chain) is for. We begin by splitting the radius variable ρ into M 1 discrete values and the field value at the 'th position into N 1 discrete values: where in the present context one might for example take a fiducial value φ 0 ≈ −a and ξ = 2a/N , with M ν = ∆ρ. Thus our Ising interaction J ij is an We must now separate those spins in the annealer that correspond to fields at different values of , effectively splitting J ij and h i into N × N sub-blocks. To do this we will utilise the Ising-chain domain wall representation introduced in [47]. That is for every position we add to the Hamiltonian (12) As shown in [47], taking Λ to be much larger than every other energy scale in the overall Hamiltonian, these terms will constrain the system to remain in the ground subspace of the Hamiltonian, where exactly one spin position, α say, is frustrated for each . These states are of the form where in the above the discretised field value is represented by the position α of the frustrated domain wall. Conversely the field value at the 'th position can be found by making the measurement which only receives a contribution from frustrated spin position with j = α . For later, it is useful to note that this is equivalent to In terms of J ij and h i , adding the full set of Ising-chain Hamiltonians given by Eq.(12) corresponds to and an h that is independent of , h (chain) This separates the system of spins into blocks of size N , each of which represents a field value.
Moving on to H (QF T ) , the potential is somewhat easier to deal with than the kinetic terms, because it can be encoded entirely in h i . This is only to be expected because the φ are independent of each other in the potential which gives entirely localised contributions to the Hamiltonian. The value of U (φ(ρ )) at each point follows directly from Eq. (14): This yields an additional contribution to the h which is also independent of : that is for all we have It can also be convenient to write this in terms of U derivatives as which correctly gives φ(ρ ) of Eq. (15) in the event that we take U (φ) = φ (because we know that σ Z N = 1). Note that in a system with arbitrary c = 0, we would need to evaluate h (U ) ≡ dρρ c U , so that h N +i would acquire a prefactor of ( ν) c .
Up to this point the M -factors have been inert and there has been no coupling between the fields at different positions in ρ . At this stage the system would simply relax to M decoupled values of φ(ρ ) that minimise U in either one of its two vacua. This changes once we include the derivatives in the kinetic terms, which contribute to the bilinear interactions, J. These terms are discretised in ρ as where ν = ∆ρ/M scales so as to keep ∆ρ constant. Inserting the discrete representation of the field values as well using Eq.(15), we find Hence the bilinear terms receive the additional contribution: Now it is the N × N indices that are inert, because every i couples to every j.
Note that the diagonal parts of Eq.23 could be embedded in the h i terms, using the fact that for valid single domain wall states we have σ Z N +i σ Z N +j = σ Z N +j − σ Z N +i + 1 for j > i. As bilinear terms may be hard to engineer on real devices, this may be desirable, but for the present study it is more convenient to keep the kinetic terms entirely.
Finally we must add terms to enforce a boundary condition. In the c = 0 case it is sufficient to fix the endpoints of the solution in the two minima (so that, at the risk of confusion, the instanton solution itself approximates a physical domain wall). This can be done by adding a term H (BC) = Λ 2 (φ(0) + a) 2 + Λ 2 (φ(ρ M ) − a) 2 with Λ being some other large parameter. This is simply an extra contribution to h which follows directly from Eq. (20), of the form Together with Eqs. (16,17,20,23), this completes the encoding of the field theory problem of Eq. (7).

IV. IMPLEMENTATION
In Sec. III we have devised a method which encodes the problem of finding a solution to a quantum field theoretical problem, i.e. of finding a solution to Eq. 7, into finding the ground state of the Hamiltonian of an Ising model. The latter can then be given an interpretation as the solution to Eq. 7 through Eq. 13, for each ρ l with l ∈ [1, ..., M ]. To show that our approach is valid and converges to the correct solution φ(ρ), we now implement the method onto various annealing samplers, as provided by D-Wave [48].
The quantum states are characterised by N M -tuples of the form |11...100...0 and the Hilbert space of the Ising model is therefore 2 N M dimensional. Sampling such a large vector space classically, with an exact sampler, while calculating the expectation value H for each state quickly becomes a computationally prohibitive task for N M 20. Conversely, a discretisation with N M 20 cannot give a reasonable approximation for the derivatives of Eq. 21.
Yet, not unlike protein-folding, in which a unique ground state is selected from an estimated number of 3 300 so-called conformations within microseconds (known as Levinthal's paradox [49]), a quantum annealer can in principle find a ground state of a Hamiltonian acting on a highly complex Hilbert space on a similar time scale, assuming there is a gap between the ground state and the other states of the system.
While the next generation of annealing processors will have approximately 5,000 qubits, they will have limited connectivity [50]. Therefore in order to accommodate the more general Ising model required for our encoding, we resorted to a hybrid asynchronous decomposition sam-pler (the Kerberos solver [51]), which can solve problems of arbitrary structure and size. To find the ground state efficiently, it applies in parallel classical tabu search algorithms, simulated annealing and D-Wave subproblem sampling on variables that have high-energy impact. Using this method we calculate the solution φ(ρ) to Eq. 7 for N = M = 50 in Fig. 2b.

V. CONCLUSION
Consequently, near-term applications of quantum devices can significantly enhance our ability to perform highly complex quantum field theoretical calculations. Focussing on the problem of calculating the classical instanton solution in Sec. II, we have proposed a method to formulate such problems as an Ising model in Sec. III, which we have solved on a D-Wave quantum annealer. Our method is highly flexible and can straightforwardly be generalised to encode any multi-dimensional differointegral equation as an Ising model.