Quantum phase transition with inhomogeneous driving in the Lechner-Hauke-Zoller model

TL;DR

This paper investigates the phase diagram of the Lechner-Hauke-Zoller model, deriving analytic expressions for free energy and critical points, and proposes an inhomogeneous driving method to avoid first-order transitions, enhancing quantum annealing efficiency.

Contribution

It introduces an inhomogeneous transverse field with a time-dependent parameter to evade first-order quantum phase transitions in the model.

Findings

Analytic expressions for free energy and critical coefficients derived.

Inhomogeneous driving improves ground-state preparation efficiency.

Numerical verification confirms theoretical predictions.

Abstract

We study the zero-temperature phase diagram of the Lechner-Hauke-Zoller model. An analytic expression for the free-energy and critical coefficients for finite-size systems and in the thermodynamic limit are derived and numerically verified. With the aim to improve standard quantum annealing, we introduce an inhomogeneously driven transverse field with an additional time-dependent parameter that allows one to evade the first-order quantum phase transition and, thus, improve the efficiency of the ground-state preparation considerably.

Equations

$${\mathcal{H}}_{\text{LHZ}}(s)$$

$${\mathcal{H}}_I(s)$$

$${\mathcal{H}}_P(s)$$

$${\mathcal{H}}_{\text{LHZ}}(s,r)$$

$$h_k(s,r)$$

$$s_k$$

$$s_{k+1}$$

$$Z_{\text{LHZ}}=\text{Tr}\left[e^{-\beta {\mathcal{H}}_{\text{LHZ}}(s,r)}\right].$$

$$Z_{\text{LHZ}}=\underset{M\to \infty }{\lim }{Z}_{\text{LHZ,M}}=\underset{M\to \infty }{\lim }\text{Tr}{\left[e^{-\beta s{\mathcal{H}}_P/M}{e}^{-\beta {\mathcal{H}}_I/M}\right]}^M$$

$$=\underset{M\to \infty }{\lim }\sum _{\{{\sigma }^z\}}⟨\{{\sigma }^z\}|\{\exp \left[\frac{s\beta }{M}\sum _{k=1}^{N_p}{J}_k{\sigma }_k^z+\frac{\beta }{M}\sum _{k=1}^{N_p}{h}_k{\sigma }_k^x\right]$$

$${\times \exp \left[\frac{s\beta }{M}\sum _{l=1}^{{\mathrm{N}}_{\mathrm{c}}}{C}_l{\sigma }_{l,n}^z{\sigma }_{l,e}^z{\sigma }_{l,s}^z{\sigma }_{l,w}^z\right]\}}^M|\{{\sigma }^z\}⟩$$

$$E_{N_l}(m)=-C\left(\frac{N_l^2}{2}-\frac{5}{2}{N}_l+3\right)m^4-C(N_l-2)m^3.$$

$$E_{N_p}(m)$$

$$-C(\sqrt{0.25+2N_p}-1.5)m^3.$$

$$Z_{\text{LHZ,M}}=\prod _{\alpha =1}^M\int dm(\alpha )d\tilde{m}(\alpha )$$

$$\exp \{\sum _{\alpha =1}^M\frac{\beta sC}{M}[\left(N_p-\sqrt{1+8N_p}+2\right)m(\alpha {)}^4$$

$$+\left(\sqrt{0.25+2N_p}-1.5\right)m(\alpha {)}^3]-N_p\tilde{m}(\alpha )m(\alpha )\}\times$$

$$\exp \left\{\sum _{k=1}^{N_p}\ln \mathrm{tr}\prod _{\alpha =1}^M\exp \left[\left(\tilde{m}(\alpha )+\frac{\beta s}{M}{J}_k\right){\hat{\sigma }}^z\right]\exp \left[\frac{\beta }{M}{h}_k{\hat{\sigma }}^x\right]\right\}$$

$$=\prod _{\alpha =1}^M\int dm(\alpha )d\tilde{m}(\alpha )\exp [-N_p\beta f_{N_p,M}(\{m(\alpha )\})],$$

$$\tilde{m}(\alpha )$$

$$+\frac{\sqrt{2.25+18N_p}-4.5}{N_p}m(\alpha {)}^2],$$

$$f_{N_p,M}(\{m(\alpha )\})=\frac{sC}{M}\sum _{\alpha =1}^M\left(3+\frac{6-\sqrt{9+72N_p}}{N_p}\right)m(\alpha {)}^4$$

$$+\left(\frac{\sqrt{1+8N_p}-3}{N_p}\right)m(\alpha {)}^3$$

$$-\frac{1}{\beta N_p}\ln \mathrm{tr}\prod _{\alpha =1}^M\exp \{\frac{\beta s}{M}[C\left(4-\frac{\sqrt{16+128N_p}+8}{N_p}\right)m(\alpha {)}^3$$

$$+C\left(\frac{\sqrt{2.25+18N_p}-4.5}{N_p}\right)m(\alpha {)}^2+J_k]{\hat{\sigma }}^z\}$$

$$\times \exp \left(\frac{\beta }{M}{h}_k{\hat{\sigma }}^x\right).$$

$$f(m,s,{\tau }^{\prime },C,J,N_p)$$

