Performance enhancement of quantum annealing under the Lechner-Hauke-Zoller scheme by non-linear driving of the constraint term

We analyze the performance of quantum annealing as formulated by Lechner, Hauke, and Zoller (LHZ), by which a Hamiltonian with all-to-all two-body interactions is reduced to a corresponding Hamiltonian with local many-body interactions. Mean-field analyses show that problematic first-order quantum phase transitions that exist in the original LHZ formulation can be avoided, and thus an exponential speedup is achieved, if we drive the coefficient of the manybody term, which represents the constraint, non-linearly as a function of time. This result applies not only to a simple ferromagnetic model but also to the spin glass problem if a parameter in the spin glass model is chosen appropriately. Numerical studies of small-size systems are consistent with the mean-field predictions.


Introduction
Quantum Annealing (QA) is a metaheuristic for solving combinatorial optimization problems using quantum fluctuations, [1][2][3][4][5][6][7] and is closely related with adiabatic quantum computation. 8,9) The goal of QA is to obtain the ground state of a classical Ising model, to which a combinatorial optimization problem can be reduced. 10) It is usually the case that the resulting Ising Hamiltonian has long-range interactions including all-to-all interactions, but the current annealing device, the D-Wave quantum annealer, does not directly implement long-range interactions. One therefore has to employ the procedure of embedding 11,12) to represent long-range interactions in terms of a combination of short-range interactions. This leads to an overhead in the qubit count and also tends to cause errors induced by imperfect realization of embedding in the real device.
Lechner, Hauke, and Zoller (LHZ) 13) proposed an ingenious scheme to partly mitigate the above problems by mapping allto-all interactions to single-qubit terms in the Hamiltonian supplemented by local four-body interactions introduced to guarantee the equivalence of two formulations. Although the issue of overhead in the qubit count still exists, the reduction of all-to-all interactions to a local representation is certainly advantageous in the device implementation as well as from the viewpoint of mitigation of errors caused by imperfections in embedding. Several proposals have been made to realize the LHZ scheme. [14][15][16][17][18] There have also been attempts to analyze the LHZ scheme theoretically. Leib,Zoller,and Lechner 14) proposed a method to reduce four-body interactions in the LHZ Hamiltonian to two-body interactions by introducing auxiliary qubits. They also showed that a proper control of the coefficient of the constraint terms is likely to improve the performance. Hartmann and Lechner 19) used non-stoquastic counter-diabatic drivers for better performance, and the same authors recently introduced a mean-field like method to analyze the effect of inhomogeneity in the transverse field for increased success probabilities. 20) * y-susa@bx.jp.nec.com We have been inspired by these developments and have analyzed the LHZ scheme within the framework of mean-field theory. The result shows that a non-linearity in the coefficient of the constraint term as a function of time leads to avoidance of first-order phase transitions that exist in the original linear time dependence of the coefficient. This implies an exponential speedup from the view point of adiabatic quantum computation because a first-order phase transition usually accompanies an exponentially-closing energy gap as a function of the system size, meaning an exponential computation time according to the adiabatic theorem of quantum mechanics. 21,22) Data from numerical computations for small-size systems are consistent with this theoretical prediction. This paper is organized as follows. After an introduction to the LZH model and the p-spin model in Sec. 2, we analyze the problem analytically and numerically in Sec. 3. Conclusion is described in Sec. 4.

