Classical Simulation and Theory of Quantum Annealing in a Thermal Environment

We study quantum annealing in the quantum Ising model coupled to a thermal environment. When the speed of quantum annealing is sufficiently slow, the system evolves following the instantaneous thermal equilibrium. This quasistatic and isothermal evolution, however, fails near the end of annealing because the relaxation time grows infinitely, therefore yielding excess energy from the thermal equilibrium. We develop a phenomenological theory based on this picture and derive a scaling relation of the excess energy after annealing. The theoretical results are numerically confirmed using a novel non-Markovian method that we recently proposed based on a path-integral representation of the reduced density matrix and the infinite time evolving block decimation. In addition, we discuss crossovers from weak to strong coupling as well as from the adiabatic to quasistatic regime, and propose experiments on the D-Wave quantum annealer.

Figure 1

(a) Energy expectation values $⟨H_S(t/t_a){⟩}_t$ and $⟨H_S(s){⟩}_{s,T_B}^{\mathrm{eq}}$ per spin of the time-dependent state ${\rho }_S(t)$ and the instantaneous thermal equilibrium state ${\rho }_S^{\mathrm{eq}}(s,T_B)$, respectively, as functions of the rescaled time $s$ for $t_a=200$ and $\eta =0.18$ at $T_B=1$. We fixed $\alpha =1$. The dashed line and the solid vertical line indicate $(1/N)⟨H_I{⟩}_{t_a}t/t_a$ and $s^{*}\equiv T_B/T^{*}$, respectively, where $T^{*}$ is determined by the minimization of the Kullback-Leibler (KL) divergence. (b) KL divergence $D_{\mathrm{KL}}(T)$ between the final state after QA and the Gibbs state of $H_I$ with temperature $T$. See the main text for a detailed definition. (c) Excess energy ${\mathcal{E}}_{\mathrm{exc}}$ per spin from the thermal expectation value after QA as a function of $t_a$ for various $\alpha$. Lines indicate the best power-law fits ${\mathcal{E}}_{\mathrm{exc}}=a{t}_a^{-b}$ to the data for $t_a>100$ with the fitting parameters $a$ and $b$. (d) Exponent $b$ as a function of $\alpha$. The numerical results (symbols) are compared to the theoretical prediction shown by the solid line. The parameters used in the numerical simulations are ${\omega }_c=5$, ${\tau }_c=10$, $\mathrm{\Delta }t=0.05$, and $N\to \infty$. The bond dimensions are up to 128.

Figure 2

Excess energy per spin after QA scaled by $t_a^{-1/3}$ as a function of $\eta$ for $T_B=1$, $\alpha =1$, and various $t_a$ ranging from 50 to 1000. With increasing $t_a$, the data collapse into a single nonmonotonic curve, which implies ${\mathcal{E}}_{\mathrm{exc}}\sim t_a^{-1/3}$ for large $t_a$. The inset shows the data rescaled by $(\eta t_a{)}^{-1/3}$. The constancy of the data with large $t_a$ near $\eta =0$ corresponds to Eq. (19). The parameters in the numerical simulations are the same as those in Fig. 1.

Figure 3

Energy expectation values per spin as functions of the rescaled time $s$. The solid and dashed lines show the energies per spin of the instantaneous equilibrium states of the closed system at $T=0$ and $T=1$, respectively. The filled symbols show the energy of the time-dependent state for $t_a=1$ of the closed system ($\eta =0$). The empty symbols denote energies of the time-dependent states of the dissipative system with $\eta =0.02$ and $T_B=1$ for $t_a=1$ and $t_a=200$. We fixed $\alpha =1$. The parameters in the numerical simulations are the same as those in Fig. 1.