Quantum Monte Carlo simulation of a particular class of non-stoquastic Hamiltonians in quantum annealing

Quantum annealing is a generic solver of the optimization problem that uses fictitious quantum fluctuation. Its simulation in classical computing is often performed using the quantum Monte Carlo simulation via the Suzuki–Trotter decomposition. However, the negative sign problem sometimes emerges in the simulation of quantum annealing with an elaborate driver Hamiltonian, since it belongs to a class of non-stoquastic Hamiltonians. In the present study, we propose an alternative way to avoid the negative sign problem involved in a particular class of the non-stoquastic Hamiltonians. To check the validity of the method, we demonstrate our method by applying it to a simple problem that includes the anti-ferromagnetic XX interaction, which is a typical instance of the non-stoquastic Hamiltonians.

in the non-diagonal elements 29,30 .However, this task is highly nontrivial and there are no generic solutions.The negative sign problem belongs to the class of the NP-hard problem 31 .
In the present study, we show a method to simulate a particular class of the non-stoquastic Hamiltonians including the anti-ferromagnetic XX interaction in the quantum Monte Carlo simulation.The numerical scheme will be demonstrated below.The present study paves the way to simulate quantum annealing with an elaborate driver Hamiltonian.The future development of hardware devices to perform quantum annealing aims at implementing the non-stoquastic Hamiltonian such as the anti-ferromagnetic XX interaction beyond the classically simulatable world.Our contribution in the present study is to establish a test bed to simulate a particular class of non-stoquastic Hamiltonians in order to verify the performance of future hardware.Moreover, our study will enable invention of many types of algorithms inspired by quantum annealing in the classical computer but with the non-stoquastic Hamiltonian.
In standard QA, we employ a system with the transverse field as where Γ represents the strength of the transverse field, and σ i x is the x component of the Pauli matrix.The classical Ising model to be solved is H 0 (σ), where σ σ σ σ =  ( , , , ) . The standard QA procedure decreases the strength of the quantum driver Hamiltonian to find the ground state of the cost function to be optimized.The adiabatic theorem ensures that we obtain the ground state after sufficient slow sweep of the quantum driver Hamiltonian.The computational time of QA can be evaluated by the energy gap between the ground state and the first excited state during the decrease in the transverse fields 8,9 .As has been extensively studied, the quantum phase transition hampers the efficient computation of QA.In particular, the first-order phase transition characterized by the exponential closure of the energy gap depending on the number of spins signals the long-time computation for optimization via QA.Seki and Nishimori successfully avoided the first order phase transition by utilizing additional fluctuations, anti-ferromagnetic XX interactions.

∑ ∑
However, their quantum fluctuation yields the negative sign problem in a naive application of the Suzuki-Trotter decomposition in the quantum Monte Carlo method 25 .Seki and Nishimori avoid the negative sign problem without rotating the basis to represent the Hamiltonian.They utilize the mean-field analysis to compute the free energy and the thermodynamic quantities such as the magnetization, energy, etc.A recent study on the performance in the optimization problem by use of the anti-ferromagnetic XX interaction is carried out by the numerical diagonalization for a limited size of the system 32 .
We then elucidate the essential part of the mean-field study to avoid the negative sign problem and apply its scheme to simulate the non-stoquastic Hamiltonian with a more generic form.The scheme involves the following form of the non-stoquastic Hamiltonian The classical Hamiltonian here is not restricted to the all-to-all connections as in the previous studies carried out by Seki and Nishimori 22,23 .The proposed method in the present study is applicable to all kinds of the classical Hamiltonian.For instance, the spin-glass Hamiltonian with all-to-all connections as and the finite-dimensional spin-glass Hamiltonian with limited connections as where, ij stands for the summation over the subset of the interactions, which are located at each bond on the specific lattice.The dimensionality is not limited because our method is based on the quantum Monte Carlo simulation rather than the mean-field analyses as in the previous studies 22,23 .In addition, the various forms of the Hamiltonian representing the optimization problem with many-body interactions can be within the scope of application of our method.On the other hand, the quantum fluctuation term is limited to the specific form, in which its argument is given by the summation over all-components of the x-component of the Pauli operators.This is a kind of extension of the quantum fluctuation used in QA.The case of the standard QA is described by f(m x ) = Γ m x .The case of the QA with the anti-ferromagnetic XX interaction and transverse field is included in this form through 2 .In our scheme, we utilize the quantum Monte Carlo simulation while avoiding the negative sign problem for a particular class of the non-stoquastic Hamiltonian through the simple calculations by replacement of the quantum fluctuation with the adaptive-changing transverse field.As discussed above, there is no generic solution for the negative sign problem.To avoid the obstacle in the numerical computation, one needs to find an adequate basis in which the negative sign problem disappears.Contrasting with the ordinary approach overcoming the difficulty involved in the negative sign problem, our method avoids the obstacle via the simple transformation and calculations rather than finding the proper basis to represent the Hamiltonian.We propose two approaches to simulate a particular class of the non-stoquastic Hamiltonian.The first one is called the adaptive quantum Monte Carlo simulation and the other is the data-analysis approach.These methods are detailed below.
The following section demonstrates the validity of our proposed method for a simple model: the (ferromagnetic) infinite-range model with the transverse field and anti-ferromagnetic XX interaction defined as 2 .The infinite model has the second-order phase transition at the critical point Γ c = 1 without any (longitudinal) magnetic field and anti-ferromagnetic XX interactions.The model we dealt with here is solvable by mean-field analysis and thus adequate for validation of our method.We emphasize, however, that the application range of our method is beyond that of the mean-field analysis.When we perform the naive application of the quantum Monte Carlo simulation through the standard approach via the Suzuki-Trotter decomposition, the numerical simulations suffer from the negative sign problem to obtain precise estimations of the physical quantities for the case away from f(m x ) = Γ m x .However, we do not encounter the negative sign problem by the use of the replacement of the non-stoquastic terms in Eq. ( 3) via the saddle-point method and integral representation of the delta function as detailed below.

