**Optimization by Use of Nature in Physics Beyond Classical Simulated Annealing**

Masayuki Ohzeki

20 Will-be-set-by-IN-TECH

[22] Martinez-Rios, F. & Frausto-Solis, J. [2008c]. Simulated annealing for sat problems using dynamic markov chains with linear regression equilibrium, *MICAI 2008, IEEE*

[23] Metropolis, N., Rosenbluth, A. W. R. M. N. & Teller, A. H. [1953]. Equation of state calculations by fast computing machines, *The journal of Chemicla Physics* 21: 1087–1092. [24] Mezard, M., Parisi, G. & Zecchina, R. [2002]. Analytic and algorithmic solution of

[25] Mezard, M. & Zecchina, R. [2002]. The random k-satisfiability problem: from an analytic

[26] Miki, M., Hiroyasu, T. & Ono, K. [2002]. Simulated annealing with advanced adaptive neighborhood, *Second international workshop on Intelligent systems design and application*,

[27] Munakata, T. & Nakamura, Y. [2001]. Temperature control for simulated annealing,

[29] PassMark [2007]. Cpu burnintest. http://www.passmark.com/download/index.htm.

[30] Sanvicente-Sanchez, H. & Frausto-Solis, J. [2004]. A method to establish the cooling scheme in simulated annealing like algorithms, *Lecture Notes in Computer Science* . [31] Schaefer, T. J. [1978]. The complexity of satisfiability problems, *Proceedings of the tenth annual ACM symposium on Theory of computing*, STOC '78, ACM, pp. 216–226.

[32] van Wandelen, C. J. [2007]. Cpubench. http://cpubench.softonic.com/eie/21660.

[33] Vorte, B. [2007]. Cpumathmark. http://files.aoaforums.com/code.php?file=931.

random satisfiability problems, *Science* 297(5582): 812–815.

solution to an efficient algorithm, *Phys. Rev.* E66(056126).

[28] Papadimitriou, C. H. [1994]. *Computational Complexity*, Addison-Wesley.

Dynamic Publishers, Inc., pp. 113–118.

URL: *http://www.passmark.com/download/index.htm*

URL: *http://doi.acm.org/10.1145/800133.804350*

URL: *http://files.aoaforums.com/code.php?file=931*

URL: *http://cpubench.softonic.com/eie/21660*

*PHYSICAL REVIEW E* 64.

pp. 182–187.

Additional information is available at the end of the chapter

http://dx.doi.org/10.5772/50636

## **1. Introduction**

We prefer to find the most appropriate choice in daily life for convenience and efficiency. When we go to a destination, we often use a searching program to find the fastest way, the minimum-length path, or most-reasonable one in cost. In such a searching problem, we mathematically design our benefit as a multivariable function (cost function) depending on many candidates and intend to maximize it. Such a mathematical issue is called the optimization problem. Simulated annealing (SA) is one of the generic solvers for the optimization problem [14]. We design the lowest-energy state in a physical system, which corresponds to the minimizer/maximizer of the cost function. The cost function to describe the instantaneous energy of the system is called as the Hamiltonian *H*0(*σ*1, *σ*2, ··· , *σN*), where *σ<sup>i</sup>* is the degrees of freedom in the system and *N* is the number of components related with the problem size. The typical instance of the Hamiltonian is a form of the spin glass, which is the disordered magnetic material, since most of the optimization problems with discrete variables can be rewritten in terms of such a physical system,

$$H\_0(\sigma\_1, \sigma\_2, \dots, \sigma\_N) = -\sum\_{\langle ij \rangle} I\_{ij} \sigma\_i \sigma\_{j \prime} \tag{1}$$

where *σ<sup>i</sup>* indicates the direction of the spin located at the site *i* in the magnetic material as *σ<sup>i</sup>* = ±1. The summation is taken over all the connected bonds (*ij*) through the interaction *Jij*. The configuration of *Jij* depends on the details of the optimization problem.

Then we introduce an artificial design of stochastic dynamics governed by the master equation.

$$\frac{d}{dt}P(\sigma;t) = \sum\_{\sigma'} M(\sigma|\sigma';t)P(\sigma';t),\tag{2}$$

where *P*(*σ*; *t*) is the probability with a specific configuration of *σ<sup>i</sup>* simply denoted as *σ* at time *t*. The transition matrix is written as *M*(*σ*� |*σ*; *t*), which obeys the conservation of probability ∑*<sup>σ</sup> M*(*σ*|*σ*� ; *t*) = 1 and the detailed balance condition

$$M(\sigma|\sigma';t)P\_{\rm eq}(\sigma';t) = M(\sigma'|\sigma;t)P\_{\rm eq}(\sigma;t). \tag{3}$$

**Figure 1.** Energy structure of a spin-glass system. We map the possible 2*<sup>N</sup>* configurations to the one-dimensional horizontal axis for simplicity. The vertical axis represents the value of the energy for each configuration.

Here we denote the instantaneous equilibrium distribution (Gibbs-Boltzmann distribution) as

$$P\_{\text{eq}}(\sigma; t) = \exp(-\beta(t)E(\sigma; t)) / Z\_{t\prime} \tag{4}$$

**Figure 2.** Behavior of the instantaneous state in SA. The left panel shows the stochastic searching in SA in the high-temperature region *β* ≈ 0, while the right one describes that in the low-temperature one. The dashed curves express the structure of the energy as in Fig. 1, while the thick ones denote the probability

processes [6]. The convergence theorem states that we reach the equilibrium distribution, if we obey the following schedule or slower rate to change the inverse temperature as

where *α* is exponentially small in *N* and *p* is a constant independence of *N*. An intuitive way to understand the performance of SA is as follows. The inverse temperature controls the range of searching, roughly speaking. The instantaneous state keeps hopping from state to state in a relatively high-temperature region as depicted in Fig. 2. By gradually decrease of the inverse temperature, we narrow the range of searching. The lower energy state means its realization with a higher possibility following the instantaneous equilibrium distribution as in Fig. 2. Demand of a sufficiently slow control of the inverse temperature implies that we need enough time to find the states with relatively lower energies by the stochastic searching

Basically, SA is based on the behavior closely to the instantaneous equilibrium state. Therefore we need to perform the change of the inverse temperature with a sufficiently slow control. In order to improve the performance, in particular to shorten the necessary time, we need to consider the protocol away from the equilibrium state, that is nonequilibrium process.

In this chapter, we show a novel direction to solve efficiently the optimization problem by use of the nature in nonequilibrium physics. In statistical physics, the interest of researchers in nonequilibrium dynamical behavior has increased. Among several remarkable developments, the Jarzynski equality (JE), which is known as a generalization of the second law of thermodynamics, might be possible to change the paradigm in optimization problem by use of the physical nature. The Jarzynski equality relates an average over all realizations during a predetermined nonequilibrium process with an expectation in an equilibrium state. As seen later, the mathematical structure of JE does not depend on the schedule and the rate of changing the external parameter. It means that, if we implement JE to solve the optimization problem, we do not need to demand slow control of the driver of the system. The challenge of the implementation of JE have been performed in several researchers. Although not yet have been studied the performance in the actual application to the optimization problem, we show

*pN* log (*α<sup>t</sup>* <sup>+</sup> <sup>1</sup>), (9)

Optimization by Use of Nature in Physics Beyond Classical Simulated Annealing 43

*<sup>β</sup>*(*t*) = <sup>1</sup>

before the barrier avoids globally searching for the lower energy state.

the possibility of the novel method from several analyses.

for realization of each configuration schematically.

where the instantaneous energy *E*(*σ*; *t*) is the value of the Hamiltonian *H*<sup>0</sup> and *Zt* denotes a normalization factor termed as the partition function. In order to satisfy these conditions, we often use the transition matrix with Metropolis rule as

$$M(\sigma|\sigma';t) = \min(1, \exp(-\beta \Delta E(\sigma|\sigma';t))),\tag{5}$$

where

$$
\Delta E(\sigma|\sigma';t) = E(\sigma;t) - E(\sigma';t), \tag{6}
$$

or heat-bath rule as

$$M(\sigma|\sigma';t) = \delta\_1(\sigma,\sigma') \frac{\exp\left(-\frac{\beta}{2}\Delta E(\sigma|\sigma';t)\right)}{2\cosh\left(\frac{\beta}{2}\Delta E(\sigma|\sigma';t)\right)},\tag{7}$$

where

$$\delta\_1(\sigma, \sigma') = \delta(2, \sum\_{i=1}^N (1 - \sigma\_i \sigma'\_i)). \tag{8}$$

The master equation simulates behavior of relaxation toward a specific distribution associated with the energy of the system. If we evolve the system for a long time with a virtual parameter *β*(*t*) being a constant *β*, which is the inverse temperature in context of physics, the probability distribution converges to the equilibrium distribution. To generate lower energy state, let us set *β* � 1. Unfortunately, the spin-glass system in the low-temperature often exhibits the extremely long time for equilibration of the system. The most relevant reason is on the complicated structure of the energy of the spin-glass system as schematically depicted in Fig. 1. There are barriers between the valleys to avoid hopping from state to state in the low-temperature, where the energy effect is dominant. Therefore it is difficult to reach the equilibrium distribution by a direct simulation with a constant temperature. Instead, by tuning a virtual parameter *β*(*t*) from zero to a large number, we perform stochastic searching while keeping to trace the instantaneous equilibrium distribution. The mathematical guarantee to converge to the instantaneous equilibrium state by gradually changing the inverse temperature has been proved by Geman and Geman based on the classical inhomogenious (time-dependent) Markov chain representing nonequilibrium

2 Will-be-set-by-IN-TECH

**Figure 1.** Energy structure of a spin-glass system. We map the possible 2*<sup>N</sup>* configurations to the one-dimensional horizontal axis for simplicity. The vertical axis represents the value of the energy for

often use the transition matrix with Metropolis rule as