The LHZ model and the p-spin model
The conventional QA has the Hamiltonian whereĤ P is the Ising Hamiltonian representing a combinatorial optimization problem, andV is the transverse field to induce quantum fluctuations, withσ z(x) i denoting the z(x) component of the Paul operator at site (qubit) i. The parameter s = t/T is the normalized time running from 0 to 1 as the time t proceeds from 0 to T, and thus T is the total computation time (annealing time).
One starts at s = 0 from the trivial ground state ofV and increases s with the expectation that the ground state of the Ising HamiltonianĤ P is reached at s = 1 (t = T). According to the adiabatic theorem of quantum mechanics, the computa- Large green circles denote physical qubits and the four qubits around each small red circle (plaquette) consist a four-body interaction. The number in a green or red circle is the index of a physical qubit k or a plaquette l, respectively. The state of auxiliary physical qubit at the bottom row (yellow) is fixed to |↑ . tion time T necessary for the system to stay close enough to the instantaneous ground state is proportional to the inverse polynomial of the minimum energy gap where E 0 (s) and E 1 (s) are the instantaneous ground and firstexcited state energy, respectively. If this energy gap closes exponentially as a function of the system size N l (the number of logical qubits) as is usually the case at a first-order quantum phase transition, the computation time T is grows exponentially e aN l (a > 0), which means that the problem is hard to solve by QA. It is therefore highly desirable to avoid or remove first-order phase transitions.
The LHZ scheme 13) reduces the all-to-all interactions implied in the Ising Hamiltonian eq. (2) to single-body terms supplemented by four-body constraint terms to enforce equivalence to the original problem, through the correspondence The number of physical qubits N in the LHZ Hamiltonian eq. (4) is the number of all-to-all interactions in the original model, and the number of constraints in the second term on the righthand side of eq. (4) is The four-body term in eq. (4) consists of four neighboring qubits (three at the bottom boundary) as depicted in Fig. 1, and is therefore local and is possibly amenable to direct experimental implementation. A recent contribution by Hartmann and Lechner 20) showed that an infinite-range (mean-field) version of the four-body term serves as a good approximation to the original nearestneighbor (short-range) interactions, which greatly facilitates analytical studies. We therefore follow their idea and introduce the following Hamiltonian, where p is four originally but we keep it arbitrary to allow for comparison with more general cases. This is the problem Hamiltonian we study in the present paper.

Mean-field and numerical analyses
It turns out to be convenient to introduce an additional parameter τ to control the time dependence of the constraint term. The Hamiltonian of the original problem and its meanfield version are then written aŝ The total Hamiltonian iŝ whereĤ p is either eq. (9) or eq. (10). The parameters s(t) and τ(t) are no longer linear in general as a function of t and change from s = τ = 0 at t = 0 to s = τ = 1 at t = T.

Mean-field analysis
It is straightforward to apply the standard procedure to derive the free energy per qubit as a function of the ferromagnetic order parameter m. 20,[23][24][25][26][27][28] We therefore just write the result for the free energy and its minimization condition, i.e. the self-consistent equation, where β is the inverse temperature and the brackets [· · · ] i stand for the average over the values of J i , 1 N i · · · . Let us first focus on the simplest case of zero temperature β → ∞ and a uniform interactions J i = J. Then eqs.(12a) and (12b) reduce to Numerical solutions to these equations reveal the phase diagram on the s-τ plane for p = 4 and J = 0.5 as shown in Fig. 2(a), where the thick blue line denotes a line of first-order phase transitions terminating at a critical point marked in orange. The precise location of this critical point can be derived following the standard prescription that the derivatives up to third order should vanish at a critical point, 29) The conventional protocol of quantum annealing corresponds to the straight line τ = s in the phase diagram, which crosses the line of first-order phase transitions. If we instead choose a trajectory τ = s r with r > 1, 56, the annealing process does not encounter a phase transition. The critical point is touched when r = 1.56. Correspondingly, the minimum energy gap between the ground state and the first excited state closes exponentially as a function of the system size for τ = s whereas it is polynomial for r = 1.56 as depicted in Fig. 2(b). When r > 1.56, the gap is expected to reach a constant in the large-N limit because there is no phase transition, but the numerical data show a slow decay. This would probably due to the proximity of the curves τ = s 2 and τ = s 3 to the critical point and the asymptotic region for N 1 is not yet reached. It is anyway the case that an exponential speedup of the computation time can be achieved by the choice of r ≥ 1.56 in comparison with the conventional annealing with r = 1 because an exponential gap closing is avoided.
To confirm that the above equilibrium statistical-mechanical analysis for the large-size, equilibrium limit is consistent with small-size, finite-time dynamics, we numerically solved the time-dependent Schrödinger equation for N = 10 and measured the probability to find ground state at t = T for τ = s and τ = s 3 , The result is given in Fig. 2(c). Although the difference in the energy gap is very small between τ = s and τ = s 3 for N = 10 as seen in Fig. 2(b), the final success probability for τ = s 3 is consistently larger by a small margin than that for τ = s. Similar results are obtained for different parameter values. Examples are shown in Fig. 3(a) for p = 3, 4, and 5 with J = 0.5, Fig. 3 A more complex case of random interactions with the distribution function is interesting because this is for random and frustrated all-toall interactions in the original problem, corresponding to the Sherrington-Kirkpatrick model of spin glasses. 30) The groundstate free energy for this problem can be derived from eq. (12a) with the result The phase diagram is drawn in Fig. 4(a)  a completely random spin-glass model, which is known as a very difficult problem to solve. 31) Nevertheless, even when the line of first-order transitions traverses the phase diagram as in the case of = 0.8, the jump in magnetization across a firstorder transition can be tuned much smaller than the naive case of τ = s by an ingenious choice of the trajectory in the phase diagram as seen in Fig. 4(b), which shows the magnetization jump along the first-order transition line. This implies that quantum tunneling probability, which strongly depends on the width of an energy barrier represented by the magnetization jump, can be tuned to be larger by an appropriate choice of a trajectory connecting s = τ = 0 and s = τ = 1. The extreme case of = 0.5 has no such properties.