Results
In this section, we propose our first method: the adaptive quantum Monte Carlo simulation.We set the spin systems with N = 4, 8, 16 and 32.In order to perform the quantum Monte Carlo simulation, we consider the replicated system at extremely low temperatures via the Suzuki-Trotter decomposition.We set the Trotter number as τ = 128 and the inverse temperature (temperature) as β = 50 (T = 0.02).To check the precision of the computation, we estimate various thermodynamic quantities: the internal energy, (longitudinal) magnetization, and the transverse magnetization including the order parameters of the system.In Fig. 1, we show estimations for three thermodynamic quantities by computing 80000 MCS time average after 20000 MCS equilibration.For comparison, we show the exact solution carried out using the mean-field analysis, because this model is tractable.The case with γ = 1 is for the anti-ferromagnetic XX interaction, while that with γ = 0 is for absence of the anti-ferromagnetic XX interaction.The results obtained by the adaptive quantum Monte Carlo simulation are close to the exact result with γ = 1 rather than that with γ = 0.This fact is evidence that the adaptive quantum Monte Carlo simulation can estimate correctly the thermodynamic quantities even in the non-stoquastic Hamiltonian.Aside from the finite-size effect, the obtained results asymptotically coincide with the exact solutions.We remark that the convergence of the estimations is relatively slow compared to the ordinary quantum Monte Carlo simulation because we have to estimate the expectation of the thermodynamic quantities as the effective transverse field f ′ (m x ) depending on the tentative value of the transverse magnetization.Although substantial slowing down of the equilibration is not observed in the simple model, the computational time to correctly estimate the thermodynamic quantities might become longer depending on the model that is related to its complexity.
Another method to estimate the thermodynamic quantities of the non-stoquastic Hamiltonian is the data analysis approach.Similar to the previous check, we test our method in the simple infinite-range model defined by Eq. ( 6).We again set the spin systems with N = 4, 8, 16 and 32, the Trotter number as τ = 128, and the inverse temperature (temperature) as β = 50 (T = 0.02).We here perform the "standard" quantum Monte Carlo simulation for the model only with the transverse field while changing the value of the transverse field from Γ = 0 to Γ = 4 by Δ Γ = 0.05.For validation of our method, we estimate the internal energy, (longitudinal) magnetization, and the transverse magnetization as in Fig. 2. The computational time is short because we perform the "standard" quantum Monte Carlo simulation.We take 80000 MCS time average after 10000 MCS equilibration.
Similar to the results of the adaptive quantum Monte Carlo simulation, we confirm the validity of our method within the finite-size effects of the simulation.The obtained results become closer to the exact solution with γ = 1, which is the non-stoquastic Hamiltonian, rather than with γ = 0.In the results obtained using the data-analysis approach, we find several adjacent points in the different value of γ.This means that the same results are obtained in several cases after estimation of the thermodynamic quantities for the non-stoquastic Hamiltonian.We can remove lack of numerical precision by increasing the step size of the transverse magnetic field in the "standard" quantum Monte Carlo simulation.We confirm that both the proposed methods based on the quantum Monte Carlo simulation can estimate the thermodynamic quantities.