$$-\left[{\int }_0^{1}d{\tau }^{\prime }\sqrt{s^2{\left(C\left(4-\frac{\sqrt{16+128N_p}+8}{N_p}\right)m^3+C\left(\frac{\sqrt{2.25+18N_p}-4.5}{N_p}\right)m^2+J\right)}^2+h({\tau }^{\prime }{)}^2}\right],$$

$$h({\tau }^{\prime })=\{\begin{array}{ll}1& \text{for}0<{\tau }^{\prime }<1-\tau ,\\ 0& \text{for}1-\tau <{\tau }^{\prime }<1.\end{array}$$

$$f(m,s,\tau ,C,J)$$

$$-(1-\tau )\sqrt{s^2(4C{m}^3+J{)}^2+1}],$$

$$\displaystyle\dfrac{\partial}{\partial m}f(m,s,\tau,C,J)\bigg{|}_{m=m_{c}}=0,$$

$$\displaystyle\dfrac{\partial^{2}}{\partial m^{2}}f(m,s,\tau,C,J)\bigg{|}_{m=m_{c}}=0,$$

$$\displaystyle\dfrac{\partial^{3}}{\partial m^{3}}f(m,s,\tau,C,J)\bigg{|}_{m=m_{c}}=0,$$

$$m_c\approx 0.679795,s_c\approx 0.219232,{\tau }_c\approx 0.38911.$$

$$F_4$$

$$F_3$$

$$f(m,N_p)$$

$$-\frac{1}{\beta N_p}\sum _{k=1}^{N_p}\ln 2\cosh \beta \sqrt{s^2{\left(C\left(4-\frac{\sqrt{16+128N_p}+8}{N_p}\right)m^3+C\left(\frac{\sqrt{2.25+18N_p}-4.5}{N_p}\right)m^2+J_k\right)}^2+h_k^2},$$

$$f(m,s,{\tau }^{\prime },J,N_p)=s\left[3m^4\sqrt{N_p}+\sqrt{8}\left(m^3-3m^4\right)+\frac{1}{\sqrt{N_p}}\left(m^4\left(6-\frac{3}{4\sqrt{2}}\right)+m^3\left(\frac{1}{4\sqrt{2}}-3\right)\right)-\frac{1}{128\sqrt{2}}\frac{1}{N_p}\left(m^3-3m^4\right)\right]$$

$$-\left[{\int }_0^{1}d{\tau }^{\prime }\sqrt{s^2{\left(4m^3\sqrt{N_p}+\sqrt{8}\left(\frac{3}{2}{m}^2-4m^3\right)-\frac{1}{\sqrt{N_p}}\left(8\frac{1}{\sqrt{2}}{m}^3-m^2\left(\frac{3}{8\sqrt{2}}-\frac{9}{2}\right)\right)-\frac{1}{N_p}\left(\frac{1.5m^2-4m^3}{128\sqrt{2}}\right)+J\right)}^2+h({\tau }^{\prime }{)}^2}\right],$$

Full text

I Introduction

With the recent experimental progress, intermediate scale quantum computers are now available in the laboratory (Preskill, 2018; Alvarez-Rodriguez et al., 2018; Blatt and Roos, 2012; Boixo et al., 2016; O’Malley et al., 2016; Bernien et al., 2017; DiCarlo et al., 2010; de Léséleuc et al., 2018; Georgescu et al., 2014; Boixo et al., 2018; Houck et al., 2012; Aspuru-Guzik and Walther, 2012; Bloch et al., 2012; Raimond et al., 2001). Although these experiments are not yet ready for scalable quantum computing with error correction, these highly developed platforms are suitable for next generation quantum simulations with full control over the individual degrees of freedom. This has motivated the concept of computation by quantum simulation, i.e. to use the simulation toolbox to solve computational problems. In particular adiabatic quantum computing (Albash and Lidar, 2018) (also known as quantum annealing) as a metaheuristic to solve combinatorial optimization problems has been studied extensively (Kadowaki and Nishimori, 1998; Finnila et al., 1994; Brooke et al., 1999; Farhi et al., 2000, 2001; Santoro et al., 2002; Santoro and Tosatti, 2006; Johnson et al., 2011; Boixo et al., 2013; Dickson et al., 2013; Boixo et al., 2014; Glaetzle et al., 2017). Despite the sizable theoretical and experimental efforts, the path towards demonstrating any quantum advantage in adiabatic quantum computing is still elusive.

From the viewpoint of statistical mechanics, solving optimization problems with adiabatic quantum computing can be understood as driving a random transverse Ising model through a zero-temperature quantum phase transition (QPT) in a one-dimensional phase diagram. In the thermodynamic limit, the QPT is associated with a minimal energy gap that closes exponentially for first-order phase transitions and polynomially for second-order transitions. The limiting factor for the efficiency of quantum annealing is the scaling of the inverse of the squared minimal energy gap between ground and first excited state (Tsuda et al., 2013; Mishra et al., 2018). For a random all-to-all model one would expect a second-order phase transition; however, it was recently shown that additional first-order transitions appear in the regime of small transverse fields Knysh (2016).

