Finding a Hadamard Matrix by Simulated Quantum Annealing

13 Pages·2017·0.75 MB·English

by  Andriyan Bayu Suksmono

entropy
Article
Finding a Hadamard Matrix by Simulated
Quantum Annealing
AndriyanBayuSuksmono ID
TelecommunicationEngineeringScientificandResearchGroup(TESRG),SchoolofElectricalEngineeringand
InformaticsandTheResearchCenteronInformationandCommunicationTechnology(PPTIK-ITB),
InstitutTeknologiBandung,Jl.GaneshaNo.10,Bandung40132,Indonesia;[email protected]
Received:2January2018;Accepted:16February2018;Published:22February2018
Abstract: Hardproblemshaverecentlybecomeanimportantissueincomputing. Variousmethods,
includingaheuristicapproachthatisinspiredbyphysicalphenomena,arebeingexplored.Inthis
paper, we propose the use of simulated quantum annealing (SQA) to find a Hadamard matrix,
which is itself a hard problem. We reformulate the problem as an energy minimization of spin
vectors connected by a complete graph. The computation is conducted based on a path-integral
Monte-Carlo(PIMC)SQAofthespinvectorsystem,withanappliedtransversemagneticfieldwhose
strengthisdecreasedovertime. Inthenumericalexperiments,theproposedmethodisemployedto
findlow-orderHadamardmatrices,includingtheonesthatcannotbeconstructedtriviallybythe
Sylvestermethod. Thescalingpropertyofthemethodandthemeasurementofresidualenergyafter
asufficientlylargenumberofiterationsshowthatSQAoutperformssimulatedannealing(SA)in
solvingthishardproblem.
Keywords: quantumannealing;adiabaticquantumcomputing;hardproblems;Hadamardmatrix;
binaryoptimization
1. Introduction
1.1. Background
Finding a solution to a hard problem is a challenging task in computing. Such a problem is
characterized by its complexity, as it grows beyond the polynomial against the size of the input.
A class of particularly important ones are NP (non-deterministic polynomial) problems, in which
verifying a solution can be conducted in polynomial time, whereas finding the solution is of
exponentialorder.Examplesofsuchproblemsare,amongothers,theTSP(travelingsalesmanproblem),
SAT(Booleansatisfiability),graphcoloring,graphisomorphism,andsubsetsums.
Aninterestingapproachtothehardproblemsisamethodinspiredbyphysicalphenomena,such
asclassicalannealing(CA)orquantumannealing(QA).BothCAandQAarephysicalprocessesthat
obtainanordered(physical)systemfromanunorderedone,whichcanbedoneeitherthermally(asis
thecaseinCA)orquantum-mechanically(asisthecaseinQA).Tosimulatethephysicalprocesseson
a(classical/non-quantum)computer,numericalmethods,suchasMC(MonteCarlo)forCAandPIMC
(path-integralMonteCarlo)forQA,havebeendeveloped. Thealgorithmorcomputationalmethod
inspiredbyclassical/thermalannealingiscalledsimulatedannealing(SA),whereastheonebased
onquantumannealingiscalledsimulatedquantumannealing(SQA).Bothofthesemethodsmake
useofthemethodsinnumericalCAornumericalQA.TheyencodetheproblemintoaHamiltonian
ofaspinsystem[1]andthenevolvethesystemfromahighenergystatedowntothegroundstate.
Theannealingprocessenablesthesystemtoavoidlocalminimatrappingandthereforeiscapableof
achievingaglobaloptimum,whichrepresentsthebestsolutionoftheproblem. Themaindifference
Entropy2018,20,141;doi:10.3390/e20020141 www.mdpi.com/journal/entropy
Entropy2018,20,141 2of13
betweenSAandSQAisintheevolutionofthesystems;whereasSAusesclassical/thermalannealing,
SQAemploysquantummechanism.
InSA[2–4],onestartsthesystemintotalrandomnesswithregardtoahightemperaturestate.
Thetemperatureisthenloweredandthesystemisevolved,whichcausestheenergytodecreaseso
thatthesystembecomesincreasinglyordered.Toavoidlocal-optimatrapping,aparticularupdating
rule,suchastheMetropolis[2],isapplied. Theruleallowsthesystemto(sometimes)movetoahigher
energystate. Uponcompletionofthealgorithm,thesystemachievesthegroundstate,atwhichpoint
asolutionisfound.
In [5], Kadowaki and Nishimori introduced quantum fluctuations to replace the thermal
fluctuationsinSAtoacceleratetheconvergence. TheyappliedthemethodonanIsingmodel,wherea
transversefieldplaystheroleoftemperatureinclassicalSA,enablingthesystemtoachievetheground
statewithgreaterprobability.Santoroetal.[6]comparedclassicalandquantumMonteCarloannealing
protocolsonatwo-dimensionalIsingmodel. TheyfoundthatthequantumMonteCarloannealingis
superiortoclassicalannealing.In[7],Boixoetal.showexperimentalresultsona 108qubitD-WaveOne,
whichisakindofhardwareimplementationofQA.AstrongcorrelationbetweenD-WaveandSQA,
comparedtothedevicewithclassicalannealing,wasfound,whichindicatesthattheD-Waveperforms
quantumannealing. ThisresultraisedtheimportantissueofwhetherQAactuallyoutperformsSA[8].
Rønnowetal. [9]showedhowquantumspeedupshouldbedefinedandmeasured. Inanexperiment
with random spin glass instances on 503 qubits of D-Wave Two, they did not find any evidence of
suchspeedup.
Regardlessoftheseissues,differentresultshavebeenachievedviaSQA.Isakov[10]performed
quantumMonteCarlo(QMC)simulationsandfoundthattheQMCtunnelingratedisplayedscaled
accordingtosystemsize. HealsofoundquadraticspeedupinQMCsimulationswhen, insteadof
periodicconditions,openboundaryconditionswereemployed. In[11],Mazzolaetal. demonstrated
thatQMCsimulationscanrecoverthescalingofground-statetunnelingrates,whichvalidatesQAin
termsofsolvingcombinatorialproblems.
Someclassesofhardproblems,includingoneswithexponentialorcombinatorialcomplexity,
havebeenasubjectofinterestinSQAresearch. Martonaketal. [12]introducedanapplicationofSQA
tosolvetheTSPproblem.TheyfoundthataPIMCalgorithmwasmoreefficientthanSAintermsof
findinganapproximatelyminimaltourinagivengraph.SQAhasalsobeenusedtosuccessfullyaddress
otherhardproblemsrelatedtographs,suchasgraphcoloring[13]andgraphisomorphism[14].
Inthispaper,weproposeSQAasameantofindaHadamardmatrix(H-matrix).Previously,in[15],
wesuccessfullyemployedSAtoperformasimilartask,inwhichlow-orderH-matriceswerefound.
ComparedtoexistingH-matrixconstructionmethods,anSA-basedmethodismoregeneralinterms
ofitscapabilityoffinding(orconstructingprobabilistically)anm =4korderH-matrix,withoutany
restrictiononthepropertyoftheorderm,whereastheSylvestermethodrequiresm =2n,wherekand
narepositiveintegers. ThispaperextendsthisclassicalSAmethodtoitsquantumversion,where
PIMCbasedonSuzuki–Trotterformulation[16,17]isemployedtosimulatethequantumprocess.
1.2. FindingAHadamardMatrix
AHadamardmatrix,orH-matrix,isanorthogonalbinary{±1}matrixofsize4k×4k,wherekis
apositiveinteger. ThismatrixwasdiscoveredbyJ.J.Sylvester[18]in1867andthenstudiedmore
extensively by J. Hadamard [19] during his investigation of the maximal determinant problem.
TheorthogonalpropertymakestheH-matrixpopularinappliedareas,suchasinformationcoding
and signal transform. In the 1960s and 1970s, Hadamard code was used in space exploration for
informationtransmission[20,21].Inarecenttechnologicalcase,CDMA(code-divisionmultipleaccess),
whichiswidelyusedincellularmobilephonesystems,employsWalsh–Hadamardsignalstoreduce
interferencebetweenitsusers[22,23].
Entropy2018,20,141 3of13
OneofthemostimportantissuesinthetheoryofH-matrixisitsexistence.Any2l orderH-matrix
withlapositiveintegercanbeconstructedusingSylvester’smethod. Furthermore,ifthereisanm
orderH-matrix,m =4kcanbeshownforapositiveintegerk. Ontheotherhand,nooneyetknows
ifthereisalwaysa4korderH-matrix[20,21]. ThelattercaseisformulatedastheHadamardmatrix
conjecture. Uptothiswriting,thesmallestunknown4korderH-matrixis668.
Variousreconstructionmethodshavebeenproposed[24–29]. Nevertheless,thesemethodsforce
theordermtofollowaparticularrule. In[15],ageneralm = 4k orderalgorithmemployingSAis
proposed. ThemethodworksonaspecialH-matrixcalledaseminormalizedHadamard(SH)matrix,
in which the first column is a 4k order unity vector(cid:126)v = (1,··· ,1)T, and the rest are 4k order SH
0
vectors(cid:126)v ∈V.
i
Abrute-forcemethodneedstoverifyall N ofthe4korderbinarymatrixtofindanH-matrix,
B
whereN (4k) =216k2 [15]. Letallmatricesconstructedwhere(cid:126)v isthefirstcolumnandacombination
B 0
of(cid:126)v ∈Vconstitutestheremaining(4k−1)columnsbecalledquasi-SH(QSH)matrices.Sincethereare
i
N =C(4k,2k)SHvectors,thereareaboutN (4k) ≈ (cid:16) 24k (cid:17)4k uniqueQSH-matrices. Althoughthe
V QU 8k3/2
numberhasbeengreatlyreducedcomparedtoN ,exhaustivecheckingstillrequiresagreatamount
B
ofcomputationalresources. TheSAmethodproposedin[15]iscapableoffindingafewlow-orderSH
matricesinamorereasonabletime.
Followingtheconventioninourpreviouspaper[15],theroleofthespin,i.e.,its±1eigenvalues,
isreplacedbySHspinvectors(cid:126)v ∈V. Tofinda4korderSH-matrix,oneneeds(4k−1)fullyconnected
i
SHspinvectors,whichinitiallyaresetrandomly. WithadefinedenergyE(Q(cid:126)),theSHspinvectorsare
randomlychangedinaccordancewithconditionswherebyatransitionintoanotherSHspinvectoris
allowedbutatransitionintoanon-SH-spin-vectorisforbidden.
2. Methods
2.1. SimulatedQuantumAnnealing
TheHamiltonianofanIsingsystemwithspinconfiguration{σˆ },wherek∈K ={1,2,··· ,i,j,···}
k
isthesetofthelattice’sindices,canbeexpressedas
Hˆ =−∑J σˆzσˆz−∑hσˆz (1)
ij i j i i
i(cid:54)=j i
where J isacouplingconstant/strengthbetweenaspinatsiteiwithaspinatsitej,h isthemagnetic
ij j
strengthatsitej,and{σˆz,σˆx}arePauli’smatricesatsitei. InSQA,quantumfluctuationiselaborated
i i
by introducing a transverse magnetic field Γ. The Hamiltonian of the system takes the following
form[5]:
Hˆ =−∑J σˆzσˆz−∑hσˆz−Γ∑σˆx. (2)
QA ij i j i i i
i(cid:54)=j i i
InEquation(2),thetransversefieldischanged(reduced)overtime,i.e.,Γ ≡ Γ(t). Ontheright
handsideoftheequation,thefirsttwotermscorrespondstopotentialenergyHˆ ,whilethethirdone
pot
istheHamiltonianintroducedbythetransversefield,whichisrelatedtokineticenergyHˆ ;i.e,we
kin
candefine
Hˆ ≡−∑J σˆzσˆz−∑hσˆz (3)
pot ij i j i i
i(cid:54)=j i
Hˆ ≡−Γ∑σˆx. (4)
kin i
i
Entropy2018,20,141 4of13
Ingeneral,Hˆ andHˆ donotcommute,so[Hˆ ,Hˆ ](cid:54)=0. DenotingtheHamiltonianofthe
pot kin pot kin
potentialasafunctionofspinconfigurationsHˆ ≡ Hˆ (cid:0){σˆz}(cid:1),wecanalsoexpressEquation(2)ina
pot i
moregeneralformasfollows:
Hˆ =−Hˆ ({σˆz})−Γ∑σˆx. (5)
QA i i
i
TosimulateaquantumsystemdescribedbyEquation(5)usingtheclassicalmethod,wehaveto
formulatePIMCbyintroducingimaginarytime. ItcanbethenapproximatedbytheSuzuki–Trotter
transformbyaddingonedimensionintheimaginarytimedirection,which,for(P×N)degreesof
freedom,takesthefollowingform[13,30]:
(cid:32) (cid:33)
1 ∑P (cid:0) (cid:1) P∑−1∑N ∑N
HST = P Hpot {Si,p} −JΓ Si,pSi,p+1+ Sj,1Sj,p (6)
p=1 p=1 i j
where N isthenumberofspinsinthelattice, PisthenumberofTrotter’sreplicas,S = ±1arethe
i
eigenvaluesofthespinmatrices,and
PT (cid:18) Γ (cid:19)
JΓ =− lntanh >0 (7)
2 PT
isthenearest-neighborcouplingofthetransversemagneticfield[30].
2.2. SQAFormulationoftheSHSpinVector
Similar to the previous paper [15], we employ a seminormalized Hadamard spin vector,
abbreviated here as an SH spin vector, instead of an ordinary spin. In a 4k order SH spin vector,
for a given positive integer k, 2k spins are −1and another2k spins are +1. Therefore, an SH spin
vectortransitionisallowedonlyifthesebalancenumbersareconserved;otherwise,suchatransition
isforbidden. WealsotreattheSHspinvectorasasingleentity,eventhoughitconsistsof4kspins,
andisdenotedas(cid:126)v ∈V,whereVisthesetofall4k-orderSHvectors. Weformulatetheenergyofa
i
particularconfigurationofspinvectors{(cid:126)v}asfollows:
i
(cid:12) (cid:12)
(cid:12) (cid:12)
E({(cid:126)v})=(cid:12)∑(cid:126)v ·(cid:126)v +∑(cid:126)1·(cid:126)v(cid:12)−16k2 (8)
i (cid:12) i j i(cid:12)
(cid:12)i(cid:54)=j i (cid:12)
where(cid:126)v ·(cid:126)v denotestheinnerproductofthevector(cid:126)v with(cid:126)v.
i j i j
Figure1showsanIsingsystemwithfourSHspinvectorswithanadditionalTrotter’sdimension.
InthelowerpartofFigure1a,eachcirclerepresentsabinaryspin,whereasthesolidlinerepresentsthe
connectionamongthespins. InteractingspiniwithbinaryvariableS andspinjwithbinaryvariable
i
S contributestheterm J SS totheHamiltonian. Fora4kordercase,every4knon-connectedspins
j ij i j
aregroupedintooneSHvector(cid:126)v,whichisillustratedasadashedline. Tosimplifythediagram,each
i
SHvectorisrepresentedbyafilledcircle;thus,weobtaintheupperpartofFigure1a,whichiscalleda
sliceorareplica. InthePIMC,thesliceisreplicatedP-times,andtheseslicesarearrangedaslayers
inimaginarytime. EachneighboringSHvectorinareplica,i.e.,(cid:126)vi,p with(cid:126)vi,p−1 and(cid:126)vi,p with(cid:126)vi,p+1,
interacts. Theextension(inimaginarytime)isillustratedinFigure1b. TheHamiltonianinEquation(6)
becomesaHamiltonianofanSHvectorspinsystemH thatcanberewrittenasfollows:
QV
(cid:32) (cid:33)
1 ∑P (cid:0) (cid:1) P∑−1∑ ∑
HQV = P Hpot {(cid:126)vi,p} −JΓ (cid:126)vi,p·(cid:126)vi,p+1+ (cid:126)vi,1·(cid:126)vi,p (9)
p=1 p=1 i i
Entropy2018,20,141 5of13
(cid:0) (cid:1)
where JΓ ≡ JΓ(t)andHpot {(cid:126)vi,p} representcomplete-graphconnectionsamongtheSHspinvectors,
similartoEquation(8),whichisgivenby
(cid:12) (cid:12)
(cid:12) (cid:12)
H (cid:0){(cid:126)v }(cid:1)=(cid:12)∑(cid:126)v ·(cid:126)v +∑(cid:126)1·(cid:126)v (cid:12)−16k2. (10)
pot i,p (cid:12) i,p j,p i,p(cid:12)
(cid:12)i(cid:54)=j i (cid:12)
TheevolutionofH inEquation(9)leadstothesolutiontotheH-matrixsearchproblem.
QV
(a) (b)
Figure1. Connectiondiagramsofthespinsandspinvectors. Weconsiderafour-orderSHvector
inthisexample: (a)fourSHspinsareconnectedbyacompletegraphK ,andeachcolumnisthen
4
groupedintoasingleSHspinvector;(b)anextensionoffullyconnectedSHspinvectorsintoaTrotter
dimension(imaginarytime)τ.
We will now formulate the SQA method for finding the H-matrix into an algorithm, which
is displayed as pseudo-code in Algorithm 1. It takes the matrix order, the number of replicas,
theinitialtemperature,theinitialvalueofΓ,andtheamountofiterationsandsub-iterationsasinputs.
ThisalgorithmyieldseitheranSH-matrixoraQSH-matrixthathasmoreorthogonalcolumnvectors
thantheinitialone. ThealgorithmstartswitharandominitializationofreplicaswithQSH-matrices,
whichare(4k−1)setsofSHvectors,andthencalculatesitsinitialenergy. Followingthescheduleofa
lineartransversefield,atrialtransitionisperformedforeachreplica. Theacceptanceandrejection
ofthetransitionisbasedontheMetropoliscriterion. Theiterationwillbestoppedwheneitherthe
numberofmaximumiterationsisreachedoranSH-matrixisfound.
Entropy2018,20,141 6of13
Algorithm1FindinganH-MatrixviaSimulatedQuantumAnnealing
1: Input: OrderofSH-matrix4k,numberofreplicasP,T0,Γ0, MaxIter,SubIter.
(cid:126) (cid:126)
2: Output: A4k-orderSH-matrixHF orapartiallyorthogonalmatrixQ.
3: InitializeT = T0,Γ = Γ0
4: Initializeall-replicasRwithrandomlygeneratedQSH-matrix: R ← {Q(cid:126)1,...,Q(cid:126)P}
5: idx ←0
6: H(cid:126)F ←(cid:126)0
7: FLAG←0
8: while(idx < MaxIter)or(FLAG==0)do
9: Calculate JΓ(idx;Γ,P,T)
10: Calculatecurrentall-replicasenergy: Erep = HQV(R,JΓ)
11: r ←0
12: whiler < Pdo
(cid:126)
13: Selectareplicaatpositionr: Qr
14: Calculatepotentialenergyofthereplica: Epot = Hpot(Q(cid:126)r)
15: ifEpot >0then
16: m ←0
17: while (m < SubIter)and(FLAG==0) do
18: FlipSHspinvectorrandomly: Q(cid:126)r → Q(cid:126)r(cid:48)
19: Calculateenergyoftheupdatedreplica: Epot1 = Hpot(Q(cid:126)r(cid:48))
20: ifEpot1==0 then
21: Epot ← Epot1
22: H(cid:126)F ← Q(cid:126)r(cid:48)
23: FLAG←1
24: r ← P
25: else
26: Updateall-replicas: R → R(cid:48)
27: Calculateenergyofupdatedall-replicasErep1 ← HQV(R(cid:48),JΓ)
28: ∆Erep ← Erep1−Erep
29: ∆Epot ← Epot1−Epot
30: Performatransitionifallowed(Metropolisupdaterule):
31: if (∆Epot <0)or(∆Erep <0)or(e−∆ETrep >rand)then
32: Acceptthetransition: R ← R(cid:48),Erep ← Erep1
33: endif
34: endif
35: m ← m+1
36: endwhile
37: else
38: H(cid:126)F ← Q(cid:126)r
39: FLAG←1
40: r ← P
41: endif
42: r ←r+1
43: endwhile
44: endwhile
Entropy2018,20,141 7of13
3. NumericalExperimentsandAnalysis
3.1. Findinga12-OrderSH-MatrixUsingSQA
Wehaveperformednumericalexperimentstofindlow-orderH-matrices. Herewepresentresults
fortheH-matrixoforder12fordetailedanalysis,sinceitisthelowest-orderH-matrixthatcannot
beconstructedbytheSylvestermethod. Initially,alloftheslices(replica)werefilledwithrandomly
generated(cid:126)v ∈V. NotethattherearetwonestediterationsinAlgorithm1.Thefirstoneisaniteration
i
ofallreplicaswiththemaximumnumbersettok·M×M,whereM=12istheH-order.Thesecondone
isaniterationofflippingwithinasliceofareplica,whosenumberisc·M,ccanbeanysmallnumber.
The energy evolution during the iteration is shown in Figure 2. The figure shows curves of
replicaenergyE ,meanpotentialenergyEp ,andminimumpotentialenergyEp . Thereplica
rep mean min
energyisdefinedsimilarlytoEquation(9), i.e., E ≡ H , whereasthepotentialenergyisgiven
rep QV
byEquation(10) E ≡ H . Themeanandminimumvalueshavebeentakenacrossthereplicas.
pot pot
Basedonthefigure,bothEp andE fluctuateovertime,buttheytendtodecrease. Theminimum
mean rep
energyofalatticeinthereplicaEp alsotendstodecrease. WhenEp =0,theH-matrixisfound.
min min
Figure2.EnergyevolutionduringtheSQAalgorithmrunstofindanSH-matrixoforder12.Fourcurves
aredrawninthegraph,whicharethemeanpotentialenergyEpmean,theminimumpotentialenergy
Epmin,thereplicaenergyErep,andthedeviationstandardofthepotentialenergyEpstd.WhenEpmin
equalszero,theiterationisstoppedsinceanSH-matrixhasbeenfound.TheEp curveindicateshigh
std
variationintheconfigurationofreplicasattheinitialstage,whichisthenreducedinlaterstages.
The degree of orthogonality of the matrix Q(cid:126) is displayed by the indicator matrix D(cid:126) ≡ Q(cid:126)TQ(cid:126).
Figure3showstheinitialQSH-matrixanditsrelatedindicatormatrix. Wealsoshowtheinitialand
final indicators for the first and last slices of the replica in Figure 4. It is expected that all of the
QSH-matricesbecomemoreorthogonal,indicatedbyalowernumberofzerosinoff-diagonalentries.
Thelastfigureshowingthelastsliceofthereplicaconditionaftertheiterationsarecompletedclearly
show this case. The found H-matrix is shown in the left part of Figure 5, with its corresponding
indicatorshownontheright,whichisadiagonalmatrix.
Entropy2018,20,141 8of13
Figure3.TheinitialstateofthefoundH-matrix:(a)TheQSH-matrix,whitesquaresindicate+1,black
squaresindicate−1.(b)Orthogonalityindicator,graysquaresshowthenon-orthogonalitycondition
ofrelatedpairofvectors.
Figure4. Indicatormatricesofthereplicacontent:(a)thefirstreplicaattheinitialstage;(b)thelast
replicaattheinitialstage;(c)thefirstreplicaatthefinalstage,and(d)thelastreplicaatthefinalstage.
Thematricesattheinitialstagesshowmostofthevectorsasnon-orthogonal,whereasthoseatthefinal
stagesshowmostofthevectorsasorthogonal.
Entropy2018,20,141 9of13
Figure5.Finalresults:(a)thefoundH-matrixand(b)itsorthogonalityindicator.Thediagonalformof
theindicatormatrixindicatesthatallofthecolumnvectorsarenoworthogonal.
3.2. TheNumberofTheReplicasandConvergenceIssue
Intheory,thenumberofTrotter’sreplicasPshouldbeaslargeaspossible. However,inpractice,
weshouldalsoconsidertheconvergenceissuewhenarunningtimerestriction(iterationnumber)
isgiven.Asexplainedin[13,30],replicasprovidediversityofsolutions;agreater Pselectsthebest
solutionwithminimumenergy. Ontheotherhand,thereplicasarenotmerelyrunningMonteCarlo
onseveralreplicas;theinteractionsbetweenreplicas JΓ(t)alsodefinetheirbehavior,i.e.,alargevalue
ofΓattheinitialstageimpliesalowvalueof JΓ,whichloosenstheconnections,andtheinteractions
thenbecomeindependent. AlowΓvalueattheendofaniterationimpliesahigh JΓ value,which
tightensthereplicaconnectionssuchthattheybecomesimilar. Tomeasurethesevariations,weuseda
simpledeviationstandardofenergyacrossthereplicas. Figure6showsthecurvesofvariationofthe
energyevolutionforP =5,10,15,and20infindinga12-orderH-matrix.
Figure6. TheeffectofthereplicanumberPinthealgorithm: althoughideallyalargePisdesired,
italsoneedstobeadjustedtotheproblem.Variationinreplicaenergy(intermsofdeviationstandards
oftheenergyacrossthereplicas)whensearchingforanH-matrixisshown.Thenumbersofreplicas
P=10,15,20yieldlargevariationsuptotheendoftheiteration,whereasP=5yieldsabetterresult
withsteadyvaluesattheend.Inallofthesecases,fortheconstructionofa12-orderH-matrix,thetotal
maximumiterationissetto20,000,consistingofaglobaliterationcountof20,000foreachP.
Sinceinitiallythereplicasweresetrandomly,theywillhavealmostidenticalenergy,sovariation
intheenergywillbeverylow. Inlateriterations,thevaluewillincreaseasanewconfigurationis
explored, and this will be followed by a decrease, which indicates that the replicas have become
homogeneous. Thiscycleofincreasing–decreasingenergyshouldbeobservedifPischosenproperly
Entropy2018,20,141 10of13
withrespecttothedimensionoftheproblem(H-order)andasufficientnumberofiterations. WhenP
istoosmall,thesystemwillperformakintoclassicalSA,whereasaPthatistoolargewillcausethe
systemtofail. Thefigureshowsthat,foragivennumberofmaximumiterations20,000,thenumberof
replicasP =5isthemostsuitable;anythinghigheristoohigh. Thisalsoshowsthatfrequentupdates
onalimitednumberofreplicas,comparedtolessfrequentupdatesonalargernumberofreplicas,
betterachieveconvergence.
3.3. PerformanceComparison: SQAvs. SA
To compare performances, in the first experiment, we measured the residual error of both
algorithms. Sincethegroundstateisachievedwhenthematrixbecomesorthogonal,inwhichcase
Equation(10)willequalzero,theresidualerror(cid:101)willbedefinedastheminimum H overallofthe
pot
replicas,i.e.,(cid:101) =min(H ). WehavechosentheorderoftheH-matrixtobesufficientlylargesothat
pot
wewillstillhavearesidualerrorattheendoftheexecutionofthealgorithm,i.e.,sothattheH-matrix
isnotfound. Weconsideredorder M =28tobesufficientforthispurpose,whereweactuallyhave
283 =21,952spins. WealsochoseaTrottersliceofP = 5andplottedthecurveforiterations50up
to5,000,000.
Following[30],theannealingschedulewaslinear;i.e.,thetemperatureTwasreducedlinearly
inSA,andwasthetransversemagneticstrengthΓ. EventhoughTisreducedlinearly,thethreshold
probabilityP willchangeexponentially. Byusingthefunction
thresh
1 −1
P (t) =1− eT(t) (11)
thresh
2
thethresholdwillstartabithigherthan0.5,whichasymptoticallyapproaches1.0attheendofiteration
timet. Figure7showsthecurveofT(t),Pthresh(t),Γ(t),and JΓ(t).
(a) (b)
Figure7.TheannealingschedulesinSAandSQA:(a)Lineartemperaturescheduleandcorresponding
thresholdscheduleinSA.(b)Lineartransverse-fieldΓ(t)andcorrespondingJΓ(t)inSQA.
Theexperimentswererepeated10timesforeachcase. Theaveragesofresidualerrorsforeach
iterationnumbersareplottedinFigure8forbothSAandSQA.
Thefigureshowsthat,althoughinitiallytheresidualerrorofSQAislargerthanSA,theslope
issteeper. Withahighernumberofiterations,whichinthiscaseisaround100,000,SQAissuperior.
ConsideringthatSQAshowstheleastamountoferroramongthereplicaslices,itseemsthatvariation
inthereplicaisanidealsolution. InSA,onceasolutionisselected,thechangeinspinconfiguration
will be less significant by the time the system reaches a lower energy state. Therefore, in terms of
findinganH-matrix,SQAissuperiortoSA.