Path optimization for U(1) gauge theory with complexified parameters
aa r X i v : . [ h e p - l a t ] S e p KUNS-2822
Path optimization for U (1) gauge theory with complexified parameters Kouji Kashiwa ∗ and Yuto Mori † Fukuoka Institute of Technology, Wajiro, Fukuoka 811-0295, Japan Department of Physics, Faculty of Science, Kyoto University, Kyoto 606-8502, Japan
In this article, we apply the path optimization method to handle the complexified parametersin the 1+1 dimensional pure U (1) gauge theory on the lattice. Complexified parameters make itpossible to explore the Lee-Yang zeros which helps us to understand the phase structure and thuswe consider the complex coupling constant with the path optimization method in the theory. Weclarify the gauge fixing issue in the path optimization method; the gauge fixing helps to optimizethe integration path effectively. With the gauge fixing, the path optimization method can treat thecomplex parameter and control the sign problem. It is the first step to directly tackle the Lee-Yangzero analysis of the gauge theory by using the path optimization method. I. INTRODUCTION
Exploring the phase structure of theories and modelswith finite external parameters such as the temperature( T ), the chemical potential ( µ ) and the external mag-netic field ( B ) is an important subject to understand ouruniverse. For example, the phase structure of quantumchromodynamics (QCD) at finite T , µ and B is directlyrelated to the early universe, current heavy ion collisionexperiments, neutron star physics and so on; see Ref. [1]as an example.One of the interesting approaches to investigate thephase structure is the Lee-Yang zero analysis [2, 3]. In theanalysis, we complexify external parameters and searchzeros of the partition function in the complex plane of theexternal parameters. Then, an approaching tendency ofzeros to the real axis indicates the existence of the phasetransition because singularities of the partition functionare the origin of the ordinary phase transition. Particu-larly, the experimental observation [4] and the quantumcomputation by using a quantum computer [5] for thezeros are currently possible and thus the analysis has at-tracted much more attention recently.There have been some attempts to perform the Lee-Yang zero analysis in the gauge theory; an interestingexample is QCD with the complexified µ [6–8]. In the cal-culation, one first gathers numerical data at finite imagi-nary chemical potential ( µ I ) and after they construct thegrand canonical partition function with the complex µ byusing the Fourier transformation and the fugacity expan-sion; see Ref. [9]. However, the Fourier transformationbecomes much more difficult with decreasing T becausethe oscillation becomes severer and thus this approachcannot tell us the phase structure at low T ; for exam-ple, see Ref. [10]. The reason why we use the imaginarychemical potential to perform the Lee-Yang zero analysisin QCD is that there is the sign problem at complexifiedexternal parameters and then the Monte-Carlo calcula-tion sometimes fails; see Ref. [11] for details of the sign ∗ kashiwa@fit.ac.jp † [email protected] problem and Ref. [12] for details of the imaginary chemi-cal potential. If we can directly perform the Monte-Carlocalculation with finite complexified parameters, there isthe possibility that we can better understand propertiesof the phase structure via the Lee-Yang zero analysis.In addition, such complexification of chemical potentialmay be related to the investigation of the confinement-deconfinement transition at finite density [13, 14].In the Lefschetz thimble and path optimization meth-ods, dynamical variables are complexified and then theintegral path and/or the configurations are generatedsuch that those obeying the sign-problem weaken mani-fold. Since these approaches can weaken the sign prob-lem and thus it is natural to expect that these approachescan be applied to explore the system with complexifiedexternal parameters. In this study, we concentrate onthe application of the path optimization method to thesystem with complexified parameters.The path optimization method is based on the stan-dard path integral formulation with the complexificationof dynamical variables [15, 16]; the actual procedure isperformed as follows:1. The cost function, which reflects the seriousness ofthe sign problem, is prepared.2. Dynamical variables are Complexified.3. The integral path in the complex domain is modi-fied to minimize the cost function.After taking the prescription, we can have a better inte-gral path (manifold) which has larger | e iθ | ; 0 ≤ | e iθ | ≤ λφ theory [16], the Polyakov-loop extended Nambu–Jona-Lasinio model [17, 18] and the 0 + 1 dimensionalQCD [19].In this article, we apply the path optimization methodto deal with the complexified parameters. We here em-ploy one plaquette system in the 1+1 dimensional U (1)gauge theory with complexified coupling constant on thelattice; some results for this theory are obtained by us-ing a modification of the integral path, see Refs. [20, 21].In Sec. II, we show the formulation of the theory onthe lattice and the explanation of the path optimizationmethod. Numerical results are shown in Sec. III. Sec-tion IV is devoted to summary. II. FORMULATION
In this section, we summarize detailed formulation ofthe 1+1 dimensional U (1) lattice gauge theory and ex-plain how we apply the path optimization method to thetheory. We here consider the one plaquette system andthus the following formulation is corresponding to thesystem. A. Action and partition function
The Wilson’s plaquette action [22] for only one pla-quette in the case of the U (1) gauge theory is writtenas S G = β n P + P − o , (1)where β = 1 /g is the lattice gauge coupling constantand P ( P − ) denotes the plaquette (its inverse). Thedefinition of P is given by P := U U U − U − , (2)where U n are the U (1) link variables defined as U n := e igA n ∈ U(1) , (3)here A n denotes the U (1) gauge field, g = 1 /β and n = 1 , , , Z = Z Y n dU n e − S G = I ( β ) , (4)where I denotes the modified Bessel function of the firstkind. It should be noted that we can obtain analyticresult of the partition function for any system volumewith periodic or open boundary conditions on the lat-tice [20, 21, 23]. For real β , I ( β ) is always positive andthus there are no zeros of Z in Eq. (4), but the gauge cou-pling constant is now complex, β ∈ C , and thus the parti-tion function can be 0 which is nothing but the Lee-Yangzeros. These zeros play an important role to understandthe phase structure.For the gauge theory, the action is invariant under the gauge transformation. In this case, to use the gaugetransformation, one can reduce the degree of freedom to n deg = 1 , · · · , U n = ( U n n = 1 , · · · , n deg I otherwise . (5) B. Path optimization method
In the path optimization method, we extend dynamicalvariables from real ( t ∈ R ) to complex ( z ∈ C ). In thepresent case, we need to extend the plaquette and thelink variable as P = U U U − U − , (6) U n = e ig A n =: U n e z n , (7)where A n ∈ C and then z n ∈ R represents the modifi-cation of the integral path. To represent z n , we employthe artificial neural network as the model to generate theintegral path. We here use the simple neural networkwhich has the mono input, hidden and output layers. Inthis network, the variables in the hidden layer nodes ( y j )and the output variables ( z n ) are given as follows. z n = ω n F ( w (2) jn y j + b j ) ,y j = F ( w (1) ij t i + b j ) , (8)where t i denotes the parametric variable which is set toRe U n ′ , Im U n ′ , i = 1 , · · · , × n deg , w , b and ω areparameters of the neural network (weight and bias) and F is so called the activation function. In this study, wechose the tangent hyperbolic function as the activationfunction.To perform the path optimization, we need the costfunction ( F ); we here use F [ z ( t )] = Z d n t | e iθ ( t ) − e iθ | × | J ( t ) e − S G ( t ) | , (9)where J ( t ) is the Jacobian, θ denotes the phase of the av-erage phase factor and θ ( t ) is the phase of J ( t ) e − S G ( t ) = e iθ ( t ) | J ( t ) e − S G ( t ) | . Of course, we do not know the ac-tual value of θ except with β ∈ R , θ = 0, and thuswe replace θ with h θ pre i EMA where h θ pre i EMA is the ex-ponential moving average (EMA) of the phase obtainedin the previous optimization steps. Minimization of theabove cost function makes phases of e − S G as a function of z take similar values on the modified integral path whenthe regions are relevant to the final result. In other words,there is no need to care for the phase of the Boltzmannweight in irrelevant regions which should be automati-cally suppressed.Since there is the sign problem in the case of the com-plex coupling constant by definition, we use the phasereweighting to perform the Monte-Carlo calculation as hOi = hO e iθ i pq h e iθ i pq , (10)where O represents any operator and h· · · i pq means thephase quenched expectation values where | Je − S G | is usedas the Boltzmann weight. Since the Boltzmann weight, | Je − S G | , is now real, we can perform the Monte-Carlocalculation exactly. It should be noted that we are notrestricted to the original integral path in the estimationof Eq. (10) unlike the ordinary reweighting calculation.The machine learning technique was first introducedto the path optimization method in Ref. [16] to representthe modified integral path with a weaker sign problem.The machine learning technique was also introduced tothe generalized Lefschetz thimble method [24] to learnthe integral manifold where the sign problem is mild inRef. [25] a few days before Ref. [16]. C. Setting of numerical calculation
Numbers of the unit in the input, hidden and out-put layers are N input = 2 × n deg , N hidden = 10 and N output = n deg , respectively. To determine the parame-ters in the neural network, we optimize these by using theADADELTA [26], one of the stochastic gradient meth-ods, as an optimizer with the Xavier initialization [27],the batch normalization [28], the mini-batch training andthe exponential moving average; see Ref. [16] for detailsof the optimization.Actual configurations are generated by using thepath optimization method with the hybrid Monte-Carlo(HMC) method [29] in the systems which includes thesingle plaquette with the open boundary condition. Itshould be noted that we here use the HMC with thereplica exchange method (exchange HMC) [30, 31] be-cause there is the global sign problem even on the mod-ified integral path [32] which means that there are someseparated regions on the modified integral path relevantto the integration. Integration over these separated re-gions is quite difficult to pick up by using ordinary HMC:We prepare the N rep replicas characterized by the param-eters in neural network as C x = xN rep C, (11)where x means the replica number, x = 1 , · · · , N rep , and C represents the parameters of the optimized neural net-work ( C = w, b, ω in Eq. (8)). We set N rep = 5 in thenumerical calculation. We use the exchange probabilityof the replicas P ( U x ↔ U x ′ ) as P ( U x ↔ U x ′ ) = min (cid:18) , P ( U x ; C x ′ ) P ( U x ′ ; C x ) P ( U x ; C x ) P ( U x ′ ; C x ′ ) (cid:19) ,P ( U x ; C x ) = | J ( U x ; C x ) e − S G ( U ( U x ; C x )) | . (12) initial con- fi gurationsHMC on J o POM withmulti-batchtrainingexchangeHMCon J m returnis (cid:1) e i θ (cid:2) su ffi cient?practicalcon fi gu-rations noyes FIG. 1. The flowchart of the algorithm to generate practicalconfigurations in this work. Symbols, J o and J m , denotethe original and the modified integral path, respectively. Theclosed loop in the flowchart is the one cycle of the optimizationprocedure of the integral path. In Ref. [16], we used the sampling of configurations basedon the symmetry of the modified integral path, but themethod is only available if we know the good modifiedintegral path which can have the symmetry. In this pa-per, we employ the exchange HMC method to generallyperform the path optimization.The expectation values are calculated with 2500 config-urations and then the corresponding errors are estimatedby using the Jackknife method with the bin-size 50. Forthe parameter of the theory, we set β = β R + iβ I with β R , β I ∈ R .Figure 1 shows the flowchart of the algorithm to gen-erate practical configurations where J o and J m are theoriginal and modified integral path, respectively. III. NUMERICAL RESULTS
In this section, we show the numerical results of the1+1 dimensional U (1) gauge theory on the lattice withthe complex coupling by using the path optimizationmethod. Here, we show the results of the 1 + 1 dimen-sional U (1) gauge theory only with single plaquette; i.e.the simplest setting of the theory. Actually, it is nothing -1-0.5 0 0.5 1 -1 -0.5 0 0.5 1 I m P ( P × e i θ ) Re P ( P × e i θ ) -1 0 1 2 3 4 5 -3 -2 -1 0 1 2 3 I m P ( P × e i θ ) Re P ( P × e i θ ) FIG. 2. The scatter plot of P and P × e iθ on the X - Y planewith β = 2 i where X = Re P and Y = Im P . The top andbottom panels show the results without and with the pathoptimization. Plus signs and crosses indicate P and P × e iθ ,respectively. but the n deg -dimensional integral.Figure 2 shows the scatter plot on the Re P -Im P planewith β = 2 i . Here, we show the results with the gaugefixing: U n = ( U n n = 1 I n = 1 . (13)We can clearly see the modification of the integral pathfrom fig. 2. In addition, we can see the bias of the dis-tribution in P × e iθ which plays a crucial role in thecalculation of hPi ; this bias makes hPi finite because thephase quenched expectation value of P becomes 0. Thehistogram of the phase for the case of Fig. 2 is shownin Fig. 3. On the original integral path, θ distributionis widely spread, but the distribution on the modifiedintegral path is well localized; we can generate configu-rations strongly localized around two separated regions.The replica exchange method well works in both cases.Figure 4 shows the average phase factor, h e iθ i EMA , as afunction of the optimization step in one epoch with β = i and 2 i ; one epoch is defined so that one sequence of themini-batch training is finished. Here, we estimate the av-erage phase factor by using EMA. From the figure, we can θ β = i β = 2i θ β = i β = 2i FIG. 3. The normalized histogram of the phase, θ , for β = i and 2 i . The top and bottom panels show the results with-out and with the path optimization. The line in the toppanel shows the probability distribution on the original path, P ( θ ) = [ πβ I sin { arccos( θ/β I ) } ] − . clearly see that the average phase factor cannot be en-hanced without the gauge fixing. With the gauge fixing,the average phase factor approaches 1 with β = i . In thecase with β = 2 i , there is the serious global sign problemas shown in Fig. 3 and thus we have the upper limit of theimprovement, but we can well enhance the average phasefactor via the path optimization. It should be noted thatthe modified integral path sometimes provides the expec-tation value with the larger error-bar compared with theoriginal one even if the average phase factor is enhanced.This indicates that there is the competition between theimprovement of the original sign problem and the mod-ification of the integral path which is responsible to thestatistical error via the path optimization method. Inaddition, we can see from Fig. 5 that the average phasefactor becomes larger with reduction of n deg by using thegauge fixing; it may be related to the fact that we havelarger degree of freedom without the gauge fixing andthen the neural network cannot show sufficient perfor-mance to optimize the suitable integral path.For the reader’s convenience, we finally show the ex-pectation value of the plaquette as a function of β I with β R = 0 where zeros exist in Fig. 6. It is clearly seen thatthe modified integral path reproduces the analytic result < e i θ > E M A Step numberwith gauge fixingwithout gauge fixingI ( β ) FIG. 4. The average phase factor, h e iθ i EMA , during the opti-mization with and without the gauge fixing in one epoch. Thedash-dot line shows the exact average phase factor on the orig-inal path. In the case with the gauge fixing, we here perfectlyfix the gauge degree of freedom; i.e. n deg = 1. Upper- andbottom-side lines are results with β = i and 2 i , respectively.In the case with β = 2 i , the serious global sign problem exists. < e i θ > E M A Step numbern deg =1n deg =2n deg =3n deg =4 FIG. 5. The average phase factor, h e iθ i EMA , during the op-timization with and without the gauge fixing in two epochswith β = i . except the region near the partition function zeros. At β I ∼ . , . , .
6, we have the zeros and then there shouldbe the strong modification of the integral path and/or theserious global sign problem.
IV. SUMMARY
In this study, we have considered the 1+1 dimensional U (1) lattice gauge theory which has the single plaquettewith the complex coupling constant as a first step to di-rectly investigate the Lee-Yang zeros in the gauge theory.We have estimated how the average phase factor is im-proved via the modification of integral variables. Since -4-2 0 2 4 0 2 4 6 8 10 I m < P > β I -1 0 1 0 2 4 6 8 10 R e < e i θ > β I FIG. 6. The top and bottom panels show the β I -dependence ofthe expectation value of the plaquette and the average phasefactor, respectively. Symbols show numerical results obtainedvia the path optimization with the gauge fixing and solid linesdenote the analytic results; the solid line in the bottom panelis corresponding to the results on the original integral path. there is the sign problem when the coupling constant iscomplex, we employ the path optimization method toperform the Monte-Carlo calculation. To represent themodified integral path in the complex domain of the in-tegral variables, we employ the artificial neural network.We have shown that the modification of the integralpath represented by using the neural network can wellenhance the average factor if we impose the gauge fixingto the theory, but we cannot without the gauge fixing;this suggests the importance of the gauge fixing to controlthe sign problem in the path optimization method on thelattice. This issue is demonstrated in the system of thesingle plaquette. Also, we have checked that the replicaexchange method well works to generate configurationlocalized in the well separated regions which are realizedvia the path optimization. From these results, we haveclarified how to use the path optimization method in thegauge theory. It should be important in the applicationof the path optimization method to the more complicatedgauge theory such as QCD.In the present study, we have restricted the system sizeto be small because we are interested in the possibilityof applying the path optimization method to the systemwith complexified external parameters. Usually, the signproblem becomes exponentially worse in terms of the sys-tem volume; there is the competition between the expo-nential suppression of it from the system volume and theenhancement of it from the optimization. In the largervolume case, we should introduce some methods to re-duce the numerical cost to calculate the Jacobian, whosenumerical cost is proportional to the square of the systemvolume. Examples of promising methods are the diago-nal ansatz of the Jacobian [33] and the nearest-neighborlattice-sites ansatz [34]. We will revisit this issue in ourfuture work.This study is a first step in the path optimizationmethod to explore the phase structure of the gauge the-ory in the complexified parameter space which is impor-tant to understand properties of the phase transition; e.g.for investigation of the distribution of the Lee-Yang zeros.In the present work, we employ the simple gauge theory, but we believe that it sheds light on the complexificationof the integral variables and also the parameters. In thefuture, we will apply the path optimization method to the SU (2) gauge theory with the complex coupling constant.It was reported that the complex Langevin method failsin some parameter regions; see Ref. [35]. ACKNOWLEDGMENTS
We are very grateful to A. Ohnishi and Y. Namekawafor their careful reading of the manuscript and usefulcomments. We appreciate S. Hawkins for reading themanuscript via the English writing correction service pro-vided from Fukuoka Institute of Technology. We ac-knowledge the Lattice Tool Kit (LTKf90) since our codeis in part based on LTKf90 [36]. This work is supportedby the Grants-in-Aid for Scientific Research from JSPS(No. 18K03618, 18J21251 and 19H01898). [1] K. Fukushima and T. Hatsuda, Rept. Prog. Phys. ,014001 (2011), arXiv:1005.4814 [hep-ph].[2] T. D. Lee and C.-N. Yang, Phys. Rev. , 410 (1952).[3] M. Biskup, C. Borgs, J. T. Chayes, L. J. Kleinwaks, andR. Koteck`y, Physical Review Letters , 4794 (2000).[4] X. Peng, H. Zhou, B.-B. Wei, J. Cui, J. Du, and R.-B.Liu, Physical review letters , 010601 (2015).[5] A. Krishnan, M. Schmitt, R. Moessner, and M. Heyl,Phys. Rev. A , 022125 (2019), arXiv:1902.07155[quant-ph].[6] A. Nakamura and K. Nagata, PTEP , 033D01(2016), arXiv:1305.0760 [hep-ph].[7] K. Nagata, K. Kashiwa, A. Nakamura, and S. M. Nishi-gaki, Phys. Rev. D91 , 094507 (2015), arXiv:1410.0783[hep-lat].[8] M. Wakayama, V. G. Bornyakov, D. L. Boyda, V. A. Goy,H. Iida, A. V. Molochkov, A. Nakamura, and V. I. Za-kharov, Phys. Lett.
B793 , 227 (2019), arXiv:1802.02014[hep-lat].[9] A. Roberge and N. Weiss, Nucl.Phys.
B275 , 734 (1986).[10] R. Fukuda, A. Nakamura, and S. Oka, Phys. Rev.
D93 ,094508 (2016), arXiv:1504.06351 [hep-lat].[11] P. de Forcrand, PoS
LAT2009 , 010 (2009),arXiv:1005.0539 [hep-lat].[12] K. Kashiwa, Symmetry , 562 (2019).[13] K. Kashiwa and A. Ohnishi, Phys. Lett. B750 , 282(2015), arXiv:1505.06799 [hep-ph].[14] K. Kashiwa and A. Ohnishi, Phys. Lett.
B772 , 669(2017), arXiv:1701.04953 [hep-ph].[15] Y. Mori, K. Kashiwa, and A. Ohnishi, Phys. Lett. B , 688 (2018), arXiv:1705.03646 [hep-lat].[16] Y. Mori, K. Kashiwa, and A. Ohnishi, Phys. Rev. D ,111501 (2017), arXiv:1705.05605 [hep-lat].[17] K. Kashiwa, Y. Mori, and A. Ohnishi, Phys. Rev. D99 ,014033 (2019), arXiv:1805.08940 [hep-ph].[18] K. Kashiwa, Y. Mori, and A. Ohnishi, Phys. Rev.
D99 ,114005 (2019), arXiv:1903.03679 [hep-lat].[19] Y. Mori, K. Kashiwa, and A. Ohnishi, PTEP ,113B01 (2019), arXiv:1904.11140 [hep-lat]. [20] J. M. Pawlowski, M. Scherzer, C. Schmidt, F. P. G.Ziegler, and F. Ziesch´e, in (2020) arXiv:2001.09767 [hep--lat].[21] W. Detmold, G. Kanwar, M. L. Wagman, and N. C.Warrington, Phys. Rev. D , 014514 (2020).[22] K. G. Wilson, Phys. Rev.
D10 , 2445 (1974), [,45(1974);,319(1974)].[23] R. Balian, J. M. Drouffe, and C. Itzykson, Phys. Rev.
D10 , 3376 (1974).[24] A. Alexandru, G. Basar, P. F. Bedaque, G. W. Ridg-way, and N. C. Warrington, JHEP , 053 (2016),arXiv:1512.08764 [hep-lat].[25] A. Alexandru, P. F. Bedaque, H. Lamm, andS. Lawrence, Phys. Rev. D96 , 094505 (2017),arXiv:1709.01971 [hep-lat].[26] M. D. Zeiler, arXiv preprint arXiv:1212.5701 (2012).[27] X. Glorot and Y. Bengio, in
Proceedings of the Thir-teenth International Conference on Artificial Intelligenceand Statistics (2010) pp. 249–256.[28] S. Ioffe and C. Szegedy, (2015), arXiv:1502.03167[cs.LG].[29] S. Duane, A. D. Kennedy, B. J. Pendleton, andD. Roweth, Phys. Lett.
B195 , 216 (1987).[30] R. H. Swendsen and J.-S. Wang, Physical review letters , 2607 (1986).[31] M. Fukuma and N. Umeda, PTEP , 073B01 (2017),arXiv:1703.00861 [hep-lat].[32] H. Fujii, D. Honda, M. Kato, Y. Kikukawa, S. Komatsu,and T. Sano, JHEP , 147 (2013), arXiv:1309.4371[hep-lat].[33] A. Alexandru, P. F. Bedaque, H. Lamm, andS. Lawrence, Phys. Rev. D97 , 094510 (2018),arXiv:1804.00697 [hep-lat].[34] F. Bursa and M. Kroyter, JHEP , 054 (2018),arXiv:1805.04941 [hep-lat].[35] H. Makino, H. Suzuki, and D. Takeda, Phys. Rev. D92 ,085020 (2015), arXiv:1503.00417 [hep-lat]. [36] S. Choe, S. Muroya, A. Nakamura, C. Nonaka, T. Saito, and F. Shoji, Nucl. Phys. B Proc. Suppl.106