The Lechner-Hauke-Zoller (LHZ) mapping (Lechner et al., 2015) is an alternative to the Ising spin glass model. LHZ consists of four-body interactions which are problem independent and random local fields. The associated quantum phase transition is thus fundamentally different compared to the all-to-all spin glass. However, LHZ is similar to an ordered p-spin model with a random local field and one expects a first-order phase transition in the limit where the four-body interactions are dominant. Recently, Susa and co-workers (Susa et al., 2018a, b) (see also Refs. Dziarmaga and Rams (2010); Adame and McMahon (2018); Rams et al. (2016)) showed for $p$-spin models, that inhomogeneneous driving of the transverse field can circumvent the first-order phase transition and improve the ground-state preparation.

In this paper, we derive the expression for the free-energy of the LHZ model and the critical coefficients of the associated quantum phase transition for the thermodynamic limit and estimate finite-size effects. We apply an inhomogeneous driving protocol for the transverse field which has been recently introduced for $p$-spin models (Susa et al., 2018a, b) to alter the phase diagram. The inhomogeneous driving introduces an additional dimension in the phase diagram and in this two-dimensional parameter space we are able to evade a particular first-order phase transition in the adiabatic protocol. We numerically demonstrate the implementation of the inhomogeneously driven transverse fields in LHZ and find an enhanced final ground-state fidelity and an enlarged minimal energy gap compared to standard quantum annealing.

II free-energy of the LHZ model

II.1 Four-body transverse Ising model

Transverse Ising models are a cornerstone of modern statistical mechanics and their quantum phase transitions have been studied extensively (for a review, see, e.g., Ref. Stinchcombe (1973)). The quantum phase transition of random transverse Ising spin models (so-called spin glasses) has recently regained considerable interest with the emergence of quantum annealing as a possible application. Quantum annealing is a metaheuristic that aims at solving combinatorial optimization problems which are encoded in Ising spin glasses (Kadowaki and Nishimori, 1998; Finnila et al., 1994). In this scheme, finding the minimum energy of the spin glass is equivalent to determining the solution of the optimization problem (Lucas, 2014). In typical examples of encoding optimization problems in the form of Ising models $\mathcal{H}_{P}=\sum_{i<j}J_{ij}\sigma_{i}^{z}\sigma_{j}^{z}$, the interaction matrix $J_{ij}$ has infinite range and is random.

In quantum annealing, the ground-state of $H_{P}$ is obtained by adiabatically connecting it to a trivial Hamiltonian, e.g. $\mathcal{H}_{I}=\sum_{k}^{N}\sigma_{k}^{x}$. The system is initially prepared in the ground-state of $H_{I}$ and further sufficiently slowly transferred to the problem Hamiltonian $\mathcal{H}_{P}$ via the protocol $\mathcal{H}(s)=[1-f(s)]\mathcal{H}_{I}+f(s)\mathcal{H}_{P}$ where $f(s)$ is a smooth function in the normalized time $s=t/t_{f}$ with $f(s=0)=0$ and $f(s=1)=1$ and $t_{f}$ is the running time of the sweep. In switching from $H_{I}$ to $H_{P}$, the system undergoes a quantum phase transition which limits its efficiency. For an Ising spin glass, one would expect a second-order phase transition at critical time $s^{*}$ with a polynomial closing gap. However, it was recently shown that additional exponentially closing gaps are present for $s>s^{*}$ Knysh (2016).

An alternative to the spin glass encoding of optimization problems has been recently introduced by LHZ (Lechner et al., 2015). In this model, physical qubits describe the relative configuration of two logical spins taking the values $+1$ for parallel (i.e., $\uparrow\uparrow$, $\downarrow\downarrow$) and $-1$ for antiparallel ($\uparrow\downarrow$, $\downarrow\uparrow$) alignment. The time-dependent Hamiltonian in LHZ reads

[TABLE]

where $\sigma_{k}^{x}$ and $\sigma_{k}^{z}$ are the $x$- and $z$-Pauli matrices for the physical qubit at site $k$ and the strengths of all local fields $h_{k}$, $J_{k}$ and constraints $C_{l}$, respectively, depend on time. Here, $\mathcal{H}_{I}(s)$ is the driver term, and $\mathcal{H}_{P}(s)$ is the encoded problem Hamiltonian to be solved.

The strengths of the controllable local magnetic fields $h_{k}$ and $J_{k}$ in Equations (2) and (3) are applied to all $N_{p}=N_{l}(N_{l}-1)/2$ physical qubits where $N_{l}$ is the number of logical spins in the original model. The third sum runs over $N_{c}=N_{p}-N_{l}+1$ four-body constraints among nearest neighbor qubits on a square lattice and $C_{l}$ is the strength of a four-body constraint at plaquette $l$. The introduction of these four-body constraints accounts for the increased number of degrees of freedom from $N_{l}$ logical to $N_{p}$ physical qubits. This notation excludes $N_{a}=N_{l}-2$ auxiliary physical qubits in the bottom row of the LHZ architecture to obtain four-body constraints on the whole square lattice. The indices $(l,n)$, $(l,w)$, $(l,s)$, and $(l,e)$ denote the northern, western, southern and eastern physical qubits of the constraint $l$, respectively (more details in Ref. (Lechner et al., 2015)). Given that the constraints are the dominant energy scale, the model is thus similar to the $p$-spin model.

II.2 Inhomogeneous transverse field