Discussion
We propose two methods to simulate the non-stoquastic Hamiltonian on the basis of the quantum Monte Carlo simulation, which is often affected by the negative sign problem due to the computational basis.Although we usually consider the rotation of the basis to remove the negative sign problem, we avoid the negative sign problem by introducing the auxiliary variable and its conjugate one m x and ∼ m x .The model which we deal with is non-stoquastic but the negative sign problem is not insurmountable.Thus, we avoid the negative sign problem though the model is in the definition of the non-stoquastic Hamiltonian.
We list the possible questions and answers on our proposed method below: • Our stand point.
In the present study, we simply check the validity of our proposed methods by comparing the numerical results and the exact solution via the mean-field analysis.Since, the proposed method considers a toy model, we do not obtain any nontrivial result.However, this means readily that our work is trivial.This is the first step to establish a systematic approach for the simulation of the non-stoquastic Hamiltonian.In the numerical simulation of the non-stoquastic Hamiltonian with a large number of components, efficient methods to compute the thermodynamic quantities are the exact diagonalization of the Hamiltonian and the renormalization group analysis.Both the methods often suffer from size limitation.We propose an alternative choice to approach the non-stoquastic Hamiltonian.• Limitation of our method.
Our method is not a generic solution to the negative sign problem.The applicable scope of our method is limited to the case with the quantum fluctuation determined by the collection of the x-component of the Pauli operators.The form of the function f(m) is very flexible.The condition of the function is that P(m x ) → 0 as m x → ± ∞ .Furthermore, when the spin operators characterizing the quantum fluctuation can be reduced to a few quantities, our proposed method can be applied.For instance, our method can be generalized to the case with , where Λ p is a subset of the indices (locations) of the spin operators.This means that the inhomogeneous XX interactions can be dealt with.The combination of the inhomogeneous XX interaction is directly related to the computational cost.We would expect that the boundary between the classical and the quantum computational capacity lies in the inhomogeneity of the XX interactions.This problem is presently not solvable.
We employ the saddle-point method for determining the value of ∼ m x as detailed below.Thus, in order to enhance the numerical precision, we take the system size to be relatively large.However, we may simulate the value of m x using the Langevin stochastic process.Then, we need to choose the direction of m x for stabilization of its behavior depending on the form of f(m x ).
• Capability of acceleration.Most of the ordinary simulations by the Langevin stochastic process and the Markov-chain Monte Carlo simulation require detailed balance condition.Various studies have ensured that the violation of the detailed balance condition exhibits acceleration of the convergence to the pre-determined steady state [33][34][35][36] .We may also utilize the exchange Monte Carlo simulation to accelerate the convergence to the steady state 3 .In the present study, we employ this method to obtain our results.
• Avoiding the first-order phase transition.
Seki and Nishimori demonstrated the possibility to remove the first-order phase transition by using the anti-ferromagnetic XX interaction.The effect of a kind of the quantum fluctuation can be recast as the modified transverse field depending on the tentative value of the transverse magnetization.This different approach provides a deep understanding of the change of the phase transition by using the elaborate quantum fluctuations.As shown above, the resulting thermodynamic quantity can be characterized by the cross points between the transverse magnetization under the effective transverse field and the inverse function determined by the form of the quantum fluctuation f(m x ).The first-order phase transition is involved in the discontinuous jump of the thermodynamic quantity.The discontinuous jump emerges from co-existing phases at the same parameter.In other words, the thermodynamic quantity is described by a multivalued function.The thermodynamically stable point is selected by comparing the values of the free energy at the multiple solutions at the same parameter.When we employ a nontrivial function as f(m x ) to avoid multiple cross points, we convert the first-order phase transition to the second-order (continuous) phase transition.The relationship between the computational complexity of QA with non-stoquastic Hamiltonian and the limitation of our approach would be a very exciting realm and will be reported elsewhere as a future study.
A recent study on the non-stoquastic Hamiltonian is carried out by numerical diagonalization for a kind of the spin glass model with long-range interaction with the anti-ferromagnetic XX interactions.However, if one implements our method, one can investigate the property of the model.In this sense, our methods open an alternate way to approach nontrivial aspects of the non-stoquastic Hamiltonian by performing the quantum Monte Carlo simulation for large-component systems as compared to the sizes considered in the numerical diagonalization.We will report such nontrivial results by use of our method in the near future.

