Mean field analysis of reverse annealing for code-division multiple-access multiuser detection
aa r X i v : . [ c ond - m a t . d i s - nn ] A p r Mean field analysis of reverse annealing for code-division multiple-access multiuser demodulator
Shunta Arai , ∗ Masayuki Ohzeki , , , and Kazuyuki Tanaka Graduate School of Information Sciences,Tohoku University, Sendai 980-8579, Japan Institute of Innovative Research,Tokyo Institute of Technology, Nagatsuta-cho 4259,Midori-ku, Yokohama, Kanagawa 226-8503, Japan Sigma-i Co., Ltd., Tokyo, Japan (Dated: April 24, 2020)In this study, we evaluate the typical reverse annealing (RA) performance of the code-division multiple-access(CDMA) multiuser demodulator by means of statistical mechanics using the replica method. If we prepare forthe proper initial states, first-order phase transition, which is troublesome in estimating the original signals, canbe avoided or mitigated. As we increase the di ffi culty to obtain the ground state, further information regardingthe ground state of the original problems is required to avoid first-order phase transition. In our theoreticalanalysis, we assume replica symmetry and static approximation. To verify our analytical results, we performquantum Monte Carlo simulations. The analytical results are consistent with the numerical results, except for theintermediate values of the annealing parameter. Moreover, we derive the Almeida–Thouless (AT) condition forthe CDMA model in RA. In our problem settings, the AT condition holds. The deviation between the analyticaland numerical results is owing to the breaking of the static approximation. Thus, to investigate the e ff ects ofthe static approximation, we perform RA without quantum fluctuations. In this case, the numerical results agreewith the analytical results. This study is the first analytical demonstration of the application of RA to practicalinference problems. I. INTRODUCTION
The code-division multiple-access (CDMA) multiuser de-modulator has been used in various communication systems[1]. Furthermore, the theoretical performance of CDMA mul-tiuser detection has been analyzed by means of statistical me-chanics [2–5]. CDMA multiuser detection is regarded as atype of signal recovery problem, similar to compressed sens-ing [6–9]. Statistical–mechanical analyses for signal recov-ery problems focus on the inference of the original informa-tion from the degraded information with noise. The noise canbe physically regarded as thermal fluctuations. By tuning thestrength of the thermal fluctuations, the original signals can beestimated from the degraded ones.In addition to thermal fluctuations, quantum fluctuationsmay be used to estimate the signals. Several studies havedemonstrated that quantum fluctuations such as the transversefield do not necessarily improve the performance of the infer-ences for image restoration, Sourlas codes, and CDMA [10–12]. The optimal decoding performance with quantum fluctu-ations is inferior to that with thermal fluctuations in Bayes op-timal cases. However, in certain non-Bayes optimal cases; forexample, where a lower temperature than the true noise levelis set, the decoding performance with finite quantum fluctu-ations and thermal fluctuations is superior to that with onlythermal fluctuations. This implies the potential of the combi-nation of quantum and thermal fluctuations for inference prob-lems.The performance of an optimization algorithm with quan-tum fluctuations, which is known as quantum annealing (QA)[13–18] or adiabatic quantum computation (AQC) [19, 20], is ∗ [email protected] equal to or better than that of an optimization algorithm withthermal fluctuations [21, 22], which is known as simulated an-nealing [23]. The physical implementation of QA is realizedby the quantum annealer [24–28]. The quantum annealer hasbeen implemented in numerous applications, such as portfoliooptimizations [29, 30], biological problems [31–33], electionforecasting [34], tra ffi c optimization [35], item listing for E-commerce [36], automated guided vehicles in factories [37],machine learning [38–42], quantum simulation [43–45], andmaterial design [46].In a closed system, QA begins from the ground state of thetransverse field term and the transverse field strength is grad-ually reduced. Following the Schrdinger equation, the triv-ial ground state evolves adiabatically into a nontrivial groundstate, which corresponds to the solution of combinatiorial op-timization problems. The quantum adiabatic theorem guaran-tees a theoretically su ffi cient condition to obtain the groundstate in QA [47]. The theorem indicates that the total compu-tational time for obtaining the ground state is characterized bythe minimum energy gap between the ground state and firstexited state. The energy gap is related to the phase transitionorder. In the case of fist-order phase transition, the computa-tional time for searching the ground state increases exponen-tially [48–51], which is the worst case of QA.Although the quantum annealer can ideally solve combina-torial optimization problems e ffi ciently, except for the worstcase, it often fails to obtain the ground state owing to thermalfluctuations and noise [52–54]. The outputs provided by thequantum annealer follow the Gibbs–Boltzmann distribution ata finite temperature [55]. In certain cases, thermal fluctuationscan aid in obtaining the ground state [56, 57]. In the quantumannealer, the annealing schedule can be changed by utilizing“pause” or “quench.” By changing the annealing schedule, wecan enhance the success probability to obtain the ground state[58, 59].Another useful implementation in the quantum annealer isreverse annealing (RA), which can mainly be classified intotwo methods: adiabatic RA (ARA) and iterated RA (IRA).ARA was proposed as a heuristic algorithm in the context ofthe AQC for the satisfiability problem [60]. The developedmethod was proposed as a local search algorithm with quan-tum fluctuations [61]. The procedure of the IRA may be fea-sible in the current quantum annealer. Recent research [62]demonstrated that the IRA is useful for open systems.In this research, we focus on ARA, the procedure of whichis outlined as follows:We start from the initial Hamiltonian, the ground state ofwhich is a candidate solution that is su ffi ciently close to theground state of the original problem we wish to solve. Next,we gradually increase the e ff ects of the quantum fluctuationsand search locally around the candidate solution. Thereafter,we gradually decrease the e ff ects of the quantum fluctuations.When the e ff ects of the quantum fluctuations disappear, theground state or lower energy state of the original problem canbe obtained. Theoretical analysis of the ARA for the p -spinmodel was carried out with respect to the static property anddynamics [63, 64]. However, to the best of the authors’ knowl-edge, it remains unknown whether or not the ARA is usefulfor certain practical problems.In this study, we investigate the e ffi ciency of the ARA forpractical problems. We apply the ARA to the CDMA mul-tiuser detection, which is a representative example in signalrecovery problems. The CDMA model can be expressed asthe quadratic unconstrained binary optimization model, with2-body interactions and the random field. The CDMA modelcan be embedded into the quantum annealer at the expenseof additional physical qubits. The CDMA model is mainlycharacterized by the pattern ratio of the number of users tothat of the measurements. In the low-temperature regions,the CDMA model has two solutions. This phenomenon re-veals the existence of the first-order phase transition, whichdegrades the demodulation e ffi cacy. We use ARA to miti-gate or avoid the demodulation di ffi culty. We set the initialHamiltonian in the ARA process. The initial Hamiltonian isinterpreted as prior information of the estimated signals in thecontext of the inference problem by using the probabilisticmodel. We expect that the prior information of the problemwill mitigate the demodulation di ffi culty.In this study, we consider the marginal posterior mode(MPM) estimation. The estimated signals correspond to theexpectation of the signals over the Gibbs–Boltzmann distribu-tion. We analyze the average demodulation performance withARA at a finite temperature using the replica method. The dif-ficulty of obtaining the original signals depends on the patternratio and noise. We consider the case in which the noise is lowand the pattern ratio is su ffi ciently small. If the ground statesof the initial Hamiltonian are su ffi ciently close to the originalground states, we can avoid the first-order phase transition. Inthe low-temperature regions, the classical CDMA model ex-hibits replica symmetry breaking (RSB) [5]. To investigate thee ff ects of the quantum fluctuations and initial Hamiltonian onthe stability of the replica symmetry (RS) solutions, we de-rive the Almeida–Thouless (AT) condition [65] in the ARA. Although the implementation of the ARA in the quantum an-nealer has not yet been realized, our results provide the firstdemonstration of the ARA as a practical technique for signalrecovery problems.The remainder of this paper is organized as follows. In Sec-tion II, we present the formulation of the CDMA model withquantum fluctuations. In Section III, we extend the formula-tion for the ARA. We calculate the partition function underthe RS ansatz and static approximation. We derive the saddle-point equations and stability condition of the RS solutions. InSection IV, we numerically solve the saddle-point equationsand illustrate the phase diagrams with and without the ARA.To verify our theoretical analysis, we perform quantum MonteCarlo simulations. Finally, we summarize our results and dis-cuss the future research directions in Section V. II. CDMA MODEL WITH QUANTUM FLUCTUATIONS
Firstly, we formulate the classical CDMA model and moveonto the quantum system. The main concept of the CDMAmodel is as follows: The digital signals of users are modulatedand transmitted to a base station through fully synchoronouschannels. By demodulating the received signals composed ofthe multiuser signals and noises, we infer the original signalsfrom the provided information.We consider that N users communicate via fully syn-chronous channels. At the base station, the receiver obtainsthe signals as follows: y k = √ N N X i = η ki ξ i + ǫ k , (1)where ξ i ∈ {± } , ( i = , . . . , N ) is the original information and η ki ∈ {± } ( i = , . . . , N , k = , . . . , K ) is the spreading codefor each user i . The length of the spreading codes for each user i is represented by K . The channel noise ǫ k is added into thereceived signals. The received signals (1) can be expressed as y = √ N A ξ + ǫ , (2)for which the following notations are used: y = (cid:16) y , . . . , y K (cid:17) T , ξ = ( ξ , . . . , ξ N ) T , ǫ = (cid:16) ǫ , . . . , ǫ K (cid:17) T , (3) A = η η · · · η N η η · · · η N ... ... . . . ...η K η K · · · η KN . (4)We assume that the spreading codes and original signals areindependently generated from the uniform distribution: P ( A ) = NK , (5) P ( ξ ) = N . (6)We consider the Gaussian channels and ǫ k is independentlygenerated from the Gaussian distribution as follows: P ( ǫ ) = P ( y | ξ ) = q πσ K exp − σ || ǫ || = r β π K exp − β (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) y − A ξ √ N (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) , (7)where β = /σ is the true noise level.In the CDMA multiuser detection, we estimate the origi-nal signals from the received output signals and the spreadingcodes that are prepared for each user in advance. Becausethe output signals fluctuate owing to noise, we formulate thisproblem as Bayesian inference. Subsequently, we introducethe posterior distribution as P ( σ | y ) = P ( y | σ ) P ( σ )Tr P ( y | σ ) P ( σ ) , (8)where the estimated signals are given by σ = ( σ , . . . , σ N ) T .The likelihood is expressed as P ( y | σ ) = r β π K exp − β (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) y − A σ √ N (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) , (9)where β = / T is the inverse temperature in statistical me-chanics and corresponds to the estimated channel noise level.If the true noise level is known, the demodulation performanceis the best and Bayes optimal. According to Eqs. (8) and(9), the posterior distribution can be written using the Gibbs–Boltzmann distribution with the Hamiltonian H ( σ ), as fol-lows: P ( σ | y ) = Z exp {− β ( H ( σ ) + H init ( σ )) } , (10) Z = Tr exp {− β ( H ( σ ) + H init ( σ )) } , (11) H ( σ ) = N X i , j K X k = η ki η kj σ i σ j − √ N N X i = K X k = η ki y k σ i , (12)where Z is the partition function and H init ( σ ) is the initialHamiltonian, which represents the prior information of theestimated signals. We generally assume that the prior of theestimated signals follows the uniform distribution P ( σ ) = N . (13)In this case, we can omit the initial Hamiltonian from Eqs.(10) and (11) . We adopt the MPM estimation to estimate theoriginal signals. The estimation performance can be evalu-ated by the overlap between the original and estimated sig-nals as m = / N P Ni = ξ i h σ i i , where h·i is the expectation overthe posterior distribution P ( σ | y ) [66]. The overlap is phys-ically interpreted as a magnetization. This quantity is ex-pected to exhibit a “self-averaging” property in the thermo-dynamics limit N → ∞ . This means that the observables, such as the overlap for a quenched realization of the data y , A , and ξ , are equivalent to the expectation of itself over thedata distribution. In this case, the overlap can be expressedas lim N →∞ m = [ ξ i h σ i i ], where the bracket [ · ] indicates theexpectation over the data distribution P ( A ) P ( ξ ) P ( y | ξ ).It is straightforward to extend the above formulation intothe quantum mechanical version:ˆ H = s ˆ H + (1 − s ) ˆ H TF , (14)ˆ H = N X i , j K X k = η ki η kj ˆ σ zi ˆ σ zj − √ N N X i = K X k = η ki y k ˆ σ zi , (15)ˆ H TF = − N X i = ˆ σ xi , (16)where ˆ σ zi and ˆ σ xi are the z and x components of the Pauli ma-trices at site i , respectively. In this case, ˆ H consists of the z components of the Pauli matrices and ˆ H TF is composed ofthe x components of the Pauli matrices. We parameterize theHamiltonian (14) with the annealing parameter s for appli-cation to the ARA. As in the classical case, we consider theMPM estimation with quantum fluctuations. The performanceof the MPM estimation with quantum fluctuations can be eval-uated by m = / N P Ni = ξ i (cid:16) Tr ˆ σ zi ˆ ρ (cid:17) , where the density matrix isˆ ρ ≡ e − β ˆ H / Tr e − β ˆ H . III. MEAN FIELD ANALYSIS FOR ARA
Following Ref. [63], we formulate the CDMA model withquantum fluctuations in the ARA. In the ARA, we have thefollowing Hamiltonian:ˆ H = s ˆ H + (1 − s )(1 − λ ) ˆ H init + (1 − s ) λ ˆ H TF , (17)ˆ H init = − N X i = τ i ˆ σ zi , (18)where λ (0 ≤ λ ≤
1) is the RA parameter. In the conven-tional QA, we monotonically increase the annealing parame-ter s from s = s = λ =
1. In the ARA, we initiallyset s = λ =
0. The ground state of Eq. (17) is the ground stateof Eq. (18), ˆ σ zi = τ i ( ∀ i ), where τ i = ± ξ i .Next, we increase the parameters s and λ towards s = λ = P ( σ z | τ ) ∝ exp (cid:16) − β ˆ H init (cid:17) .The typical behaviors of the order parameters such as theoverlap can be obtained via the free energy. We calculate thepartition function and derive the RS free energy in the limit of N , K → ∞ , while maintaining the pattern ratio α ≡ K / N = O (1). We assume the “self-averaging” property in the thermo-dynamic limit.The free energy per user f can be evaluated as − β f = lim N →∞ / N [ln Z ]. Firstly, we employ the Suzuki–Trotter decomposition[67] in the partition function: Z = Tr exp (cid:16) − β ˆ H (cid:17) = lim M →∞ Tr (cid:26) exp (cid:18) − β M (cid:16) ˆ H + ˆ H init (cid:17)(cid:19) exp (cid:18) − β M ˆ H TF (cid:19)(cid:27) M = lim M →∞ Z M , (19)where Z M = Tr M Y t = exp − β s N M X i , j K X k = η ki η kj σ zit σ zjt + β sM √ N N X i = K X k = η ki y k σ zit + β (1 − s )(1 − λ ) M N X i = τ i σ zit + β (1 − s ) λ M N X i = σ xit , (20)in which the symbol t is the index of the Trotter slice and M is the Trotter number. To evaluate [ln Z M ], we use the replica method[68]: [ln Z M ] = lim n → [ Z nM ] − n . (21)We introduce the order parameters as follows: m t ( a ) = N N X i = ξ i σ zit ( a ) , (22) q tt ′ ( a , b ) = N N X i = σ zit ( a ) σ zit ′ ( b ) , (23) R tt ′ ( a ) = N N X i = σ zit ( a ) σ zit ′ ( a ) , (24) m xt ( a ) = N N X i = σ xit ( a ) . (25)Moreover, we introduce the auxiliary parameters ˜ m t ( a ) , ˜ q tt ′ ( a , b ) , ˜ R tt ′ ( a ) , ˜ m xt ( a ) of the order parameters. Under the RS ansatz andstatic approximation: m t ( a ) = m , q tt ′ ( a , b ) = q , R tt ′ ( a ) = R , m xt ( a ) = m x , ˜ m t ( a ) = ˜ m , ˜ q tt ′ ( a , b ) = ˜ q , ˜ R tt ′ ( a ) = ˜ R , ˜ m xt ( a ) = ˜ m x , we canfinally obtain the RS free energy: − β f RS = α −
12 ln(1 − β s ( q − R )) + β s R − + β s + β β + m − q − (1 + β − )1 − β s ( q − R ) + β (1 − s ) λ m x − m ˜ m − m x ˜ m x − R ˜ R + q ˜ q + c Z Dz ln Z Dy q g + + ( ˜ m x ) + (1 − c ) Z Dz ln Z Dy q g − + ( ˜ m x ) , (26)where g ± = ( ˜ m ± β (1 − λ )(1 − s )) + p ˜ qz + q R − ˜ qy , (27)in which Dz means that the Gaussian measure Dz : = / √ π dze − z / and Dy is the same as Dz . The number c (0 ≤ c ≤
1) denotesthe fraction of the ground state τ i = ξ i in the initial state as follows: c = N N X i = δ τ i ,ξ i . (28)The detailed calculations for the derivation the RS free energy in Eq. (26) are presented in Appendix A. The extremization ofEq. (26) yields the following saddle-point equations: m = c Z DzY − + Z Dy g + u + ! sinh u + + (1 − c ) Z DzY − − Z Dy g − u − ! sinh u − , (29) q = c Z Dz ( Y − + Z Dy g + u + ! sinh u + ) + (1 − c ) Z Dz ( Y − − Z Dy g − u − ! sinh u − ) , (30) R = c Z DzY − + Z Dy ( ˜ m x ) u + ! sinh u + + g + u + ! cosh u + + (1 − c ) Z DzY − − Z Dy ( ˜ m x ) u − ! sinh u − + g − u − ! cosh u − , (31) m x = c Z DzY − + Z Dy ˜ m x u + ! sinh u + + (1 − c ) Z DzY − − Z Dy ˜ m x u − ! sinh u − , (32)˜ m = αβ s + β s ( R − q ) , (33)˜ q = αβ s (cid:16) q − m + + β − (cid:17) (1 + β s ( R − q )) , (34)2 ˜ R − ˜ q = αβ s ( R − q )1 + β s ( R − q ) , (35)˜ m x = β (1 − s ) λ, (36) Y ± ≡ Z Dy cosh u ± , (37) u ± ≡ q g ± + ( ˜ m x ) . (38)Next, we consider the stability of the RS solutions. Two instabilities exist in the RS solutions: the local instability and globalinstability of the RS solutions. The local stability condition of the RS solutions under the static approximation is expressed as α c β s (1 + β s ( R − q )) Z Dz Y − + Z Dy g + u + ! sinh u + ! − Y − + Z Dy ( ˜ m x ) u + ! sinh u + + Z Dy g + u + ! cosh u + + α (1 − c ) β s (1 + β s ( R − q )) Z Dz Y − − Z Dy g − u − ! sinh u − ! − Y − − Z Dy ( ˜ m x ) u − ! sinh u − + Z Dy g − u − ! cosh u − < . (39)This condition corresponds to the AT condition in the ARA. We can achieve this condition by taking into account the perturba-tions to the RS solutions. This result is consistent with the previous result in Ref. [5] for the classical limit s = λ =
1. Thedetailed calculations for deriving the AT condition in Eq. (39) are presented in Appendix B. The global instability condition ofthe RS solutions is related to the negative entropy. The existence of the global instability corresponds to the freezing behavior[69]. To detect the freezing behavior, we calculate the RS entropy as follows: S = − ∂∂ T f RS = − α { ln (1 + β s ( R − q )) } + R − q m − ˜ q ) ++ ˜ RR − q ˜ q + c Z Dz ln Y + + (1 − c ) Z Dz ln Y − − β ( c Z DzY − + Z Dyu + sinh u + + (1 − c ) Z DzY − − Z Dyu − sinh u − ) . (40)In the case of s = λ =
1, this result is also consistent with the classical one.
IV. EXPERIMENTAL RESULTS
We numerically solve the saddle-point equations in Eqs.(29) to (36) with T = . λ =
1. The phase diagrams for T = , .
05 and 0 . α s spinodal 1spinodal 2critical (a) α s spinodal 1spinodal 2critical (b) α s spinodal 1spinodal 2 (c) FIG. 1. Phase diagram of CDMA model with quantum fluctuations. The horizontal axis denotes the pattern ratio. The vertical axis denotes theannealing parameter. The experimental settings are (a) T =
0, (b) T = .
05, and (c) T = .
1. The “spinodal 1” and “spinodal 2” lines denotethe solutions from the two di ff erent branches. The “critical” line denotes the point at which the RS free energy takes the same value. α m branch 1branch 2QMC (a) α m x branch 1branch 2QMC (b) α R branch 1branch 2QMC (c) FIG. 2. Dependence of the order parameters on the pattern ratio for the fixed annealing parameter s = .
9. The vertical axes denote theseorder parameters: (a) magnetization, (b) transverse magnetization, and (c) correlation between Trotter slices. The blue solid line and orangedashed line denote the two di ff erent branches that are obtained from the saddle-point equations. The circles represent the results obtained bythe quantum Monte Carlo simulations. green dotted line denotes the critical point at which the RSfree energy takes the same value. As we cannot distinguishthe critical point from the spinodal points in this scale, we donot write down the line in Fig. 1(c). Higher noise results ina narrower region in which the two solutions coexist. Withthese problem settings, the AT condition is not broken and thefreezing behavior does not occur.To verify the RS ansatz and static approximation, we per-form quantum Monte Carlo simulations using the CDMAmodel without RA, which means that we fix the RA parame-ter as λ =
1. We set the system size as N = M =
50, the temperature as T = .
1, and the truenoise scale as T =
0. We use a 100000 Monte Carlo step(MCS) average after 50000 MCS equilibrations for each in-stance. We take the configuration average over the spreadingcodes and the original signals by randomly generating 50 in-stances. The error bar is given by the standard deviation. Weplot the behavior of the magnetization with respect to the pat-tern ratio for the fixed annealing parameter s = . α = . ff erences between the numerical andanalytical results in the low pattern ratio probably result fromthe breaking of the static approximation. The e ff ects of thequantum fluctuations are underestimated with the static ap-proximation because we neglect the dependence of the orderparameters on the Trotter slices. According to Fig. 3, thenumerical results for the magnetization and transverse mag-netization are consistent with the analytical results, except forthe intermediate values of the annealing parameter. As in Fig.3(c), the correlation between the Trotter slices obtained bythe saddle-point equations is overestimated owing to the staticapproximation [70] . According to Figs. 2 and 3, when the an-nealing parameter is large and close to 1 (the transverse fieldstrength is weak), the numerical results are consistent with theanalytical results.Moreover, we investigate the AT stability against highertrue noise than the temperature T under the transverse field.In Fig. 4, the blue solid line denotes the AT line, above whichthe AT condition in Eq. (39) is broken. The orange dash-dotted line is the zero-entropy line S =
0, above which the s m branch 1branch 2QMC (a) s m x branch 1branch 2QMC (b) s R branch 1branch 2QMC (c) FIG. 3. Dependence of order parameters on annealing parameter for fixed pattern ratio α = .
6. The same symbols as those in Fig. 2 are used. α s AT linezero-entropy linem max m = m classical (a) α s AT linezero-entropy linem max m = m classical (b) α s AT linezero-entropy linem max m = m classical (c)
FIG. 4. The blue solid line denotes the AT line; the orange dash-dotted line denotes the zero-entropy line; the green dashed line denotes themaximum of the magnetization for each pattern ratio; and the red dotted line represents the point at which the magnetization value is equal tothe classical one s =
1. The experimental settings are (a) T = . T = .
2, (b) T = . T = .
5, and (c) T = .
05 and T = .
2. Bothaxes are the same as those in Fig. 1.
RS entropy is negative. We can regard these solutions as non-physical solutions. The green dashed line is the maximum ofthe magnetization. The red dotted line denotes the point atwhich the magnetization is equal to the classical value s = T = .
1, and the true noise level as T = . T = . T = .
05 whilemaintaining the true noise level at T = .
2. Firstly, the globalstability of the RS solutions is broken at the low pattern ratioaround α ≃ .
566 when the annealing parameter is increased.According to Fig. 4, the maximum values of the magnetiza-tion are above the AT line. When the AT condition is broken,the maximum values of the magnetization are at s =
1. In Fig.4(c), the red dotted line is under the AT line around α = . ff ects of the quantumfluctuations.Subsequently, we apply the ARA to the CDMA model withquantum fluctuations. The experimental settings are the sameas those in Fig.1(a). Figure 5 presents the phase diagram of theCDMA model in the ARA for α = . , .
5, and 0 .
4. We con- sider four initial conditions: c = . , . , .
9, and 0 .
95. Eachline represents a point of the first-order phase transition. Thedi ffi culty of estimating the original signals is related to the freeenergy landscape. It is easy to estimate the original signals onthe right side of spinodal line 2 in Fig. 1(a), because the freeenergy exhibits a minimum, which is a good estimator. Whenwe set the pattern ratio as α = .
6, first-order phase transitionexists. The free energy landscape has two valleys. At spin-odal line 2, the free energy landscape is transformed into thesimple valley. In this case, it is comparatively easy to estimatethe original signals. For α = .
5, spinodal line 2 does not ex-ist. The free energy landscape maintains two valleys. In thiscase, the time for searching the original signals increases ex-ponentially. If the proper initial conditions are known, whichmeans that the initial states are close to the original signals,it is easy to search them. The critical point does not exist for α = .
4. The minima of the free energy do not provide uswith an e ff ective estimation. In this case, it is most di ffi cultto obtain the original signals. We can observe from Figs. 5(a)and 5(b) that the first-order phase transition can be avoidedif the proper initial conditions are provided. As the informa-tion regarding the original signals is increased, the region foravoiding the first-order phase transition is broadened. Thismeans that the di ffi culty of estimating the original signals ismitigated by introducing prior information regarding the orig-inal signals. In Fig. 5(b), the region in which the first-orderphase transition can be avoided is narrower than that in Fig. λ s c=0.95c=0.9c=0.8c=0.7 s (a) λ s c=0.95c=0.9c=0.8c=0.7 s (b) λ s c=0.7 spinodalc=0.8 criticalc=0.8 spinodal s (c) λ s spinodalcritical s (d) λ s spinodalcritical s (e) FIG. 5. Phase diagrams of CDMA model in ARA for four di ff erent values of c . The horizontal axis denotes the RA parameter. The vertical axisdenotes the annealing parameter. These lines represent the points at which the first-order phase transitions occur. The experimental settingsare (a) α = .
6, (b) α = .
5, (c) α = . c = .
7, and c = .
8, (d) α = . c = .
9, and (e) α = . c = . ffi culty in estimating the original signals. In this case, fur-ther information regarding the original signals is required toavoid the first-order phase transition. From Figs. 5(c) to 5(e),we can observe di ff erent behaviors to those of Figs.5(a) and5(b). In Fig. 5(c), we cannot avoid the first-order phase tran-sition. The spinodal lines are the same as spinodal 1. For c = .
8, the dashed line represents the point at which the first-order phase transition occurs. Between the dashed lines, theminima of the free energy are e ff ective estimators. In otherregions, the minima of the free energy are ine ff ective estima-tors. In Figs. 5(d) and 5(e), the qualitative behaviors of thesystems can be separated by the dotted line. On the right sideof the dotted line, the dash-dotted lines represent the spinodal1 line. In this region, e ff ective estimators only exist betweenthe dashed critical lines. The dashed lines are the same as inFig. 5(c). On the left side of the dotted line, there exists a re-gion in which the first-order phase transition can be avoided.We also investigate the stability of the RS solutions and findthat the AT stability holds in these problem settings.To analyze the extent to which the di ffi culty in obtaining theoriginal signals is mitigated by the ARA, we plot the di ff er-ences in the magnetization ∆ m between the two local minimaat the first-order phase transition in the case of α = . . ff erences in the magnetization result in thetwo local minima of the free energy being further. Figure 6 in-dicates that ∆ m decreases as c increases. Even in the case ofthe low pattern ratio α = . ∆ m is smaller than that of the vanilla QA λ =
1. By using the ARA, the two local minima ofthe free energy become closer. This result demonstrates thatthe ARA enhances the e ff ects of the quantum tunneling if weprepare for su ffi cient proper initial conditions.We consider the noise e ff ects for the CDMA model in theARA. The experimental settings are the same as those illus-trated in Fig. 1(b). Figure 7 displays the phase diagrams ofthe CDMA in the ARA for α = . , .
57, and 0 .
5. Thequalitative behaviors of the systems are approximately thesame as those in the noiseless case. The regions in whichthe first-order phase transition can be avoided are larger thanthose of the noiseless cases because the first-order phase tran-sition is weakened owing to the noise e ff ects. Furthermore,the RS solutions are stable in these problem settings. Figure8 presents the di ff erences in the magnetization between thetwo local minima at the first-order phase transition in the caseof α = .
62 and 0 .
57. Although the qualitative behaviors of ∆ m are the same as those in Fig. 6, ∆ m is smaller than in thenoiseless case.To validate the replica analysis results, we perform quan-tum Monte Carlo simulations. The experimental settings arethe same as those in Fig. 3. We set the RA parameter as λ = .
8, and the initial conditions as c = . .
9. Fig-ure 9 indicates that the numerical results for the magnetiza-tion and transverse magnetization are consistent with the an-alytical results, except for the intermediate values of the an-nealing parameter. The numerical results for the correlationbetween the Trotter slices do not agree with the analytical re-sults, except for the region in which the e ff ects of the quan- λ Δ m c=0.95c=0.9c=0.8c=0.7 (a) λ Δ m c=0.95c=0.9c=0.8c=0.7 (b) FIG. 6. Di ff erences in magnetization between two local minima at first-order phase transition in Figs. 5(a) and 5(b). The vertical axis denotesthe di ff erences in the magnetization between the two local minima at the first-order phase transition. The horizontal axis denotes the RAparameter. The experimental settings are (a) α = . α = . tum fluctuations are weak. The order parameters for c = . c = .
9, it can be observed that the first-order phasetransition can be avoided, as in Figs. 9(d) to 9(f). Althoughthe order parameter behaviors from the numerical simulationsdo not completely match the analytical results owing to thebreaking of the static approximation, the behaviors for avoid-ing the first-order phase transition do not change. Therefore,the analytical results obtained by the replica method underthe static approximation provide us with the variable predic-tion concerning avoiding the first-order phase transition in ourproblem settings.Finally, we consider the ARA without quantum fluctuationsto verify whether or not the di ff erences between the analyti-cal and numerical results arise from the breaking of the staticapproximation. When we set the RA parameter as λ =
0, thequantum part in Eq. (17) disappears. We do not have to usethe static approximation in this case. Figures 5 and 7 indi-cate that the first-order phase transition can be avoided in theRA parameter λ =
0. The experimental settings are approx-imately the same as those in Figs. 2 and 3. We set α = . T = α = .
62 and T = .
05 in Fig.10(b). We consider three initial conditions: c = .
7, 0 .
8, and0 .
95. The error bar is given by the standard deviation. Eachline represents the analytical results and each symbol denotesthe numerical results obtained by the Markov chain MonteCarlo simulations. It can be observed that the numerical re-sults are consistent with the analytical results. By applyingthe ARA without quantum fluctuations, we can also avoid thefirst-order phase transition with quantum fluctuations.
V. SUMMARY
We performed mean field analysis of the ARA for CDMAmultiuser detection. In CDMA multiuser detection, first-orderphase transition is encountered at intermediate pattern ratios.This first-order phase transition degrades the estimation per-formance. To avoid first-order phase transition, we applied theARA to the CDMA multiuser detection. If we prepare for theproper initial conditions, we can avoid the first-order phasetransition. The di ffi culty of obtaining the original signals isdistinguished by the pattern ratio. We found that, as the prob- lems become more di ffi cult, the number of ground states inthe initial Hamiltonian increases.Firstly, we presented the phase diagrams of the CDMAmodel without the ARA to determine the existence of the first-order phase transition. At the intermediate pattern ratios andunder a weak transverse field strength, this model exhibitedfirst-order phase transition. The di ffi culty of obtaining theoriginal signals is characterized by the pattern ratio and truenoise level. A lower pattern ratio results in less informationregarding the signals. A higher true noise level results in anarrower region in which the first-order phase transition ex-ists. We performed quantum Monte Carlo simulations to val-idate our results. We considered two cases. Firstly, we fixedthe annealing parameter at s = . α = . ff ects of the quantum fluctuations because we neglected thedependence on the Trotter slices.Moreover, we derived the AT condition for the CDMAmodel in the ARA. We first investigated the AT stability of theRS solutions in our numerical simulations without the ARA.0 λ s c=0.95c=0.9c=0.8c=0.7 s (a) λ s c=0.95c=0.9c=0.8c=0.7 s (b) λ s c=0.7 spinodalc=0.8 spinodalc=0.8 critical s (c) λ s spinodalcritical s (d) λ s spinodalcritical s (e) FIG. 7. Phase diagrams of CDMA model in ARA for four di ff erent values of c . Both axes are the same as those in Fig. 5. The experimentalsettings are (a) α = .
62, (b) α = .
57, (c) α = . c = .
7, and c = .
8, (d) α = . c = .
9, and (e) α = . c = . The AT condition held in our problem settings. To examinethe noise e ff ects, we considered higher noise than the temper-ature. In this case, the stability of the RS solutions was broken.In certain regions below the AT line, the magnetization underthe transverse field was larger than that in the classical case.This indicated that the instability of the RS solutions was re-covered by the e ff ects of the quantum fluctuations. In a futurestudy, we will compare the behaviors of the order parametersobtained from the replica method with those obtained fromthe numerical simulations when RSB occurs and the static ap-proximation is broken.Next, we presented the phase diagrams of the CDMAmodel in the ARA with and without noise. The first-orderphase transition could be avoided by applying the ARA. Thismeans that the free energy landscape was transformed into thesimple valley. Even if the first-order phase transition could notbe avoided, its di ffi culty was mitigated. In the ARA, the dif-ferences in the magnetization between the two local minimaat the first-order phase transition were smaller than those inthe vanilla QA. We performed quantum Monte Carlo simula-tions to validate our analysis. The numerical results were con-sistent with the analytical results, except for the intermediatevalues of the annealing parameter, owing to the breaking ofthe static approximation. The behaviors of the order parame-ters attained from the numerical simulations to avoid the first-order phase transition did not change as a result of the e ff ectsof the static approximation. The analytical results under thestatic approximation are useful for understanding the behav-iors of quantum spin-glass systems such as the CDMA model.Finally, we considered the ARA without quantum fluctuations to neglect the e ff ects of the static approximation. The numer-ical results were consistent with the analytical results. Theanalytical results were valid in our problem settings. The re-gions avoiding the first-order phase transition when using theARA with quantum fluctuations were broader than those with-out quantum fluctuations. These results reflected the e ff ects ofthe quantum tunneling. By using the ARA, the two minima ofthe free energy became closer if we prepared for the properinitial states. The ARA enhanced the probability of escap-ing the local minimum with quantum tunneling. Our resultsdemonstrate that the iterative algorithm can overcome the dif-ficulty of inference problems.In this study, we analyzed ideal cases in which we fixedthe number of ground states in the initial Hamiltonian. In realsituations, whether or not we can prepare for the proper initialconditions should be investigated. Moreover, we can considerthe sparsity and the orthogonality with respect to the spreadingcodes [72–74]. In a future study, we will analyze the e ff ects ofthe quantum fluctuations for the CDMA with such spreadingcodes. ACKNOWLEDGMENTS
S.A. was partially supported by Grants-in-Aid for JSPS Fel-lows for the Promotion of Science (No. 19J21790). M.O.was supported by KAKENHI (No. 19H01095), and theNext Generation High- Performance Computing Infrastruc-tures and Applications R & D Program by MEXT. K.T. wassupported by JSPS KAKENHI (No. 18H03303). This work1 λ Δ m c=0.95c=0.9c=0.8c=0.7 (a) λ Δ m c=0.95c=0.9c=0.8c=0.7 (b) FIG. 8. Di ff erences in magnetization between two local minima at first-order phase transition in Figs. 7(a) and 7(b). Both axes are the same asthose in Fig. 6. The experimental settings are (a) α = .
62 and (b) α = . was partly supported by JST-CREST (No. JPMJCR1402). Appendix A: DERIVATION OF FREE ENERGY
We derive the free energy under the RS ansatz and static approximation. We introduce the following terms: u k = √ N N X i = η ki ξ i , (A1) u kt ( a ) = √ N N X i = η ki σ zit ( a ) . (A2)Furthermore, we introduce the delta function and its Fourier integral representation for Eqs. (A1) and (A2). The partitionfunction can be rewritten as[ Z nM ] = NK N X A Tr ξ Tr { σ z } Z K Y k = dy k exp (cid:26) − β y k − u k ) (cid:27) K Y k = Z du k d ˜ u k π exp i ˜ u k u k − √ N N X i = η ki ξ i × Y k , a , t Z du kt ( a ) d ˜ u kt ( a )2 π exp i ˜ u kt ( a ) u kt ( a ) − √ N N X i = η ki σ zit ( a ) exp − β s M X k , a , t u kt ( a ) + β sM X k , a , t y k u kt ( a ) × exp β (1 − s )(1 − λ ) M X a , t , i τ i σ zit ( a ) + β (1 − s ) λ M X a , t , i σ xit ( a ) . (A3)We implement the expectation for the spreading codes A in the above expression as follows: L η i = NK X H K Y k = N Y i = exp − i √ N ˜ u k ξ i − X a , t ˜ u kt ( a ) σ zit ( a ) η ki ≈ exp − N X k , i ˜ u k ξ i − X a , t ˜ u kt ( a ) σ zit ( a ) = Y k , i exp − N ( ˜ u k ) − N ˜ u k ξ i X a , t ˜ u kt ( a ) σ zit ( a ) − N n X a = M X t = ˜ u kt ( a ) σ zit ( a ) − N X a < b X t , t ′ ˜ u kt ( a ) ˜ u kt ′ ( b ) σ zit ( a ) σ zit ′ ( b ) . (A4)2 s m ReplicaQMC (a) s m x ReplicaQMC (b) s R ReplicaQMC (c) s m ReplicaQMC (d) s m x ReplicaQMC (e) s R ReplicaQMC (f)
FIG. 9. Dependence of order parameters on annealing parameter for fixed RA parameter λ = .
8. The experimental settings are as follows: (a)to (c) are c = . c = .
9. The vertical axes are the same as those in Fig. 3.
We introduce the δ function and its Fourier integral representation for Eqs. (25)-(28) as follows: Y a , t Z dm t ( a ) δ m t ( a ) − N N X i = ξ i σ zit ( a ) = Y a , t Z Ndm t ( a ) d ˜ m t ( a )2 π iM exp − ˜ m t ( a ) M Nm t ( a ) − N X i = ξ i σ zit ( a ) , (A5) Y a , t , t ′ Z dR tt ′ ( a ) δ R tt ′ ( a ) − N N X i = σ zit ( a ) σ zit ′ ( a ) = Y a , t , t ′ Z NdR tt ′ ( a ) d ˜ R tt ′ ( a )2 π iM exp − ˜ R tt ′ ( a ) M NR tt ′ ( a ) − N X i = σ zit ( a ) σ zit ′ ( a ) , (A6) Y a < b , t , t ′ Z dq tt ′ ( a , b ) δ q tt ′ ( a , b ) − N N X i = σ zit ( a ) σ zit ′ ( b ) = Y a < b , t , t ′ Z Ndq tt ′ ( a , b ) d ˜ q tt ′ ( a , b )2 π iM exp − ˜ q tt ′ ( a , b ) M Nq tt ′ ( a , b ) − N X i = σ zit ( a ) σ zit ′ ( b ) , (A7) Y a , t Z dm xt ( a ) δ m xt ( a ) − N N X i = σ xit ( a ) = Y a , t Z Ndm xt ( a ) d ˜ m xt ( a )2 π iM exp − ˜ m xt ( a ) M Nm xt ( a ) − N X i = σ xit ( a ) . (A8)Finally, the partition function can be rewritten as3 s m c = 0.7c = 0.8c = 0.95 (a) s m c = 0.7c = 0.8c = 0.95 (b) FIG. 10. Dependence of magnetization on annealing parameter for fixed RA parameter λ =
0. The experimental settings are as follows: (a) is α = . T =
0, and (b) is α = .
62 and T = .
05. Both axes are the same as those in Figs. 9(a) and 9(d). [ Z nM ] = Y a , t Z Ndm t ( a ) d ˜ m t ( a )2 π iM Y a , t , t ′ Z NdR tt ′ ( a ) d ˜ R tt ′ ( a )2 π iM Y a < b , t , t ′ Z Ndq tt ′ ( a , b ) d ˜ q tt ′ ( a , b )2 π iM Y a , t Z Ndm xt ( a ) d ˜ m xt ( a )2 π iM e G + G + G , (A9) e G ≡ K Y k = Z du k d ˜ u k π Y k , a , t du kt ( a ) d ˜ u kt ( a )2 π K Y k = dy k (cid:18) β π (cid:19) exp − β (cid:16) y k − u k (cid:17) + i X k , a , t u kt ( a ) ˜ u kt ( a ) −
12 ( ˜ u k ) − ˜ u k X a , t ˜ u kt ( a ) m t ( a ) − X a < b X t , t ′ ˜ u kt ( a ) ˜ u kt ′ ( b ) q tt ′ ( a , b ) − n X a = X t , t ′ ˜ u kt ( a ) ˜ u kt ′ ( a ) R tt ′ ( a ) − β s M X k , a , t (cid:18)(cid:16) u kt ( a ) (cid:17) − y k u kt ( a ) (cid:19) , (A10) e G ≡ N Tr ξ Tr { σ z } exp M X a , t ˜ m t ( a ) N X i = ξ i σ zit ( a ) + M X a , t , t ′ ˜ R tt ′ ( a ) N X i = σ zit ( a ) σ zit ′ ( a ) + M X a < b X t , t ′ ˜ q tt ′ ( a , b ) N X i = σ zit ( a ) σ zit ′ ( b ) + β (1 − s )(1 − λ ) M X a , t , i τ i σ zit ( a ) + M X a , t ˜ m xt ( a ) N X i = σ xit ( a ) , (A11) e G ≡ exp − NM X a , t ˜ m t ( a ) m t ( a ) − NM X a , t ˜ m xt ( a ) m xt ( a ) − X t , t ′ , a NM ˜ R tt ′ ( a ) R tt ′ ( a ) − X a < b X t , t ′ NM ˜ q tt ′ ( a , b ) q tt ′ ( a , b ) + β (1 − s ) λ NM X a , t m xt ( a ) . (A12)We calculate each free energy part e G , e G , and e G . Firstly, we perform integration for u k and ˜ u k in Eq. (A10), and we canobtain e G as follows: e G = Y k , a , t Z du kt ( a ) d ˜ u kt ( a )2 π K Y k = π s β π + β Z dy k exp β + β ) iy k − X a , t ˜ u kt ( a ) m t ( a ) − X a < b X t , t ′ ˜ u kt ( a ) ˜ u kt ′ ( b ) q tt ′ ( a , b ) − n X a = X t , t ′ ˜ u kt ( a ) ˜ u kt ′ ( a ) R tt ′ ( a ) + i X a , t u kt ( a ) ˜ u kt ( a ) − β s M X a , t (cid:18)(cid:16) u kt ( a ) (cid:17) − y k u kt ( a ) (cid:19) . (A13)The integration over y k and u kt ( a ) can be performed as follows: e G = Y k Z Dv Y a , t s M πβ s exp β sv M + β β − M β s ( ˜ u kt ( a )) − β s M X a , t ˜ u ′ kt ( a ) m t ( a ) + iv s + β β ˜ u kt ( a ) + ˜ u kt ( a ) X a ′ , t ′ ˜ u kt ′ ( a ′ ) m t ′ ( a ′ ) − i β sv M s + β β X a ′ , t ′ ˜ u kt ′ ( a ) m t ′ ( a ′ ) × exp − X a < b X t , t ′ ˜ u kt ( a ) ˜ u kt ′ ( b ) q tt ′ ( a , b ) − n X a = X t , t ′ ˜ u kt ( a ) ˜ u kt ′ ( a ) R tt ′ ( a ) , (A14)4where we introduce the Hubbard–Stratonovich transformation,exp x ! = Z Dv exp ( xv ) . (A15)We assume that the RS ansatz and static approximation: m t ( a ) = m , q tt ′ ( a , b ) = q , R tt ′ ( a ) = R ( t , t ′ ) , m xt ( a ) = m x , ˜ m t ( a ) = ˜ m , ˜ q tt ′ ( a , b ) = ˜ q , ˜ R tt ′ ( a ) = ˜ R ( t , t ′ ) , ˜ m xt ( a ) = m x . (A16)Under the RS ansatz and static approximation, Eq. (A14) can be expressed as e G = Y k Z Dv Z Dv Y a Z Dv Y t Z d ˜ u kt ( a )2 π s M πβ s exp β sv M + β β − M − β s ( R − β s ( ˜ u kt ( a )) + iv s + β β (1 − n β sm ) + v q m − q − n β sm + v p q − R ˜ u kt ( a ) , (A17)where we use the following relationships: X a X t , t ′ ˜ u kt ( a ) ˜ u kt ′ ( a ) R tt ′ ( a ) = R X a X t ˜ u kt ( a ) − ( R − X a , t (cid:16) ˜ u kt ( a ) (cid:17) , (A18) X a < b X t , t ′ ˜ u kt ( a ) ˜ u kt ′ ( a ) q tt ′ ( a , b ) = q X a , t ˜ u kt ( a ) − X a X t ˜ u kt ( a ) , (A19)as well as the Hubbard–Stratonovich transformation on (cid:16)P a , t (cid:17) and P a (cid:0)P t (cid:1) . We perform integration over ˜ u kt ( a ) as follows: e G = Y k Z Dv Z Dv Y a Z Dv Y t Z d ˜ u kt ( a )2 π s M πβ s exp β sv M + β β exp ( − M − β s ( R − β s (cid:16) ˜ u kt ( a ) − β sM − β s ( R − iv s + β β (1 − n β sm ) + v q m − q − n β sm + v p q − R × exp β s M − β s ( R − iv s + β β (1 − n β sm ) + v q m − q − n β sm + v p q − R = Y k Z Dv Z Dv Y a Z Dv Y t s MM − β s ( R −
1) exp β sv M + β β × exp β s M − β s ( R − iv s + β β (1 − n β sm ) + v q m − q − n β sm + v p q − R . (A20)In the limit of M → ∞ , we note that the coe ffi cient term in Eq. (A20) is reduced to Y t s MM − β s ( R − ≈ exp (cid:26) β s R − (cid:27) . (A21)After carrying out integration over the Gaussian variables v , v , and v , we can obtain e G as follows: e G = exp α nN −
12 ln(1 − β s ( q − R )) + β s R − + β s + β β + m − q − (1 + β − )1 − β s ( q − R ) . (A22)5We calculate e G under the RS ansatz and static approximation as follows: e G = N Y i =
12 Tr ξ Tr σ Z Dz exp ˜ mM X a , t ξ i σ zit ( a ) + R − ˜ q M X a M X t = σ zit ( a ) + √ ˜ qM X a , t σ zit ( a ) z + β (1 − s )(1 − λ ) M X a , t τ i σ zit ( a ) + ˜ m x M X a , t σ xit ( a ) , = N Y i =
12 Tr ξ Z Dz n Y a = Z Dy M Y t = Tr σ exp ˜ mM ξ i σ zit ( a ) + p R − ˜ qM σ zit ( a ) y + √ ˜ qM σ zit ( a ) z + β (1 − s )(1 − λ ) M τ i σ it ( a ) + ˜ m x M σ xit ( a ) ) = N Y i =
12 Tr ξ Z Dz Z Dy q g ( ξ i , τ i ) + ( ˜ m x ) ! n = N Y i =
12 Tr ξ exp ( n Z Dz " ln Z Dy q g ( ξ i , τ i ) + ( ˜ m x ) i ) = N Y i =
12 Tr ξ exp ( nc Z Dz ln Z Dy q g + ( ξ i ) + ( ˜ m x ) + n (1 − c ) Z Dz ln Z Dy q g − ( ξ i ) + ( ˜ m x ) ) ≈ exp ( nN c Z Dz ln Z Dy q g + + ( ˜ m x ) + (1 − c ) Z Dz ln Z Dy q g − + ( ˜ m x ) !) , (A23)where g ( ξ i , τ i ) = ˜ m ξ i + β (1 − λ )(1 − s ) τ i + p ˜ qz + q R − ˜ qy , (A24) g ± ( ξ i ) = ( ˜ m ± β (1 − λ )(1 − s )) ξ i + p ˜ qz + q R − ˜ qy . (A25)Here, the brackets [ · · · ] i represent the average over the sites:[ · · · ] i ≡ N N X i = ( · · · ) . (A26)Under the RS ansatz and static approximation, e G is expressed as e G = exp (cid:26) − NM nMm ˜ m − NM nMm x ˜ m x − NnM M R ˜ R − N M M ( n − n ) q ˜ q (cid:27) = exp ( Nn − m ˜ m − m x ˜ m x − R ˜ R + q ˜ q + β (1 − s ) λ m x + O ( n ) !) . (A27)In the thermodynamics limit N → ∞ , the saddle-point method can be used and the RS free energy is expressed as − β f RS = lim n → [ Z n ] − Nn = α − ln(1 + β s ( R − + β s ( R − + β s + β β − + β − + q − m + β s ( R − q ) + β (1 − s ) λ m x − m ˜ m − m x ˜ m x − R ˜ R + q ˜ q + c Z Dz ln Z Dy q g + + ( ˜ m x ) + (1 − c ) Z Dz ln Z Dy q g − + ( ˜ m x ) . (A28) Appendix B: DERIVATION OF AT CONDITION
In this Appendix, we derive the AT condition for the CDMA model with quantum fluctuations. The local stability of theRS solutions against the RSB perturbation is analyzed by constructing the 1-step (1RSB) solutions. The detailed derivation of6the 1RSB solutions is as follows. We assume the RS ansatz and static approximation for m t ( a ) , R tt ′ ( a ) , m xt ( a ) , ˜ m t ( a ) , ˜ R tt ′ ( a ), and˜ m xt ( a ). We divide the replicas into two blocks and introduce two order parameters, as follows: q tt ′ ( a l , b l ) = q ( l ∈ block ) q ( l < block ) , ˜ q tt ′ ( a l , b l ) = ˜ q ( l ∈ block )˜ q ( l < block ) , (B1)where l = , , . . . , n / m is the block number, m is Parisi’s breaking parameter, and a l , b l = , , . . . , m is the index inside ablock.By using the 1RSB scheme, we can divide these terms in Eq. (A14) as follows: X a < b X t , t ′ q tt ′ ( a , b ) ˜ u kt ( a ) ˜ u kt ′ ( b ) = q X l , a l , t ˜ u kt ( a l ) + q − q X l X a l , t ˜ u kt ( a l ) − q X l , a l X t ˜ u kt ( a l ) , (B2) X a , t , t ′ R tt ′ ( a ) ˜ u kt ( a l ) ˜ u kt ′ ( a l ) = R X l , a l X t ˜ u kt ( a l ) − ( R − X l , a l , t (cid:16) ˜ u kt ( a l ) (cid:17) . (B3)Using Eqs. (B2) and (B3), we can rewrite e G as e G = Y k Z Dv Y l , a l , t d ˜ u kt ( a l )2 π s M πβ s exp β sv M + β β − M β s (cid:16) ˜ u kt ( a l ) (cid:17) − β sm M X l , a l , t ˜ u kt ( a l ) + iv s + β β ˜ u kt ( a l ) + m ˜ u kt ( a l ) X l ′ , a ′ l , t ′ ˜ u kt ′ ( a ′ l ) − iv β smM s + β β X l ′ , a ′ l , t ′ ˜ u kt ′ ( a ′ l ) × exp − q X l , a l , t ˜ u kt ( a l ) − q − q X l X a l , t ˜ u kt ( a l ) + q X l , a l X t ˜ u kt ( a l ) − R X l , a l X t ˜ u kt ( a l ) + ( R − X l , a l , t (cid:16) ˜ u kt ( a l ) (cid:17) = Y k Z Dv Z Dv Y l Z Dv Y a l Z Dv Y t d ˜ u kt ( a l )2 π s M πβ s exp β sv M + β β exp ( − M − β s ( R − β s (cid:16) ˜ u kt ( a l ) (cid:17) + iv (1 − n β sm ) s + β β + v q m − q − n β sm + v √ q − q + v p q − R ˜ u kt ( a l ) = Y k Z Dv Z Dv Y l Z Dv Y a l Z Dv Y t s MM − β s ( R − × exp β s M − β s ( R − iv (1 − n β sm ) s + β β + v q m − q − n β sm + v √ q − q + v p q − R , (B4)where the Hubbard–Stratonovich transformation is used on (cid:16)P l , a l , t (cid:17) , P l (cid:16)P a l , t (cid:17) , and P l , a l (cid:0)P t (cid:1) . The coe ffi cient term in Eq. (B4)converges to the same as Eq. (A25) in the limit of M → ∞ . By performing the Gaussian integration over v , v , v , and v , e G can be expressed as e G = α nN − ln (1 + β s ( R − q )) − m ln + β sm ( q − q )1 + β s ( R − q ) ! + + β β β s − β s (1 + β − + q − m )1 + β s ( R − q + m ( q − q )) . (B5)7We calculate e G under the 1RSB scheme and static approximation as follows: e G = N Y i =
12 Tr ξ Tr σ Z Dz exp ˜ mM X l , a l , t ξ i σ zit ( a l ) + R − ˜ q M X l , a l M X t = σ zit ( a l ) + ˜ q − ˜ q M X l X a l , t σ zit ( a l ) + √ ˜ q M X l , a l , t σ zit ( a l ) z + β (1 − s )(1 − λ ) M X l , a l , t τ i σ zit ( a l ) + ˜ m x M X l , a l , t σ xit ( a l ) = N Y i =
12 Tr ξ Z Dz nm Y l = Z Dy m Y a l = Z Dx M Y t = Tr σ exp ˜ mM ξ i σ zit ( a l ) + p R − ˜ qM σ zit ( a l ) x + √ ˜ q − ˜ q M σ zit ( a l ) y + √ ˜ q M σ zit ( a l ) z + β (1 − s )(1 − λ ) M τ i σ it ( a l ) + ˜ m x M σ xit ( a l ) ) = N Y i =
12 Tr ξ Z Dz (Z Dy Z Dx q g ( ξ i , τ i ) + ( ˜ m x ) ! m ) nm = N Y i =
12 Tr ξ exp ( nm Z Dz " ln Z Dy Z Dx q g ( ξ i , τ i ) + ( ˜ m x ) ! m i ) = N Y i =
12 Tr ξ exp ( ncm Z Dz ln Z Dy Z Dx q g + ( ξ i ) + ( ˜ m x ) ! m + n (1 − c ) m Z Dz ln Z Dy Z Dx q g − ( ξ i ) + ( ˜ m x ) ! m ) ≈ exp ( nNm c Z Dz ln Z Dy Z Dx q g + ( ξ i = + ( ˜ m x ) ! m + (1 − c ) Z Dz ln Z Dy Z Dx q g − ( ξ i = + ( ˜ m x ) ! m !) , (B6)where g ( ξ i , τ i ) = ˜ m ξ i + β (1 − λ )(1 − s ) τ i + p ˜ q z + p ˜ q − ˜ q y + q R − ˜ q x , (B7) g ± ( ξ i ) = ( ˜ m ± β (1 − λ )(1 − s )) ξ i + p ˜ q z + p ˜ q − ˜ q y + q R − ˜ q x . (B8)Under the 1RSB scheme and static approximation, e G is expressed as e G = exp ( − NM nm m Mm ˜ m − NM nm m Mm x ˜ m x − NM nm m ( M − M ) R ˜ R − NM nm m M − Nn M M (( m − q ˜ q + ( n − m ) q ˜ q ) (cid:27) = exp ( Nn − m ˜ m − m x ˜ m x − R ˜ R − m − q ˜ q + m q ˜ q − M + O ( n ) !) . (B9)Finally, we can obtain the 1RSB free energy of the CDMA model with quantum fluctuations, as follows: − β f = α − ln (1 + β s ( R − q )) − m ln + β sm ( q − q )1 + β s ( R − q ) ! + + β β β s − β s (1 + β − + q − m )1 + β s ( R − q + m ( q − q )) − m ˜ m − m x ˜ m x − R ˜ R − m − q ˜ q + m q ˜ q + cm Z Dz ln Z Dy Z Dx q g + ( ξ i = + ( ˜ m x ) ! m + (1 − c ) m Z Dz ln Z Dy Z Dx q g − ( ξ i = + ( ˜ m x ) ! m . (B10)8The extrimization condition of the 1RSB free energy provides us with the saddle-point equations, as follows: m = c Z DzX − + Z DyY m − + Z Dx g + u + ! sin u + + (1 − c ) Z DzX − − Z DyY m − − Z Dx g − u − ! sin u − , (B11) q = c Z DzX − + Z DyY m − + (Z Dx g + u + ! sin u + ) + (1 − c ) Z DzX − − Z DyY m − − (Z Dx g − u − ! sin u − ) , (B12) q = c Z Dz ( X − + Z DyY m − + Z Dx g + u + ! sin u + ) + (1 − c ) Z Dz ( X − − Z DyY m − − Z Dx g − u − ! sin u − ) , (B13) R = c Z DzX − + Z DyY m − + Z Dx ( ˜ m x ) u + ! sinh u + + g + u + ! cosh u + + (1 − c ) Z DzX − − Z DyY m − − Z Dx ( ˜ m x ) u − ! sinh u − + g − u − ! cosh u − , (B14) m x = c Z DzX − + Z DyY m − + Z Dx ˜ m x u + ! sin u + + (1 − c ) Z DzX − − Z DyY m − − Z Dx ˜ m x u − ! sin u − , (B15)˜ m = αβ s + β s ( R − q ) , (B16)˜ q = ˜ q + αβ s ( q − q )(1 + β s ( R − q + m ( q − q )))(1 + β s ( R − q )) , (B17)˜ q = αβ s (1 + β − + q − m )(1 + β s ( R − q + m ( q − q ))) , (B18)2 ˜ R − ˜ q = αβ s ( R − q )1 + β s ( R − q ) , (B19)˜ m x = β (1 − s ) λ, (B20) X ± = Z Dy ( Y ± ) m , (B21) Y ± = Z Dx cosh u ± . (B22)When we set q = q = q and ˜ q = ˜ q = ˜ q , the 1RSB solutions are reduced to the RS solutions. The stability of the RS solutionsis evaluated by the stability analysis of the 1RSB solutions of q = q = q and ˜ q = ˜ q = ˜ q . We define ∆ = q − q and˜ ∆ = ˜ q − ˜ q . We apply Taylor expansion to these as follows: ∆ = c Z Dz Y − + Z Dy g + u + ! sinh u + ! − Y − + Z Dy ( ˜ m x ) u + ! sinh u + + g + u + ! cosh u + + (1 − c ) Z Dz Y − − Z Dy g − u − ! sinh u − ! − Y − − Z Dy ( ˜ m x ) u − ! sinh u − + g − u − ! cosh u − ˜ ∆ + O ( ˜ ∆ ) (B23)˜ ∆ ≈ ∂ ˜ ∆ ∂ ∆ (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) q = q = q ∆ + O ( ∆ ) ≈ αβ s (1 + β s ( R − q )) ∆ . (B24)We substitute Eq. (B23) for Eq. (B24). Finally, we can obtain the stability condition: α c β s (1 + β s ( R − q )) Z Dz Y − + Z Dy g + u + ! sinh u + ! − Y − + Z Dy ( ˜ m x ) u + ! sinh u + + Z Dy g + u + ! cosh u + + α (1 − c ) β s (1 + β s ( R − q )) Z Dz Y − − Z Dy g − u − ! sinh u − ! − Y − − Z Dy ( ˜ m x ) u − ! sinh u − + Z Dy g − u − ! cosh u − < . (B25)9 [1] S. Verdu, Multiuser Detection , 1st ed. (Cambridge UniversityPress, USA, 1998).[2] H. Nishimori,
Statistical Physics of Spin Glasses and Informa-tion Processing: an Introduction (Oxford University Press, Ox-ford; New York, 2001).[3] T. Tanaka, EPL (Europhysics Letters) , 540 (2001).[4] T. Tanaka, IEEE Transactions on Information Theory , 2888 (2002).[5] M. Yoshida, T. Uezu, T. Tanaka, and M. Okada,Journal of the Physical Society of Japan , 054003 (2007).[6] D. L. Donoho, IEEE Transactions on Information Theory , 1289 (2006).[7] D. Donoho and J. Tanner,Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences , 4273 (2009).[8] Y. Kabashima, T. Wadayama, and T. Tanaka,Journal of Statistical Mechanics: Theory and Experiment , L09003 (2009).[9] S. Ganguli and H. Sompolinsky,Phys. Rev. Lett. , 188701 (2010).[10] J.-i. Inoue, Phys. Rev. E , 046114 (2001).[11] Y. Otsubo, J.-i. Inoue, K. Nagata, and M. Okada,Phys. Rev. E , 051138 (2012).[12] Y. Otsubo, J.-i. Inoue, K. Nagata, and M. Okada,Phys. Rev. E , 012126 (2014).[13] T. Kadowaki and H. Nishimori, Phys. Rev. E , 5355 (1998).[14] G. E. Santoro, R. Martoˇn´ak, E. Tosatti, and R. Car,Science , 2427 (2002).[15] G. E. Santoro and E. Tosatti,Journal of Physics A: Mathematical and General , R393 (2006).[16] A. Das and B. K. Chakrabarti,Rev. Mod. Phys. , 1061 (2008).[17] S. Morita and H. Nishimori,Journal of Mathematical Physics , 125210 (2008).[18] R. D. Somma, D. Nagaj, and M. Kieferov´a,Phys. Rev. Lett. , 050501 (2012).[19] E. Farhi, J. Goldstone, S. Gutmann, J. Lapan, A. Lundgren, andD. Preda, Science , 472 (2001).[20] T. Albash and D. A. Lidar, Rev. Mod. Phys. , 015002 (2018).[21] V. S. Denchev, S. Boixo, S. V. Isakov, N. Ding, R. Bab-bush, V. Smelyanskiy, J. Martinis, and H. Neven,Phys. Rev. X , 031015 (2016).[22] T. Albash and D. A. Lidar, Phys. Rev. X , 031016 (2018).[23] S. Kirkpatrick, C. D. Gelatt, and M. P. Vecchi,Science , 671 (1983).[24] M. W. Johnson, P. Bunyk, F. Maibaum, E. Tolkacheva, A. J.Berkley, E. M. Chapple, R. Harris, J. Johansson, T. Lanting,I. Perminov, E. Ladizinsky, T. Oh, and G. Rose, SuperconductorScience and Technology , 065004 (2010).[25] A. J. Berkley, M. W. Johnson, P. Bunyk, R. Harris, J. Johansson,T. Lanting, E. Ladizinsky, E. Tolkacheva, M. H. S. Amin, andG. Rose, Superconductor Science and Technology , 105014(2010).[26] R. Harris, M. W. Johnson, T. Lanting, A. J. Berkley, J. Johans-son, P. Bunyk, E. Tolkacheva, E. Ladizinsky, N. Ladizinsky,T. Oh, F. Cioata, I. Perminov, P. Spear, C. Enderud, C. Rich,S. Uchaikin, M. C. Thom, E. M. Chapple, J. Wang, B. Wilson,M. H. S. Amin, N. Dickson, K. Karimi, B. Macready, C. J. S.Truncik, and G. Rose, Phys. Rev. B , 024511 (2010).[27] S. Boixo, T. F. Rønnow, S. V. Isakov, Z. Wang,D. Wecker, D. A. Lidar, J. M. Martinis, and M. Troyer,Nature Physics , 218 EP (2014).[28] H. G. Katzgraber, F. Hamze, and R. S. Andrist,Phys. Rev. X , 021008 (2014).[29] G. Rosenberg, P. Haghnegahdar, P. God- dard, P. Carr, K. Wu, and M. L. de Prado,IEEE Journal of Selected Topics in Signal Processing , 1053 (2016).[30] D. Venturelli and A. Kondratyev,Quantum Machine Intelligence , 17 (2019).[31] A. Perdomo-Ortiz, N. Dickson, M. Drew-Brook, G. Rose, andA. Aspuru-Guzik, Scientific Reports , 571 EP (2012).[32] M. Hernandez, A. Zaribafiyan, M. Aramon, and M. Naghibi,ArXiv e-prints (2016).[33] R. Y. Li, R. Di Felice, R. Rohs, and D. A. Lidar,npj Quantum Information , 14 (2018).[34] M. Henderson, J. Novak, and T. Cook,Journal of the Physical Society of Japan , 061009 (2019).[35] F. Neukart, G. Compostella, C. Seidel, D. von Dollen,S. Yarkoni, and B. Parney, Frontiers in ICT , 29 (2017).[36] N. Nishimura, K. Tanahashi, K. Suganuma, M. J. Miyama, andM. Ohzeki, Frontiers in Computer Science , 2 (2019).[37] M. Ohzeki, A. Miki, M. J. Miyama, and M. Terabe,Frontiers in Computer Science , 9 (2019).[38] D. Crawford, A. Levit, N. Ghadermarzy, J. S.Oberoi, and P. Ronagh, ArXiv e-prints (2016),arXiv:1612.05695 [quant-ph].[39] M. Benedetti, J. Realpe-G´omez, R. Biswas, and A. Perdomo-Ortiz, Phys. Rev. A , 022308 (2016).[40] F. Neukart, D. Von Dollen, C. Seidel, and G. Compostella,Frontiers in Physics , 71 (2018).[41] M. H. Amin, E. Andriyash, J. Rolfe, B. Kulchytskyy, andR. Melko, Phys. Rev. X , 021050 (2018).[42] A. Khoshaman, W. Vinci, B. Denis, E. An-driyash, H. Sadeghi, and M. H. Amin,Quantum Science and Technology , 014001 (2018).[43] A. D. King, J. Carrasquilla, J. Raymond, I. Ozfidan, E. An-driyash, A. Berkley, M. Reis, T. Lanting, R. Harris, F. Altomare,K. Boothby, P. I. Bunyk, C. Enderud, A. Fr´echette, E. Hoskin-son, N. Ladizinsky, T. Oh, G. Poulin-Lamarre, C. Rich, Y. Sato,A. Y. Smirnov, L. J. Swenson, M. H. Volkmann, J. Whittaker,J. Yao, E. Ladizinsky, M. W. Johnson, J. Hilton, and M. H.Amin, Nature , 456 (2018).[44] R. Harris, Y. Sato, A. J. Berkley, M. Reis, F. Altomare, M. H.Amin, K. Boothby, P. Bunyk, C. Deng, C. Enderud, S. Huang,E. Hoskinson, M. W. Johnson, E. Ladizinsky, N. Ladizinsky,T. Lanting, R. Li, T. Medina, R. Molavi, R. Neufeld, T. Oh,I. Pavlov, I. Perminov, G. Poulin-Lamarre, C. Rich, A. Smirnov,L. Swenson, N. Tsai, M. Volkmann, J. Whittaker, and J. Yao,Science , 162 (2018).[45] P. Weinberg, M. Tylutki, J. M. R¨onkk¨o, J. Westerholm,J. A. Åstr¨om, P. Manninen, P. T¨orm¨a, and A. W. Sandvik,Phys. Rev. Lett. , 090502 (2020).[46] K. Kitai, J. Guo, S. Ju, S. Tanaka, K. Tsuda, J. Shiomi, andR. Tamura, Phys. Rev. Research , 013319 (2020).[47] S. Suzuki and M. Okada, Journal of the Physical Society ofJapan , 1649 (2005).[48] T. J¨org, F. Krzakala, J. Kurchan, and A. C. Maggs,Phys. Rev. Lett. , 147204 (2008).[49] A. P. Young, S. Knysh, and V. N. Smelyanskiy,Phys. Rev. Lett. , 020502 (2010).[50] T. J¨org, F. Krzakala, J. Kurchan, A. C. Maggs, and J. Pujos,EPL (Europhysics Letters) , 40004 (2010).[51] T. J¨org, F. Krzakala, G. Semerjian, and F. Zamponi,Phys. Rev. Lett. , 207206 (2010).[52] M. Keck, S. Montangero, G. E. Santoro, R. Fazio, andD. Rossini, New Journal of Physics , 113029 (2017). [53] L. Arceci, S. Barbarino, D. Rossini, and G. E. Santoro,Phys. Rev. B , 064307 (2018).[54] A. Mishra, T. Albash, and D. A. Lidar,Nature Communications , 2917 (2018).[55] M. H. Amin, Phys. Rev. A , 052323 (2015).[56] N. G. Dickson, M. W. Johnson, M. H. Amin, R. Harris,F. Altomare, A. J. Berkley, P. Bunyk, J. Cai, E. M. Chapple,P. Chavez, F. Cioata, T. Cirip, P. deBuen, M. Drew-Brook,C. Enderud, S. Gildert, F. Hamze, J. P. Hilton, E. Hoskinson,K. Karimi, E. Ladizinsky, N. Ladizinsky, T. Lanting, T. Ma-hon, R. Neufeld, T. Oh, I. Perminov, C. Petro ff , A. Przybysz,C. Rich, P. Spear, A. Tcaciuc, M. C. Thom, E. Tolkacheva,S. Uchaikin, J. Wang, A. B. Wilson, Z. Merali, and G. Rose,Nature Communications , 1903 (2013).[57] G. Passarelli, G. De Filippis, V. Cataudella, and P. Lucignano,Phys. Rev. A , 022319 (2018).[58] G. Passarelli, V. Cataudella, and P. Lucignano,Phys. Rev. B , 024302 (2019).[59] J. Marshall, D. Venturelli, I. Hen, and E. G. Rie ff el,Phys. Rev. Applied , 044083 (2019).[60] A. Perdomo-Ortiz, S. E. Venegas-Andraca, and A. Aspuru-Guzik, Quantum Information Processing , 33 (2011).[61] N. Chancellor, New Journal of Physics , 023024 (2017).[62] G. Passarelli, K.-W. Yip, D. A. Lidar, H. Nishimori, and P. Lu-cignano, Phys. Rev. A , 022331 (2020).[63] M. Ohkuwa, H. Nishimori, and D. A. Lidar,Phys. Rev. A , 022314 (2018).[64] Y. Yamashiro, M. Ohkuwa, H. Nishimori, and D. A. Lidar,Phys. Rev. A , 052321 (2019). [65] J. R. L. de Almeida and D. J. Thouless,Journal of Physics A: Mathematical and General , 983 (1978).[66] In usual, the estimated signal is described as ˜ σ i = sgn h σ i i where sgn( · ) is the signum function. The average perfor-mance of the estimation can be evaluated by M ( β ) = P A Tr ξ R d y P ( A ) P ( ξ ) P ( y | ξ ) ξ i sgn h σ i i . In this paper, we adopt m to evaluate the performance of the MPM estimation for simplic-ity. If the estimated signals correspond to the original signals,we can attain m =
1. Therefore we may the overlap m as a mea-surement of the average performance of the MPM estimation.[67] M. Suzuki, Communications in Mathematical Physics , 183 (1976).[68] D. Sherrington and S. Kirkpatrick,Phys. Rev. Lett. , 1792 (1975).[69] Krauth, Werner and M´ezard, Marc,J. Phys. France , 3057 (1989).[70] We check the dependence of the order parameters on the Trotternumber in our simulations. As we increase the Trotter number,the deviation of the correlation between the Trotter slices getssmall. To check the validity of the static approximation is outof scope in our paper.[71] T. Obuchi, H. Nishimori, and D. Sherrington,Journal of the Physical Society of Japan , 054002 (2007).[72] M. Yoshida and T. Tanaka, in (2006) pp. 2378–2382.[73] J. Raymond and D. Saad,Journal of Physics A: Mathematical and Theoretical , 12315 (2007).[74] K. Takeda, S. Uda, and Y. Kabashima,Europhysics Letters (EPL)76