The $p$-spin model (Nishimori and Takada, 2017) with a standard homogeneously driven tranverse field undergoes a first-order quantum phase transition in the zero-temperature phase diagram. Thus, in LHZ–with its similarity to the $p$-spin model for $p=4$–a first-order QPT is also expected during a quantum annealing sweep. As the minimal energy gap of Hamiltonian (1) between the ground-state and the first excited state decreases exponentially with increasing system size $N$ at the critical point (i.e., $\propto e^{-aN},\,a>0$), the computation time $t_{f}$ grows exponentially (i.e., $t_{f}\propto|\langle 1|d\mathcal{H}/dt|0\rangle|/\Delta^{2}$ with $\Delta$ as the minimal energy gap and $\ket{0}$ and $\ket{1}$ as the instantaneous ground-state and first excited state, respectively) with the system size according to the adiabatic theorem and Landau-Zener’s formula.

Spatiotemporal inhomogeneous driving (Mohseni et al., 2018; Susa et al., 2018a, b) of the transverse field introduces an additional parameter and, thus, an additional dimension in the phase diagram. This allows one to avoid first-order phase transitions by connecting $H_{I}$ and $H_{P}$ via a continuous path around the critical point. This is achieved by switching off the strength of the transverse field inhomogeneously. With the goal to apply this to LHZ in mind, we modify our Hamiltonian (1) as

[TABLE]

where $h_{k}(s,r)$ is the strength of the inhomogeneously driven transverse field. In this paper, we choose a protocol for the strength of the transverse field that reads

[TABLE]

This protocol first switches off the transverse field of the qubits in the first row and the auxiliary qubits in the last row (see Fig. S1 in the Appendix).

The protocol $h_{k}(s,r)$ is chosen as a continuous piecewise-differentiable function to avoid diverging derivatives of the Hamiltonian (4). Here, we have included a new parameter $r$ which enters in an additional time-dependent function $\tau=s^{r}$ with $0\leq\tau\leq 1$. In this spatiotemporal formulation, $s$ and $\tau$ are both controlled as a function of time with $s=\tau=0$ at time $t=0$ and $s=\tau=1$ at time $t=t_{f}$, the total sweep time.

II.3 Inhomogeneously driven LHZ

In the following, we derive the free-energy of Hamiltonian Eq. (4) with an inhomogeneous driving field. This derivation follows the Suzuki-Trotter decomposition used in Refs. (Jörg et al., 2010; Seki and Nishimori, 2012). The partition function of Hamiltonian Eq. (4) reads

[TABLE]

Using the Suzuki-Trotter decomposition $e^{(A+B)}=\lim_{n\to\infty}\left(e^{A/n}e^{B/n}\right)^{n}$ with $A,B$ being quantum operators, the partition function reads

[TABLE]

where $\sum_{\{\sigma^{z}\}}$ refers to the summation over all $2^{N_{p}}$ possible spin configurations in the $z$-basis and with $\{\sigma^{z}\}=\bigotimes_{k=1}^{N_{p}}\ket{\sigma_{k}^{z}}$. We introduce $M$ replicas of the quantum state $\ket{\sigma(\alpha)}$, each labeled $\alpha(=1,\dots,M)$ such that $\hat{\mathbb{1}}(\alpha)=\sum_{\{\sigma^{z}(\alpha)\}}\ket{\{\sigma^{z}(\alpha)\}}\bra{\{\sigma^{z}(\alpha)\}}\times\sum_{\{\sigma^{x}(\alpha)\}}\ket{\{\sigma^{x}(\alpha)\}}\bra{\{\sigma^{x}(\alpha)\}}$. In these replicas, $\alpha$ can be understood as an imaginary time in the dynamic evolution through all these replicas.

Next we derive the expression for the energy of the four-body term in Eq. (4), i.e. $E$ as a function of the magnetization $m$ for a given number of physical qubits $N_{p}$.

LHZ consists of $N_{l}-2$ three-body plaquettes and $(N_{l}-1)(N_{l}-2)/2-N_{l}+2$ four-body plaquettes. Thus, the model resembles a mixture of a $p$-spin model with $p=3$ and $p=4$ and sparse connectivity. Counting the numbers of constraint terms, the energy as a function of magnetization expressed as a function of the number of logical qubits $N_{l}$ reads as

[TABLE]

The same equation expressed as a function of the physical qubits $N_{p}$ reads as

[TABLE]

Here, we assumed that the constraints are the dominant energy in the system and we neglected the random field terms. Note, that local fields that are drawn randomly form a distribution with mean at [math] and their contribution to the energy averages out for all magnetizations.

In order to verify the energy expression (9), we calculate numerically the energy of LHZ averaged over MS samples of each magnetization and fit the result to a function $f(m)=am^{4}+bm^{3}$. We repeat this for system sizes between $N_{p}=6$ up to $N_{p}=5886$. In order to sample the magnetization, we randomly shuffle the configurations of spins being up or down whereas keeping the total number of spins up constant. For example for the case of the magnetization value $m=0.0$, $N_{p}/2$ physical qubits are spin-up and $N_{p}/2$ spin-down. We calculate the energy for randomly shuffled configurations in LHZ with $N_{p}/2$ qubits being spin-up and $N_{p}/2$ spin-down, and compute the mean energy and standard deviation of the four-body term. Similarly we proceed for all other possible magnetization values for a chosen number of physical qubits in LHZ.