Methods
In the standard approach of the quantum Monte Carlo simulation, we construct the replicated partition function via the Suzuki-Trotter decomposition to reduce the transverse field to the interactions between the replicated systems.For a class of the non-stoquastic Hamiltonian written as in Eq. ( 3), we replace the quantum fluctuation with the adaptively-changing transverse field via the integral representation of the delta function and saddle-point method.Thus, the method we detail below can be more precise when the number of degrees of freedom is large.
The partition function of the present system is written as follows: We employ the Suzuki-Trotter decomposition to divide the exponentiated Hamiltonian into the diagonal and non-diagonal parts as Here, we utilize an identity through the delta function as The partition function is written as This is a kind of micro-canonical form of the present system.We further employ the Fourier integral form of the delta function as The resulting partition function is the same as the Ising model with the transverse field.In the thermodynamic limit, we may take the saddle point in the integral.The saddle point is evaluated by , where f ′ (m) is the derivative of the function f(m).In other words, the transverse field is determined by the auxiliary variable m xt , which corresponds to the transverse magnetization.The transverse field can be rewritten as the interaction between the different Trotter slices t and t + 1, according to the prescription of the quantum Monte Carlo simulation as below.We define the joint probability distribution for ({σ t }, m x ) as x t Simultaneously, we define the conditional probability as ( ) ( 14) where Z(m x ) is the normalization constant and the marginal distribution as In the standard procedure of the Markov-chain Monte Carlo simulation, we generate a series of realizations of the Ising variable following the probability distribution conditioned on some fixed parameters such as the temperature and transverse field.However, in our case, the effective transverse field can change the tentative value of the auxiliary variable as x .The argument of the function f′ (m x ) also fluctuates following the marginal distribution.We may then simulate the stochastic process governed by the Langevin equation to generate the auxiliary variable m x following the marginal distribution as where dW is the Weiner process.Here, we define where the bracket denotes the expectation with respect to the weight of the conditional probability P(σ|m x ).In the case of the number of the spin variables N → ∞ , the effect of the Weiner process is negligible in the Langevin equation.In other words, the fluctuation around the saddle point σ = m x i x can be ignored.We can then generate two strategies to employ our approach.The first strategy is to directly simulate the non-stoquastic Hamiltonian by using the adaptive transverse field: 1. Perform the quantum Monte Carlo simulation following Eq.( 14) and estimate the value of the expectation of the transverse magnetization.To estimate the expectation, we generally require long-time equilibration.In the case of actual application, we compute the approximate value of the expectation by the empirical average over the whole simulation after the first relaxation.2. Change the transverse field ∼ m x following the saddle-point solution as x x 3. Repeat until the physical quantities converge.
The second strategy is performed by the data analysis for the results of the "standard" quantum Monte Carlo simulation of the Ising model, only with the transverse field.We generate the results obtained by the wide range of the values of the transverse field ∼ m x .Then, we plot the transverse magnetization vs. the transverse field, and we read out its cross point to the curve determined by is the inverse function of f(m).Then, the cross point identifies the realization of the transverse magnetization of the non-stoquastic Hamiltonian.In the case of the transverse field and the anti-ferromagnetic XX interaction, the inverse function is simply given by the linear function f ′ −1 (m x ) = (Γ − m x )/γ.In this example, for the given Γ and γ, we find the cross point and re-plot the realized transverse magnetization and other corresponding thermodynamic quantities vs. the values of Γ and γ, as shown in Fig. 3.

Figure 1 .
Figure 1.Results of the adaptive quantum Monte Carlo simulation.The dashed curves denote the exact solution for the case without antiferromagnetic interactions and solid ones represent those with antiferromagnetic interactions.The gradating red dots represent the results of the adaptive quantum Monte Carlo simulation (denser red point stands for larger size of system).

Figure 2 .
Figure 2. Results of data analysis.The same symbols are used as those in Fig. 1.

Figure 3 .
Figure 3. Data analysis for the non-stoquastic Hamiltonian.The red dots are the results of the Ising model, only with the transverse field and the blue dots represent the inverse function as = ∼ − m f m ( ) x x 1 .