Quantum Monte Carlo tunneling from quantum chemistry to quantum annealing

Quantum tunneling is ubiquitous across different fields, from quantum chemical reactions and magnetic materials to quantum simulators and quantum computers. While simulating the real-time quantum dynamics of tunneling is infeasible for high-dimensional systems, quantum tunneling also shows up in quantum Monte Carlo (QMC) simulations, which aim to simulate quantum statistics with resources growing only polynomially with the system size. Here we extend the recent results obtained for quantum spin models [Phys. Rev. Lett. 117, 180402 (2016)], and we study continuous-variable models for proton transfer reactions. We demonstrate that QMC simulations efficiently recover the scaling of ground-state tunneling rates due to the existence of an instanton path, which always connects the reactant state with the product. We discuss the implications of our results in the context of quantum chemical reactions and quantum annealing, where quantum tunneling is expected to be a valuable resource for solving combinatorial optimization problems.

Figure 1

Illustration of a double-well potential energy $V\left(x\right)$ and energy levels. The degenerate levels $E_L$ and $E_R$ correspond to the localized states $|{\psi }_L\rangle$ (blue) and $|{\psi }_R\rangle$ (green). The degeneracy is lifted by the linear combination of localized states that produce the true eigenstates ${\psi }_0=1/\sqrt{2}\left(\right|{\psi }_L\rangle +|{\psi }_R\rangle )$ (red curve) and ${\psi }_1=1/\sqrt{2}\left(\right|{\psi }_L\rangle -|{\psi }_R\rangle )$. The tunneling splitting $\mathrm{\Delta }$ can be calculated from the overlap of the localized states.

Figure 2

Left: Illustration of the typical instantonic paths in configuration space, with PBC, ${\mathbf{x}}^{**}\left(\tau \right)$ (cyan), and OBC in imaginary time, ${\mathbf{x}}^{*}\left(\tau \right)$ (pink). These paths are transition states of the PIMC and path integral ground-state (PIGS) pseudodynamics, respectively (in the space of imaginary-time trajectories) in double-well models (sketched in the gray-scale heat map; see Fig. 6 for a more realistic example). Right: Instantonic trajectories (projected on the reaction coordinate $x$ axis) as a function of the imaginary time $\tau$. Notice that PIMC instantons have to cross twice the barrier to fulfill the PBC constraint.

Figure 3

Average mean first tunneling time (MFTT) with PIMC (for PBCs and OBCs) as a function of $x_0$ for different values of $\lambda$, at $\beta =20$, corresponding to a temperature always much lower than the barrier height. We use a dimensionless mass parameter $m=1/2$. The inset shows the shape of the double-well potential $V\left(x\right)$, whose barrier width (at the top) is $2x_0$. In the OBC case, deviations from the expected $1/\mathrm{\Delta }$ behavior occur when the tunneling rate becomes larger.

Figure 4

Top panel: position distributions (histograms) obtained considering the center (blue) or the tail (orange) of the OBC path. The distributions are area-normalized, respectively, with the exact $\rho \left(x\right)\approx |{\psi }_0{|}^2$ distribution (red) and the exact ground state ${\psi }_0\left(x\right)$ (green). We plot only for $x>0$ and we use $x_0=3$ and $\lambda =0.14$ in Eq. (5). The difference between the sampled distributions and the reference ones are negligible. We perform simulations at low temperatures, $\beta =20\gg \mathrm{\Delta }V$. Middle and lower panel: the position of the particle as the simulation progresses for PBCs and OBCs (both for the center and the tail). As expected, the tunneling rate is much larger for OBCs.

Figure 5

Average MFTT tunneling time with PIMC (with PBC) as a function of ${\omega }_y$ for different values of $g=0.2,0.3,0.4$ (pink, blue, and red data series, respectively) and two values of $y_0=0$ (empty symbols) and 1 (full symbols). The potential considered is $V_A$ as in Eq. (8). Lines represent a fit to the exact ${\mathrm{\Delta }}^{-2}$ gap values, obtained with the DVR method. The proportionality constant $\alpha \left(g\right)$ that multiplies the inverse gap squared is different for each $g$ value, and it is fitted using only the $y_0=0$ data series. Notice the logarithmic scale on both axes.

Figure 6

Heat-map plot of the two-dimensional model potential $V_S(x,y)$ of malonaldehyde of Eq. (9). Red lines represents a collection of 20 instanton paths sampled with PIMC. These samples are uncorrelated as they correspond to independent tunneling events after full reinitialization of the starting path in the reactant well. That is, when the MFTT identification criteria are met (see the text), we stop the PIMC simulation and collect the last path that has been generated. We use a sufficiently large inverse temperature ($\beta =400$, $P=512$), Trotter slices, and OBCs in imaginary time. We single out, in orange, one of these instances in order to appreciate their instantonic character. Indeed, most of the Trotter slices are located at the bottom of the two wells, and only a very few of them are located on the barrier (cf. Fig. 2); these correspond to the middle of the imaginary-time trajectory ($\tau \approx \beta /2$). The PIMC instanton paths are not smooth, given the large number of Trotter slices, and they represent fluctuations around an average transition path that is qualitatively very close to the RPI solution taken from Ref. [50] (green; see the text). We also plot the MEP (blue) for comparison. The proton tunneling paths typically take place on a region quite far from the saddle point $(0,\gamma /{\omega }_y^2\approx 1.7)$ of the potential.

Figure 7

(a) Unstable manifold with a cusp singularity and three typical imaginary-time paths emanating from one of the minima of the potential. (b) Projection of the unstable manifold onto the coordinate plane ($X_1,X_2$). The two folds project onto caustics, while projections of the three trajectories intersect at the point $X_0$ that lies on the switching line, showed as a dashed line. (c) The action $S\left(X_0\right)$ is a three-valued function in between the caustics. Two of its lower branches intersect along the switching line. After Refs. [67, 69].