Figure 1 depicts the energy of the constraints in Eq. (3) and its fluctuations using the sampling method described above. The origin of the fluctuations in the energy as a function of magnetization is twofold. One is the result of finite sampling MS, and the other is the result of entropy, i.e., configurations with the same $m$ can have different energies. Let us first consider the error from sampling.

Figure 1(a) depicts the energy as a function of the magnetization for various system sizes $N_{p}$. For small systems, the energy is asymmetric and resembles a cubic function. For large system sizes, the energy approaches a quartic and symmetric function in $m$. The standard deviation of the energies is depicted in Fig. 1(b) for the same system sizes. Figure 1(c) depicts the relative errors as a function of the shuffling parameter MS. The standard deviation scales, as expected with the system size as $1/\sqrt{N_{p}}$. For small and large values of shuffling parameter MS, we plotted the relative error of the energy for magnetization $m=0.0$ for different system sizes $N_{p}$. Figure 1(d) shows the relative error for magnetizations $m=0.0$ and $m=\pm 0.5$ and different system sizes $N_{p}$. The inset in Fig. 1(d) shows the relative error multiplied by the inverse of the scaling, i.e. $\sqrt{N_{p}}$, for magnetization $m=0.0$ and increasing system sizes $N_{p}$. Thus for small systems, the entropic energy fluctuations persist for all system sizes.

With the numerical data for the energy, we can now verify the individual terms in the analytical expression of Eq. (8). Figure 2 depicts the comparison of Eq. (9) and the numerical data with a fit $f(m)=am^{4}+bm^{3}$ with parameters $a$ and $b$. The analytical expressions from Eq. (9) are in excellent agreement with the data for both, the cubic and the quartic terms.

As expected from the $N_{p}$ dependence in Eq. (9), the quartic $m^{4}$ term dominates in the thermodynamic limit as shown in Fig. S2(b) in the Appendix for the case of $N_{l}=109$ logical and, thus, $N_{p}=5886$ physical qubits in LHZ.

Plugging expression Eq. (9) into Eq. (7), the decomposed partition function $Z_{\textrm{LHZ,M}}$ reads with the help of the integral definition of the $\delta$ function, i.e., $\delta(N_{p}m(\alpha)-\sum_{k=1}^{N_{p}}\sigma_{k}^{z}(\alpha))=\int d\tilde{m}(\alpha)\exp[-\tilde{m}(\alpha)(N_{p}m(\alpha)-\sum_{k=1}^{N_{p}}\sigma_{k}^{z}(\alpha))]$, as

[TABLE]

where $f_{N_{p},M}(\{m(\alpha)\})$ is the free-energy of the system consisting of $N_{p}$ physical qubits as a function of the magnetization $m$. The parameter $M$ denotes the number of imaginary time slices and $\beta$ is the reciprocal temperature.

We are interested in an expression for $\tilde{m}$ that minimizes the integrand of the partition function and thus the free-energy of LHZ. This saddlepoint condition for $\tilde{m}(\alpha)$, i.e. $\partial Z_{\textrm{LHZ,M}}/\partial m=0$ and solving for $\tilde{m}$, reads

[TABLE]

and the free-energy, thus, becomes

[TABLE]

We now apply the static approximation $m=m(\alpha)$ for all $\alpha$ and take the reverse operation of the Suzuki-Trotter decomposition for $M\to\infty$. We can further rewrite the strength $h_{k}$ of the inhomogeneously driven transverse field as a continuous function $h(\tau^{\prime})$ with $\tau^{\prime}\in[0,1]$ as denoted in Eq. (5) to obtain the integral form of the free-energy in the zero-temperature limit $\beta\to\infty$, i.e., $1/\beta\ln\,2\,\cosh\,\beta\to 1$, as

[TABLE]

where the square brackets $[\cdots]$ over the integral denote the average value over the distribution of the strengths $J_{k}$ of the longitudinal magnetic field denoted as $J$ according to the law of large numbers $\lim_{N_{p}\to\infty}1/N_{p}\sum_{k=1}^{N_{p}}(\cdots)=[(\cdots)]$. This means the integral is evaluated for uniformly distributed values of $J_{k}$. Equation (13) together with the protocol $h(s,r)$ in Eq. (5) for the inhomogeneously driven transverse field [i.e., $h(\tau^{\prime})$ in Eq. (13)] describes the free-energy of LHZ for finite sizes with respect to the number $N_{p}$ of physical qubits (see the Appendix for further information).

For finite-size systems with $N_{p}$ physical qubits in LHZ, we use the expression $h(s,\tau)$ of Eq. (5) for $h(\tau^{\prime})$ in Eq. (13). In the thermodynamic limit $N_{p}\gg 1$, the interval $[s_{k+1},s_{k}]$ for $s$ of Eq. (5) becomes infinitesimally small, and thus, the protocol for the strength of the inhomogeneously driven transverse field can be written as

[TABLE]

For this choice of inhomogeneous transverse field strength $h$, the free-energy of LHZ in the thermodynamic limit can be written as

[TABLE]