*M*(*σ*|*σ*�

*M*(*σ*|*σ*�

Δ*E*(*σ*|*σ*�

*δ*1(*σ*, *σ*�

; *t*) = *δ*1(*σ*, *σ*�

Here we denote the instantaneous equilibrium distribution (Gibbs-Boltzmann distribution) as

where the instantaneous energy *E*(*σ*; *t*) is the value of the Hamiltonian *H*<sup>0</sup> and *Zt* denotes a normalization factor termed as the partition function. In order to satisfy these conditions, we

; *t*) = min(1, exp(−*β*Δ*E*(*σ*|*σ*�

; *t*) = *E*(*σ*; *t*) − *E*(*σ*�

2 cosh

*N* ∑ *i*=1

The master equation simulates behavior of relaxation toward a specific distribution associated with the energy of the system. If we evolve the system for a long time with a virtual parameter *β*(*t*) being a constant *β*, which is the inverse temperature in context of physics, the probability distribution converges to the equilibrium distribution. To generate lower energy state, let us set *β* � 1. Unfortunately, the spin-glass system in the low-temperature often exhibits the extremely long time for equilibration of the system. The most relevant reason is on the complicated structure of the energy of the spin-glass system as schematically depicted in Fig. 1. There are barriers between the valleys to avoid hopping from state to state in the low-temperature, where the energy effect is dominant. Therefore it is difficult to reach the equilibrium distribution by a direct simulation with a constant temperature. Instead, by tuning a virtual parameter *β*(*t*) from zero to a large number, we perform stochastic searching while keeping to trace the instantaneous equilibrium distribution. The mathematical guarantee to converge to the instantaneous equilibrium state by gradually changing the inverse temperature has been proved by Geman and Geman based on the classical inhomogenious (time-dependent) Markov chain representing nonequilibrium

*β*

(1 − *σiσ*�

) exp − *β*

) = *δ*(2,

*P*eq(*σ*; *t*) = exp(−*β*(*t*)*E*(*σ*;*t*))/*Zt*, (4)

<sup>2</sup> Δ*E*(*σ*|*σ*�

<sup>2</sup> Δ*E*(*σ*|*σ*�

; *t*) 

; *t*)

; *t*))), (5)

, (7)

; *t*), (6)

*<sup>i</sup>*)). (8)

each configuration.

where

where

or heat-bath rule as

**Figure 2.** Behavior of the instantaneous state in SA. The left panel shows the stochastic searching in SA in the high-temperature region *β* ≈ 0, while the right one describes that in the low-temperature one. The dashed curves express the structure of the energy as in Fig. 1, while the thick ones denote the probability for realization of each configuration schematically.

processes [6]. The convergence theorem states that we reach the equilibrium distribution, if we obey the following schedule or slower rate to change the inverse temperature as

$$\beta(t) = \frac{1}{pN} \log\left(\alpha t + 1\right),\tag{9}$$

where *α* is exponentially small in *N* and *p* is a constant independence of *N*. An intuitive way to understand the performance of SA is as follows. The inverse temperature controls the range of searching, roughly speaking. The instantaneous state keeps hopping from state to state in a relatively high-temperature region as depicted in Fig. 2. By gradually decrease of the inverse temperature, we narrow the range of searching. The lower energy state means its realization with a higher possibility following the instantaneous equilibrium distribution as in Fig. 2. Demand of a sufficiently slow control of the inverse temperature implies that we need enough time to find the states with relatively lower energies by the stochastic searching before the barrier avoids globally searching for the lower energy state.

Basically, SA is based on the behavior closely to the instantaneous equilibrium state. Therefore we need to perform the change of the inverse temperature with a sufficiently slow control. In order to improve the performance, in particular to shorten the necessary time, we need to consider the protocol away from the equilibrium state, that is nonequilibrium process.

In this chapter, we show a novel direction to solve efficiently the optimization problem by use of the nature in nonequilibrium physics. In statistical physics, the interest of researchers in nonequilibrium dynamical behavior has increased. Among several remarkable developments, the Jarzynski equality (JE), which is known as a generalization of the second law of thermodynamics, might be possible to change the paradigm in optimization problem by use of the physical nature. The Jarzynski equality relates an average over all realizations during a predetermined nonequilibrium process with an expectation in an equilibrium state. As seen later, the mathematical structure of JE does not depend on the schedule and the rate of changing the external parameter. It means that, if we implement JE to solve the optimization problem, we do not need to demand slow control of the driver of the system. The challenge of the implementation of JE have been performed in several researchers. Although not yet have been studied the performance in the actual application to the optimization problem, we show the possibility of the novel method from several analyses.

#### **2. Population annealing**

We introduce a couple of theories in nonequilibrium statistical physics in short, before we show the actual application. They provide the supplement to make the protocol of SA faster. The Jarzynski equality is the most important key.

#### **2.1. Jarzynski equality**

Among several recent developments in nonequilibrium statistical mechanics, we take JE as an attempt to improve the performance of SA. The Jarzynski equality relates quantities at two different thermal equilibrium states with those of nonequilibrium processes from *t* = 0 to *t* = *T* connecting these two states as [10, 11]

$$
\left< \mathbf{e}^{-\beta \mathcal{W}} \right>\_{0 \to T} = \frac{Z\_T}{Z\_0},
\tag{10}
$$

where we use the formal solution of the master equation by the exponentiated transition

e*δtM*0(*σ*1|*σ*0;*t*0)

*Zt*<sup>1</sup> *Z*0

Repetition of the above manipulation in Eq. (14) yields the quantity in the right-hand side of

*n*−1 ∏ *k*=0

The Jarzynski equality is a consequence from the fluctuation theorem [2–4], which relates the

This leads to the more generic result, for an observable *O*({*σt*}) depending on the

where *O*<sup>r</sup> denotes the observable which depends on the backward process *T* → 0. The brackets with subscript 0 → *T* express the average with the weight *P*0→*T*({*σt*}) over possible

If we choose an observable depending only on the final state, which is denoted as *OT*, instead of *O*({*σt*}), *O*<sup>r</sup> reads an observable at the initial state in the backward process. Then �*O*r�*T*→<sup>0</sup> equals to the ordinary thermal average at the initial equilibrium state with *β<sup>T</sup>* represented by

<sup>0</sup>→*<sup>T</sup>* <sup>=</sup> �*O*�*β<sup>T</sup>*

By looking over the above calculations, we can understand the roll of the exponentiated pseudo work. The resultant distribution after SA is given by *P*0→*T*({*σt*}). The biased distribution with the exponentiated pseudo work, *P*0→*T*({*σt*}) exp(−*Y*) always gives the equilibrium one. In this sense, the exponentiated pseudo work plays a roll to fill the gap between the equilibrium distribution and the resultant one after SA. Therefore, if we skillfully use the exponentiated pseudo work to keep the instantaneous distribution close to

*ZT Z*0

*Ztk*<sup>+</sup><sup>1</sup> *Ztk*

<sup>e</sup>−*<sup>Y</sup>* <sup>=</sup> *ZT Z*0

<sup>0</sup>→*<sup>T</sup>* <sup>=</sup> �*O*r({*σt*})�*T*→<sup>0</sup>

*ZT Z*0

<sup>=</sup> *ZT Z*0

*P*eq(*σ*0; *t*0)

Optimization by Use of Nature in Physics Beyond Classical Simulated Annealing 45

. (15)

. (16)

. (17)

, (18)

. (19)

matrix. We take the first product of the above equation as,

∑*σ*0 

JE as,

**2.2. Fluctuation theorem**

instantaneous spin configurations,

�· · · �*β<sup>T</sup>* , and Eq. (18) leads to

realizations {*σt*}. For *O* = *O*<sup>r</sup> = 1, Eq. (18) reduces to JE.

*<sup>O</sup>*({*σt*})e−*<sup>Y</sup>*

 *OT*e−*<sup>Y</sup>* 

the equilibrium one, we can invent the improved version of SA.

e−*Y*(*σ*1;*t*0)

= *P*eq(*σ*1; *t*1)

*P*eq(*σn*; *tn*)

*P*0→*T*({*σt*}) *PT*→0({*σt*})

∑*σn*

probability *P*0→*T*({*σt*}) with that of the inverse process *PT*→0({*σt*}) as

where the partition functions appearing in the ratio on the right-hand side are for the initial (*t* = 0) and final Hamiltonians (*t* = *T*). The quantity on the right-hand side can be represented by the exponentiated difference of the free energy exp(−*β*Δ*F*) between the initial and final conditions. The brackets on the left-hand side express the nonequilibrium average over all the instantaneous realizations of the degrees of freedom, for instance spin configurations, during the nonequilibrium process with the following path probability defined as

$$P\_{0 \to T}(\{\sigma\_l\}) = \prod\_{k=0}^{n-1} \left\{ \mathbf{e}^{\delta t M(\sigma\_{k+1}|\sigma\_k; t\_k)} \right\} P\_{\mathbf{eq}}(\sigma\_0; t\_0). \tag{11}$$

It implies that the observations of the nonequilibrium behavior can estimate the equilibrium quantity represented by the partition functions, that is the free energy. This equality is regarded as a generalization of the well-known inequality, the second law of thermodynamics �*W*�0→*<sup>T</sup>* ≥ Δ*F*, which can be reproduced by the Jensen inequality. One of the important features is that JE holds independently of the pre-determined schedule of the nonequilibrium process.

In order to consider the improvement of SA, let us apply the nonequilibrium process with change of the temperature. We then have to employ the pseudo work instead of the ordinary performed work due to the energy difference as

$$Y(\sigma; t\_k) = \left(\beta\_{k+1} - \beta\_k\right) E(\sigma),\tag{12}$$

where we use discrete time expressions as *t*<sup>0</sup> = 0 and *tn* = *T* for simplicity and we assume that the instantaneous energy does not depend on the time as *E*(*σ*). The Jarzynski equality holds formally in the case with change of temperature,

$$
\left<\mathbf{e}^{-Y}\right>\_{0\to T} = \frac{Z\_T}{Z\_0}.\tag{13}
$$

We show a simple proof of JE for the particular dynamics in SA. Let us consider a nonequilibrium process in a finite-time schedule governed by the master equation. The left-hand side of JE is written as

$$\left\langle \mathbf{e}^{-Y} \right\rangle\_{0 \to T} = \sum\_{\{\sigma\_k\}} \prod\_{k=0}^{n-1} \left\{ \mathbf{e}^{-Y(\sigma\_{k+1}t\_k)} \mathbf{e}^{\delta t M(\sigma\_{k+1}|\sigma\_k; t\_k)} \right\} P\_{\mathbf{eq}}(\sigma\_0; t\_0). \tag{14}$$

where we use the formal solution of the master equation by the exponentiated transition matrix. We take the first product of the above equation as,

$$\sum\_{\sigma\_0} \left\{ \mathbf{e}^{-Y(\sigma\_0; t\_0)} \mathbf{e}^{\delta t M\_0(\sigma\_1 | \sigma\_0; t\_0)} \right\} P\_{\text{eq}}(\sigma\_0; t\_0)$$

$$= P\_{\text{eq}}(\sigma\_1; t\_1) \frac{Z\_{t\_1}}{Z\_0}. \tag{15}$$

Repetition of the above manipulation in Eq. (14) yields the quantity in the right-hand side of JE as,

$$\sum\_{\sigma\_{\mathfrak{u}}} P\_{\mathbf{eq}}(\sigma\_{\mathfrak{u}}; t\_{\mathfrak{u}}) \prod\_{k=0}^{n-1} \frac{Z\_{t\_{k+1}}}{Z\_{t\_k}} = \frac{Z\_T}{Z\_0}. \tag{16}$$

#### **2.2. Fluctuation theorem**

4 Will-be-set-by-IN-TECH

We introduce a couple of theories in nonequilibrium statistical physics in short, before we show the actual application. They provide the supplement to make the protocol of SA faster.

Among several recent developments in nonequilibrium statistical mechanics, we take JE as an attempt to improve the performance of SA. The Jarzynski equality relates quantities at two different thermal equilibrium states with those of nonequilibrium processes from *t* = 0 to

where the partition functions appearing in the ratio on the right-hand side are for the initial (*t* = 0) and final Hamiltonians (*t* = *T*). The quantity on the right-hand side can be represented by the exponentiated difference of the free energy exp(−*β*Δ*F*) between the initial and final conditions. The brackets on the left-hand side express the nonequilibrium average over all the instantaneous realizations of the degrees of freedom, for instance spin configurations, during

<sup>0</sup>→*<sup>T</sup>* <sup>=</sup> *ZT Z*0

e*δtM*(*σk*+1|*σk*;*tk*)

It implies that the observations of the nonequilibrium behavior can estimate the equilibrium quantity represented by the partition functions, that is the free energy. This equality is regarded as a generalization of the well-known inequality, the second law of thermodynamics �*W*�0→*<sup>T</sup>* ≥ Δ*F*, which can be reproduced by the Jensen inequality. One of the important features is that JE holds independently of the pre-determined schedule of the nonequilibrium

In order to consider the improvement of SA, let us apply the nonequilibrium process with change of the temperature. We then have to employ the pseudo work instead of the ordinary

where we use discrete time expressions as *t*<sup>0</sup> = 0 and *tn* = *T* for simplicity and we assume that the instantaneous energy does not depend on the time as *E*(*σ*). The Jarzynski equality

We show a simple proof of JE for the particular dynamics in SA. Let us consider a nonequilibrium process in a finite-time schedule governed by the master equation. The

e−*Y*(*σk*<sup>+</sup>1;*tk*)

<sup>0</sup>→*<sup>T</sup>* <sup>=</sup> *ZT Z*0

e*δtM*(*σk*+1|*σk*;*tk*)

 e−*<sup>Y</sup>* 

*n*−1 ∏ *k*=0

*Y*(*σ*; *tk*) = (*βk*+<sup>1</sup> − *βk*) *E*(*σ*), (12)

, (10)

*P*eq(*σ*0; *t*0). (11)

. (13)

*P*eq(*σ*0; *t*0). (14)

 e−*β<sup>W</sup>*

the nonequilibrium process with the following path probability defined as

*n*−1 ∏ *k*=0

*P*0→*T*({*σt*}) =

performed work due to the energy difference as

left-hand side of JE is written as

 e−*<sup>Y</sup>* 

holds formally in the case with change of temperature,

<sup>0</sup>→*<sup>T</sup>* <sup>=</sup> <sup>∑</sup> {*σk*}

**2. Population annealing**

**2.1. Jarzynski equality**

process.

The Jarzynski equality is the most important key.

*t* = *T* connecting these two states as [10, 11]

The Jarzynski equality is a consequence from the fluctuation theorem [2–4], which relates the probability *P*0→*T*({*σt*}) with that of the inverse process *PT*→0({*σt*}) as

$$\frac{P\_{0 \to T}(\{\sigma\_l\})}{P\_{T \to 0}(\{\sigma\_l\})} \mathbf{e}^{-Y} = \frac{Z\_T}{Z\_0}.\tag{17}$$

This leads to the more generic result, for an observable *O*({*σt*}) depending on the instantaneous spin configurations,

$$\left\langle O(\{\sigma\_{\rm t}\}) \mathbf{e}^{-Y} \right\rangle\_{0 \to T} = \langle O\_{\rm r}(\{\sigma\_{\rm t}\}) \rangle\_{T \to 0} \frac{\mathbf{Z}\_{T}}{\mathbf{Z}\_{0}} \,. \tag{18}$$

where *O*<sup>r</sup> denotes the observable which depends on the backward process *T* → 0. The brackets with subscript 0 → *T* express the average with the weight *P*0→*T*({*σt*}) over possible realizations {*σt*}. For *O* = *O*<sup>r</sup> = 1, Eq. (18) reduces to JE.

If we choose an observable depending only on the final state, which is denoted as *OT*, instead of *O*({*σt*}), *O*<sup>r</sup> reads an observable at the initial state in the backward process. Then �*O*r�*T*→<sup>0</sup> equals to the ordinary thermal average at the initial equilibrium state with *β<sup>T</sup>* represented by �· · · �*β<sup>T</sup>* , and Eq. (18) leads to

$$
\left< \mathcal{O}\_T \mathbf{e}^{-Y} \right>\_{0 \to T} = \left< \mathcal{O} \right>\_{\beta \tau} \frac{Z\_T}{Z\_0}.\tag{19}
$$

By looking over the above calculations, we can understand the roll of the exponentiated pseudo work. The resultant distribution after SA is given by *P*0→*T*({*σt*}). The biased distribution with the exponentiated pseudo work, *P*0→*T*({*σt*}) exp(−*Y*) always gives the equilibrium one. In this sense, the exponentiated pseudo work plays a roll to fill the gap between the equilibrium distribution and the resultant one after SA. Therefore, if we skillfully use the exponentiated pseudo work to keep the instantaneous distribution close to the equilibrium one, we can invent the improved version of SA.

#### **2.3. Population annealing**

We introduce an improvement of SA by use of the property of JE. Let us consider to implement JE in numerical simulation. We parallelize the instantaneous spin configurations {*σt*} as {*σt*}*i*=1,···,*C*. Each of the configuration independently will be evolved by the master equation. We regard the exponentiated pseudo work on the left-hand side of JE, Eq. (14), as the weight for each realization of the configuration. By computing the pseudo work for each realization and multiplying the weight given by the exponentiated pseudo work, we simultaneously perform the stochastic dynamics governed by the master equation. At the last stage of repetition of the above procedure, we obtain the ratio of the partition function as in the right-hand side of JE. In order to calculate its value, we estimate the empirical average as, after parallel computing of the master equation,

$$\frac{1}{C} \sum\_{i=1}^{C} \exp\left(-\sum\_{k=1}^{n} Y(t\_k; \sigma)\right). \tag{20}$$

optimization problem. The advantage of PA is cutting the computational time compared to SA, since PA follows the property of JE. It means that we find the optimal solution by use of

Below, we propose an ambitious use of PA to evaluate the equilibrium property in the

At first we briefly review a useful analysis in several spin-glass systems, which provides a powerful technique to discuss the possibility of PA. Let us consider a simple model of spin glasses, the ±*J* Ising model, on an arbitrary lattice. The Hamiltonian is the same form as Eq.

The partition function, which is the most important quantity through the free energy, is

{*σi*} ∏ �*ij*�

where the product *βJ* is rewritten as *K*. Both of the above quantities depend on the specific configuration {*τij*}. In order to evaluate the physical property of spin glasses in equilibrium, we strive the difficult task to deal with the free energy depending on the non-uniform interactions. Instead of the direct manipulation, the averaged quantity over all the possible configurations of {*τij*} may be considered based on the self-averaging property as, in the

> 1 *N*

where the square bracket denotes the average over all the combinations of {*τij*} (configurational average). The self-averaging property is valid for other observables, which

Here let us define a local transformation by the simultaneous change of the interactions and

*F*(*K*; {*τij*})

*Jτijσiσj*, (21)

exp(*Kτijσiσj*). (23)

, (25)

*P*(*τij*) = *pδ*(*τij* − 1)+(1 − *p*)*δ*(*τij* + 1). (22)

Optimization by Use of Nature in Physics Beyond Classical Simulated Annealing 47

− *βF*(*K*; {*τij*}) = log *Z*(*K*; {*τij*}), (24)

*τij* → *�i�jτij* (26) *σ<sup>i</sup>* → *�iσi*. (27)

*<sup>H</sup>* = − ∑ �*ij*�

where we extract the sign of the interaction *τij* following the distribution function as

*<sup>Z</sup>*(*K*; {*τij*}) = ∑

*<sup>N</sup> <sup>F</sup>*(*K*; {*τij*}) <sup>→</sup>

can be obtained from the free energy per site like the internal energy.

1

spin variables as, by the binary variables *�<sup>i</sup>* = ±1 [18, 19]

low-temperature in spin glasses by use of the special symmetry [21].

PA faster than SA.

**2.4. Spin glass**

(1) as

defined as

large-limit *N*,

The free energy is then given by

**2.5. Gauge transformation**

While estimating the ratio of the partition functions by JE, implementation of Eq. (19) gives the thermal average of the observable through their ratio. This is the typical implementation of JE in a numerical simulation, which is called as population annealing (PA) [7, 9, 17] as depicted in Fig. 3,

**Figure 3.** Schematic picture of the process of PA by *C* = 4. The size of the circles denotes the weight given by the multiplication of the exponentiated pseudo work during PA.

We remark that, as proposed in the literatures [7, 9], we have to employ a skillful technique, resampling, to efficiently generate the relevant copies to estimate the nonequilibrium average and maintain the stability of the method. The population annealing with resampling method indeed shows outstanding performance comparable to a standard technique to equilibrate the spin-glass system known as the exchange Monte Carlo method [8]. If we successfully generate the equilibrium distribution in the low-temperature region, we efficiently find the lowest energy state, which corresponds to the optimal solution in context of the optimization problem. Therefore PA is also relevant for the improvement of SA as a solver for the optimization problem. The advantage of PA is cutting the computational time compared to SA, since PA follows the property of JE. It means that we find the optimal solution by use of PA faster than SA.

Below, we propose an ambitious use of PA to evaluate the equilibrium property in the low-temperature in spin glasses by use of the special symmetry [21].

#### **2.4. Spin glass**

6 Will-be-set-by-IN-TECH

We introduce an improvement of SA by use of the property of JE. Let us consider to implement JE in numerical simulation. We parallelize the instantaneous spin configurations {*σt*} as {*σt*}*i*=1,···,*C*. Each of the configuration independently will be evolved by the master equation. We regard the exponentiated pseudo work on the left-hand side of JE, Eq. (14), as the weight for each realization of the configuration. By computing the pseudo work for each realization and multiplying the weight given by the exponentiated pseudo work, we simultaneously perform the stochastic dynamics governed by the master equation. At the last stage of repetition of the above procedure, we obtain the ratio of the partition function as in the right-hand side of JE. In order to calculate its value, we estimate the empirical average as,

While estimating the ratio of the partition functions by JE, implementation of Eq. (19) gives the thermal average of the observable through their ratio. This is the typical implementation of JE in a numerical simulation, which is called as population annealing (PA) [7, 9, 17] as depicted

**Figure 3.** Schematic picture of the process of PA by *C* = 4. The size of the circles denotes the weight

We remark that, as proposed in the literatures [7, 9], we have to employ a skillful technique, resampling, to efficiently generate the relevant copies to estimate the nonequilibrium average and maintain the stability of the method. The population annealing with resampling method indeed shows outstanding performance comparable to a standard technique to equilibrate the spin-glass system known as the exchange Monte Carlo method [8]. If we successfully generate the equilibrium distribution in the low-temperature region, we efficiently find the lowest energy state, which corresponds to the optimal solution in context of the optimization problem. Therefore PA is also relevant for the improvement of SA as a solver for the

given by the multiplication of the exponentiated pseudo work during PA.

*Y*(*tk*; *σ*)

. (20)

**2.3. Population annealing**

in Fig. 3,

after parallel computing of the master equation,

1 *C*

*C* ∑ *i*=1 exp − *n* ∑ *k*=1 At first we briefly review a useful analysis in several spin-glass systems, which provides a powerful technique to discuss the possibility of PA. Let us consider a simple model of spin glasses, the ±*J* Ising model, on an arbitrary lattice. The Hamiltonian is the same form as Eq. (1) as

$$H = -\sum\_{\langle ij \rangle} f \pi\_{ij} \sigma\_i \sigma\_{j\prime} \tag{21}$$

where we extract the sign of the interaction *τij* following the distribution function as

$$P(\tau\_{ij}) = p\delta(\tau\_{ij} - 1) + (1 - p)\delta(\tau\_{ij} + 1). \tag{22}$$

The partition function, which is the most important quantity through the free energy, is defined as

$$Z(K; \{\tau\_{i\bar{j}}\}) = \sum\_{\{\sigma\_l\}} \prod\_{\langle i\bar{j} \rangle} \exp(K \tau\_{i\bar{j}} \sigma\_l \sigma\_{\bar{j}}).\tag{23}$$

The free energy is then given by

$$-\beta F(\mathbf{K}; \{\tau\_{i\bar{j}}\}) = \log Z(\mathbf{K}; \{\tau\_{i\bar{j}}\}),\tag{24}$$

where the product *βJ* is rewritten as *K*. Both of the above quantities depend on the specific configuration {*τij*}. In order to evaluate the physical property of spin glasses in equilibrium, we strive the difficult task to deal with the free energy depending on the non-uniform interactions. Instead of the direct manipulation, the averaged quantity over all the possible configurations of {*τij*} may be considered based on the self-averaging property as, in the large-limit *N*,

$$\frac{1}{N}F(K;\{\tau\_{i\bar{j}}\}) \to \frac{1}{N}\left[F(K;\{\tau\_{i\bar{j}}\})\right],\tag{25}$$

where the square bracket denotes the average over all the combinations of {*τij*} (configurational average). The self-averaging property is valid for other observables, which can be obtained from the free energy per site like the internal energy.

#### **2.5. Gauge transformation**

Here let us define a local transformation by the simultaneous change of the interactions and spin variables as, by the binary variables *�<sup>i</sup>* = ±1 [18, 19]

$$
\pi\_{\vec{l}\vec{j}} \to \mathfrak{e}\_{\vec{l}} \mathfrak{e}\_{\vec{l}} \pi\_{\vec{l}\vec{j}} \tag{26}
$$

$$
\sigma\_{\dot{l}} \to \mathfrak{e}\_{\dot{l}} \sigma\_{\dot{l}}.\tag{27}
$$

#### 8 Will-be-set-by-IN-TECH 48 Simulated Annealing – Advances, Applications and Hybridizations Optimization by Use of Nature in Physics Beyond Classical Simulated Annealing <sup>9</sup>

This is called as the gauge transformation. Notice that the gauge transformation does not alter the value of the physical quantity given by the double average over *τij* and *σ<sup>i</sup>* since it changes only the order of the summations. The Hamiltonian can not change its form after the gauge transformation since the right-hand side is evaluated as

$$-\sum\_{\langle ij\rangle} f \mathfrak{e}\_{\hat{l}\hat{l}} \mathfrak{e}\_{\hat{l}} \mathfrak{e}\_{\hat{l}} \mathfrak{e}\_{\hat{l}} \mathfrak{e}\_{\hat{l}} \mathfrak{e}\_{\hat{j}} \mathfrak{v}\_{\hat{j}} = H\_{\prime} \tag{28}$$

We take the summation over {*�i*} in advance of that over {*τij*} and then find the partition

where *NB* is the number of bonds. Going back to the definition (31), both of the partition

Similarly, we can evaluate the rigorous upper bound on the specific heat. The condition *Kp* = *K* confines the special subspace in which we can perform the exact analysis for spin glasses.

By use of the gauge transformation as above introduced briefly, let us consider the application of the relations (18) and (19) to spin glasses, namely PA for such a complicated system in a tricky way. We analyze JE for the spin-glass model for several interesting quantities below.

We apply Eq. (19) to a gauge-invariant quantity *G*({*τij*}) for the purpose of evaluation of the

The quantity on the left-hand side is the configurational and nonequilibrium averages of the observable *G*({*τij*}) at final stage of PA, that is after the protocol *K*<sup>0</sup> → *KT* with the factor *<sup>e</sup>*−*βW*. On the other hand, �*G*({*τij*})�*KT* on the right-hand side expresses the configurational

All the quantities in this equation are invariant under the gauge transformation. The

�*G*({*τij*})�*KT*

�*G*({*τij*})�*KT* <sup>∏</sup>�*ij*� <sup>e</sup>*Kpτij�i�<sup>j</sup>* (2 cosh *Kp*)*NB*

�*G*({*τij*})�*KTZ*(*Kp*; {*τij*}) 2*N*(2 cosh *Kp*)*NB*

equilibrium quantity in spin glasses. The configurational average for Eq. (19) yields

*Kp* =  {*σi*}

∑ {*τij*}

*d dK* exp

= −*NBJ* tanh *K*. (37)

Optimization by Use of Nature in Physics Beyond Classical Simulated Annealing 49

*Z*(*Kp*; {*τij*})

(2 cosh *Kp*)*NB* × �*H*�*K*. (36)

*Kτijσiσ<sup>j</sup>*

*Z*(*KT*; {*τij*}) *Z*(*K*0; {*τij*})

*Kp*

*Z*(*KT*; {*τij*}) *Z*(*K*0; {*τij*})

*Z*(*KT*; {*τij*}) *Z*(*K*0; {*τij*})

. (38)

. (39)

. (40)

<sup>2</sup>*<sup>N</sup>* ∑ {*τij*}

functions on the denominator and numerator can be cancelled when *Kp* = *K* as

<sup>2</sup>*N*(2 cosh *Kp*)*NB* ∑

[�*H*�*K*]*Kp* <sup>=</sup> <sup>1</sup>

[�*H*�*K*]*<sup>K</sup>* <sup>=</sup> <sup>−</sup>*<sup>J</sup>*

This subspace is called as the Nishimori line (NL) [18, 19].

**2.6. Jarzynski equality for spin glasses**

*2.6.1. Gauge-invariant quantities like internal energy*

*GT*({*τij*})e−*<sup>Y</sup>*

*K*0→*KT*

and thermal averages of the equilibrium state for the final Hamiltonian.

The gauge transformation *σ<sup>i</sup>* → *σi�i*, *τij* → *τij�i�<sup>j</sup>* (∀*i*, *j*) leads us to

*Kp*

*Kp*

= ∑ {*τij*}

> = ∑ {*τij*}

*Z*(*KT*; {*τij*}) *Z*(*K*0; {*τij*})

*Z*(*KT*; {*τij*}) *Z*(*K*0; {*τij*})

summation over {*�i*} and division by 2*<sup>N</sup>* gives

�*G*({*τij*})�*KT*

�*G*({*τij*})�*KT*

function with *Kp* instead of *K*.

As this case, if the physical quantity is invariant under the gauge transformation (gauge invariant), we can evaluate its exact value even for finite-dimensional spin glasses. The key point of the analyses by the gauge transformation is on the form of the distribution function. Before performing the gauge transformation, the distribution function can take the following form as

$$P(\tau\_{lj}) = \frac{\mathbf{e}^{K\_p \tau\_{lj}}}{2 \cosh K\_p} \,\tag{29}$$

where exp(−2*Kp*)=(1 − *p*)/*p*. The gauge transformation changes this into

$$P(\tau\_{l\bar{j}}) = \frac{\mathbf{e}^{K\_p \tau\_{l\bar{j}} \varepsilon\_l \varepsilon\_{l\bar{j}}}}{2 \cosh K\_p}.\tag{30}$$

Let us evaluate the internal energy by aid of the gauge transformation here. The thermal average of the Hamiltonian is given by

$$\langle H \rangle\_{\mathcal{K}} = \sum\_{\{\sigma\_i\}} \frac{1}{\mathbb{Z}(\mathbf{K}; \{\tau\_{ij}\})} H \prod\_{\langle ij \rangle} \exp \left( \mathbf{K} \tau\_{ij} \sigma\_i \sigma\_j \right) \tag{31}$$

$$\dot{\lambda} = -J \frac{d}{dK} \log Z(K; \{\tau\_{\dot{l}\dot{l}}\}). \tag{32}$$

We can use the self-averaging property here and thus take the configurational average as

$$\mathbb{E}\left[\langle H\rangle\_{\mathbb{K}}\right]\_{\mathbb{K}\_p} = \sum\_{\{\tau\_{ij}\}} \prod\_{\langle ij\rangle} \frac{\exp(K\_p \tau\_{ij})}{2 \cosh K\_p} \times \langle H \rangle\_{\mathbb{K}}.\tag{33}$$

where [··· ]*Kp* denotes the configurational average by the distribution function (29) with *Kp*. Then we perform the gauge transformation, which does not change the value of the internal energy due to gauge invariance,

$$\mathbb{E}\left[\langle H\rangle\_{K}\right]\_{K\_{p}} = \sum\_{\{\pi\_{\emptyset}\}} \prod\_{\langle ij\rangle} \frac{\exp(K\_{p}\pi\_{\emptyset\rangle}\varepsilon\_{i}\varepsilon\_{j})}{2\cosh K\_{p}} \times \langle H\rangle\_{K}.\tag{34}$$

Therefore we here take the summation over all the possible configurations of {*σi*} and divide it by 2*<sup>N</sup>* (the number of configurations) as

$$\left[\left\_{K\_{p}}\right]\_{K\_{p}} = \frac{1}{2^{N}} \sum\_{\{\varepsilon\_{l}\}} \sum\_{\{\pi\_{l}\}} \prod\_{\{ij\}} \frac{\exp(K\_{p}\tau\_{ij}\varepsilon\_{i}\varepsilon\_{j})}{2\cosh K\_{p}} \times \langle H \rangle\_{K}.\tag{35}$$

We take the summation over {*�i*} in advance of that over {*τij*} and then find the partition function with *Kp* instead of *K*.

$$\mathbb{E}\left[\langle H\rangle\_{\mathbb{K}}\right]\_{\mathbb{K}\_p} = \frac{1}{\mathfrak{D}^N} \sum\_{\{\tau\_{\vec{l}}\}} \frac{Z(K\_p; \{\tau\_{\vec{l}j}\})}{(\mathfrak{D}\cosh K\_p)^{N\_{\mathfrak{B}}}} \times \langle H\rangle\_{\mathbb{K}}.\tag{36}$$

where *NB* is the number of bonds. Going back to the definition (31), both of the partition functions on the denominator and numerator can be cancelled when *Kp* = *K* as

$$\begin{split} \left[ \langle H \rangle\_{\mathcal{K}} \right]\_{\mathcal{K}} &= \frac{-J}{2^{N} (2 \cosh \mathcal{K}\_{p})^{N\_{B}}} \sum\_{\{\sigma\_{i}\}} \sum\_{\{\tau\_{ij}\}} \frac{d}{d\mathcal{K}} \exp \left( \mathcal{K} \tau\_{ij} \sigma\_{i} \sigma\_{j} \right) \\ &= -N\_{B} J \tanh \mathcal{K}. \end{split} \tag{37}$$

Similarly, we can evaluate the rigorous upper bound on the specific heat. The condition *Kp* = *K* confines the special subspace in which we can perform the exact analysis for spin glasses. This subspace is called as the Nishimori line (NL) [18, 19].

#### **2.6. Jarzynski equality for spin glasses**

8 Will-be-set-by-IN-TECH

This is called as the gauge transformation. Notice that the gauge transformation does not alter the value of the physical quantity given by the double average over *τij* and *σ<sup>i</sup>* since it changes only the order of the summations. The Hamiltonian can not change its form after the gauge

As this case, if the physical quantity is invariant under the gauge transformation (gauge invariant), we can evaluate its exact value even for finite-dimensional spin glasses. The key point of the analyses by the gauge transformation is on the form of the distribution function. Before performing the gauge transformation, the distribution function can take the following

*<sup>P</sup>*(*τij*) = <sup>e</sup>*Kpτij*

*<sup>P</sup>*(*τij*) = <sup>e</sup>*Kpτij�i�<sup>j</sup>*

Let us evaluate the internal energy by aid of the gauge transformation here. The thermal

1 *Z*(*K*; {*τij*})

We can use the self-averaging property here and thus take the configurational average as

where [··· ]*Kp* denotes the configurational average by the distribution function (29) with *Kp*. Then we perform the gauge transformation, which does not change the value of the internal

Therefore we here take the summation over all the possible configurations of {*σi*} and divide

∑ {*τij*} ∏ �*ij*�

{*τij*} ∏ �*ij*�

{*τij*} ∏ �*ij*�

<sup>2</sup>*<sup>N</sup>* ∑ {*�i*} 2 cosh *Kp*

2 cosh *Kp*

*H* ∏ �*ij*� exp *Kτijσiσ<sup>j</sup>*

exp(*Kpτij*) 2 cosh *Kp*

exp(*Kpτij�i�j*) 2 cosh *Kp*

> exp(*Kpτij�i�j*) 2 cosh *Kp*

*Jτij�i�j�iσi�jσ<sup>j</sup>* = *H*, (28)

, (29)

. (30)

× �*H*�*K*. (33)

× �*H*�*K*. (34)

× �*H*�*K*. (35)

(31)

*dK* log *<sup>Z</sup>*(*K*; {*τij*}). (32)

transformation since the right-hand side is evaluated as

average of the Hamiltonian is given by

energy due to gauge invariance,

it by 2*<sup>N</sup>* (the number of configurations) as

�*H*�*<sup>K</sup>* = ∑

{*σi*}

<sup>=</sup> <sup>−</sup>*<sup>J</sup> <sup>d</sup>*

[�*H*�*K*]*Kp* = ∑

[�*H*�*K*]*Kp* = ∑

[�*H*�*K*]*Kp* <sup>=</sup> <sup>1</sup>

form as

− ∑ �*ij*�

where exp(−2*Kp*)=(1 − *p*)/*p*. The gauge transformation changes this into

By use of the gauge transformation as above introduced briefly, let us consider the application of the relations (18) and (19) to spin glasses, namely PA for such a complicated system in a tricky way. We analyze JE for the spin-glass model for several interesting quantities below.

#### *2.6.1. Gauge-invariant quantities like internal energy*

We apply Eq. (19) to a gauge-invariant quantity *G*({*τij*}) for the purpose of evaluation of the equilibrium quantity in spin glasses. The configurational average for Eq. (19) yields

$$\left[ \left< G\_T(\{\tau\_{i\bar{l}}\}) \mathbf{e}^{-Y} \right>\_{K\_0 \to K\_T} \right]\_{K\_p} = \left[ \langle G(\{\tau\_{i\bar{l}}\}) \rangle\_{K\_T} \frac{Z(K\_T; \{\tau\_{i\bar{l}}\})}{Z(K\_0; \{\tau\_{i\bar{l}}\})} \right]\_{K\_p} \tag{38}$$

The quantity on the left-hand side is the configurational and nonequilibrium averages of the observable *G*({*τij*}) at final stage of PA, that is after the protocol *K*<sup>0</sup> → *KT* with the factor *<sup>e</sup>*−*βW*. On the other hand, �*G*({*τij*})�*KT* on the right-hand side expresses the configurational and thermal averages of the equilibrium state for the final Hamiltonian.

The gauge transformation *σ<sup>i</sup>* → *σi�i*, *τij* → *τij�i�<sup>j</sup>* (∀*i*, *j*) leads us to

$$\left[ \langle G(\{\tau\_{\vec{\imath}j}\}) \rangle\_{\mathcal{K}\tau} \frac{Z(\mathcal{K}\_{\mathcal{T}}; \{\tau\_{\vec{\imath}j}\})}{Z(\mathcal{K}\_{0}; \{\tau\_{\vec{\imath}j}\})} \right]\_{\mathcal{K}\_{p}} = \sum\_{\{\tau\_{\vec{\imath}j}\}} \frac{\langle G(\{\tau\_{\vec{\imath}j}\}) \rangle\_{\mathcal{K}\tau} \prod\_{\{\vec{\imath}j\}} \mathbf{e}^{\mathbf{K}\_{p}\tau\_{\vec{\imath}}\varepsilon\_{\vec{\imath}}} \frac{Z(\mathcal{K}\_{\mathcal{T}}; \{\tau\_{\vec{\imath}j}\})}{Z(\mathcal{K}\_{0}; \{\tau\_{\vec{\imath}j}\})} . \tag{39}$$

All the quantities in this equation are invariant under the gauge transformation. The summation over {*�i*} and division by 2*<sup>N</sup>* gives

$$\left[ \langle G(\{\tau\_{\vec{l}}\}) \rangle\_{\mathbb{K}\_{\mathsf{T}}} \frac{Z(\mathsf{K}\_{\mathsf{T}}; \{\tau\_{\vec{l}j}\})}{Z(\mathsf{K}\_{\mathsf{0}}; \{\tau\_{\vec{l}j}\})} \right]\_{\mathbb{K}\_{p}} = \sum\_{\{\tau\_{\vec{l}}\}} \frac{\langle G(\{\tau\_{\vec{l}j}\}) \rangle\_{\mathbb{K}\_{\mathsf{T}}} Z(\mathsf{K}\_{p}; \{\tau\_{\vec{l}j}\}) }{2^{N} (2 \cosh \mathsf{K}\_{p})^{N\_{\mathsf{B}}}} \frac{Z(\mathsf{K}\_{\mathsf{T}}; \{\tau\_{\vec{l}j}\})}{Z(\mathsf{K}\_{\mathsf{0}}; \{\tau\_{\vec{l}j}\})}. \tag{40}$$

#### 10 Will-be-set-by-IN-TECH 50 Simulated Annealing – Advances, Applications and Hybridizations Optimization by Use of Nature in Physics Beyond Classical Simulated Annealing <sup>11</sup>

On the other hand, let us evaluate the quantity �*G*({*τij*})�*KT Kp* . Similarly to the above calculation, the following identity can be obtained by the gauge transformation,

$$\mathbb{E}\left[\langle G(\{\tau\_{\vec{l}\vec{j}}\})\rangle\_{\mathbb{K}\_{\mathbb{T}}}\right]\_{\mathbb{K}\_{\mathbb{P}}} = \sum\_{\{\tau\_{\vec{l}\vec{j}}\}} \frac{\langle G(\{\tau\_{\vec{l}\vec{j}}\})\rangle\_{\mathbb{K}\_{\mathbb{T}}} Z(K\_{p};\{\tau\_{\vec{l}\vec{j}}\})}{2^{N} (2 \cosh K\_{p})^{N\_{\mathbb{B}}}}.\tag{41}$$

By Setting *Kp* = *K*<sup>0</sup> in Eq. (40) and *Kp* = *KT* in the above equation, we reach the following nonequilibrium relation,

$$\left[ \langle \mathbf{G}\_T(\{\tau\_{i\bar{j}}\}) \mathbf{e}^{-Y} \rangle\_{\mathbf{K}\_0 \rightarrow \mathbf{K}\_\tau} \right]\_{\mathbf{K}\_0} = \left[ \langle G(\{\tau\_{i\bar{j}}\}) \rangle\_{\mathbf{K}\_\tau} \right]\_{\mathbf{K}\_\tau} \left( \frac{\cosh K\_T}{\cosh K\_0} \right)^{N\_\mathcal{B}}.\tag{42}$$

If we set *GT*({*τij*}) = 1 in the resultant equation, the Jarzynski equality for spin glass is obtained,

$$
\left[ \left< \mathbf{e}^{-Y} \right>\_{K\_{\mathbf{0}} \to K\_{T}} \right]\_{K\_{\mathbf{0}}} = \left( \frac{\cosh K\_{T}}{\cosh K\_{\mathbf{0}}} \right)^{N\_{\mathbf{0}}}.\tag{43}
$$

O

*Z*(*KT*; {*τij*}) *Z*(*K*0; {*τij*})

�*σi*�*KT �<sup>i</sup>* ∏

�*ij*�

<sup>2</sup>*N*(2 cosh *Kp*)*NB* �*σi*�*K*<sup>0</sup> �*�i*�*Kp* . (50)

*NB*

cosh *KT* cosh *K*<sup>0</sup> *Kp*

e*Kpτij�i�<sup>j</sup>* 2 cosh *Kp*

*Z*(*KT*; {*τij*}) *Z*(*K*0; {*τij*})

. (47)

. (48)

. (49)

. (51)

*K*

*T*

*K*

0

*K*

Optimization by Use of Nature in Physics Beyond Classical Simulated Annealing 51

Para

**Figure 4.** The nonequilibrium process as in SA/PA on the left-hand side as in Eq. (45) drawn as an arrow on the phase diagram. The corresponding equilibrium state as on the right-hand side of Eq. (45) is at the point on NL. The solid curves are the phase boundaries and the dashed line of 45o represents NL.

Therefore it is important to observe the behavior of the first momentum of spin variable in equilibrium. For this purpose, we choose *σi*(*T*) for *O* in Eq. (19) and consider the possibility

�*σi*�*KT*

*Z*(*KT*; {*τij*}) *Z*(*K*0; {*τij*})

We again sum both sides of this equation over all the possible configurations of {*�i*} and

*Z*(*Kp*; {*τij*})

*Z*(*Kp*; {*τij*})

= [�*σi*�*K*<sup>0</sup> ]

*KT*

<sup>2</sup>*N*(2 cosh *Kp*)*NB* �*σi*�*KT* �*�i*�*Kp*

Ferro

SG

*Kp*

�*σi*�*KT*

�*σi*�*KT*

divide the obtained quantity by 2*<sup>N</sup>* to find

�*σi*�*KT*

*Z*(*KT*; {*τij*}) *Z*(*K*0; {*τij*})

[�*σi*�*K*<sup>0</sup> ]

*σi*(*T*)e−*<sup>Y</sup>*

*Z*(*KT*; {*τij*}) *Z*(*K*0; {*τij*})

*Kp*

of the application of PA. After the configurational average, we obtain

Gauge transformation for the right-hand side in this equation yields

*Kp*

= ∑ {*τij*}

The following relation can also be obtained in a similar manipulation,

*Kp* = ∑ {*τij*}

*Z*(*KT*; {*τij*}) *Z*(*K*0; {*τij*})

By setting *Kp* = *K*<sup>0</sup> in Eq. (49) and *Kp* = *KT* in Eq. (50), we reach a relation

*K*0

*K*0→*KT*

= ∑ {*τij*}

*Kp* = 

Equation (43) leads to the lower bound on the pseudo work, using Jensen's inequality for the average of e−*Y*,

$$\left[ \langle Y \rangle\_{\mathcal{K}\_0 \to \mathcal{K}\_\Gamma} \right]\_{\mathcal{K}\_0} \ge -N\_B \log \left( \frac{\cosh K\_T}{\cosh K\_0} \right). \tag{44}$$

By substituting *GT*({*τij*}) = *H* into Eq. (42), we obtain

$$\left[ \left< H\mathbf{e}^{-\beta W} \right>\_{K\_{\mathbb{B}} \to K\_{\mathbb{T}}} \right]\_{K\_{\mathbb{B}}} = \left[ \langle H \rangle\_{K\_{\mathbb{T}}} \right]\_{K\_{\mathbb{T}}} \left( \frac{\cosh K\_{\mathbb{T}}}{\cosh K\_{0}} \right)^{N\_{\mathbb{B}}}.\tag{45}$$

This equation shows that the internal energy after the cooling as in SA or heating process starting from a temperature on NL, which is in the present case *K*<sup>0</sup> = *Kp*, is proportional to the internal energy in equilibrium on NL corresponding to the final temperature *KT* = *Kp* as in Fig. 4. We here assume to perform PA and take an average over all results after many repetitions. The nonequilibrium process starts from NL (1/*K*0, 1/*K*0) and ends at the point away from NL (1/*K*0, 1/*KT*). Notice that the ordinary procedure in PA gives the estimation of the equilibrium quantity at the last condition as (1/*K*0, 1/*KT*). However the corresponding internal energy is at the different point but on NL as (1/*KT*, 1/*KT*). It means that we can obtain the equilibrium quantities in the different amount of the randomness from the initial condition through the configurational average of the results from PA.

#### *2.6.2. Gauge-non-invariant quantities*

In statistical physics, it is important to detect the order of the instantaneous spin configuration in the system. For instance, as in Fig. 4, there are several phases, ferromagnetic, paramagnetic and spin-glass ones, involved in the spin-glass model. They have the characteristic quantities to distinguish themselves, termed as the order parameter. The order parameter to identify the phase boundary between the ferromagnetic and paramagnetic phases is the magnetization defined as

$$m = \frac{1}{N} \sum\_{i=1}^{N} \sigma\_i \tag{46}$$

**Figure 4.** The nonequilibrium process as in SA/PA on the left-hand side as in Eq. (45) drawn as an arrow on the phase diagram. The corresponding equilibrium state as on the right-hand side of Eq. (45) is at the point on NL. The solid curves are the phase boundaries and the dashed line of 45o represents NL.

Therefore it is important to observe the behavior of the first momentum of spin variable in equilibrium. For this purpose, we choose *σi*(*T*) for *O* in Eq. (19) and consider the possibility of the application of PA. After the configurational average, we obtain

$$\left[ \left< \sigma\_{\dot{\imath}}(T) \mathbf{e}^{-Y} \right>\_{\mathcal{K}\_{0} \to \mathcal{K}\_{T}} \right]\_{\mathcal{K}\_{p}} = \left[ \left< \sigma\_{\dot{\imath}} \right>\_{\mathcal{K}\_{T}} \frac{Z(\mathcal{K}\_{T}; \{\tau\_{\dot{\imath}j}\})}{Z(\mathcal{K}\_{0}; \{\tau\_{\dot{\imath}j}\})} \right]\_{\mathcal{K}\_{p}}.\tag{47}$$

Gauge transformation for the right-hand side in this equation yields

10 Will-be-set-by-IN-TECH

By Setting *Kp* = *K*<sup>0</sup> in Eq. (40) and *Kp* = *KT* in the above equation, we reach the following

If we set *GT*({*τij*}) = 1 in the resultant equation, the Jarzynski equality for spin glass is

Equation (43) leads to the lower bound on the pseudo work, using Jensen's inequality for the

*<sup>K</sup>*<sup>0</sup> ≥ −*NB* log

= [�*H*�*KT* ]

This equation shows that the internal energy after the cooling as in SA or heating process starting from a temperature on NL, which is in the present case *K*<sup>0</sup> = *Kp*, is proportional to the internal energy in equilibrium on NL corresponding to the final temperature *KT* = *Kp* as in Fig. 4. We here assume to perform PA and take an average over all results after many repetitions. The nonequilibrium process starts from NL (1/*K*0, 1/*K*0) and ends at the point away from NL (1/*K*0, 1/*KT*). Notice that the ordinary procedure in PA gives the estimation of the equilibrium quantity at the last condition as (1/*K*0, 1/*KT*). However the corresponding internal energy is at the different point but on NL as (1/*KT*, 1/*KT*). It means that we can obtain the equilibrium quantities in the different amount of the randomness from the initial

In statistical physics, it is important to detect the order of the instantaneous spin configuration in the system. For instance, as in Fig. 4, there are several phases, ferromagnetic, paramagnetic and spin-glass ones, involved in the spin-glass model. They have the characteristic quantities to distinguish themselves, termed as the order parameter. The order parameter to identify the phase boundary between the ferromagnetic and paramagnetic phases is the magnetization

> *N* ∑ *i*=1

*<sup>m</sup>* <sup>=</sup> <sup>1</sup> *N*

 *K*0 =

calculation, the following identity can be obtained by the gauge transformation,

= ∑ {*τij*}

 *K*0 = 

*K*0→*KT*

 *Kp* �*G*({*τij*})�*KT*

cosh *KT* cosh *K*<sup>0</sup>

*KT*

�*G*({*τij*})�*KT*

�*G*({*τij*})�*KTZ*(*Kp*; {*τij*})

 *KT*

cosh *KT* cosh *K*<sup>0</sup>

> cosh *KT* cosh *K*<sup>0</sup>

*NB*

*NB*

*σ<sup>i</sup>* (46)

 *Kp*

cosh *KT* cosh *K*<sup>0</sup>

<sup>2</sup>*N*(2 cosh *Kp*)*NB* . (41)

*NB*

. (43)

. (44)

. (45)

. Similarly to the above

. (42)

On the other hand, let us evaluate the quantity

�*G*({*τij*})�*KT*

�*GT*({*τij*})e−*Y*�*K*0→*KT*

 e−*<sup>Y</sup>* 

By substituting *GT*({*τij*}) = *H* into Eq. (42), we obtain

*H*e−*βW*

*2.6.2. Gauge-non-invariant quantities*

defined as

[�*Y*�*K*0→*KT* ]

*K*0→*KT*

condition through the configurational average of the results from PA.

 *K*0

nonequilibrium relation,

obtained,

average of e−*Y*,

$$
\left[ \langle \sigma\_{\overline{l}} \rangle\_{\mathbb{K}\_{\mathbb{T}}} \frac{Z(\mathsf{K}\_{\mathcal{T}}; \{\tau\_{\overline{l}j}\})}{Z(\mathsf{K}\_{\mathcal{U}}; \{\tau\_{\overline{l}j}\})} \right]\_{\mathbb{K}\_{p}} = \sum\_{\{\mathsf{r}\_{\overline{l}}\}} \frac{Z(\mathsf{K}\_{\mathcal{T}}; \{\tau\_{\overline{l}j}\})}{Z(\mathsf{K}\_{\mathcal{U}}; \{\tau\_{\overline{l}j}\})} \langle \sigma\_{\overline{l}} \rangle\_{\mathbb{K}\_{\mathbb{T}}} \mathbf{e}\_{\overline{l}} \prod\_{\langle \overline{l}j \rangle} \frac{\mathbf{e}^{\mathsf{K}\_{p} \tau\_{\overline{l}} \varepsilon\_{\overline{l}} \varepsilon\_{\overline{l}}}}{2 \cosh \mathsf{K}\_{p}}.\tag{48}
$$

We again sum both sides of this equation over all the possible configurations of {*�i*} and divide the obtained quantity by 2*<sup>N</sup>* to find

$$\left[ \langle \sigma\_{\mathrm{i}} \rangle\_{\mathrm{K}\_{\mathrm{T}}} \frac{Z(\mathrm{K}\_{\mathrm{T}}; \{\tau\_{\mathrm{i}\bar{\jmath}}\})}{Z(\mathrm{K}\_{0}; \{\tau\_{\mathrm{i}\bar{\jmath}}\})} \right]\_{\mathrm{K}\_{p}} = \sum\_{\{\tau\_{\mathrm{i}}\}} \frac{Z(\mathrm{K}\_{p}; \{\tau\_{\mathrm{i}\bar{\jmath}}\})}{2^{N} (2 \cosh \mathrm{K}\_{p})^{N\_{\mathrm{B}}}} \langle \sigma\_{\mathrm{i}} \rangle\_{\mathrm{K}\_{\mathrm{T}}} \langle \varepsilon\_{\mathrm{i}} \rangle\_{\mathrm{K}\_{p}} \frac{Z(\mathrm{K}\_{\mathrm{T}}; \{\tau\_{\mathrm{i}\bar{\jmath}}\})}{Z(\mathrm{K}\_{0}; \{\tau\_{\mathrm{i}\bar{\jmath}}\})}. \tag{49}$$

The following relation can also be obtained in a similar manipulation,

$$\left[ \left< \sigma\_{l} \right>\_{\mathcal{K}\_{p}} \right]\_{\mathcal{K}\_{p}} = \sum\_{\{\tau\_{l}\}} \frac{Z(\mathcal{K}\_{p}; \{\tau\_{l\bar{l}}\})}{2^{N} (2 \cosh \mathcal{K}\_{p})^{N\_{\mathcal{B}}}} \langle \sigma\_{l} \rangle\_{\mathcal{K}\_{l}} \langle \epsilon\_{l} \rangle\_{\mathcal{K}\_{p}}.\tag{50}$$

By setting *Kp* = *K*<sup>0</sup> in Eq. (49) and *Kp* = *KT* in Eq. (50), we reach a relation

$$\left[ \langle \sigma\_{\dot{i}} \rangle\_{\mathcal{K}\_{\mathcal{T}}} \frac{Z(\mathcal{K}\_{\mathcal{T}}; \{\tau\_{\dot{i}\dot{j}}\})}{Z(\mathcal{K}\_{0}; \{\tau\_{\dot{i}j}\})} \right]\_{\mathcal{K}\_{0}} = \left[ \langle \sigma\_{\dot{i}} \rangle\_{\mathcal{K}\_{0}} \right]\_{\mathcal{K}\_{\mathcal{T}}} \left( \frac{\cosh \mathcal{K}\_{\mathcal{T}}}{\cosh \mathcal{K}\_{0}} \right)^{N\_{\mathcal{B}}}.\tag{51}$$

**Figure 5.** Equations (52) and (53) relate equilibrium quantities at the point (1/*K*0, 1/*KT*) (the lower-left dot) with other physical quantities estimated during the nonequilibrium process shown in the arrow. The lower-left dot is in the spin glass phase whereas the corresponding arrow is in the ferromagnetic phase.

As a result, we obtain a nonequilibrium relation ,

$$\left[ \langle \sigma\_i(T) \mathbf{e}^{-Y} \rangle\_{K\_0 \to K\_T} \right]\_{K\_0} = \left[ \langle \sigma\_i \rangle\_{K\_0} \right]\_{K\_T} \left( \frac{\cosh K\_T}{\cosh K\_0} \right)^{N\_B}.\tag{52}$$



that of the empirical average of the estimation given by PA.

associated with this procedure with PA for spin glasses.

diff

diff

0.5 1 1.5 2 2.5 3 3.5 4 4.5 5

Optimization by Use of Nature in Physics Beyond Classical Simulated Annealing 53

T

0.5 1 1.5 2 2.5 3 3.5 4 4.5 5

T

**Figure 7.** Difference between the exact values and the empirical averages over 100 realizations. The cross points give the deviation of the empirical averages of the exact values of ratios of the partition functions from the exact value on the right-hand side of Eq. (43), while the tilted-cross ones represent

do not correctly reproduce the quantity of the right-hand side. We describe the test results for estimation of the ratio of the partition functions as in Eq. (13) by PA for 100 realizations of {*τij*} in Fig. 6. We perform PA for the ±*J* Ising model on the square lattice with a linear size *L* = 6, and Monte-Carlo step *R* = 1000. The number of copies is *C* = 100. The population annealing can correctly estimate the ratio of the partition functions for each realization, but their simple average does not coincide with the quantity on the right-hand side of Eq. (43) as in Fig. 7. Both of the results are away from the exact solutions due to the sample-to-sample fluctuation and show nontrivial behavior depending on the linear size. These facts imply lack of self-averaging property. Therefore, if we exploit all the above results given by the configurational average of the exponentiated pseudo work, we have to overcome this violation due to lack of the self-averaging property. This is one of the remaining problem

**Figure 6.** Superimposed plots of Differences between the exact value and the estimation given by PA. The horizontal axis denotes the instantaneous temperature during PA. The vertical axis represents the

difference from the exact values given by the transfer matrix method in advance.

The same method yields another relation for the correlation functions to similarly measure the magnitude of order in the system

$$\left[ \langle \sigma\_0(T)\sigma\_r(T)\mathbf{e}^{-Y} \rangle\_{K\_0 \to K\_T} \right]\_{K\_0} = \left[ \langle \sigma\_0 \sigma\_r \rangle\_{K\_0} \right]\_{K\_T} \left( \frac{\cosh K\_T}{\cosh K\_0} \right)^{N\_B}.\tag{53}$$

The obtained relations (52) and (53) relate the equilibrium physical quantities away from NL (the right-hand sides) with other quantities measured during the nonequilibrium process from a point on NL to another point away from NL (the left-hand sides) as depicted in Fig. 5. The spin-glass system in the low-temperature region exhibits the extremely slow relaxation toward equilibrium. This feature hampers to observe the equilibrium behavior of spin glasses. However our results imply that the configurational average of PA would overcome the difficulty. One may attempts the heating process from NL in order to evaluate the low-temperature property through Eqs. (52) and (53) as depicted in Fig. 5. The Jarzynski equality holds irrespectively of the schedule to control the external field. It means that we can investigate the low-temperature behavior for spin glasses without suffering from critical slowing down. Unfortunately, however, the exponentiated pseudo work does not hold the self-averaging property. It means that the sample-to-sample fluctuation between different configurations of {*τij*} remains to be relevant even in a large-*N* system. Therefore, if we estimate the empirical average of realizations of {*τij*} following the obtained equalities, we

12 Will-be-set-by-IN-TECH

Ferro

**Figure 5.** Equations (52) and (53) relate equilibrium quantities at the point (1/*K*0, 1/*KT*) (the lower-left dot) with other physical quantities estimated during the nonequilibrium process shown in the arrow. The lower-left dot is in the spin glass phase whereas the corresponding arrow is in the ferromagnetic phase.

= [�*σi*�*K*<sup>0</sup> ]

= [�*σ*0*σr*�*K*<sup>0</sup> ]

The same method yields another relation for the correlation functions to similarly measure the

The obtained relations (52) and (53) relate the equilibrium physical quantities away from NL (the right-hand sides) with other quantities measured during the nonequilibrium process from a point on NL to another point away from NL (the left-hand sides) as depicted in Fig. 5. The spin-glass system in the low-temperature region exhibits the extremely slow relaxation toward equilibrium. This feature hampers to observe the equilibrium behavior of spin glasses. However our results imply that the configurational average of PA would overcome the difficulty. One may attempts the heating process from NL in order to evaluate the low-temperature property through Eqs. (52) and (53) as depicted in Fig. 5. The Jarzynski equality holds irrespectively of the schedule to control the external field. It means that we can investigate the low-temperature behavior for spin glasses without suffering from critical slowing down. Unfortunately, however, the exponentiated pseudo work does not hold the self-averaging property. It means that the sample-to-sample fluctuation between different configurations of {*τij*} remains to be relevant even in a large-*N* system. Therefore, if we estimate the empirical average of realizations of {*τij*} following the obtained equalities, we

*KT*

*KT*

cosh *KT* cosh *K*<sup>0</sup>

> cosh *KT* cosh *K*<sup>0</sup>

 *K*0

> *K*0

Para

SG

*Kp*

As a result, we obtain a nonequilibrium relation ,

�*σi*(*T*)e−*Y*�*K*0→*KT*

�*σ*0(*T*)*σr*(*T*)e−*Y*�*K*0→*KT*

magnitude of order in the system

O

*K*

*NB*

*NB*

. (52)

. (53)

*T*

*K*

0

*K*

**Figure 6.** Superimposed plots of Differences between the exact value and the estimation given by PA. The horizontal axis denotes the instantaneous temperature during PA. The vertical axis represents the difference from the exact values given by the transfer matrix method in advance.

**Figure 7.** Difference between the exact values and the empirical averages over 100 realizations. The cross points give the deviation of the empirical averages of the exact values of ratios of the partition functions from the exact value on the right-hand side of Eq. (43), while the tilted-cross ones represent that of the empirical average of the estimation given by PA.

do not correctly reproduce the quantity of the right-hand side. We describe the test results for estimation of the ratio of the partition functions as in Eq. (13) by PA for 100 realizations of {*τij*} in Fig. 6. We perform PA for the ±*J* Ising model on the square lattice with a linear size *L* = 6, and Monte-Carlo step *R* = 1000. The number of copies is *C* = 100. The population annealing can correctly estimate the ratio of the partition functions for each realization, but their simple average does not coincide with the quantity on the right-hand side of Eq. (43) as in Fig. 7. Both of the results are away from the exact solutions due to the sample-to-sample fluctuation and show nontrivial behavior depending on the linear size. These facts imply lack of self-averaging property. Therefore, if we exploit all the above results given by the configurational average of the exponentiated pseudo work, we have to overcome this violation due to lack of the self-averaging property. This is one of the remaining problem associated with this procedure with PA for spin glasses.

#### **3. Quantum annealing**

Observant readers may begin to recognize the possibility to use physical nature to drive the system in searching the lowest energy (ground state) instead of thermal fluctuation controlled by the inverse temperature. We show another strategy to find the ground state recently studied in a field of physics, quantum annealing (QA) [13].

#### **3.1. Quantum adiabatic computation**

In quantum-mechanical system, we can use a parallel way to drive all the candidates of the desired solution in optimization problem by use of superposition. Quantum annealing uses quantum fluctuation between superposed states to search for the ground state. One of the successful strategies is to use the adiabatic evolution known as quantum adiabatic computation (QAC) [5]. In QAC, as the procedure of SA, we control to gradually decrease the strength of quantum fluctuations to drive the system. Similarly to SA, the convergence into the optimal solution of QAC (the ground state) is also guaranteed by a mathematical proof [15].

In QAC, we introduce a non-commutative operator to drive the system by quantum nature in addition to the original Hamiltonian *H*0, which is designed to represent the optimization problem to be solved, as

$$H(t) = f(t)H\_0 + \left(1 - f(t)\right)H\_{1\prime} \tag{54}$$

where *f*(*t*) is assumed to be a monotonically increasing function satisfying *f*(0) = 0 and *f*(*T*) = 1. For instance, *f*(*t*) = *t*/*T*, where *T* denotes the computation time for QAC. In order to exemplify the explicit instance, we again assume that *H*<sup>0</sup> is the spin glass Hamiltonian as

*<sup>H</sup>*<sup>0</sup> = − ∑ �*ij*� *Jijσ<sup>z</sup> i σz <sup>j</sup>* , (55)

The adiabatic theorem guarantees that the instantaneous state at time *t*, |Ψ(*t*)�, is very close to the instantaneous ground state for a sufficiently large *T* (implying slow control) as |0(*t*)�, if the instantaneous ground state |0(*t*)� is non-degenerate. The condition for |0(*t*)�, �0(*t*)|Ψ(*t*)� ≈

> *dH*(*t*) *dt* |0(*t*)�

where |1(*t*)� is the instantaneous first excited state, and Δ(*t*) is the energy gap between the ground state and first excited one. The maximum and minimum should be evaluated between 0 and *T*. In the present case, since *dH*(*t*)/*dt* ∝ 1/*T*, the above adiabatic condition is reduced

*δ*minΔ2(*t*)

It means that if we desire to solve the optimization problems by use of QAC, which one of the specialized version of QA, we take the computational time proportional to the inverse square of the energy gap. If the problem involved with the exponential closure of the energy gap for increasing of *N*, QAC must take extremely long time to obtain the ground state with a high probability [12, 25]. Interestingly, we can reproduce the convergence theorem of SA (9) from the above adiabatic condition with recourse to a mathematical mapping of the procedure of SA into the quantum dynamics [15, 23]. It implies that the nature of QAC can be understood

Below, we would provide a new paradigm to solve faster than the ordinary scheme of SA. A fast sweep of the system yields nonequilibrium behavior. Although we have not yet understood deeply the nonequilibrium phenomena, there are a few well-established theories which rises to applications to the optimization problem. One possibility is PA for the quantum system by use of JE and its alternatives [20]. Here we again employ JE to give another scheme

In order to consider the nonequilibrium behavior away from the adiabatic dynamics of QAC,

To directly use JE in the protocol to find the ground state of the spin-glass Hamiltonian *H*<sup>0</sup> as in Eq. (55), we prepare the dynamical quantum system following the time-dependent Hamiltonian (58). In addition, we pick up a state from the canonical ensemble for *H*(0) =

time-dependent Schrödinger equation. We measure the performed work in the isolated quantum system as *W* = *Em*(*T*) − *En*(0), which is simply given by the difference between the outputs of projective measurements of the initial and final energies. Here *m* and *n* denote the indices of the instantaneous eigenstates as *H*(*T*)|*m*(*T*)� = *Em*(*T*)|*m*(*T*)� and *H*(0)|*n*(0)� = *En*(0)|*n*(0)�, respectively. The time-evolution operator is given by the following

> i *T* 0

*UT* <sup>=</sup> <sup>T</sup> exp

*<sup>i</sup>* at the initial stage of the procedure and then let it evolve following the

*dtH*(*t*) 

*<sup>T</sup>* <sup>∝</sup> <sup>1</sup>

 

min <sup>Δ</sup>2(*t*) <sup>=</sup> *<sup>δ</sup>*, (59)

Optimization by Use of Nature in Physics Beyond Classical Simulated Annealing 55

. (60)

, (61)

max �1(*t*)<sup>|</sup>

through that of SA by the mathematical mapping technique.

of QA while considering the nonequilibrium behavior.

**3.2. Jarzynski equality for isolated quantum system**

we shortly review JE for an isolated quantum system [1, 24].

<sup>1</sup> <sup>−</sup> *<sup>δ</sup>*2(*<sup>δ</sup>* � <sup>1</sup>) to hold can be written as [15]

into

*<sup>H</sup>*<sup>1</sup> <sup>=</sup> <sup>−</sup>Γ<sup>0</sup> <sup>∑</sup>*<sup>i</sup> <sup>σ</sup><sup>x</sup>*

unitary operator as

where *σ<sup>z</sup> <sup>i</sup>* is the *z* component of the Pauli operators defined as

$$
\sigma^x = \begin{pmatrix} 0 \ 1 \\ 1 \ 0 \end{pmatrix}, \quad \sigma^y = \begin{pmatrix} 0 \ \mathrm{i} \\ -\mathrm{i} \ 0 \end{pmatrix}, \quad \sigma^z = \begin{pmatrix} 1 & 0 \\ 0 & -1 \end{pmatrix}. \tag{56}
$$

We take the computational basis of the eigenstates of the *z*-component of the Pauli matrix (Ising variables) to represent the instantaneous state as <sup>|</sup>Ψ(*t*)� <sup>=</sup> <sup>|</sup>*σ<sup>z</sup>* <sup>1</sup> , *<sup>σ</sup><sup>z</sup>* <sup>2</sup> , ··· , *<sup>σ</sup><sup>z</sup> <sup>N</sup>*�. The transverse-field operator is often used as quantum fluctuations for implementing QAC for the spin-glass model

$$H\_1 = -\Gamma\_0 \sum\_{i=1}^N \sigma\_i^{\chi}.\tag{57}$$

where Γ<sup>0</sup> is the strength of the transverse field. The whole Hamiltonian of QAC (although widely used also for QA) thus becomes

$$H(t) = f(t) \sum\_{\langle i \rangle} f\_{\overline{i}\overline{j}} \sigma\_i^z \sigma\_{\overline{j}}^z + (1 - f(t)) \,\Gamma\_0 \sum\_{i=1}^N \sigma\_i^x. \tag{58}$$

The quantum adiabatic computation starts from a trivial ground state of *H*1. In the present case, the ground state of the transverse-field operator *H*<sup>1</sup> is simply written by a uniform linear combination as <sup>∑</sup>{*σ*} <sup>|</sup>*σ<sup>z</sup>* <sup>1</sup> , *<sup>σ</sup><sup>z</sup>* <sup>2</sup> , ··· , *<sup>σ</sup><sup>z</sup> <sup>N</sup>*�/ √2 *N*.

The adiabatic theorem guarantees that the instantaneous state at time *t*, |Ψ(*t*)�, is very close to the instantaneous ground state for a sufficiently large *T* (implying slow control) as |0(*t*)�, if the instantaneous ground state |0(*t*)� is non-degenerate. The condition for |0(*t*)�, �0(*t*)|Ψ(*t*)� ≈ <sup>1</sup> <sup>−</sup> *<sup>δ</sup>*2(*<sup>δ</sup>* � <sup>1</sup>) to hold can be written as [15]

$$\frac{\max\left|\langle 1(t) | \frac{dH(t)}{dt} | 0(t) \rangle\right|}{\min \Delta^2(t)} = \delta\_\prime \tag{59}$$

where |1(*t*)� is the instantaneous first excited state, and Δ(*t*) is the energy gap between the ground state and first excited one. The maximum and minimum should be evaluated between 0 and *T*. In the present case, since *dH*(*t*)/*dt* ∝ 1/*T*, the above adiabatic condition is reduced into

$$T \propto \frac{1}{\delta \text{min} \Delta^2(t)}.\tag{60}$$

It means that if we desire to solve the optimization problems by use of QAC, which one of the specialized version of QA, we take the computational time proportional to the inverse square of the energy gap. If the problem involved with the exponential closure of the energy gap for increasing of *N*, QAC must take extremely long time to obtain the ground state with a high probability [12, 25]. Interestingly, we can reproduce the convergence theorem of SA (9) from the above adiabatic condition with recourse to a mathematical mapping of the procedure of SA into the quantum dynamics [15, 23]. It implies that the nature of QAC can be understood through that of SA by the mathematical mapping technique.

Below, we would provide a new paradigm to solve faster than the ordinary scheme of SA. A fast sweep of the system yields nonequilibrium behavior. Although we have not yet understood deeply the nonequilibrium phenomena, there are a few well-established theories which rises to applications to the optimization problem. One possibility is PA for the quantum system by use of JE and its alternatives [20]. Here we again employ JE to give another scheme of QA while considering the nonequilibrium behavior.

#### **3.2. Jarzynski equality for isolated quantum system**

14 Will-be-set-by-IN-TECH

Observant readers may begin to recognize the possibility to use physical nature to drive the system in searching the lowest energy (ground state) instead of thermal fluctuation controlled by the inverse temperature. We show another strategy to find the ground state recently

In quantum-mechanical system, we can use a parallel way to drive all the candidates of the desired solution in optimization problem by use of superposition. Quantum annealing uses quantum fluctuation between superposed states to search for the ground state. One of the successful strategies is to use the adiabatic evolution known as quantum adiabatic computation (QAC) [5]. In QAC, as the procedure of SA, we control to gradually decrease the strength of quantum fluctuations to drive the system. Similarly to SA, the convergence into the optimal solution of QAC (the ground state) is also guaranteed by a mathematical

In QAC, we introduce a non-commutative operator to drive the system by quantum nature in addition to the original Hamiltonian *H*0, which is designed to represent the optimization

where *f*(*t*) is assumed to be a monotonically increasing function satisfying *f*(0) = 0 and *f*(*T*) = 1. For instance, *f*(*t*) = *t*/*T*, where *T* denotes the computation time for QAC. In order to exemplify the explicit instance, we again assume that *H*<sup>0</sup> is the spin glass Hamiltonian as

�*ij*�

 0 i −i 0

We take the computational basis of the eigenstates of the *z*-component of the Pauli matrix

transverse-field operator is often used as quantum fluctuations for implementing QAC for

where Γ<sup>0</sup> is the strength of the transverse field. The whole Hamiltonian of QAC (although

The quantum adiabatic computation starts from a trivial ground state of *H*1. In the present case, the ground state of the transverse-field operator *H*<sup>1</sup> is simply written by a uniform linear

*N* ∑ *i*=1 *σx*

*<sup>j</sup>* + (1 − *f*(*t*)) Γ<sup>0</sup>

*H*<sup>1</sup> = −Γ<sup>0</sup>

*Jijσ<sup>z</sup> i σz* *Jijσ<sup>z</sup> i σz*

, *σ<sup>z</sup>* =

*<sup>H</sup>*<sup>0</sup> = − ∑

*<sup>i</sup>* is the *z* component of the Pauli operators defined as

, *σ<sup>y</sup>* =

(Ising variables) to represent the instantaneous state as <sup>|</sup>Ψ(*t*)� <sup>=</sup> <sup>|</sup>*σ<sup>z</sup>*

 0 1 1 0 

*H*(*t*) = *f*(*t*)∑

<sup>2</sup> , ··· , *<sup>σ</sup><sup>z</sup>*

�*ij*�

*<sup>N</sup>*�/ √2 *N*.

*σ<sup>x</sup>* =

widely used also for QA) thus becomes

<sup>1</sup> , *<sup>σ</sup><sup>z</sup>*

*H*(*t*) = *f*(*t*)*H*<sup>0</sup> + (1 − *f*(*t*)) *H*1, (54)

 1 0 0 −1

> *N* ∑ *i*=1 *σx*

*<sup>j</sup>* , (55)

<sup>1</sup> , *<sup>σ</sup><sup>z</sup>*

*<sup>i</sup>* . (57)

. (56)

<sup>2</sup> , ··· , *<sup>σ</sup><sup>z</sup>*

*<sup>i</sup>* . (58)

*<sup>N</sup>*�. The

**3. Quantum annealing**

proof [15].

where *σ<sup>z</sup>*

the spin-glass model

combination as <sup>∑</sup>{*σ*} <sup>|</sup>*σ<sup>z</sup>*

problem to be solved, as

studied in a field of physics, quantum annealing (QA) [13].

**3.1. Quantum adiabatic computation**

In order to consider the nonequilibrium behavior away from the adiabatic dynamics of QAC, we shortly review JE for an isolated quantum system [1, 24].

To directly use JE in the protocol to find the ground state of the spin-glass Hamiltonian *H*<sup>0</sup> as in Eq. (55), we prepare the dynamical quantum system following the time-dependent Hamiltonian (58). In addition, we pick up a state from the canonical ensemble for *H*(0) = *<sup>H</sup>*<sup>1</sup> <sup>=</sup> <sup>−</sup>Γ<sup>0</sup> <sup>∑</sup>*<sup>i</sup> <sup>σ</sup><sup>x</sup> <sup>i</sup>* at the initial stage of the procedure and then let it evolve following the time-dependent Schrödinger equation. We measure the performed work in the isolated quantum system as *W* = *Em*(*T*) − *En*(0), which is simply given by the difference between the outputs of projective measurements of the initial and final energies. Here *m* and *n* denote the indices of the instantaneous eigenstates as *H*(*T*)|*m*(*T*)� = *Em*(*T*)|*m*(*T*)� and *H*(0)|*n*(0)� = *En*(0)|*n*(0)�, respectively. The time-evolution operator is given by the following unitary operator as

$$\mathcal{U}\_T = \mathcal{T} \exp\left(\mathrm{i} \int\_0^T dt H(t)\right),\tag{61}$$

where T denotes the time ordered product. Therefore the transition probability between the initial and final stages is given as

$$P\_{m,n}(0 \to T) = |\langle \Psi\_m(T) | \mathcal{U}\_T | \Psi\_n(0) \rangle|^2. \tag{62}$$

The resultant equation suggests that we can estimate the thermal average under the Hamiltonian *H*<sup>0</sup> on the right-hand side through the left-hand side of JE. This fact may be useful in the evaluation of equilibrium average, since the left-hand side is evaluated without slow adiabatic processes. Differently from QAC, we must sweep the quantum system repeatedly to correctly estimate the nonequilibrium average as in JE, but in short time. We propose such a procedure as the alternative of QAC, the non-adiabatic quantum annealing (NQA) based on the property of JE as above established. First we discuss the possibility as a solver of the

Optimization by Use of Nature in Physics Beyond Classical Simulated Annealing 57

In order to investigate the property of the ground state, we tune the inverse temperature into a very large value *β* � 1. The nonequilibrium average on the left-hand side of JE involves a non-extensive quantity, the exponentiated work, whose value fluctuates significantly from process to process. Therefore the average on the left-hand side must be calculated by many trials of annealing processes. Thus, rare events with large values of the exponentiated work (i.e. *β*|*W*| � Γ0) would contribute to the average significantly, and we have to repeat the annealing process very many, typically exponentially, times in order to reach the correct value of the average. This property would be a bottleneck of the simple implementation of NQA, instead of long time involved by the closure of the energy gap in QAC. What about PA in the classical counter part? In order to generate the relevant contributions in PA, we use the biased distribution with the exponentiated pseudo work through resampling technique [7, 9]. Without resampling technique, we cannot efficiently reproduce the prediction given by JE. In this fact, if we implement the biased distribution in the quantum system, we would use NQA without suffering from rare events. It means that we can perform NQA in order to find the ground state in a short time with several repetitions. So far, it is the very difficult task to realize the biased distribution in the quantum system. However it is worthwhile to consider

Unfortunately, we have not reached any positive answers on the performance of NQA. Instead let us here evaluate several properties in nonequilibrium process as in NQA for the particular spin glasses. We can exactly analyze nonequilibrium behavior by combination of JE with the gauge transformation, although there are few exact results in nonequilibrium quantum

Following the prescription of JE, let us consider a repetition of NQA starting from the equilibrium ensemble. The initial Hamiltonian of NQA is given only by the transverse field *H*(0) = *H*1. It turns out that the starting point of our analyses is the specialized JE to the case

�e−*βW*�QA <sup>=</sup> *ZT*(*β*, {*Jij*})

We assume that the interactions {*Jij*} follow the distribution function for the ±*J* Ising model

*<sup>P</sup>*(*Jij*) = exp(*β<sup>p</sup> Jij*)

2 cosh *βp J*

(2 cosh *<sup>β</sup>*Γ0)*<sup>N</sup>* . (68)

, (69)

optimization problem below.

its possibility in the future.

for NQA as

**4.1. Several analyses of non-adiabatic quantum annealing**

dynamical system with many components [22].

(22), which is better to be rewritten as

**4. Non-adiabatic quantum annealing**

The path probability for the nonequilibrium process starting from the equilibrium ensemble is then evaluated as

$$P\_{\rm m,n}(0 \to T) \frac{\exp(-\beta E\_{\rm n}(0))}{Z\_0(\beta; \{I\_{\rm ij}\})},\tag{63}$$

where we express the instantaneous partition function for the instantaneous Hamiltonian at each time *t* as *Zt*(*β*; {*Jij*}).

By directly evaluating the nonequilibrium average of the exponentiated work but for the isolated system, we reach JE applicable to non-adiabatic version of QA. We define the nonequilibrium average of the exponentiated work as

$$\left\langle \mathbf{e}^{-\beta W} \right\rangle\_{\mathrm{QA}} = \sum\_{m,n} \mathbf{e}^{-\beta W} P\_{m,n}(0 \to T) \frac{\exp(-\beta E\_{\mathrm{il}}(0))}{Z\_0(\beta; \{I\_{\mathrm{ij}}\})},$$

which becomes the left-hand side of JE. The quantity defined here is evaluated as

$$
\left\langle \mathbf{e}^{-\beta W} \right\rangle\_{\mathrm{QA}} = \sum\_{m,n} \frac{\mathbf{e}^{-\beta E\_m(T)}}{Z\_0(\beta; \{I\_{lj}\})} P\_{m,n}(0 \to T)
$$

$$
= \sum\_m \frac{\mathbf{e}^{-\beta E\_m(T)}}{Z\_0(\beta; \{I\_{lj}\})}
$$

$$
= \frac{Z\_T(\beta; \{I\_{lj}\})}{Z\_0(\beta; \{I\_{lj}\})} \,\,\,\,\tag{64}
$$

where we used the fact that the performed work *W* was a classical number and

$$\sum\_{n} P\_{m,n}(0 \to T) = \sum\_{n} \langle \Psi\_{m}(T) | \mathcal{U}\_{T} | \Psi\_{n}(0) \rangle \langle \Psi\_{n}(0) | \mathcal{U}\_{T}^{\dagger} | \Psi\_{i}(T) \rangle$$

$$= \sum\_{m} \langle \Psi\_{m}(T) | \mathcal{U}\_{T} \mathcal{U}\_{T}^{\dagger} | \Psi\_{m}(T) \rangle = 1. \tag{65}$$

If we measure the physical observable *O*ˆ *<sup>T</sup>* at the last of the nonequilibrium process, we obtain another equation as, similarly to the classical version,

$$
\langle \bullet\_T \text{e}^{-\beta W} \rangle\_{\text{QA}} = \langle \bullet \rangle\_{\beta} \frac{Z\_T(\nexists, \{I\_{\text{ij}}\})}{Z\_0(\beta; \{I\_{\text{ij}}\})} \,\text{}\,\text{}\tag{66}
$$

where the subscript on the square brackets in the right-hand side denotes the thermal average in the last equilibrium state with the inverse temperature *β*. The ratio of Eqs. (64) and (66) gives

$$\frac{\langle \bullet\_{T} \mathbf{e}^{-\beta W} \rangle\_{\mathbf{QA}}}{\langle \mathbf{e}^{-\beta W} \rangle\_{\mathbf{QA}}} = \langle \bullet \rangle\_{\beta}. \tag{67}$$

The resultant equation suggests that we can estimate the thermal average under the Hamiltonian *H*<sup>0</sup> on the right-hand side through the left-hand side of JE. This fact may be useful in the evaluation of equilibrium average, since the left-hand side is evaluated without slow adiabatic processes. Differently from QAC, we must sweep the quantum system repeatedly to correctly estimate the nonequilibrium average as in JE, but in short time. We propose such a procedure as the alternative of QAC, the non-adiabatic quantum annealing (NQA) based on the property of JE as above established. First we discuss the possibility as a solver of the optimization problem below.

#### **4. Non-adiabatic quantum annealing**

16 Will-be-set-by-IN-TECH

where T denotes the time ordered product. Therefore the transition probability between the

The path probability for the nonequilibrium process starting from the equilibrium ensemble

where we express the instantaneous partition function for the instantaneous Hamiltonian at

By directly evaluating the nonequilibrium average of the exponentiated work but for the isolated system, we reach JE applicable to non-adiabatic version of QA. We define the

<sup>e</sup>−*<sup>β</sup>WPm*,*n*(<sup>0</sup> <sup>→</sup> *<sup>T</sup>*)

e−*βEm*(*T*) *Z*0(*β*; {*Jij*})

e−*βEm*(*T*) *Z*0(*β*; {*Jij*})

�Ψ*m*(*T*)|*UT*|Ψ*n*(0)��Ψ*n*(0)|*U*†

*ZT*(*β*; {*Jij*}) *Z*0(*β*; {*Jij*})

<sup>=</sup> *ZT*(*β*; {*Jij*}) *Z*0(*β*; {*Jij*})

�Ψ*m*(*T*)|*UTU*†

If we measure the physical observable *O*ˆ *<sup>T</sup>* at the last of the nonequilibrium process, we obtain

where the subscript on the square brackets in the right-hand side denotes the thermal average in the last equilibrium state with the inverse temperature *β*. The ratio of Eqs. (64) and (66)

exp(−*βEn*(0))

*Pm*,*n*(0 → *T*)

which becomes the left-hand side of JE. The quantity defined here is evaluated as

QA <sup>=</sup> ∑ *m*,*n*

> = ∑ *m*

where we used the fact that the performed work *W* was a classical number and

�*O*<sup>ˆ</sup> *<sup>T</sup>*e−*βW*�QA <sup>=</sup> �*O*ˆ�*<sup>β</sup>*

�*O*<sup>ˆ</sup> *<sup>T</sup>*e−*βW*�QA �e−*βW*�QA

= ∑ *m*

nonequilibrium average of the exponentiated work as

 e−*β<sup>W</sup>*

*Pm*,*n*(<sup>0</sup> <sup>→</sup> *<sup>T</sup>*) = <sup>∑</sup>*<sup>n</sup>*

another equation as, similarly to the classical version,

QA <sup>=</sup> ∑ *m*,*n*

 e−*β<sup>W</sup>*

∑*n*

gives

*Pm*,*n*(<sup>0</sup> <sup>→</sup> *<sup>T</sup>*) = |�Ψ*m*(*T*)|*UT*|Ψ*n*(0)�|2. (62)

*<sup>Z</sup>*0(*β*; {*Jij*}) , (63)

, (64)

*<sup>T</sup>*|Ψ*i*(*T*)�

*<sup>T</sup>*|Ψ*m*(*T*)� = 1. (65)

<sup>=</sup> �*O*ˆ�*β*. (67)

, (66)

exp(−*βEn*(0)) *<sup>Z</sup>*0(*β*; {*Jij*}) ,

*Pm*,*n*(0 → *T*)

initial and final stages is given as

is then evaluated as

each time *t* as *Zt*(*β*; {*Jij*}).

In order to investigate the property of the ground state, we tune the inverse temperature into a very large value *β* � 1. The nonequilibrium average on the left-hand side of JE involves a non-extensive quantity, the exponentiated work, whose value fluctuates significantly from process to process. Therefore the average on the left-hand side must be calculated by many trials of annealing processes. Thus, rare events with large values of the exponentiated work (i.e. *β*|*W*| � Γ0) would contribute to the average significantly, and we have to repeat the annealing process very many, typically exponentially, times in order to reach the correct value of the average. This property would be a bottleneck of the simple implementation of NQA, instead of long time involved by the closure of the energy gap in QAC. What about PA in the classical counter part? In order to generate the relevant contributions in PA, we use the biased distribution with the exponentiated pseudo work through resampling technique [7, 9]. Without resampling technique, we cannot efficiently reproduce the prediction given by JE. In this fact, if we implement the biased distribution in the quantum system, we would use NQA without suffering from rare events. It means that we can perform NQA in order to find the ground state in a short time with several repetitions. So far, it is the very difficult task to realize the biased distribution in the quantum system. However it is worthwhile to consider its possibility in the future.

#### **4.1. Several analyses of non-adiabatic quantum annealing**

Unfortunately, we have not reached any positive answers on the performance of NQA. Instead let us here evaluate several properties in nonequilibrium process as in NQA for the particular spin glasses. We can exactly analyze nonequilibrium behavior by combination of JE with the gauge transformation, although there are few exact results in nonequilibrium quantum dynamical system with many components [22].

Following the prescription of JE, let us consider a repetition of NQA starting from the equilibrium ensemble. The initial Hamiltonian of NQA is given only by the transverse field *H*(0) = *H*1. It turns out that the starting point of our analyses is the specialized JE to the case for NQA as

$$
\langle \mathbf{e}^{-\beta W} \rangle\_{\mathrm{QA}} = \frac{Z\_T(\boldsymbol{\beta}, \{I\_{l\bar{l}}\})}{(2 \cosh \beta \Gamma\_0)^N}. \tag{68}
$$

We assume that the interactions {*Jij*} follow the distribution function for the ±*J* Ising model (22), which is better to be rewritten as

$$P(f\_{lj}) = \frac{\exp(\beta\_p f\_{lj})}{2 \cosh \beta\_p f} \tag{69}$$

#### 18 Will-be-set-by-IN-TECH 58 Simulated Annealing – Advances, Applications and Hybridizations Optimization by Use of Nature in Physics Beyond Classical Simulated Annealing <sup>19</sup>

where we do not use *K* = *βJ* for transparency, and exp(−2*β<sup>p</sup> J*)=(1 − *p*)/*p*.

#### **4.2. Gauge transformation for quantum spin systems**

For several special spin glasses as the ±*J* Ising model, the gauge transformation is available for analyses on the dynamical property even under quantum fluctuations. The time-dependent Hamiltonian as in Eq. (58) is invariant under the following local transformation,

$$
\sigma\_i^{\mathbf{x}} \to \sigma\_i^{\mathbf{x}}, \sigma\_i^{\mathbf{y}} \to \mathfrak{xi}\_i \\
\sigma\_i^{\mathbf{y}}, \sigma\_i^{\mathbf{z}} \to \mathfrak{ξ}\_i \\
\sigma\_i^{\mathbf{z}}, \mathfrak{f}\_{\mathbf{ij}} \to \mathfrak{f}\_{\mathbf{ij}} \mathfrak{f}\_i \\
\mathfrak{f}\_{\mathbf{j}} \quad (\forall i, j), \tag{70}
$$

P

SG

Ǫ

the performed work during NQA

NL (*β<sup>p</sup>* = *β*).

work as

distribution,

Ǫ

 �*W*�QA

> 0 ≥ 1

 �*W*�QA

�e−*β*1*W*�QA

<sup>0</sup> ≥ − *<sup>N</sup>*

F

Ǫ

**Figure 8.** Two different processes of NQA in Eq. (75). The left-hand side of Eq. (75) represents the process toward the upper-right red dot and the right-hand side ends at the lower-left dot. Three phases expressed by the same symbols as in Fig. 4 are separated by solid curves and the dotted line describes

in Eq. (75), (implying *p* = 1/2, the symmetric distribution), we find a nontrivial relation on

The symmetric distribution (*β*<sup>2</sup> = 0 on the left-hand side) makes it possible to reduce the right-hand side to the above trivial expression. It is remarkable that NQA, which involves very complex dynamics, satisfies such a simple identity irrespective of the speed of annealing *T*. The Jensen inequality for the above equality leads us to the lower bound for the performed

where *z* is the coordination number as *z* = *NB*/*N*. Here we generalize the inverse temperature to *β* from the specific choice *β*1. This lower bound is loose, since the direct application of the Jensen inequality to JE for NQA yields, after the configurational average with the symmetric

*<sup>β</sup> <sup>D</sup>*(0|*β*) <sup>−</sup> *<sup>N</sup>*

<sup>0</sup> <sup>=</sup> (cosh *<sup>β</sup>*<sup>1</sup> *<sup>J</sup>*)*NB*

*<sup>β</sup>* log (cosh *<sup>β</sup>J*)*<sup>z</sup>* cosh *β*Γ<sup>0</sup>

> *<sup>β</sup>* log (cosh *<sup>β</sup>J*)*<sup>z</sup>* cosh *β*Γ<sup>0</sup>

(cosh *<sup>β</sup>*1Γ0)*<sup>N</sup>* . (76)

. (77)

. (78)

Ǫ

Optimization by Use of Nature in Physics Beyond Classical Simulated Annealing 59

Ǫ

Ǫ

where *ξi*(= ±1) is called as a gauge variable. This transformation is designed to preserve the commutation relations between different components of Pauli matrices [16].

#### **4.3. Relationship between two different paths of NQA**

Below, we reveal several properties inherent in NQA by the gauge transformation. Let us take the configurational average of Eq. (68) over all the realizations of {*Jij*} for the special case with *β* = *β*<sup>1</sup> and *β<sup>p</sup>* = *β*<sup>2</sup> as

$$\left[ \langle \mathbf{e}^{-\beta\_1 W} \rangle\_{\rm QA} \right]\_{\beta\_2} = \left[ \frac{Z\_T(\beta\_1; \{l\_{l\bar{j}}\})}{(2 \cosh \beta\_1 \Gamma\_0)^N} \right]\_{\beta\_2}.\tag{71}$$

The right-hand side is written explicitly as

$$\left[ \langle \mathbf{e}^{-\beta\_1 W} \rangle\_{\mathrm{QA}} \right]\_{\beta\_2} = \sum\_{\{l\_{\vec{l}}\}} \frac{\exp \left( \beta\_2 \sum\_{\langle \vec{l} \rangle} l\_{\vec{l} \vec{j}} \right)}{\left( 2 \cosh \beta\_2 I \right)^{N\_B}} \frac{Z\_T(\beta\_1; \{l\_{\vec{l} \vec{j}}\})}{\left( 2 \cosh \beta\_1 \Gamma\_0 \right)^N}. \tag{72}$$

Let us here apply the gauge transformation as introduced above. Since the time-dependent Hamiltonian is invariant, we may sum over all possible configurations of the gauge variables {*ξi*} and divide the result by 2*<sup>N</sup>* in order to obtain the quantity on the left-hand side,

$$\left[ \langle \mathbf{e}^{-\beta\_1 W} \rangle\_{\mathrm{QA}} \right]\_{\beta\_2} = \sum\_{\{I\_{\vec{l}}\}} \frac{Z\_T(\beta\_2; \{I\_{\vec{l}\vec{j}}\}) Z\_T(\beta\_1; \{I\_{\vec{l}\vec{j}}\})}{2^N (2 \cosh \beta\_2 I)^{N\_B} (2 \cosh \beta\_1 \Gamma\_0)^N}. \tag{73}$$

A similar quantity of the average of the exponentiated work for spin glass with the inverse temperature *β*<sup>2</sup> and the parameter for the quenched randomness *β*<sup>1</sup> gives

$$\left[ \langle \mathbf{e}^{-\beta\_2 W} \rangle\_{\mathrm{QA}} \right]\_{\beta\_1} = \sum\_{\{I\_{\vec{l}}\}} \frac{Z\_T(\beta\_2; \{I\_{\vec{l}\vec{j}}\}) Z\_T(\beta\_1; \{I\_{\vec{l}\vec{j}}\})}{2^N (2 \cosh \beta\_1 I)^{N\_B} \left( 2 \cosh \beta\_2 \Gamma\_0 \right)^N}. \tag{74}$$

Comparing Eqs. (73) and (74), we find the following relation between two different non-adiabatic processes,

$$\left[ \langle \mathbf{e}^{-\beta\_1 W} \rangle\_{\mathrm{QA}} \right]\_{\beta\_2} = \left[ \langle \mathbf{e}^{-\beta\_2 W} \rangle\_{\mathrm{QA}} \right]\_{\beta\_1} \left( \frac{\cosh \beta\_1 I}{\cosh \beta\_2 I} \right)^{N\_\mathrm{B}} \left( \frac{\cosh \beta\_2 \Gamma\_0}{\cosh \beta\_1 \Gamma\_0} \right)^N. \tag{75}$$

We describe the two different paths of NQA related by this equality in Fig. 8. Setting *β*<sup>2</sup> = 0

18 Will-be-set-by-IN-TECH

For several special spin glasses as the ±*J* Ising model, the gauge transformation is available for analyses on the dynamical property even under quantum fluctuations. The time-dependent

*<sup>i</sup>* <sup>→</sup> *<sup>ξ</sup>iσ<sup>z</sup>*

where *ξi*(= ±1) is called as a gauge variable. This transformation is designed to preserve the

Below, we reveal several properties inherent in NQA by the gauge transformation. Let us take the configurational average of Eq. (68) over all the realizations of {*Jij*} for the special case with

exp

{*ξi*} and divide the result by 2*<sup>N</sup>* in order to obtain the quantity on the left-hand side,

Let us here apply the gauge transformation as introduced above. Since the time-dependent Hamiltonian is invariant, we may sum over all possible configurations of the gauge variables

A similar quantity of the average of the exponentiated work for spin glass with the inverse

Comparing Eqs. (73) and (74), we find the following relation between two different

 *β*1

We describe the two different paths of NQA related by this equality in Fig. 8. Setting *β*<sup>2</sup> = 0

 *ZT*(*β*1; {*Jij*}) (2 cosh *β*1Γ0)

*<sup>β</sup>*<sup>2</sup> <sup>∑</sup>�*ij*� *Jij*

(2 cosh *β*<sup>2</sup> *J*)*NB*

*ZT*(*β*2; {*Jij*})*ZT*(*β*1; {*Jij*}) 2*N*(2 cosh *β*<sup>2</sup> *J*)*NB* (2 cosh *β*1Γ0)

*ZT*(*β*2; {*Jij*})*ZT*(*β*1; {*Jij*}) 2*N*(2 cosh *β*<sup>1</sup> *J*)*NB* (2 cosh *β*2Γ0)

> cosh *β*<sup>1</sup> *J* cosh *β*<sup>2</sup> *J*

*N* 

*β*2

*ZT*(*β*1; {*Jij*}) (2 cosh *β*1Γ0)

*NB* cosh *<sup>β</sup>*2Γ<sup>0</sup>

cosh *β*1Γ<sup>0</sup>

*<sup>N</sup>*

*<sup>i</sup>* , *Jij* → *Jijξiξ<sup>j</sup>* (∀*i*, *j*), (70)

. (71)

*<sup>N</sup>* . (72)

*<sup>N</sup>* . (73)

*<sup>N</sup>* . (74)

. (75)

where we do not use *K* = *βJ* for transparency, and exp(−2*β<sup>p</sup> J*)=(1 − *p*)/*p*.

Hamiltonian as in Eq. (58) is invariant under the following local transformation,

*<sup>i</sup>* , *<sup>σ</sup><sup>z</sup>*

*<sup>i</sup>* <sup>→</sup> *<sup>ξ</sup>iσ<sup>y</sup>*

commutation relations between different components of Pauli matrices [16].

 *β*2 =

**4.2. Gauge transformation for quantum spin systems**

*<sup>i</sup>* , *<sup>σ</sup><sup>y</sup>*

**4.3. Relationship between two different paths of NQA**

�e−*β*1*W*�QA

 *β*2 = ∑ {*Jij*}

 *β*2 = ∑ {*Jij*}

 *β*1 = ∑ {*Jij*}

temperature *β*<sup>2</sup> and the parameter for the quenched randomness *β*<sup>1</sup> gives

�e−*β*2*W*�QA

�e−*β*1*W*�QA

�e−*β*1*W*�QA

�e−*β*2*W*�QA

 *β*2 = 

*σx <sup>i</sup>* <sup>→</sup> *<sup>σ</sup><sup>x</sup>*

The right-hand side is written explicitly as

�e−*β*1*W*�QA

non-adiabatic processes,

*β* = *β*<sup>1</sup> and *β<sup>p</sup>* = *β*<sup>2</sup> as

**Figure 8.** Two different processes of NQA in Eq. (75). The left-hand side of Eq. (75) represents the process toward the upper-right red dot and the right-hand side ends at the lower-left dot. Three phases expressed by the same symbols as in Fig. 4 are separated by solid curves and the dotted line describes NL (*β<sup>p</sup>* = *β*).

in Eq. (75), (implying *p* = 1/2, the symmetric distribution), we find a nontrivial relation on the performed work during NQA

$$\left[ \langle \mathbf{e}^{-\beta\_1 W} \rangle\_{\mathrm{QA}} \right]\_0 = \frac{(\cosh \beta\_1 I)^{N\_B}}{(\cosh \beta\_1 \Gamma\_0)^N}. \tag{76}$$

The symmetric distribution (*β*<sup>2</sup> = 0 on the left-hand side) makes it possible to reduce the right-hand side to the above trivial expression. It is remarkable that NQA, which involves very complex dynamics, satisfies such a simple identity irrespective of the speed of annealing *T*. The Jensen inequality for the above equality leads us to the lower bound for the performed work as

$$\mathbb{E}\left[\langle W\rangle\_{\text{QA}}\right]\_0 \geq -\frac{N}{\beta} \log\left(\frac{(\cosh\beta I)^z}{\cosh\beta \Gamma\_0}\right). \tag{77}$$

where *z* is the coordination number as *z* = *NB*/*N*. Here we generalize the inverse temperature to *β* from the specific choice *β*1. This lower bound is loose, since the direct application of the Jensen inequality to JE for NQA yields, after the configurational average with the symmetric distribution,

$$\mathbb{E}\left[\langle W \rangle\_{\rm QA}\right]\_0 \geq \frac{1}{\beta} D(0|\beta) - \frac{N}{\beta} \log\left(\frac{(\cosh \beta I)^z}{\cosh \beta \Gamma\_0}\right). \tag{78}$$

#### 20 Will-be-set-by-IN-TECH 60 Simulated Annealing – Advances, Applications and Hybridizations Optimization by Use of Nature in Physics Beyond Classical Simulated Annealing <sup>21</sup>

where *D*(*β*|*β*� ) is the Kullback-Leibler divergence defined as

$$D(\beta|\beta') = \sum\_{\{I\_{\vec{l}}\}} \tilde{P}\_{\beta}(\{I\_{\vec{l}\vec{j}}\}) \log \frac{\tilde{P}\_{\beta'}(\{I\_{\vec{l}\vec{j}}\})}{\tilde{P}\_{\beta}(\{I\_{\vec{l}\vec{j}}\})} \tag{79}$$

P

SG

Ǫ

Eq. (85), whereas the red one is for the left-hand side.

relation

�*σz i σz* F

**Figure 9.** Two different nonequilibrium processes in NQA through Eq. (85). The same symbols are depicted as in Fig. 4. The blue circle represents the determination of the process on the right-hand side of

Nishimori line, while the other expresses for the symmetric distribution.

further consider the case for the two-point correlation *OT* = *σ<sup>z</sup>*

1

*<sup>j</sup>* <sup>e</sup>−*βW*�QA

*β*

*<sup>j</sup>* <sup>e</sup>−*βW*�QA

1

average of both sides under the condition *β<sup>p</sup>* = *β*, we find

which is another exact identity for processes of NQA.

�*σz i σz*

As shown in Fig. 9, the resultant equation leads us to a fascinating relationship of the two completely different processes through the inverse statistics. One denotes NQA toward the

We can find the exact results through the inverse statics of the inverse statics of Eq. (66). Let us

The quantity on the right-hand side becomes unity by the analysis with the gauge transformation as has been shown in the literatures [18, 19]. We thus reach a simple exact

*β*

We obtain several exact nontrivial relations between completely different paths in NQA as shown above by use of the gauge transformation, which is a specialized tool to analyze spin glasses. We should notice that such results are very rare for the nonequilibrium behavior in disordered quantum system. The importance of the above equalities is still not clear. We emphasize that, when we realize the quantum spin systems in experiments, the above

<sup>=</sup> (cosh *<sup>β</sup>*Γ0)*<sup>N</sup>* (cosh *βJ*)*NB*

*i σz*

1 �*σz i σz j* �*β* *β*

(cosh *<sup>β</sup>J*)*NB* , (87)

<sup>=</sup> (cosh *<sup>β</sup>*Γ0)*<sup>N</sup>*

*<sup>j</sup>* . Taking the configurational

. (86)

Ǫ

Optimization by Use of Nature in Physics Beyond Classical Simulated Annealing 61

Here we defined the marginal distribution for the specific configuration {*Jij*} summed over all the possible gauge transformations,

$$\mathcal{P}\_{\hat{\boldsymbol{\beta}}}(\{I\_{\hat{\boldsymbol{ij}}}\}) = \frac{1}{2^{N}} \sum\_{\{\tilde{\xi}\_{i}\}} \prod\_{\langle i\hat{\boldsymbol{j}}\rangle} P(I\_{\hat{\boldsymbol{ij}}}) = \frac{Z\_{\mathcal{T}}(\boldsymbol{\beta}; \{I\_{\hat{\boldsymbol{ij}}}\})}{2^{N} (2 \cosh \boldsymbol{\beta} I)^{N\_{\mathcal{B}}}}.\tag{80}$$

Since the Kullback-Leibler divergence does not become non-negative, the work performed by the transverse field during a nonequilibrium process in the symmetric distribution (i.e. the left-hand side of Eq. (78)) does not lower below the second quantity on the right-hand side of Eq. (78), namely Eq. (77). This fact means that Eq. (77) is a loose lower bound.

#### **4.4. Exact relations involving inverse statistics**

Beyond the above results, we can perform further non-trivial analyses for the nonequilibrium process in the special conditions. Let us next take the configurational average of the inverse of JE, Eq. (68), as

$$
\left[\frac{1}{\langle \mathbf{e}^{-\beta W} \rangle\_{\mathrm{QA}}}\right]\_{\beta\_p} = \left[\frac{(2 \cosh \beta \Gamma\_0)^N}{Z\_T(\beta; \{I\_{ij}\})}\right]\_{\beta\_p}.\tag{81}
$$

The application of the gauge transformation to the right-hand side yields

$$
\left[\frac{1}{\langle \mathbf{e}^{-\beta W} \rangle\_{\mathrm{QA}}}\right]\_{\mathfrak{P}\_{\mathbb{P}}} = \sum\_{\{I\_{\vec{l}}\}} \frac{\exp\left(\beta\_{\mathcal{P}} \sum\_{\{\vec{l}\}} J\_{\vec{l}\vec{l}} \mathfrak{T}\_{\vec{l}} \mathfrak{T}\_{\vec{l}}\right)}{\left(2 \cosh \beta\_{\mathcal{P}} I\right)^{\mathrm{N}\_{\mathbb{B}}}} \frac{\left(2 \cosh \beta \Gamma\_{\mathbb{Q}}\right)^{\mathrm{N}}}{Z\_{\mathsf{T}}(\beta; \{I\_{\vec{l}\vec{j}}\})}.\tag{82}
$$

Summation of the right-hand side over all the possible configurations of {*ξi*} and division of the result by 2*<sup>N</sup>* give

$$
\left[\frac{1}{\langle \mathbf{e}^{-\beta \mathcal{W}} \rangle\_{\mathrm{QA}}}\right]\_{\beta\_{\mathrm{p}}} = \sum\_{\{I\_{\overline{l}}\}} \frac{Z\_{T}(\beta\_{\overline{p}}; \{I\_{\overline{l}\overline{j}}\})}{2^{N}(2\cosh\beta\_{\overline{p}}I)^{N\_{\overline{p}}}} \frac{(2\cosh\beta\Gamma\_{0})^{N}}{Z\_{T}(\beta; \{I\_{\overline{l}\overline{j}}\})}.\tag{83}
$$

This equation reduces to, by setting *β<sup>p</sup>* = *β*, namely on the Nishimori line,

$$
\left[\frac{1}{\langle \mathbf{e}^{-\beta \mathcal{W}} \rangle\_{\rm QA}}\right]\_{\beta} = \frac{(\cosh \beta \Gamma\_0)^N}{(\cosh \beta f)^{N\_B}}.\tag{84}
$$

Comparison of Eqs. (76) and (84) gives

$$\left[ \langle \mathbf{e}^{-\beta \mathcal{W}} \rangle\_{\mathrm{QA}} \right]\_0 = \left( \left[ \frac{1}{\langle \mathbf{e}^{-\beta \mathcal{W}} \rangle\_{\mathrm{QA}}} \right]\_\beta \right)^{-1} . \tag{85}$$

**Figure 9.** Two different nonequilibrium processes in NQA through Eq. (85). The same symbols are depicted as in Fig. 4. The blue circle represents the determination of the process on the right-hand side of Eq. (85), whereas the red one is for the left-hand side.

As shown in Fig. 9, the resultant equation leads us to a fascinating relationship of the two completely different processes through the inverse statistics. One denotes NQA toward the Nishimori line, while the other expresses for the symmetric distribution.

We can find the exact results through the inverse statics of the inverse statics of Eq. (66). Let us further consider the case for the two-point correlation *OT* = *σ<sup>z</sup> i σz <sup>j</sup>* . Taking the configurational average of both sides under the condition *β<sup>p</sup>* = *β*, we find

$$
\left[\frac{1}{\langle \sigma\_{\!i}^{z} \sigma\_{\!j}^{z} \mathbf{e}^{-\beta W} \rangle\_{\text{QA}}}\right]\_{\beta} = \frac{(\cosh \beta \Gamma\_{0})^{N}}{(\cosh \beta I)^{N\_{\text{B}}}} \left[\frac{1}{\langle \sigma\_{\!i}^{z} \sigma\_{\!j}^{z} \rangle\_{\beta}}\right]\_{\beta}.\tag{86}
$$

The quantity on the right-hand side becomes unity by the analysis with the gauge transformation as has been shown in the literatures [18, 19]. We thus reach a simple exact relation

$$
\left[\frac{1}{\langle \sigma\_i^z \sigma\_j^z \mathbf{e}^{-\beta \mathcal{W}} \rangle\_{\text{QA}}}\right]\_{\beta} = \frac{(\cosh \beta \Gamma\_0)^N}{(\cosh \beta f)^{N\_\mathcal{B}}}\,\tag{87}
$$

which is another exact identity for processes of NQA.

20 Will-be-set-by-IN-TECH

Here we defined the marginal distribution for the specific configuration {*Jij*} summed over

Since the Kullback-Leibler divergence does not become non-negative, the work performed by the transverse field during a nonequilibrium process in the symmetric distribution (i.e. the left-hand side of Eq. (78)) does not lower below the second quantity on the right-hand side of

Beyond the above results, we can perform further non-trivial analyses for the nonequilibrium process in the special conditions. Let us next take the configurational average of the inverse

(2 cosh *β*Γ0)

*<sup>β</sup><sup>p</sup>* <sup>∑</sup>�*ij*� *Jijξiξ<sup>j</sup>*

(2 cosh *β<sup>p</sup> J*)*NB*

*ZT*(*βp*; {*Jij*}) 2*N*(2 cosh *β<sup>p</sup> J*)*NB*

<sup>=</sup> (cosh *<sup>β</sup>*Γ0)*<sup>N</sup>*

1 �e−*βW*�QA �

*β*

⎞ ⎠ −1

*ZT*(*β*; {*Jij*})

*N*

�

�

*βp*

(2 cosh *β*Γ0)

(2 cosh *β*Γ0)

*N*

*N*

(cosh *<sup>β</sup>J*)*NB* . (84)

*ZT*(*β*; {*Jij*}) . (82)

*ZT*(*β*; {*Jij*}) . (83)

. (85)

*<sup>β</sup>*({*Jij*})log

*P*(*Jij*) =

*P*˜ *<sup>β</sup>*�({*Jij*})

*P*˜

*ZT*(*β*; {*Jij*})

*<sup>β</sup>*({*Jij*}) (79)

<sup>2</sup>*N*(2 cosh *<sup>β</sup>J*)*NB* . (80)

. (81)

) is the Kullback-Leibler divergence defined as

) = ∑ {*Jij*} *P*˜

<sup>2</sup>*<sup>N</sup>* ∑ {*ξi*} ∏ �*ij*�

Eq. (78), namely Eq. (77). This fact means that Eq. (77) is a loose lower bound.

�

*βp* = �

exp �

Summation of the right-hand side over all the possible configurations of {*ξi*} and division of

*D*(*β*|*β*�

*<sup>β</sup>*({*Jij*}) = <sup>1</sup>

all the possible gauge transformations,

�

�

Comparison of Eqs. (76) and (84) gives

1 �e−*βW*�QA

> 1 �e−*βW*�QA

*P*˜

**4.4. Exact relations involving inverse statistics**

�

�

*βp*

�

�

�e−*βW*�QA

�

*βp*

1 �e−*βW*�QA

The application of the gauge transformation to the right-hand side yields

= ∑ {*Jij*}

> = ∑ {*Jij*}

This equation reduces to, by setting *β<sup>p</sup>* = *β*, namely on the Nishimori line,

1 �e−*βW*�QA

> � 0 =

�

*β*

⎛ ⎝ �

where *D*(*β*|*β*�

of JE, Eq. (68), as

the result by 2*<sup>N</sup>* give

We obtain several exact nontrivial relations between completely different paths in NQA as shown above by use of the gauge transformation, which is a specialized tool to analyze spin glasses. We should notice that such results are very rare for the nonequilibrium behavior in disordered quantum system. The importance of the above equalities is still not clear. We emphasize that, when we realize the quantum spin systems in experiments, the above

#### 22 Will-be-set-by-IN-TECH 62 Simulated Annealing – Advances, Applications and Hybridizations Optimization by Use of Nature in Physics Beyond Classical Simulated Annealing <sup>23</sup>

results would be a valuable platform to confirm their precisions and conditions. The quantum annealing is originally intended to be a tool implemented in so-called quantum computer. Therefore we need theoretical studies not only on the performance in the application but also how to make good condition to implement the protocol. In this regard, our studies shed light on the future indicator to open the way to realize the generic solver in the quantum computer. We must continue the active studies in this direction.

**Acknowledgement**

**Author details** Masayuki Ohzeki

**6. References**

*Sakyo-ku, Kyoto, 606-8501, Japan*

*Japan*

*Phys. Rev. E* 61: 2361–2366.

*on* PAMI-6(6): 721 –741.

16(2): 279–286.

78: 2690–2693.

quantum systems, *Phys. Rev. Lett.* 102: 210401.

URL: *http://link.aps.org/doi/10.1103/PhysRevLett.102.210401*

URL: *http://link.aps.org/doi/10.1103/PhysRevE.60.2721*

URL: *http://link.aps.org/doi/10.1103/PhysRevE.61.2361*

Adiabatic Evolution, *arXiv:quant-ph/0001106v1* .

*AIP Conference Proceedings* 690(1): 200–206.

URL: *http://jpsj.ipap.jp/link?JPSJ/65/1604/*

URL: *http://link.aps.org/doi/10.1103/PhysRevE.56.5018*

URL: *http://link.aps.org/doi/10.1103/PhysRevLett.78.2690*

URL: *http://link.aps.org/doi/10.1103/PhysRevLett.101.147204*

quantum annealing, *Phys. Rev. Lett.* 101: 147204.

The author thanks the fruitful discussions with Prof. Koji Hukushima of University of Tokyo and Prof. Hidetoshi Nishimori of Tokyo Institute of Technology. This work was partially

Optimization by Use of Nature in Physics Beyond Classical Simulated Annealing 63

*Department of Systems Science, Graduate School of Informatics, Kyoto University, Yoshida-Honmachi,*

[1] Campisi, M., Talkner, P. & Hänggi, P. [2009]. Fluctuation theorem for arbitrary open

[2] Crooks, G. E. [1998]. Nonequilibrium measurements of free energy differences for

[3] Crooks, G. E. [1999]. Entropy production fluctuation theorem and the nonequilibrium

[4] Crooks, G. E. [2000]. Path-ensemble averages in systems driven far from equilibrium,

[5] Farhi, E., Goldstone, J., Gutomann, S. & Sipser, M. [2000]. Quantum Computation by

[6] Geman, S. & Geman, D. [1984]. Stochastic relaxation, gibbs distributions, and the bayesian restoration of images, *Pattern Analysis and Machine Intelligence, IEEE Transactions*

[7] Hukushima, K. & Iba, Y. [2003]. Population annealing and its application to a spin glass,

[8] Hukushima, K. & Nemoto, K. [1996]. Exchange monte carlo method and application to spin glass simulations, *Journal of the Physical Society of Japan* 65(6): 1604–1608.

[9] Iba, Y. [2001]. Population monte carlo algorithms, *Trans. Jpn. Soc. Artif. Intel.*

[10] Jarzynski, C. [1997a]. Equilibrium free-energy differences from nonequilibrium

[11] Jarzynski, C. [1997b]. Nonequilibrium equality for free energy differences, *Phys. Rev. Lett.*

[12] Jörg, T., Krzakala, F., Kurchan, J. & Maggs, A. C. [2008]. Simple glass models and their

measurements: A master-equation approach, *Phys. Rev. E* 56: 5018–5035.

microscopically reversible markovian systems, *J. Stat. Phys.* 90: 1481–1487.

work relation for free energy differences, *Phys. Rev. E* 60: 2721–2726.

supported by MEXT in Japan, Grant-in-Aid for Young Scientists (B) No.24740263.

## **5. Conclusion**

We reviewed several recent active studies on the alternatives of SA with use of the novel substantial progress in statistical physics. The key point was to exploit the nonequilibrium behavior during performing the active control on the system. The Jarzynski equality states the possibility to estimate the equilibrium quantities by the average quantity through the nonequilibrium behavior. It means that we can invent several new strategies by use of JE away from the paradigm of the simulated annealing, which sticks to the quasi-static dynamics. The population annealing is a starting point of the studies in this direction. It is certain that population annealing find out the desired solution of the optimization problem from the property of JE faster than SA. Roughly speaking, the population annealing cuts the computational time (CPU) by use of the parallel dynamics (memory). The remaining problem is to evaluate its qualitative performance of PA for the optimization problem.

Not only the direct use of PA, we propose another type of its application in this chapter. Regarding on this type, we show skillful analyses by use of the special symmetry hidden in spin glasses to give several nontrivial exact relations. The resultant relations are useful to investigate the low-temperature region for spin glasses if we implement them by aid of PA, since we do not suffer from the critical slowing down peculiar in spin glass.

Meanwhile, if we employ a different rule to drive the system, we would be able to find the way to solve the optimization problem as SA. We reviewed QA, which was by use of quantum fluctuations as a driver. The specialized version of QA, QAC, is found to has a crucial bottleneck to solve a part of the optimization problem. Therefore we need to remove this problem while keeping its generality as a solver of the optimization problem. We again considered the application of JE to propose an alternative method, NQA. Although we do not assess its quantitative performance in the application to the optimization problem, our proposal gives a new paradigm to solve the optimization problem through the physical process like SA. We have to emphasize that QA was invented to solve the optimization problem in quantum computer. Therefore we must prepare the quantitative results to verify the precision and conditions in the actual experience on quantum computers. Along this line, we gave several results for the nonequilibrium behavior in the quantum system with gauge symmetry. These studies would be significant in the future development to realize the quantum computation.

Beyond the original version of SA, in order to find the desired solution as fast as possible, we must be away from the quasi-static procedure. The key point is to deal with nonequilibrium behavior. The further understanding of its peculiar behavior in statistical physics would be helpful to invent a genius and generic solver as PA and NQA.

## **Acknowledgement**

22 Will-be-set-by-IN-TECH

results would be a valuable platform to confirm their precisions and conditions. The quantum annealing is originally intended to be a tool implemented in so-called quantum computer. Therefore we need theoretical studies not only on the performance in the application but also how to make good condition to implement the protocol. In this regard, our studies shed light on the future indicator to open the way to realize the generic solver in the quantum computer.

We reviewed several recent active studies on the alternatives of SA with use of the novel substantial progress in statistical physics. The key point was to exploit the nonequilibrium behavior during performing the active control on the system. The Jarzynski equality states the possibility to estimate the equilibrium quantities by the average quantity through the nonequilibrium behavior. It means that we can invent several new strategies by use of JE away from the paradigm of the simulated annealing, which sticks to the quasi-static dynamics. The population annealing is a starting point of the studies in this direction. It is certain that population annealing find out the desired solution of the optimization problem from the property of JE faster than SA. Roughly speaking, the population annealing cuts the computational time (CPU) by use of the parallel dynamics (memory). The remaining problem

Not only the direct use of PA, we propose another type of its application in this chapter. Regarding on this type, we show skillful analyses by use of the special symmetry hidden in spin glasses to give several nontrivial exact relations. The resultant relations are useful to investigate the low-temperature region for spin glasses if we implement them by aid of PA,

Meanwhile, if we employ a different rule to drive the system, we would be able to find the way to solve the optimization problem as SA. We reviewed QA, which was by use of quantum fluctuations as a driver. The specialized version of QA, QAC, is found to has a crucial bottleneck to solve a part of the optimization problem. Therefore we need to remove this problem while keeping its generality as a solver of the optimization problem. We again considered the application of JE to propose an alternative method, NQA. Although we do not assess its quantitative performance in the application to the optimization problem, our proposal gives a new paradigm to solve the optimization problem through the physical process like SA. We have to emphasize that QA was invented to solve the optimization problem in quantum computer. Therefore we must prepare the quantitative results to verify the precision and conditions in the actual experience on quantum computers. Along this line, we gave several results for the nonequilibrium behavior in the quantum system with gauge symmetry. These studies would be significant in the future development to realize the

Beyond the original version of SA, in order to find the desired solution as fast as possible, we must be away from the quasi-static procedure. The key point is to deal with nonequilibrium behavior. The further understanding of its peculiar behavior in statistical physics would be

helpful to invent a genius and generic solver as PA and NQA.

is to evaluate its qualitative performance of PA for the optimization problem.

since we do not suffer from the critical slowing down peculiar in spin glass.

We must continue the active studies in this direction.

**5. Conclusion**

quantum computation.

The author thanks the fruitful discussions with Prof. Koji Hukushima of University of Tokyo and Prof. Hidetoshi Nishimori of Tokyo Institute of Technology. This work was partially supported by MEXT in Japan, Grant-in-Aid for Young Scientists (B) No.24740263.

## **Author details**

Masayuki Ohzeki

*Department of Systems Science, Graduate School of Informatics, Kyoto University, Yoshida-Honmachi, Sakyo-ku, Kyoto, 606-8501, Japan Japan*

## **6. References**


[12] Jörg, T., Krzakala, F., Kurchan, J. & Maggs, A. C. [2008]. Simple glass models and their quantum annealing, *Phys. Rev. Lett.* 101: 147204. URL: *http://link.aps.org/doi/10.1103/PhysRevLett.101.147204*

	- [13] Kadowaki, T. & Nishimori, H. [1998]. Quantum annealing in the transverse ising model, *Phys. Rev. E* 58: 5355–5363. URL: *http://link.aps.org/doi/10.1103/PhysRevE.58.5355*

**SA Applications** 


URL: *http://link.aps.org/doi/10.1103/PhysRevLett.105.050401*