Numerical analysis of finite-size systems
To verify that the results of mean-field analyses are applicable to the finite-size short-range problem of eq. (9), we evaluated the minimum energy gap as a function of the parameter r in τ = s r for N l = 4 and 5 with J k = 0.5 As seen in Fig. 5(a), the minimum energy gap is larger for r > 1 than for r = 1. Correspondingly, the success probability as a function of the computation time T is slightly larger for τ = s 4 than for τ = s as observed in Fig. 5(b). Although the difference in success probabilities is small due probably partly to the smallness of the difference in the energy gaps (about 0.8 for r = 1 and 0.9 for r = 4 according to Fig. 5(a)), the overall tendency is in favor of the τ = s 4 case and is consistent with Fig. 2(c).
We have also studied the distribution of energy gap for uniformly random instances of interaction with |J k | ≤ 0.5. The result is summarized in Fig. 6 as the difference of the minimum energy gap between the cases of τ = s 2 and τ = s as a function of the minimum energy gap for τ = s. It is seen that the gap is larger for τ = s 2 than for τ = s for a majority of samples.

Conclusion
We have studied the scheme of Lechner, Hauke and Zoller 13) to express long-range (all-to-all), two-body interactions by short-range, many-body interactions. Mean-field and numerical methods were used to show that non-linear driving of the four-body constraint term as a function of time is advantageous for improved performance. In particular, increasing the amplitude of the constraint term more slowly than for the intrinsic problem term can lead to an exponential speedup according to the mean-field prediction for the phase diagram although it would be difficult in practice to observe such a drastic effect because of non-ideal environmental effects as well as due to the limited applicability of mean-field theory. It is nevertheless encouraging that numerical results for small-size short-range cases also show qualitatively similar behavior though the gain is not very large.
We may learn a generic lesson that constraints are better introduced later in the process of quantum annealing compared to the main problem Hamiltonian. In other words, one may first search for good solutions without constraints and then gradually select among candidate solutions those that satisfy the constraint. It is an interesting future topic to test this idea for various problems.  6. (a) Each dot represents the difference of the minimum energy gaps ∆ τ=s 2 − ∆ τ=s for 1,000 instances with uniform random |J k | ≤ 0.5 as a function ∆ τ=s . The system size is N l = 4.