where $C$ and $J$ are the strengths of the constraints and longitudinal magnetic fields, respectively.

The critical coefficients $m_{c}$, $s_{c}$, and $\tau_{c}$ of the free-energy term Eq. (15) are obtained by a Landau-type expansion of the free-energy term with the condition that the first three derivatives vanish (Nishimori and Ortiz, 2011), i.e., we have to solve the system of equations

[TABLE]

with respect to its critical coefficients $m_{c}$, $s_{c}$ and $\tau_{c}$. Note, that the critical coefficients are a function of the constraint strength $C$ and the distribution of $J$. The values of the critical coefficients for a uniform distribution of the strength of the longitudinal magnetic field $J$ with values between $-1$ and $1$, and constraint strength $C=2$ are

[TABLE]

We can obtain these critical coefficients with our thermodynamical free-energy term Eq. (15) [or with finite-size free-energy term Eq. (13) by increasing the number of physical qubits $N_{p}$] and which can be seen in Fig. 3. Here, we have plotted the free-energy term of Eq. (15) with respect to the magnetization $m$ for different points $(s,\tau)$ in the two-dimensional phase diagram. On the first-order transition line [points (a), (b) and (d)], we see a degenerated minimum of the free-energy. At the crossing of the first-order transition line starting from (c) and going to (e), we see that the value for the magnetization that minimizes the free-energy changes discontinuously from the paramagnetic solution $m=m_{\mathrm{p}}=0.0$ to the ferromagnetic solution $m=m_{\mathrm{f}}=1$ and which depicts a quantum phase transition of first order.

III Numerical Results

Let us now apply the results to a quantum annealing protocol. For quantum annealing, an important measure of the efficiency is the ground-state fidelity $F(t_{f})=\langle\psi(t_{f})|\phi_{0}(t_{f})\rangle$ with $\ket{\psi(t_{f})}$ as the state of our system and $\ket{\phi_{0}(t_{f})}$ as the ground-state of our final Hamiltonian at time $t=t_{f}$. Another important measure is the minimal energy gap $\Delta E_{\textrm{min}}$ of the corresponding energy eigenspectra.

Figure 4 shows the statistics of the squared final ground-state fidelities $F^{2}(t_{f})$ for sweeps with different running times $t_{f}$ for an ensemble of 100 randomly uniformly distributed instances of interactions $J_{k}$ for homogeneous (1) and inhomogeneous Eq. (4) drivings, respectively. The system size is $N_{p}=6$ plus two auxiliary physical qubits, and strengths of the constraints are each $C=2$ for all three constraints and $10$ is the value of the strength of the auxiliary local fields in the bottom row of LHZ. The free parameter in the inhomogeneously driven transverse field is $r=0.5$.

One can see that inhomogeneous driving of the transverse field can enhance the performance of traditional quantum annealing considerably. Furthermore, the ratio of the squared final ground-state fidelities of inhomogeneous to homogeneous driving improves with increasing system size as shown in Fig. 4(b).

A free parameter in the protocol is the choice of the control parameter $r$ in the inhomogeneous sweep, i.e., the path we take in the two-dimensional $s-\tau$ diagram. Figure S4 in the Appendix shows that the control parameter $r$ can considerably enhance the performance of inhomogeneous driving.

The excess energy can also be lowered considerably by inhomogeneous driving of the transverse field. For the same ensemble of 100 randomly chosen instances for system size $N_{p}=6$ plus two auxiliary physical qubits, we have plotted the excess energy $E=\sum_{n}E_{n}-E_{0}$ (where $E_{n}$ denotes the energy of the $n$-th lowest eigenstate) for homogeneous (1) and inhomogenous (4) driving of the transverse field in Fig. S5 in the Appendix.

An interesting question arises whether the particular order at which the transverse fields of the qubits are switched off is relevant for the performance. Statistics of the ground-state fidelities and excess energies of homogeneous (1) and inhomogeneous Hamiltonian (4) with value $r=0.5$ and descending order are included in Fig. S6 in the Appendix. The ground-state fidelities and excess energies are the same as for the ascending order case. The method is insensitive to the precise order in which the transverse fields are switched off.

The minimal energy gap is considered one of the main limiting factors in quantum annealing. The inset of Fig. 4a depicts the statistics of the minimal energy gaps of homogeneous Hamiltonian (1) and inhomogeneous Hamiltonian (4) with value $r=0.5$ for different system sizes $N_{p}$ over an ensemble of 100 randomly chosen instances of interaction strength $J_{k}$.

Due to inhomogeneous driving of the transverse field we can considerably enlarge the minimal energy gap for all instances compared to standard quantum annealing. Also, the ratio of the minimal energy gap of inhomogeneous driving to homogeneous driving increases with system size.

IV Conclusion and Outlook

In this paper, we have introduced an inhomogeneously driven transverse field of the Hamiltonian in the LHZ lattice gauge model architecture. We find that by using inhomogeneous driving of the transverse field in LHZ the ground-state fidelities are increased considerably compared to standard quantum annealing. The method is insensitive to the order in which the fields are switched off.

As an important step, we analytically derived an energy expression of the four-body constraint term $\sigma^{z}\sigma^{z}\sigma^{z}\sigma^{z}$ in LHZ. The term in front of the $m^{4}$ term stems from the four-body constraints of LHZ whereas the term in front of $m^{3}$ stems from the three-body constraints in the lower row of LHZ. As the ratio of three-body to four-body constraints converges towards the value 0 for increasing system sizes $N_{p}$, i.e., the finite-size effect of the three-body constraints vanish, the $m^{4}$ term dominates as can be seen in Figs. 1(a) and S2 in the Appendix. In our derivation, we followed the Suzuki-Trotter decomposition, the saddlepoint approximation as well as the static approximation. We believe that these approximations are valid for LHZ as it was shown recently that it reproduces the exact free-energy term under some valid constraints (Okuyama and Ohzeki, 2018).

From our free-energy expression (13) of LHZ, we calculated the critical coefficients for different system sizes and for the thermodynamic limit, respectively, and where we further computed the line of first-order quantum phase transitions. We note here that, for small system sizes, we can always avoid first-order phase transitions for different values of the control parameter $r$, i.e., taking different paths through the two-dimensional $s-\tau$ diagram. We further note that, for quantum annealing in LHZ for larger system sizes, the strength of the constraints obey a scaling behavior which changes the values of its critical coefficients considerably. We have included an analytical expression of the free-energy expression for a scaling behavior of $C\propto\sqrt{N_{p}}$ and the calculation of the critical coefficients for different system sizes in the Appendix (see also Fig. S3).

Furthermore, we numerically demonstrated an increase in the minimal energy gap and final ground-state fidelity due to avoiding first-order quantum phase transitions by inhomogeneously driving the strength of the transverse field. We note that the ratio of the final ground-state fidelities of inhomogeneous driving Eq. (4) to homogeneous driving Eq. (1) increases with the system size $N_{p}$ of physical qubits in LHZ which we expect due to the exponential closing of the minimal energy gap with system size. This is an encouraging result which we will further study numerically for larger system sizes by using path integral Monte Carlo methods (Mazzola and Troyer, 2017; Isakov et al., 2016; Inack et al., 2018) in future work.

We note, that the minimal gap of random Ising models is not directly at the critical point but rather in the spin-glass phase close to the end of the adiabatic protocol Knysh (2016). Whether inhomogeneous driving also affects these additional exponential closing gaps is an interesting question for future studies.

As a future direction, the inhomogeneous driving scheme in the LHZ model may be applied to the counter-diabatic driving of LHZ as described in Ref. (Hartmann and Lechner, 2019) where quantum phase transitions during the sweep may decrease the efficiency of the approximate counterdiabatic term added to speedup quantum annealing. This may open a new branch of developing fast near-term quantum annealer devices.

V Acknowledgments

We thank T. Hatomura, H. Nishimori, L. Sieberer and K. Ender for valuable discussions. The research was funded by the Austrian Science Fund (FWF) through a START Grant under Project No. Y1067-N27 and the SFB BeyondC Project No. F7108-N38, the Hauser-Raspe Foundation, and the European Union’s Horizon 2020 Eesearch and Innovation Program under Grant Agreement No. 817482 PasQuanS.

Appendix

Inhomogeneous driving of the transverse field

The inhomogeneous driving of the transverse field $h_{k}(s,r)$ in the main text is assumed to be a linear function with different slopes and delays. Figure S1 shows the function $h_{k}(s,r)$ of Eq. (5) in the main text for $N_{l}=4$ logical and, thus, $N_{p}=6$ plus two auxiliary physical qubits in the bottom row in LHZ with parameter value $r=0.5$ in the inhomogeneous driving scheme. Here, $s=t/t_{f}$ is the normalized time, and $k=1$ denotes the physical qubit on the lower left in LHZ and $k=8$ denotes the additional physical qubit fixed to a value of 1 on the bottom right in LHZ, i.e., we first switch off the transverse field of the physical qubit in the lower left and at last the additional qubit in the bottom right.

Derivation of the free-energy in LHZ

The free-energy of LHZ in the thermodynamic limit $N_{p}\to\infty$ [Eq. (15) in the main text] is similar to the free-energy expression of the $p$-spin model for $p=4$ (Susa et al., 2018b, a). Compared to the $p$-spin model, the LHZ model contains $m^{4}$ and $m^{3}$ terms with a ratio that depends on the system size. In LHZ, the ratios $N_{4}/N_{c}$ of the number of four-body constraints to all constraints as well as the ratios $N_{3}/N_{c}$ the number of three-body constraints to all constraints read

[TABLE]

These scalings translate to the terms in front of $m^{4}$ and $m^{3}$ [i.e., Equations (8) and (9)] compared to the $p$-spin model with $p=4$. The ratio $F_{4}$ of four-body constraints to all constraints converges towards 1 and the ratio $F_{3}=1-F_{4}$ towards 0, meaning that the finite-size effect of three-body constraints are negligible in the thermodynamic limit and vice versa the term of the four-body constraints dominates.

With this fact in mind, Fig. S2 shows the energy which approaches a symmetric function due to the fact that the $m^{4}$ term is dominant for the case of $N_{p}=5886$ physical qubits.

For the derivation of the finite-size free-energy expression Eq. (13) in LHZ, we need to apply the static approximation $m=m(\alpha)$ for all $\alpha$ and take the reverse operation of the Suzuki-Trotter decomposition for $M\to\infty$ for the expression Eq. (12). This gives us the expression,

[TABLE]

where we further use the zero-temperature limit $T\to 0$, i.e., $\beta\to\infty$, and rewrite the sum into an integral for large $N_{p}$. In this thermodynamic limit $N_{p}\to\infty$, the actually stepwise function $h_{k}$ becomes continuous, i.e., $h_{k}(\tau^{\prime})$, and, thus, we obtain the free-energy expression for a finite-size system in LHZ as in Eq. (4) in the main text.

A scaling of the constraint strengths in the form of $C_{l}\propto\sqrt{N_{p}}$ may have to be applied in LHZ to suppress any unreasonable solutions in the emerging unreachable sub Hilbert space due to the increase in $N_{l}$ logical to $N_{p}\approx N_{l}^{2}$ qubits and, thus, increasing Hilbert space. To account for this, the finite-size free-energy term Eq. (13) in LHZ can be rewritten as

[TABLE]

Here, it can be clearly seen that the four-body term dominates with increasing system size. In the second term is a constant offset of the free-energy which stems from the finite-size three-body constraints in LHZ. The critical coefficients of the finite-size free-energy expression (A3) in LHZ can, thus, be calculated according to Eq. (16) in the main text and are plotted in Fig. S3.

We can see that the critical coefficients wander from the right of the $s_{c}-\tau_{c}$ diagram to the left. The critical coefficients for a system size $N_{l}=7$ logical spins and, thus, $N_{p}=21$ physical qubits in LHZ, read $s_{c}\approx 0.505$ and $\tau_{c}\approx 0.242$ (blue plus); whereas for a system size $N_{l}=25$ and, thus, $N_{p}=300$, the critical coefficients read $s_{c}\approx 0.029$ and $\tau_{c}\approx 0.371$ (beige plus). For the thermodynamic limits $N_{p}\to\infty$ and $C\to\infty$, the critical coefficients will be on the $\tau_{c}$ axis (i.e., $s_{c}=0$) so that first-order quantum phase transitions cannot be avoided for any choice of the control parameter $r$.

Additional numerical results

Choice of control parameter $r$

As mentioned in the main text, the free control parameter $r$ sets the path in the $s-\tau$ diagram one chooses in order to avoid first-order quantum phase transitions. Its choice can also further enhance the performance of quantum annealing by not only avoiding first-order QPTs, but also enlarging the probability of finding the ground-state of our problem Hamiltonian to be solved.

Figure S4 depicts the statistics of the squared final ground-state fidelity $F^{2}(t_{f})=|\langle\psi(t_{f})|\phi_{0}(t_{f})\rangle|^{2}$ for Hamiltonian Eq. (4) with system size $N_{p}=6$ plus two auxiliary physical qubits in LHZ for an ensemble of 100 randomly chosen instances of $J_{k}$ and different values of the control parameter $r$.

As can be seen here, an appropriate choice of the control parameter $r$ can enlarge the squared final ground-state fidelity of the inhomogeneously driven Hamiltonian considerably. For a value of $r=0.5$, the final ground-state fidelity squared reaches its maximum.

The excess energy is another measure for the performance of quantum annealing as it gives rise to the amount of transitions to higher excited eigenstates that have occurred during quantum annealing sweeps.

Figure S5 depicts the excess energies of the homogeneously driven Hamiltonian (1) and inhomogeneously driven Hamiltonian Eq. (4) with parameter $r=0.5$ and where the parameters of LHZ are as described in the main text. The excess energies of the inhomogeneously driven Hamiltonian (4) are smaller than for the homogeneously driven Hamiltonian (1) in LHZ, meaning that less transitions to higher excited states have occurred during these sweeps.

Different orders of inhomogeneous driving

An interesting question arises whether the order in which we switch off the transverse field of each qubit is of matter for our theory.

Figure S6 depicts the statistics of final ground-state fidelity and excess energies of the descending order of inhomogeneously driven transverse fields of the qubits over an ensemble of 100 uniformly distributed interaction strengths $J_{k}$ for system sizes $N_{p}=3$ plus one auxiliary and $N_{p}=6$ plus two auxiliary physical qubits, respectively. Here, we first switch off the transverse field of the qubit at the top of the triangular LHZ structure and the qubits 1 in the left lower row and auxiliary qubits in the last row. We see the same results as in the case of ascending order from the main text with function $h_{k}(s,r)$ as in Figure S1, meaning that the order in which we switch off the transverse field of the qubits does not matter for the efficiency of our method.

Energy spectrum

As we can enlarge the minimal energy gap by inhomogeneous driving of the transverse field in LHZ, we are interested in the energy spectrum during the whole sweep. Figure S7 shows the energy spectra of the homogeneous (1) and inhomogenous Hamiltonian Eq. (4) with parameter $r=0.5$ for system size $N_{p}=6$ plus two auxiliary physical qubits and a randomly chosen instance $J_{k}$ with constraint strength $C_{l}=2$ for all constraints. We can see that the minimal energy gap is enlarged and shifted in time from around $t/t_{f}\approx 0.5$ for the homogeneous Hamiltonian to $t/t_{f}\approx 0.8$ for the inhomogeneous Hamiltonian.