Chiral properties of (2+1)-flavor QCD in strong magnetic fields at zero temperature
aa r X i v : . [ h e p - l a t ] A ug Chiral properties of (2+1)-flavor QCD in strong magnetic fields at zero temperature
H.-T. Ding a , S.-T. Li b , a , A. Tomiya c , X.-D. Wang a , Y. Zhang a a Key Laboratory of Quark and Lepton Physics (MOE) and Institute of Particle Physics,Central China Normal University, Wuhan 430079, China b Institute of Modern Physics, Chinese Academy of Sciences, Lanzhou 730000, China c RIKEN BNL Research center, Brookhaven National Laboratory, Upton, NY, 11973, USA
We present lattice QCD results for masses and magnetic polarizabilities of light and strangepseudo-scalar mesons, chiral condensates, decay constants of neutral pion and neutral kaon in thepresence of background magnetic fields with eB ranging up to around 3.35 GeV ( ∼ M π ) in thevacuum. The computations were carried out in (2+1)-flavor QCD on 32 ×
96 lattices using theHighly Improved Staggered Quarks (HISQ) action with M π ≈
220 MeV at zero temperature. Wefind that the masses of neutral pseudo-scalar mesons monotonously decrease as the magnetic fieldstrength grows and then saturate at a nonzero value, while there exists a non-monotonous behaviorof charged pion and kaon masses in the magnetic field. We observe a qB scaling of the up and downquark flavor components of neutral pion mass, neutral pion decay constant as well as the quark chiralcondensates at 0.05 . eB . . We show that the next-to-leading order chiral correctionto the Gell-Mann-Oakes-Renner relation involving neutral pion is less than 6% and the correctionfor the relation involving neutral kaon is in the range of (56-72)% at eB . . We alsoderive the Ward-Takahashi identities for QCD in the magnetic field in the continuum formulationincluding the relation between neutral pseudo-scalar meson correlators and chiral condensates. I. INTRODUCTION
The properties of strong interacting matter in the external magnetic field have attracted a lot of studies in the recentyears as strong magnetic fields appear in heavy-ion collisions [1–3], and early universe [4] and magnet stars [5]. TheQCD thermodynamics in the presence of a background magnetic field is of particular interest. At zero temperature itis found from lattice QCD studies using standard staggered fermions that the order parameter of the transition, chiralcondensate, increases with the magnetic field strength eB (so-called magnetic catalysis) [6, 7]. This suggests that themagnitude of the chiral symmetry breaking becomes larger in the vacuum which leads to an expectation that the chiralcrossover transition temperature T pc increases with the magnetic field strength eB , which is also seen in the latticeQCD studies in 2-flavor and 3-flavor QCD using standard staggered fermions at finite lattice cutoffs [6, 8]. However,surprise came later that T pc actually decreases with eB as found from continuum extrapolated results in N f = 2 + 1lattice QCD studies using improved staggered (stout) fermions [9]. The discrepancy of results in Ref. [6] from those inRef. [9] is most likely due to the large discretization errors in the standard staggered fermions [8]. Accompanying thereduction of T pc a decreasing behavior of chiral condensate in eB in the proximity of transition temperature, i.e. theso-called inverse magnetic catalysis (IMC) is also found in Ref. [9, 10] and further lattice QCD studies using improveddiscretization schemes [11–14]. Many model and theoretical studies have been performed to understand the (inverse)magnetic catalysis, reduction of T pc as well as their relations [7, 15–26]. Recently it has been suggested from latticeQCD studies that the IMC is not necessarily associated with the reduction of T pc as a function of eB by performingsimulations with heavier-than-physical pions [27, 28] and it is more like a deconfinement catalysis [29].Although the increase (reduction) of T pc is often connected with the increase (reduction) of chiral condensates, whichsuggests the increase (reduction) of the magnitude of the chiral symmetry breaking, the breaking of chiral symmetryin QCD is also related to the Goldstone pions. At vanishing magnetic field the square of the Goldstone pion mass isproportional to the product of a sum of quark masses and quark condensates, where the former explicitly breaks thechiral symmetry of the QCD Lagrangian, the latter measures the strength of spontaneous symmetry breaking. Thisis the well-known Gell-Mann-Oakes-Renner (GMOR) relation [30], and the validity of GMOR the relation has beenconfirmed on lattice QCD simulations at vanishing magnetic field and in the vacuum [31]. The GMOR relation hasalso been extended to the 3-flavor case including a strange quark [32], and the next-to-leading order chiral correctionsto the GMOR relation have also been studied at zero temperature and vanishing magnetic field [33–35]. Extendingto the case of low temperature and vanishing magnetic field [36], weak magnetic field at zero temperature [37] andin both low temperature and weak magnetic field [38] it is found that the GMOR relation for neutral pions holdstrue in the leading order chiral perturbation theory ( χ PT) in the chiral limit of quark masses as well. As it is knownthat at vanishing magnetic field the transition temperature decreases with lighter pions [39, 40] one may expect thatthe reduction of the transition temperature in the magnetic field could be associated with a lighter Goldstone bosonwhich can be a neutral pion.Meson spectrum of QCD in the external magnetic field is of important interest by itself [16, 26, 41–47], and inparticular the mass of a neutral pion in the external magnetic field may be helpful to understand the reduction of T pc given that neutral pion is a Goldstone boson in the nonzero magnetic field. Besides that as implied from the QCDinequality [48] the mass sum of the up and down quark flavor components of the neutral pion is a lower bound of themass of a charged ρ meson whose condensation could signal the transition of the QCD vacuum into a superconductingstate in a sufficiently strong magnetic field [49, 50]. Most of lattice studies on the meson spectrum in external magneticfields have been performed in the quenched approximation, e.g. in the quenched two color QCD with overlap fermionsas valence quarks [51], quenched QCD with Wilson fermions [48, 52] and overlap fermions [53, 54] as valence quarksin the computation of meson correlation functions. In the earlier studies there exists a discrepancy in the behavior eB dependence of the neural pion mass from lattice QCD studies in the quenched approximation, i.e. for quenchedunimproved Wilson fermions the neutral meson mass firstly decreases and then increases as eB grows [48] whilefor quenched overlap fermions it monotonously decreases [53]. Later the discrepancy is resolved as pointed out inRef. [52] that the eB -dependence of hopping parameter κ has to be taken into account in the discretization schemeusing Wilson fermions, and a monotonous reduction of neutral pion mass with increasing eB is finally established inthe quenched QCD [52]. While for the case of a charged pion a monotonous increasing of its mass M π − with eB isfound in all the quenched QCD studies [48, 52–54]. On the other hand, lattice studies in full QCD on eB dependenceof mass spectrum focus only on the light pseudo-scalar mesons ( π , ± ) [9, 52]. It is shown in Ref. [9, 52] that behaviorof masses of π and π ± in external magnetic fields obtained from (2+1)-flavor QCD using stout fermions have similartrends as those from quenched QCD [52–54]. And it is worth mentioning that there exists a possible decreasing trendof charged pion mass in the very strong magnetic field in both quenched [54] and full QCD [9].Similar decreasing trend of the neutral pion mass and T pc in the nonzero magnetic field shown in [52] mightindicate the connection between these two quantities given that the GMOR relation holds in the nonzero magneticfield. However, it is not known yet whether the GMOR relation holds in a strong magnetic field although χ PTsuggests its validity in the weak magnetic field. Model studies on the other hand suggest that the GMOR relationholds for neutral chiral pions at nonzero magnetic field while they are violated for the charged ones at eB & [55]. Studies on the GMOR relation in the nonzero magnetic field are intricate due to the explicit breakingof rational invariance caused by the magnetic field [56]. Investigations on the charge pion decay constant have beenperformed on the lattice [57] and it was found that a further pion decay constant exists due to the possibility of anonzero pion-to-vacuum transition via the vector piece of the electroweak current. Both charged and neutral piondecay constants have been parameterized for the one-pion-to-vacuum matrix elements of the vector and axial vectorhadronic currents in the background magnetic fields [58] and studied in the NJL model [44, 59]. It is worth notingthat early studies on pion decay constants in e.g. Refs. [37, 38, 55, 60] involves only the one for neutral pion fieldsrelated to the axial vector current that is parallel to the magnetic field.In this paper we focus on the chiral properties of QCD vacuum by studying the light and strange meson massesin the pseudo-scalar channel, chiral condensates as well as the neutral pion and kaon decay constants related tothe axial vector current in a wide range of magnetic field strength from 0 to ∼ .
35 GeV in N f = 2 + 1 QCD.We will present a first lattice QCD study on the GMOR relation in the external magnetic field, and show a firstobservation of qB scaling of up and down quark flavor components of neutral pion mass, neutral pion decay constantand chiral condensates, and discuss the magnetic polarizabilities of light and strange pseudo-scalar mesons. Ourlattice simulations are performed on 32 ×
96 lattices at a single lattice cutoff a = 0 .
117 fm using Highly ImprovedStaggered fermions with a heavier-than-physical pion mass of ∼
220 MeV at zero temperature.The paper is organized as follows. In Section II we introduce basic quantities we are going to study, in Section IIIwe will describe our simulation parameters in lattice QCD, our methodology to extract the meson masses and theamplitude to compute the decay constants as well as the discussion on the π − ρ mixing in the magnetic field via thegeneralized Ward-Takahashi identities, and in Section IV we will show our results on masses of pseudo-scalar mesons,magnetic polarizabilities, chiral condensates, neutral pion and kaon decay constants as well as corrections to theGMOR relation as a function of the magnetic field strength. Finally we will summarize in Section V. In Appendix Awe will show the extended Ward-Takahashi identities in the noznero magnetic field and in Appendix B we will showthe details of the implementation of the magnetic field in the lattice QCD simulations using the Highly ImprovedStaggered fermions. Some of previous results on the masses of pseudo-scalar mesons have been reported in conferenceproceedings [61]. II. TEMPORAL CORRELATORS, MASSES OF PSEUDO-SCALAR MESONS, CHIRALCONDENSATES AND NEUTRAL PION AND KAON DECAY CONSTANTS
Hadron spectrum in the vacuum can be extracted from two-point temporal correlation functions in the Euclideanspace G ( τ ) = Z d ~x (cid:28) M ( ~x, τ ) (cid:16) M ( ~ , (cid:17) † (cid:29) , (1)where M = ¯ ψ ( τ, ~x ) Γ ψ ( τ, ~x ) is a meson operator that projects to a certain quantum channel Γ = Γ D N t a with Diracmatrices Γ D and a flavor matrix t a , and for instance Γ = γ and γ µ correspond to the pseudo-scalar and vectorchannel, respectively. The angular brackets h· · · i stand for the expectation value over the gauge field ensembles. Thetemporal correlator decays exponentially at large distance τ lim τ →∞ G ( τ ) ∼ e − m Γ τ , (2)which defines the mass m Γ of the corresponding ground state. In the case of staggered fermions the correspondingmeson operators are of the form ψ ( ~x )(Γ D N Γ ∗ T ) ψ ( x ) with ψ ( x ) a 16-component hypercubic spinor and Γ D and Γ T Dirac matrices in spin and taste space, respectively. Here we consider local operator only and the meson operatorreduces to M = ζ ( ~x ) ¯ χ ( ~x ) χ ( ~x ), where ζ ( x ) is the phase factor depending on the choice of Γ = Γ D = Γ T and χ ( ~x ) isthe staggered fermion field.The connected part of correlation function of staggered bilinear can thus be written asG( τ ) = − X x,y,z ζ ( ~n )Tr (cid:20)(cid:16) M − ( ~x, τ ; ~ , (cid:17) † M − ( ~x, τ ; ~ , (cid:21) , (3)where M − ( ~x, τ ; ~ ,
0) is the staggered propagator from ( ~ ,
0) to ( ~x, τ ). In this work we focus on the mesons in thepseudo-scalar channel built from ¯ q i q j flavor combinations, here i, j = u, d, s , and the phase factor for the pseudo-scalarchannel is ζ ( ~n ) = 1. Due to the presence of the magnetic field the iso-spin symmetry of up and down quarks is brokenby their different electric charges. The neutral pion thus is not an iso-vector state any more and a full computationof neutral pion correlation function thus involves both connected and disconnected parts of the correlation functionas well as the mixing factor between u ¯ u ( π u ) and d ¯ d ( π d ) components . It has been shown in Ref. [54] that thequark-line disconnected part is negligible in the nonzero magnetic fields, and in our current study of neutral pions wethus neglect the disconnected contributions.A typical staggered meson correlator, for a fixed separation (in lattice unit) between the source and sink, is anoscillating correlator that simultaneously couples to two sets of mesons with the same spin but with opposite parities,and it thus can be parameterized as [62] G ( n τ ) = N nosc X i =1 A nosc,i exp ( − M nosc,i n τ ) − ( − n τ N osc X i =0 A osc,i exp ( − M osc,i n τ ) , (4)where N nosc ( N osc ) is the number of non-oscillating (oscillating) meson states whose masses are denoted by M nosc,i and M osc,i , and n τ = τ /a ∈ Z and a is the lattice spacing. Note that both amplitudes, A nosc,i and A osc,i are positive.In the current study the mass M nosc is the mass of pseudo-scalar meson we are interested in.It is well known that the energy of a point-like charged particle in the nonzero magnetic field and at zero temperatureare the Landau levels, E n = M + (2 n + 1) | eB | − g s z qB + p z , n ∈ Z +0 , (5)where M is the mass of a charged particle at zero magnetic field, q is the electrical charge of the particle, s z is the spinpolarization in the z direction and the magnetic field is assumed to go along with the z direction. For pseudo-scalarmesons the gyromagnetic ratio g = 0 while for vector mesons g = 2. In the case of the lowest Landau level and zeromomentum along the z direction, i.e. n =0 and p z =0, the mass of a charged point-like pseudo-scalar meson in theexternal magnetic fields can be expressed as M ± ps ( B ) = q(cid:0) M ± ps ( B = 0) (cid:1) + | eB | . (6)While the charged particles become heavier with increasing eB , the neutral particle is supposed to remain inde-pendent of magnetic field if it remains as a point-like particle. The lightest pseudo-scalar mesons, i.e. pions are ofparticular interests as which are Goldstone bosons at vanishing magnetic field, and their masses are connected to thequark chiral condensate in the vacuum known as the Gell-Mann-Oakes-Renner relation expressed as [30]( m u + m d ) (cid:0) h ¯ ψψ i u + h ¯ ψψ i d (cid:1) = 2 f π M π (1 − δ π ) , (7) In the presence of magnetic field the operator for π should be α ¯ uγ u − β ¯ dγ d with α + β = 1 [57]. Here in our computation for theneutral pion correlation function to extract the mass of π and neutral pion decay constant f π we have α = β = 1 / √ and the above GMOR relation in the two-flavor theory has been extended to the 3-flavor case including a strangequark as follows [32] ( m s + m d ) (cid:0) h ¯ ψψ i s + h ¯ ψψ i d (cid:1) = 2 f K M K (1 − δ K ) , (8)where m π and m K are the masses of pion and kaon respectively, f π ( f K ) is the pion (kaon) decay constant, m f = u,d,s isthe mass of a up, down and strange quark and the corresponding quark chiral condensate is denoted by h ¯ ψψ i f = u,d,s .Here in our study of N f = 2 + 1 QCD, m u = m d and m s /m l = 10, and the single flavor quark chiral condensate h ¯ ψψ i f is obtained as follows h ¯ ψψ i f = 14 1 V ∂ ln Z∂m f = 14 1 V Tr D − f , (9)where Z is the partition function of QCD, V is the full volume of space-time, and the factor 1 / D f in the staggered theory. δ π and δ K are the next-to-leading order chiral corrections, andboth of them are related with some low-energy constants and have a relation of δ K = M K /M π δ π [32]. The estimateson these two quantities using χ PT combining with QCD sum rules are δ π = (6 . ± . δ K = (55 ± √ f π M π = ( m u + m d ) (cid:28) (cid:12)(cid:12)(cid:12)(cid:12) √ (cid:0) ¯ uγ u − ¯ dγ d (cid:1)(cid:12)(cid:12)(cid:12)(cid:12) π ( p = 0) (cid:29) , (10) √ f K M K = ( m d + m s ) (cid:10) (cid:12)(cid:12) ¯ dγ s (cid:12)(cid:12) K ( p = 0) (cid:11) . (11)We will also look into the up and down quark flavor components of the neutral pion decay constant ( f π u and f π d )and neutral pion mass ( m π u and m π d ), whose relation is expressed as follows √ f π u M π u = 2 m u (cid:10) | ¯ uγ u | π u ( p = 0) (cid:11) , (12) √ f π d M π d = 2 m d (cid:10) (cid:12)(cid:12) ¯ dγ d (cid:12)(cid:12) π d ( p = 0) (cid:11) . (13)The corresponding GMOR relation is thus4 m u h ¯ ψψ i u = 2 f π u M π u (1 − δ π u ) , (14)4 m d h ¯ ψψ i d = 2 f π d M π d (1 − δ π d ) . (15)The matrix elements in Eqs. 10, 11, 12 and 13 can be extracted from the correlation function of pseudo-scalar mesonat zero spatial momentum h O S ( τ )P W (0) i [64]. The amplitude of the correlation function can be written as C O S P W = h | O S | P( ~p = 0) i h P( ~p = 0) | P W | i M P V s , (16)where V s is the spatial lattice volume, P denotes the operator for π , K and π u,d with M P the corresponding mesonmasses. By inserting the corresponding operator ( O S = P W ) into Eq. 16 and combining Eqs. 10, 11, 12 and 13, thepion and kaon decay constants can be obtained as follows [64] f P π = 2 m l r V s s C P πW P πW M π , (17) f K = ( m d + m s ) r V s s C K W K W M K , (18)where a factor of 1 / √ m l = m u = m d , and P π denotes the cases for π , π u and π d . The amplitudes C P πW P πW and C K W K W as well as themasses M P π and M K will be obtained from the mass fit (cf. Eq. 4) to the correlators in the pseudo-scalar channel. III. NUMERICAL SETUPSA. Lattice setup
Most of previous lattice QCD studies for (2 + 1)-flavor QCD in the external magnetic field were performed usingthe stout staggered fermions. In our current simulations we use N f = 2 + 1 Highly Improved Staggered Quarks(HISQ) [65]. At a given value of the lattice spacing the HISQ action achieves better taste symmetry than 2stoutactions [66]. The HISQ action is constructed by the Kogut-Susskind 1-link action and Naik improvement term withsmeared links. The smeared links are obtained in the following way. Firstly level one smeared links V µ are constructedby fat-7 from thin SU(3) links U µ . Next, re-unitarized links W µ are constructed by projecting V µ on U (3). Finally,level two smeared links X µ are constructed by fat-7 from thin SU(3) links W µ with the Lepage term. The HISQDirac operator are built by the Kogut-Suskind term with X µ and the Naik term with W µ . The magnetic field onlycouples directly to quarks thus the implementation of magnetic field is done just by replacing X µ → u µ X µ in theKogut-Susskind term and W µ → u µ W µ in the Naik term [14]. Details about the implementation of magnetic fields inthe lattice QCD simulations using the HISQ action are summarized in Appendix B.The external magnetic field pointing along the z -direction ~B = (0 , , B ) is described by a fixed link variable u µ ( n )of the U (1) field, and u µ ( n ) is expressed as follows in the Landau gauge [9, 67], u x ( n x , n y , n z , n τ ) = ( exp[ − iqa BN x n y ] ( n x = N x − u y ( n x , n y , n z , n τ ) = exp[ iqa Bn x ] ,u z ( n x , n y , n z , n τ ) = u t ( n x , n y , n z , n τ ) = 1 . (19)Here the lattice size is denoted as ( N x , N y , N z , N τ ) and coordinates as n µ = 0 , · · · , N µ − µ = x, y, z, τ ). Theexternal magnetic field applied along the z -direction B = (0 , , B ) is quantized in the following way qB = 2 πN b N x N y a − , (20)where q is the electric charge of a quark, and N x ( N y ) is the number of points in the x ( y ) direction on the lattice.Since the quantization has to be satisfied for all the quarks in the system, a greatest common divisor of the electriccharge for all the quarks, i.e. | q d | = | q s | = e/ e the elementary electric charge, is chosen in our simulation. Thusthe strength of the magnetic field eB practically is expressed as follows eB = 6 πN b N x N y a − , (21)where N b ∈ Z is the number of magnetic flux through unit area in the x-y plane. The periodic boundary conditionfor U (1) links is applied for all directions except for the x -direction. As limited by the boundary condition N b isconstrained in the range of 0 ≤ N b < N x N y .For the gauge part we use a tree-level improved Symanzik gauge action. The simulations have been performed on32 ×
96 lattices at zero temperature. The inverse lattice spacing is a − ≃ .
685 GeV, and strange quark mass m s istuned to its physical value by tuning the mass of the η s meson, M η s ≃
684 MeV which is based on the leading orderchiral perturbation theory relation M η s = p M K − M π between masses of η s , π and K . The light quark mass isthen set to be m l = m s /
10 corresponding to a pion mass M π ≃
220 MeV in the vacuum. The details on the scalesetting can be found in Ref. [62, 68].In our simulations 16 different values of N b have been chosen, i.e. 0, 1, 2, 3, 4, 6, 8, 10, 12, 16, 20, 24, 32, 40, 48and 64 whose corresponding magnetic field strength eB ranges from 0 to around 3.35 GeV . All configurations havebeen produced using Rational Hybrid Monte Carlo (RHMC) algorithm and are saved by every 5 time units, and thestatistics for each N b is about O ( ∼ ). Details of eB and statistics are listed in Table I. B. Meson correlation functions
We show our results of correlation functions for π u ( u ¯ u ) and π d ( d ¯ d ) as an example in the left and right plot ofFig. 1, respectively. It can be seen clearly that with increasing magnetic field strength eB both correlation functionsbecome larger, and the u ¯ u component of the two-point correlation function increases faster. This may be understood eB [GeV ] 0 0.052 0.104 0.157 0.209 0.314 0.418 0.523 N b eB [GeV ] 0.627 0.836 1.045 1.255 1.673 2.09 2.510 3.345 N b
12 16 20 24 32 40 48 64 ×
96 lattices with β = 6 .
68 and a ≃ .
117 fm ( a − ≃ .
685 GeV) produced usingthe HISQ fermion action and a tree-level improved Symanzik gauge action in N f = 2 + 1 QCD. The pion mass used in thesimulation is tuned to be 220 MeV, obtained mass of ρ is M ρ =779(35) MeV, M K =507.03(7) MeV and M η s =684.44(6) MeV at eB = 0. π (n τ ) n τ N b =0N b =4N b =8N b =16N b =32N b =48N b =64 π (n τ ) n τ N b =0N b =4N b =8N b =16N b =32N b =48N b =64 FIG. 1. Left: The u ¯ u component of the neutral pion correlation function G π u ( n τ ) at different values of N b calculated using 12corner wall sources. Right: Same as left one but for the d ¯ d component G π d ( n τ ). G shown here are rescaled with lattice spacing a such that they are dimensionless quantities. that the internal structure of pion is probed and the u ¯ u component is more affected due to the larger absolute valueof electrical charge of the u quark.We then extract the mass of neutral as well as the charged pseudo-scalar mesons by fitting to the correspondingtemporal correlation functions with the ansatz (cf. Eq. 4). The non-oscillating states are the physical states we need,and oscillating states are also necessary to be included in the fit in particular in the case to correlation functionsfor charged particles. In the fits we used various numbers of non-oscillating as well as oscillating states, i.e. ( N nosc , N osc ) has been set to (1,0), (1,1), (2,1). All the fits have been performed in a given interval [ τ min , N τ / τ min ranges from 0 to N τ / −
1. The final results are chosen as the best fit from all different fit modes through the AICc(corrected Akaike Information criterion) [69, 70]AICc = 2 k − ln( ˆ L ) + 2 k + 2 kn − k − , (22)where k is the number of parameters, ˆ L is the likelihood function, and the last term is needed to correct over-fittingif the number of data points n is not much larger than k . Then we choose a plateau in the AICc selected results andobtain the final mass and its uncertainty from the plateau by using a Gaussian bootstrapping method [62].Since correlation functions computed using a single point source do not suppress the excited states and make theextraction of ground states difficult we thus also compute correlation functions using a corner wall source which meansputting an unit source at the origin of each 2 on a chosen z -slice [62, 71–73]. We further improve the signal by putting12 corner wall sources at (0,0,0,0), (0,0,0,8),...,(0,0,0,88). The signal-to-noise ratio δG ( τ ) / h G ( τ ) i obtained using asingle corner wall source is reduced by a factor of 6 compared to the results obtained using a single point source.The signal-to-noise ratio obtained using multiple corner wall sources is found √ u ¯ u component of a neutral pion correlation function at eB =0.84 GeV measured using a single point source(top two plots) and 12 corner wall sources (bottom two plots) in Fig. 2. We find that the usage of corner wall sources M[GeV] τ min plateau(1,0)(1,1)(2,1)AICc selected C PP V s τ min plateau(1,0)(1,1)(2,1)AICc selected M[GeV] τ min plateau(1,0)(1,1)(2,1)AICc selected C PP V s τ min plateau(1,0)(1,1)(2,1)AICc selected FIG. 2. Top: Mass (left) and its associated matrix element C PP multiplied by the spatial volume V s (right) of π u at N b = 16( eB = 0.84 GeV ) with there different fit modes: ( yields a much better signal in the ground states as compared to that of the point source and has a longer plateau withmuch smaller uncertainty. The value of ground state mass extracted using the corner wall source is consistent withinthat extracted using the point source. As seen from Fig. 2 the corner wall source also works better for the extractionof the amplitude in the same way, i.e. the amplitude obtained with corner wall source is about 3% larger than theamplitude obtained with point source while the relative error is reduced by ∼
20 times. The amplitude obtained fromthe corner wall sources is then used in the calculation of decay constants in Section IV.
C. Mixing of pion states and Ward-Takahashi identities at eB > Since the states of vector meson ρ with spin polarization s z = 0, i.e. ρ s z =0 has the same quantum number ofthe states of π in the presence of magnetic field it is supposed that the mixing between the pion and ρ s z =0 states isenabled [52], i.e. neutral π mixes with neutral ρ s z =0 while charged π ± mixes with charged ρ ± s z =0 . This is to saythat once the mixing is enabled correlation functions in the vector channel and pseudo-scalar channels receive mutualcontributions and their corresponding ground state masses should be the same. The mixing has been investigated indetails in quenched QCD in Ref. [52], and it was found that the influence to the ground state of π is negligible. Herewe show another evidence that the influence of mixing to neutral pseudo-scalar meson states is mild by looking into Note that the mixing with s z = ± the following Ward-Takahashi identities (cf. A.31, A.32 and A.33 in Appendix A 2) h ¯ ψψ i f = m f χ ps , (23) h ¯ ψψ i d + h ¯ ψψ i s = ( m d + m s ) χ K , (24)where χ ps is the space-time sum of the correlation function of neutral mesons consisting only f ¯ f in the pseudo-scalarchannel, χ ps = P Nτ − τ =0 G ps ( τ ) with f = u, d and s quark flavors. This is to say that χ ps is evaluated by usingthe u ¯ u , d ¯ d and s ¯ s components of the correlation functions of neutral pseudo-scalar mesons, i.e. π u , π d and η s . Thechiral condensate h ¯ ψψ i f with the corresponding flavor is computed via Eq. 9. And χ K is the space-time sum of theneutral kaon correlators. As shown in Appendix A.32 and A.31 the above Eq. 23 for f = s should hold at nonzeromagnetic field, while for neutral pion we only look the connected part and break it into the up and down quark flavorcomponents, i.e. for f = u and d , respectively. The ratio of h ¯ ψψ i f / ( m f χ ps f ) ( f = u, d, s ) as a function of eB is shownin the left plot of Fig. 3. We found that all the ratios agree with unity with devations less than 1.2% at all the valuesof magnetic field strength we studied and this indicates that the influence of the mixing to neutral pseudo-scalar statesis negligible. The other indication is that for the neutral pion the mixing between ¯ uγ and ¯ dγ d stays almost the sameas they were at zero magnetic field. Results for the ratio of left hand side to right hand side of Eq. 24 as derived inAppendix A 2 are shown in the right plot of Fig. 3, and the deviation from unity is less than 0.6% in eB ∈ [0 , . . ψ - ψ > f / (m f χ ps ) eB [GeV ] f=uf=df=s ψ - ψ > d + < ψ - ψ > s )/((m d +m s ) χ K )eB [GeV ] FIG. 3. Left: Ratio of h ¯ ψψ i f / ( m f χ ps f ) as a function of eB for f as up ( u ), down ( d ) and strange ( s ) quark flavors. Right:Ratio of ( h ¯ ψψ i d + h ¯ ψψ i s ) / (( m d + m s ) χ K ) as a function of eB . D. UV divergence of quark chiral condensates
To investigate the GMOR relation we need to take care of the UV divergence in the light and strange quark chiralcondensates at zero and nonzero magnetic field. Since it has been shown in Ref. [9] that the UV divergence part ofthe chiral condensate is independent of the magnetic field, we can obtain the UV-free chiral condensate at nonzeromagnetic field by subtracting the UV-divergence in the chiral condensate, i.e. h ¯ ψψ i UV f = l,s obtained at the zero magneticfield. To obtain h ¯ ψψ i UV l,s we thus look into the Dirac spectrum representation of the subtracted chiral condensate atzero magnetic field, h ¯ ψψ i sub ≡ h ¯ ψψ i l − m l m s h ¯ ψψ i s = Z ∞ m l ( m s − m l ) ρ ( λ )( λ + m l )( λ + m s ) d λ, (25)where ρ ( λ ) is the eigenvalue spectral density of fermion matrix D f (cf. Eq. 9), and the light (up or down) quark andstrange quark chiral condensate, h ¯ ψψ i l and h ¯ ψψ i s , are connected to ρ ( λ ) through the following relation h ¯ ψψ i l,s = Z ∞ m l,s ρ ( λ ) λ + m l,s d λ. (26)The UV divergence part of the quark chiral condensate that is linear in quark mass is thus absent in h ¯ ψψ i sub ,while a logarithm UV divergence in the light quark chiral condensate should be negligible. Thanks to the Chebyshevfiltering technique combined with the stochastic estimate method [14, 74–77] we are able to compute the completeDirac eigenvalue spectrum ρ ( λ ) and then reproduce h ¯ ψψ i sub as well as h ¯ ψψ i l and h ¯ ψψ i s through Eqs. 25 and 26,respectively . In the Dirac eigenvalue spectrum the UV divergence part should be represented by ρ ( λ ) with λ ≥ λ UVcut .Thus the UV divergence part of light quark condensate can be expressed as h ¯ ψψ i UV l,s = Z ∞ λ UV cut m l,s ρ ( λ ) λ + m l,s d λ. (27)Given that the logarithm divergence in quark mass is negligible the UV-divergence part in the strange quark chiral < ψ - ψ > UVs /< ψ - ψ > UVl λ cut < ψ - ψ > UVs /< ψ - ψ > UVl < ψ - ψ > UVs / < ψ - ψ > s (full) < ψ - ψ > UVl / < ψ - ψ > l (full) ψ - ψ > sub / < ψ - ψ > sub (full) λ cut FIG. 4. Left: Ratio of h ¯ ψψ i UV s / h ¯ ψψ i UV l as a function of λ cut . The inset shows the ratios of h ¯ ψψ i UV l,s to the corresponding fullquark condensates as a function of λ cut . Right: Ratio of the subtracted chiral condensates with an upper cutoff λ cut in theintegration of λ in Eq. 25 to that with a full spectrum of ρ ( λ ) as a function of λ cut . The inset shows a blowup in the y axisregion close to 1. condensate should be m s /m l times as that in the light quark chiral condensate, i.e. h ¯ ψψ i UV s = 10 h ¯ ψψ i UV l in our case.Thus h ¯ ψψ i UV l,s can be determined with the smallest value of λ UV cut which makes h ¯ ψψ i UV s / h ¯ ψψ i UV l =10.To determine the value of λ UVcut we show the ratio h ¯ ψψ i UV s / h ¯ ψψ i UV l as a function of λ cut in the left plot of Fig. 4. Wesee that the ratio approaches 10 rapidly at λ cut . . λ cut . We thus pickup a value of 0.24 to be λ UV cut which makes the ratio start to approach 10 by less than 0.5%. As seen from the inset inthe left plot of Fig. 4 h ¯ ψψ i UVf ( λ cut ) itself is a rapidly decreasing function of λ cut , we also check the uncertainty in thedetermination of λ UV cut by looking to the subtracted chiral condensate. We compute h ¯ ψψ i sub ( λ cut ) as a function of λ cut being the upper bound of the integration variable λ in Eq. 25, and show its ratio to h ¯ ψψ i sub ( λ cut = ∞ ) as a function λ cut in the right plot of Fig. 4. Thus λ UV cut should be the smallest value of λ cut at which h ¯ ψψ i sub ( λ cut ) / h ¯ ψψ i sub ( λ cut = ∞ ) ≃
1. As seen from the plot the ratio approaches unity at very small values of λ cut , e.g. at λ cut ≃ .
12 the ratio is0.5% deviation from unity.Since it is not very clear how large the logarithm divergence in quark mass contributes to the UV-divergence in thequark chiral condensates, we adopt a rather wide window for the values of λ UV cut from 0.12 to 0.36 with the centralvalue 0.24 which gives h ¯ ψψ i UV s / h ¯ ψψ i UV l ≃
10. Then the obtained h ¯ ψψ i UV l ranges from about 32% to 27% of h ¯ ψψ i l at eB = 0, while h ¯ ψψ i UV s ranges from about 83% to 71% of h ¯ ψψ i s at eB = 0. IV. RESULTSA. Masses and magnetic dipole polarizabilities of light and strange pseudo-scalar mesons
We now present our results for masses of pseudo-scalar mesons calculated at 16 different values of eB ranging from0 to ∼ .
35 GeV . In the left plot of Fig. 5 we show the ratio of masses of neutral pseudo-scalar mesons at nonzero With the order of Chebyshev polynomials being 24000 and the binwidth of Dirac eigenvalue spectrum in lattice spacing being 0.002 h ¯ ψψ i sub , h ¯ ψψ i l and h ¯ ψψ i s obtained from the stochastic noise method (cf. Eq. 9) are reproduced within an accuracy of 1% via Eqs. 25and 26. eB . We found that the masses of all neutral mesonsdecrease with increasing eB and has a tendency to saturate at eB & . . By comparing the normalized massesof π u , π d , K , η s , it is obvious that the lighter hadrons are more affected by magnetic field, i.e. in the strongestmagnetic field ( eB ≃ ) we have it can be seen that M η s and M π u ( M π d ) are about 70% and 60% of theirvalues at B =0, respectively. The amount of reduction in ¯ uu and ¯ dd components of pion mass is roughly consistentwith results presented in SU(2) gauge theory [51] and SU(3) quenched QCD [52] as well as in N f = 2 + 1 QCD withstout fermions and physical pion mass in the vacuum [9], while in the former case M π u [52] decreases faster while thelatter M π u [9] decrease slower with eB compared to our current study. This could be partly due to the fact that thehadrons with larger masses are less affected by the magnetic field, as in [52] the pion mass in the vacuum is about415 MeV, while in [9] it is 135 MeV. Due to the presence of a nonzero magnetic field the SU V (2) symmetry is brokenand mixture of the u ¯ u and d ¯ d flavor contents in the neutral pion could depend on eB [52]. To determine the mixturecoefficient is beyond the scope of our current paper, and for demonstration we nevertheless show in the left plot ofFig. 5 the ground state mass of π extracted from the averaged correlation functions of u ¯ u and d ¯ d in the pseudo-scalarchannel, i.e. G π = ( G π u + G π d ) / B = 0 case [54]. As seen from the plot the ratio for π is in between thosefor π u and π d as expected. π M(B)/M(B=0) eB [GeV ] π π π K η π (|q u B u |) / M π (|q d B d |)|qB| [GeV ] FIG. 5. Left: Masses of π u , π d , K , η s normalized by their corresponding masses at eB =0 as a function of eB . Right: Ratio M π u ( | q u B u | ) /M π d ( | q d B d | ) as a function of | qB | . Here q u and q d stand for the electrical charges of u and d quarks, and B u and B d are different values of B which makes | qB | ≡ | q u B u | = | q d B d | . Since M ( B ) /M ( B = 0) deviates from unity at all the values of eB we simulated the pseudo-scalar mesons studiedhere thus cannot be considered as neutral point-like particles whose masses should be independent of eB . Also thedifferent magnitudes of the mass reduction between π u and π d may come from the different electric charges of up anddown quarks which again indicates that the inner structure of meson has been revealed. Since the internal structureof the neutral pion is probed within our current window of magnetic field, we intend to investigate the influence ofthe electrical charge of quarks on the mass of neutral pion. We thus show the ratio of M π u to M π d as a function of qB instead of eB in the right plot of Fig. 5. We found for the first time to our knowledge that after rescaling x axisfrom eB into qB , M π u ( | q u B u | ) is almost the same as M π d ( | q d B u | ) at | qB | = | q u B u | = | q d B d | and differs at most by2%. Here q u and q d are the electric charges of u and d quarks, respectively, and B u,d stands for different values of B the quark feels to make | qB | the same for up and down quarks. We call this behavior the qB scaling. This againsupports that the internal structure of pions is probed, and may reveal that the dominant degree is represented bythe single quark already starting at the smallest magnetic field we simulated, i.e. eB ≈ .
05 GeV . Note that thissmallest value of eB is about the value of M π ( B = 0) in our simulation.We now turn to the case of charged pseudo-scalar mesons, i.e. π − and K − , and show the differences of their squaredmasses from the case of zero magnetic field, i.e. M ( eB ) − M ( eB = 0) in the left plot of Fig. 6 . We see that forboth π − and K − the differences show a non-monotonous behavior in the magnetic field, i.e. first increase and thendecrease with the magnetic field strength eB and seem to saturate at eB & . . In the small magnetic field, i.e. Due to the parity in eB , the masses of their anti-particles should be the same. π M (B)-M (B=0) [GeV ]eB [GeV ] π - K - LLL π M π - /M K - eB [GeV ] LLL2*M π - /M K - (B=0) FIG. 6. Left: Differences of squared masses between the case at B = 0 and B = 0, M ( eB ) − M ( eB = 0) for π − and K − as a function of eB . Right: Ratio of the charged pion mass to charged kaon mass at nonzero magnetic field. Dashed lines inboth plots show the results of the lowest Landau level approximation, while the grey horizontal line in the right plot gives the2 times the value of M π − /M K − at B = 0. eB . .
31 GeV ( N b ≤ π − and K − respectively, canbe well described by the lowest Landau level (LLL) approximation (cf. Eq. 6) shown as the dashed line in the plot,while at larger eB , the masses start to deviate from the results of the LLL approximation and then decrease with eB . The deviation of M π − and M K − from the LLL approximation suggests that π − and K − cannot be considered aspoint-particles any more at eB & . The decreasing behavior of M ( B ) − M ( B = 0) in eB starts at a muchweaker magnetic field, i.e. eB ≃ .
63 GeV ( N b =12) compared to those obtained from quenched QCD [54], where atendency of decreasing behavior sets in only till eB ≃ . One can also observe that at large magnetic field themass of K − is less affected than π − by eB , which is probably due to the fact that mass of K − is larger than π − inthe vacuum as the case for neutral mesons. We further show the ratio M K − ( B ) /M π − ( B ) in the right plot of Fig. 6.At vanishing magnetic field, the mass of π − is about 40% of that of K − , while as the magnetic field strength eB increases M π − / M K − firstly increases as reaching up to ∼ eB ∼ , and then slightly decreases becomingflat at a value of ∼ . eB & . It is worth noting that this value is about 2 times the value of M π − /M K − at B = 0 as indicated by the grey horizontal line. π π π N M(B)/M(B=0) (eB) [GeV ] K η π π π K η β m [GeV -3 ] obtained from fits:linear in (eB) up to O((eB) )up to O((eB) ) FIG. 7. Left: Ratio of masses of various neutral pseudo-scalar mesons to their values at zero magnetic field as a functionof ( eB ) in a small magnetic field region. The solid lines represent the fits to data points with an ansatz of M ( B ) /M ( B =0) = 1 − πβ m ( eB ) /M ( B = 0) with a fit range of ( eB ) ∈ [0 , .
03] GeV except that for π u the fit range is [0, 0.02] GeV ,while dashed lines and the dash-dotted line denote fits with an ansatz of including higher order terms up to ( eB ) and ( eB ) ,respectively, with a fit range of ( eB ) ∈ [0 , .
1] GeV . Right: The obtained magnetic polarizability for the u ¯ u and d ¯ d flavorcomponents of neutral pion, π , neutral kaon and η s from the fits shown in the left plot. The fit results including fit ranges, fitansatz and χ /d.o.f. are summarized in Table II. M ( B ) = M ( B = 0) + | qB | − πM ( B = 0) β m ( eB ) − πM ( B = 0) β hm ( eB ) + O (( eB ) ) , (28)where q is the electric charge of the particle, β m is the magnetic dipole polarizability and β hm is the first order magnetichyper-polarizability. In the weak magnetic field we thus fit the ratio of neutral pseudo-scalar meson mass at nonzeromagnetic field to its value at zero magnetic field using the following ansatz [54] M ( B ) M ( B = 0) = 1 − πM ( B = 0) β m ( eB ) − πM ( B = 0) β hm ( eB ) + O (( eB ) ) . (29)We show M ( B ) /M ( B = 0) for the case of neutral pseudo-scalar mesons in a small magnetic field range in the leftplot of Fig. 7. A clear linear behavior in ( eB ) can be observed for π d , K and η s at ( eB ) . ( N b ≤ π u in a smaller window, i.e. at ( eB ) . ( N b ≤ eB ) and the fits are denoted by solid lines in the plot.The obtained magnetic dipole polarizabilities β m are shown as red points in the right plot of Fig. 7. To check theuncertainties of β m we also performed the fits to the data in a wider range of ( eB ) < . ( N b ≤
6) by addinghigher order terms in the fit ansatz. It can be seen from the left plot of Fig. 7 that the fit ansatz including termsup to ( eB ) (denoted as dashed lines) can describe the data for π d , K and η s fairly well, while an even higher orderterm, i.e. ( eB ) is needed to describe the data for π u (denoted as a dashed-dotted line). The corresponding resultsof β m from fits including terms up to ( eB ) and ( eB ) are shown as blue and black points in the right plot of Fig. 7,respectively. It can be seen that the uncertainties of β m for π d , K and η s are mild while for π u are relatively large. β m for the d ¯ d component of neutral pion is close to the half of that for η s and K , which could be due to differentconstituent quark masses and electric-magnetic interactions. In case of π u β m is about 0 .
167 obtained from the linearfit in ( eB ) and drops (increases) by about 15% (10%) obtained from fits including terms up to ( eB ) (( eB ) ). Thevalue of β m for π u is thus in the ballpark of 4 times the value of β m for π d , which is consistent with the qB scalingshown in the right plot of Fig. 5 due to ( q u ) = 4( q d ) . The qB scaling can also been seen from the values of thefirst order hyper-polarizabilities for π u and π d , i.e. β hm,π u ≃ β hm,π d as seen from Table II with the fit range N b ≤ π , which are obtained using the same fit policy as that for π d in Fig. 7 and Table II.The quality of the fit for π is similar to that for π d , and the obtained β m ( β hm ) of π from the best fit with smallest χ /d.o.f. is about 1.5 (5) times as those of both K and η s . (B)[GeV ] LLLfits with N maxb =12N maxb =16N maxb =20N maxb =24 b M (B)[GeV ] eB [GeV ] K - π - FIG. 8. Fits to the squared masses of charged pion and kaon using the ansatz Eq. 28 to extract magnetic dipole polarizability.Solid lines represent fit results with different fit ranges of N maxb =12, 16, 20 and 24 while black dashed lines represent the resultsfrom the LLL approximation. The fit results including fit ranges, fit ansatz and χ /d.o.f. are summarized in Table II. We move forward to show the results of magnetic polarizabilities for charged pseudo-scalar mesons, i.e. π − and K − based on a fit ansatz of Eq. 28 including terms up to ( eB ) . We performed fits to M ( B ) using 4 different fit ranges,3i.e. N b ≤
12, 16, 20 and 24 and the corresponding fit results are denoted by lines with N maxb =12, 16, 20 and 24 inFig. 8, respectively . The best fit is obtained with the narrowest fit range of N maxb =12 whose χ /d.o.f. is closest tounity as listed in Table II. We find that the resulting β m is consistent with zero for both π − and K − , while valuesof β hm for π − and K − are comparable to be around 0.3 GeV − . As the fit range becomes wider the quality of thefit becomes worse, and this indicates that higher order hyper-polarizabilities are needed to describe the data which isbeyond the scope of the current paper. channel β m [GeV − ] β hm [GeV − ] χ /d.o.f. fit ansatz & range π u O (( eB ) ), N b ∈ [0 , O (( eB ) ), N b ∈ [0 , O (( eB ) ), N b ∈ [0 , π d O (( eB ) ), N b ∈ [0 , O (( eB ) ), N b ∈ [0 , π O (( eB ) ), N b ∈ [0 , O (( eB ) ), N b ∈ [0 , O (( eB ) ), N b ∈ [0 , K O (( eB ) ), N b ∈ [0 , O (( eB ) ), N b ∈ [0 , η s O (( eB ) ), N b ∈ [0 , O (( eB ) ), N b ∈ [0 , π − -0.00(5) 0.4(1) 1.3 up to O (( eB ) ), N b ∈ [0 , − , N b ∈ [0 , − , N b ∈ [0 , − , N b ∈ [0 , K − -0.02(2) 0.30(9) 1.9 up to O (( eB ) ), N b ∈ [0 , − , N b ∈ [0 , − , N b ∈ [0 , − , N b ∈ [0 , β m and the first order hyper-polarizability β hm of neutral mesons π u , π d , π , K and η s as well as charged mesons π − and K − , fit ranges and fit ansatz as well as χ /d.o.f. obtained from fits shown in Figs. 7and 8. We close this subsection by comparing our results of magnetic polarizabilities (cf. Figs. 7 and 8 and Table II)to previous computations of β m for light pseudo-scalar mesons in lattice QCD which were all done in the quenchedapproximation [52–54, 78]. In these quenched QCD studies the corresponding pion mass from the tuned valence quarkmass at zero magnetic field is generally large, e.g. M π ( eB = 0) &
512 MeV in [78] , M π ( eB = 0) &
320 MeV in [53, 54]and M π ( eB = 0) &
400 MeV in [52]. One difference to be noted is that in Ref. [52] β m for π u is about twice ofthat for π d , while in our case β m is about four times as that for π d as consistent with the | qB | scaling. This couldbe due to the fact that sufficiently small magnetic field needs to be used to extract β m and relatively larger weakestmagnetic fields are applied in [52]. As pointed in Ref. [52] the hopping parameter κ used in the Wilson propagatorneeds to be along the line of constant physics in the magnetic field, i.e. κ is eB dependent, β m obtained from [78]with κ as a constant in eB might become smaller when the eB dependence of κ was taken into account as indicatedfrom study in [52]. On the other hand, studies performed using Overlap valence quarks in quenched QCD on 18 and20 lattices give a nonzero value of β m for charged pion which are close to the experiment results obtained from theCOMPASS collaboration, i.e. β m = ( − . ± . stat ± . syst ) × − fm [79]. Note that the experiment value of β m is obtained under the assumption that electrical polarizability α π = − β m [79]. While in our study the value of β m forcharged pion is consistent with zero, thus aiming for a comparison to experiment results of magnetic polarizability astudy with continuum extrapolations at physical pion mass is crucially needed. The reason we choose the smallest value of N maxb to be 12 is that at N b ≤ B. Light quark chiral condensates
As has been pointed out in Ref. [9] the presence of the external magnetic field does not introduce any new eB − dependent divergences. To take care of the additive divergences as well as the multiplicative renormalizationin the chiral condensate, we investigate the following dimensionless quantity [10],Σ l ( B ) = 2 m l M π f π (cid:0) h ¯ ψψ i l ( B = 0) − h ¯ ψψ i l ( B = 0) (cid:1) + 1 , (30)where m l ≡ m u ≡ m d is the bare quark mass for up and down quarks, and M π and f π is the mass of pion and piondecay constant respectively at zero magnetic field. In our study M π is found to be 220.61(6) MeV while f π is 96.93(2)MeV whose determination will be shown in Section IV C. Σ l =2m l /(M π f π ) * [ ψ - ψ l (B)- ψ - ψ l (0)]+1 eB/M π eB [GeV ] l=ul=d Σ l =2m l /(M π f π ) * [ ψ - ψ l (B)- ψ - ψ l (0)]+1 eB/M π eB [GeV ] l=ul=d FIG. 9. Left: Renormalized up and down quark chiral condensates as a function of magnetic field strength eB . The left plotshows chiral condensates in the whole eB region while the right one shows those in a narrower eB region. The dashed linesin both plots denote two-parameter linear fits of chiral condensates at eB ≥ . while the sold lines in the right plotrepresent power-law fits with ansatz of h ( eB ) γ +1 at eB ≤ . . We show Σ u and Σ d as a function of eB in Fig. 9. As expected that due to different electric charges of up anddown quarks up quark and down quark chiral condensates becomes non-degenerate in the nonzero magnetic field.And the up quark chiral condensate is more affected compared to that of the down quark one by the magnetic fieldwhich is probably due to | q u | > | q d | . It is obvious to see that at large magnetic field both up and down quark chiralcondensates show linear behavior in eB . We perform linear fits for these two condensates at ( eB ) & . , andthe corresponding fit results, which are shown as dashed lines in Fig. 9, describe the data fairly well. We find that the Σ u - Σ d eB/M π eB [GeV ] 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.90 1 2 3 4 5 6 7 8 9 10111213141516 Σ u - Σ d eB/M π eB [GeV ] FIG. 10. Similar as Fig. 9 but for the difference of renormalized quark chiral condensates Σ u − Σ d . The solid line represent atwo-parameter power-law fit with an ansatz of h | eB | γ while the dashed lines denote two-parameter linear fits in eB . eB . as seen from a blowup plot of the left plot of Fig. 9, i.e. the right plot of Fig. 9. At eB . both chiralcondensates increase faster than a linear behavior in eB . We thus adopt a two-parameter power-law fit ansatz of h | eB | γ + 1 to fit the chiral condensates. We found that this fit ansatz can describe the data at eB . welland the exponent γ obtained for these two condensates are almost the same, i.e. γ = 1 . γ = 1 . u − Σ d in Fig. 10. Similar as the fits we showed in Fig. 9 we performed two-parameterlinear fits (shown as the dashed line) at eB & . and two-parameter power-law fits (shown as the solid line)at eB . . . It is expected that Σ u − Σ d posses a linear behavior in large eB , and a power-law behavior withthe same exponent as that in Fig. 9 at eB . . Σ u (eB, λ UVcut ) / Σ d (eB, λ UVcut ) eB/M π eB [GeV ] λ UVcut =0 λ UVcut =0.12 λ UVcut =0.36 λ UVcut = ∞ Σ u (|q u B u |)/ Σ d (|q d B d |) |qB| [GeV ] λ UVcut =0 λ UVcut =0.12 λ UVcut =0.36 λ UVcut = ∞ FIG. 11. Left: Ratio of renormalized up quark condensate to the down quark one Σ u ( eB ) / Σ d ( eB ) as a function of eB . Right:Similar as the left plot but as a function of | qB | = | q u B u | = | q d B d | . In both plots λ UV cut = 0 corresponds to the case that thesubtrahend in the renormalized chiral condensate is h ¯ ψψ i l ( B = 0) (c.f. Eq. 30), and λ UV cut = ∞ to the unsubtracted chiralcondensate, i.e. the subtrahend is 0, and λ UV cut = 0 .
12 and 0.36 gives the uncertainty on the estimate of UV divergence part ofthe chiral condensates as the subtrahend in Eq. 31.
We further show the ratio of Σ u / Σ d in Fig. 11. To better understand the influence of the UV divergence part ofthe chiral condensate in the ratio we also investigate the following quantity,Σ l ( B, λ UV cut ) = 2 m l M π f π (cid:0) h ¯ ψψ i l ( B ) − h ¯ ψψ i UV l ( B = 0 , λ UV cut ) (cid:1) + 1 , (31)where λ UV cut is the lower limit of λ in the integration in Eq. 27 which gives the UV divergence part of the chiralcondensate, i.e. h ¯ ψψ i UV l ( B = 0 , λ UV cut ). It is obvious that when λ UV cut = 0 Eq. 31 is the same as Eq. 30, i.e. Σ l ( B, λ UV cut =0) ≡ Σ l ( B ), and when λ UV cut = ∞ h ¯ ψψ i UV l ( B = 0 , λ UV cut = ∞ ) goes to zero and the left hand side of Eq. 31 Σ l ( B, λ UV cut = ∞ ) equals to m l M π f π h ¯ ψψ i l ( B ) + 1. We see in the left plot of Fig. 11 that the ratio Σ u / Σ d obtained using differentvalues of λ UV cut all increases with increasing strength of the magnetic field eB , which is a consequence that Σ u − Σ d increases faster in eB than Σ d . Results obtained using λ UV cut = 0 and ∞ show two extreme cases that Σ l is obtained bysubtracting the full chiral chiral condensate at zero magnetic field h ¯ ψψ i l ( B = 0), and by subtracting zero, respectively.The obtained ratio Σ u ( eB ) / Σ d ( eB ) in these two extreme cases as shown in the left plot of Fig. 11 differs at most about10%. To avoid the subtraction of the infrared contribution to h ¯ ψψ i l , we estimated the UV-divergence contribution to h ¯ ψψ i l with λ UV cut ∈ [0 . , .
36] as discussed in Section III D. As shown in the plot the obtained Σ u ( eB ) / Σ d ( eB ) lies inbetween the two extreme cases and the difference between results obtained using λ UV cut = 0 .
12 and 0 .
36 is mild.Keeping in mind the fact of qB scaling for the u and d flavor components of the neutral pion mass, we also show theratio of Σ u and Σ d at same values of | qB | = | q u B u | = | q d B d | in the right plot of Fig. 11. Results obtained using variousvalues of λ UV cut are also shown. We find that the ratio Σ u ( | q u B u | ) / Σ d ( | q d B d | ) is very close to unity with deviation atmost by 3% with qB up to about 1.1 GeV . We thus conclude that the qB scaling of light quark chiral condensatesholds within an accuracy of 3% in our current window of magnetic field.To compare with the results from χ PT we then show the average of chiral condensates, i.e. Σ avg = (Σ u + Σ d ) / χ PT. Although our results are not yet continuum extrapolated and are obtainedat a single lattice cut-off with lattice spacing a ≃ .
117 fm, as compared with results from Ref. [10] where the finest6 Σ u + Σ d )/2 eB/M π eB [GeV ] χ PT 1-loop χ PT 2-loop Σ u + Σ d )/2 eB/M π eB [GeV ] Σ u + Σ d )/2 eB/M π eB [GeV ] χ PT 1-loop χ PT 2-loop
FIG. 12. Left: Average of the up and down quark chiral condensates as a function of magnetic field strength eB . Right: Sameas the left plot but a blow up in the small values of eB . The red dashed line stands for the one-loop χ PT results [16, 80] whilethe purple solid line with a grey band represents the two loop χ PT results [81]. lattice spacing is 0.1 fm our results should be close to the results in the continuum limits besides that we are using aHighly Improved Staggered fermions which has a smaller taste system breaking effects [66]. By comparing to the oneloop χ PT results extended to nonzero pion mass (dashed lines in the plot) we find that the results from χ PT describeour lattice data only at the weakest magnetic field of eB = 0.052 GeV , which is already at the scale of M π ( eB = 0).While the two-loop χ PT results are slightly larger than those from the one-loop, it has large uncertainties (denotedas the grey band) arising from the undetermined low energy constants. It is worth noting that our pion mass at eB = 0 is heavier than the physical one, and both one-loop and two-loop chiral perturbation theory give slightlysmaller results for (Σ u + Σ d ) / C. Decay constants of neutral pion and kaon and the GMOR relation
100 150 200 250 300 350 400 0 0.5 1 1.5 2 2.5 3 3.5f P [MeV] eB [GeV ] P=K P= π P= π P= π f K / f P eB [GeV ] P= π P= π P= π FIG. 13. Left: eB dependence of decay constants of π u , π d , π and K . Right: Ratio of neutral kaon decay constant to otherthree decay constants as a function of eB . We start by showing decay constants of neutral pion and kaon in Fig. 13. At eB = 0 we obtained pion decayconstant f π = 96.93(2) MeV and kaon decay constant f K = 112 . f K /f π = 1 . f π =92.1(6) MeV, f K =110.1(5) MeV and7 f K /f π =1.1917(37) [82] . As seen from the left plot of Fig. 13 all the decay constants increase with eB , where thedecay constants from the up quark flavor component of neutral pion, f π u increases most rapidly with respect to eB while f K is the least. It is interesting to see that f K becomes degenerate with the decay constant of the downquark flavor component of neutral pion f π d in the large magnetic field. The value of neutral pion decay constant f π , as extracted from an average on the correlation function of G π u and G π d in the same way as we did for M π inSection IV A, lies in between the results of f π u and f π d . We also show the ratio of f K to other three decay constantsin the right plot of Fig. 13. It is interesting to see that the ratio f K /f π d firstly decreases with increasing eB andsaturates to be unity at eB & . , while both f K /f π u and f K /f π decreases faster in eB as compared to f K /f π d , and goes below unity at eB & . . In the large magnetic field f K /f π u and f K /f π seems to saturateto be ∼ . ∼ .
85, respectively, in the range of eB ∈ (1 . , .
5) GeV , and then slightly increase by less than 5%at our largest value of eB .As seen from Fig. 5 and Fig. 11 there exists qB scaling behavior in M π u ( M π d ) and Σ u (Σ d ) we also wonderabout the case for the up and down quark flavor components of neutral pion decay constant. We show the ratio f π u ( | q u B u | ) /f π d ( | q d B d | ) as a function of | qB | in the left plot of Fig. 14. It can be clearly seen that the ratio is veryclose to 1 and the deviation is always less than 2% in our current window of magnetic field studied. Thus the qB scaling behavior is also found in the case of neutral pion decay constant. As M π u ( M π d ), Σ u (Σ d ) and f π u ( f π d ) are allrelated to the pion correlation functions, i.e. G π u ( τ ) and G π d ( τ ), we show the ratio of G π u ( τ, q u B u )/ G π d ( τ, q d B d ) asa function of temporal distance n τ at 13 different values of | qB | in the right plot of Fig. 14. We found that at largedistance, i.e. n τ close to 48 ( N τ /
2) which is most relevant to the quantities we discussed earlier, the ratio deviatesfrom unity at most by 2%. This naturally explains the consistency of the qB scaling behavior shown in M π u ( M π d ),Σ u (Σ d ) and f π u ( f π d ). f π (|q u B u |) / f π (|q d B d |) |qB| [GeV ] 0.94 0.96 0.98 1 1.02 1.04 1.06 1.08 1.1 0 10 20 30 40 50 60 70 80 90 G π (|q u B u |)/G π (|q d B d |) n τ |qB| = 0 GeV |qB| = 0.035 GeV |qB| = 0.070 GeV |qB| = 0.104 GeV |qB| = 0.139 GeV |qB| = 0.209 GeV |qB| = 0.278 GeV |qB| = 0.348 GeV |qB| = 0.418 GeV |qB| = 0.557 GeV |qB| = 0.697 GeV |qB| = 0.836 GeV |qB| = 1.115 GeV FIG. 14. Left: Ratio of f π u to f π d at same value of | qB | = | q u B u | = | q d B d | as a function of | qB | . Right: G π u ( n τ , | q u B u | )/ G π d ( n τ , | q d B d | ) as a function of temporal distance n τ at various values of | qB | = | q u B u | = | q d B d | . Since we have obtained the masses and decay constants of neutral pseudo-scalar mesons, light and strange quarkchiral condensates, we are now ready to check the validity of the 2-flavor as well as the 3-flavor GMOR relations. Weshow the corrections to the 2-flavor and 3-flavor GMOR relations in a wide window of eB from 0 to ∼ .
35 GeV inthe left and right plots of Fig. 15, respectively. In the left plot of Fig. 15 the next-to-leading order chiral correction δ π u,d (cf. Eq. 14 and Eq. 15) for u and d quark flavor components separately. To get a UV-free chiral condensate,as discussed in Section III D, we subtract the UV-divergence part h ¯ ψψ i UVl obtained using λ UVcut = 0 .
12 and 0.36 fromthe chiral condensates h ¯ ψψ i l . The results obtained using λ UVcut = 0 .
12 and 0.36 are shown as open and filled points,respectively. At zero magnetic field the correction to the GMOR relation at most is about 6%. As the GMOR relationstrictly holds in the chiral limit of quarks at zero magnetic field from the leading order chiral perturbation theory,the next to leading order chiral corrections to the 2-flavor GMOR relation at the physical pion mass M π is (6.2 ± eB = 0 is 220 MeV, the corrections to the GMOR relation at Note that there is a factor of √ √ -0.06-0.04-0.02 0 0.02 0.04 0.06 0.08 0 0.5 1 1.5 2 2.5 3 3.50 10 20 30 40 50 60 δ i λ UVcut =0.12 λ UVcut =0.36 eB/M π eB [GeV ] i = π i = π i = π -0.06-0.04-0.02 0 0.02 0.04 0.06 0.08 0 0.5 1 1.5 2 2.5 3 3.50 10 20 30 40 50 60 δ i λ UVcut =0.12 λ UVcut =0.36 eB/M π eB [GeV ] δ K eB/M π eB [GeV ] λ UVcut = 0.12 λ UVcut = 0.24 λ UVcut = 0.36
FIG. 15. The next-to-leading order chiral corrections δ i = π ,π u ,π d for the 2-flavor GMOR relation (left) and δ K for the 3-flavorGMOR relation (right). eB =0 shown in Fig. 15 is in the same ballpark. At nonzero magnetic field we see that the correction either decreasesor increases towards unity and compared to the case at eB = 0 | δ π u,d | becomes smaller at eB & . . This is tosay that the deviation of δ π u,d from zero is at most 6% in the whole range of magnetic field strength studied. We thusconclude that in our study the GMOR relation for the up and down quark flavor components of π holds with anaccuracy of ∼
6% in eB ∈ [0 , .
35) GeV . While for π the correction δ π is obtained by averaging up and down quarkflavor components of the correlation function with the disconnected part of the correlation function being ignored inour study. We also see that the GMOR relation for neutral pion π (cf. Eq. 7) holds quite well in the range of eB from 0 to 3.35 GeV .We then have a look at the next-to-leading order chiral correction δ K to the 3-flavor GMOR relation (cf. Eq. 8)shown in the right plot of Fig. 15. The UV-divergence of the strange quark chiral condensate is taken care of inthe same way as for the light quark chiral condensate. At zero magnetic field the chiral correction δ K is muchlarger than δ π due to an enhancement factor M K /M π [32]. At eB = 0 δ K is roughly in the range of (60-72)% fromour computation while chiral perturbation theory combined with the QCD sum rules give δ K =(55 ± eB increases δ K decreases and seems to saturate at about 56% at eB & . Thus inthe nonzero magnetic field the correction to the 3-flavor GMOR relation becomes smaller, and it is reduced by about15% at our largest value of eB compared to the case at eB = 0. V. CONCLUSIONS
In this paper we investigated the masses and magnetic polarizabilities of light and strange pseudo-scalar mesons,chiral condensates, as well as neutral pion and kaon decay constants in the presence of background magnetic fields with eB . .
35 GeV by performing lattice simulations of N f = 2 + 1 QCD at zero temperature using HISQ fermions with m π ≈
220 MeV on 32 ×
96 lattices at a single lattice cutoff a = 0 .
117 fm. Our main results include the qB scaling ofvarious chiral observables and the eB -dependence of the corrections to 2- and 3-flavor GMOR relations. We find thatthe qB scaling of M π u ( M π d ), Σ u (Σ d ) and f π u ( f π d ) as well as G π u ( τ )( G π d ( τ )) in the range of eB ∈ [0 . , .
35) GeV .With complete Dirac eigenvalue spectrum we are able to estimate the UV-divergence of the quark chiral condensate.This render us to study the 2-flavor and 3-flavor GMOR relations. We found that the chiral corrections to the 2-flavorand 3-flavor GMOR relations are about 6% and 65% at eB = 0 in our setup, respectively. The correction to the2-flavor GMOR relation for the neutral pion is less than 6% in the whole window of magnetic field we studied andbecomes less than 2% at our strongest magnetic field. Thus the 2-flavor GMOR relation holds true with an accuracyof 6% at eB ∈ [0 , .
35) GeV . This makes the connection of the reduction of T pc and reduction of neutral pion massin the nonzero magnetic field more clear. On a different note that the in the case of vanishing magnetic field theboth sides of the GMOR relation has similar behavior in terms of breaking fields, e.g. quark mass at eB = 0. I.e. allthe chiral observables, i.e. chiral condensates, pion decay constants as well as the pion mass decreases with lighterquark mass. However, this is not the case at nonzero magnetic fields, as the chiral condensate and neutral pion decayconstants increase with larger eB , while the neutral pion mass decreases with larger eB as restricted by the 2-flavorGMOR relation. While for the second lightest Goldstone boson, i.e. neutral kaon in the magnetic field the correction9to the 3-flavor GMOR relation decreases monotonously as eB grows and saturates at ∼ eB up to ∼ . . The nonzeroness of (connected) neutral pion mass thus disfavors an occurrence of a superconducting phase inthe current window of magnetic field. For the charged pion and kaon, their masses show a non-monotonous behaviorin eB , i.e. they firstly increase in eB following the behavior of the LLL approximation, and then increase slower thanthe LLL approximation and finally show a turning point with subsequent decreasing behavior in eB . The decreasingbehavior of the charged pion mass is different from previous computations and this may be due to the fact that themagnetic field in terms of M π ( eB = 0) is much larger compared to previous cases. While the neutral pion cannot beconsidered as a point-like particle already at our smallest value of eB ∼ .
05 GeV ∼ M π ( B = 0), the charged pionremains point-like till eB at about 6 M π ( B = 0). The obtained magnetic polarizabilities of charged and neutral pionsfrom the eB dependence of their masses are at odd with the results from quenched QCD as well as the experimentmeasurements where the magnetic polarizability is assumed to the additive inverse of the electric polarizability. Thisand together with the different behavior of M π − at large eB leaves a plenty big room for studies in the future on thepossible effects from dynamic quarks as well as discretization effects on the lattice in the weak and strong magneticfields.On the other hand the pion and kaon decay constants and their ratio are obtained as f π = 96.93(2) MeV, f K =112 . f K /f π = 1 . eB = 0 in our study. These results deviate by 5% from the state-of-the-art lattice QCD results obtained at the physical mass point and in the continuum limit. As eB increases both theneutral and kaon decay constants increase and the ratios between f K and other three decay constants ( f π , f π u , f π d )monotonously decrease and then saturate at nonzero values at large eB . Since we only deal with the decay constantsof neutral pseudo-scalar mesons that are related to the axial vector current which is parallel to the magnetic field atzero momentum, it would be interesting to study the new decay constants at nonzero momentum and those relatedto the vector current as well in the future. ACKNOWLEDGEMENTS
We thank Swagato Mukherjee for the early involvement of the work and enlightening discussions, and thank GunnarBali, Toru Kojo, Jinfeng Liao, Zhaofeng Liu and Pengfei Zhuang for interesting discussions. This work was supportedby the National Natural Science Foundation of China under grant numbers 11535012, 11775096 and 11947237, andthe RIKEN Special Postdoctoral Researcher program. The numerical simulations have been performed on the GPUcluster in the Nuclear Science Computing Center at Central China Normal University (NSC ), Wuhan, China. Appendix A: Generalized Ward-Takahashi identities in a background magnetic field
The general Ward-Takahashi identity in the continuum QCD in the nonzero magnetic field has already been derived,see e.g. Ref. [52]. Here In this appendix we will further derive the non-anomalous Ward-Takahashi identity related tothe pseudo-scalar operator and vector operator correlation functions with chiral transformation in the presence of amagnetic field interacting with two degenerate light quark flavors. The main results are shown in Eqs. A.31, A.32, A.33and A.37, where the first three identities are discussed in Fig. 3 in Section III C.We start by showing the Ward-Takahashi identity at zero magnetic field with our conventions. The expectationvalue of an observable O in QCD is given by, h O i = 1 Z Z D G µ D ¯ ψ D ψ e − S QCD O, (A.1)where O = O ( y , · · · , y n ) is an position depending operator, h i = 1 and S QCD is the Euclidean action, S QCD = Z d x (cid:16)
12 tr G µν G µν + ¯ ψ ( / D + M ) ψ (cid:17) , (A.2)and G µν = ∂ µ G ν − ∂ ν G µ − i g [ G µ , G ν ] and G µ represents the gluon field. M = diag ( m u , m d , · · · ) is a mass matrix. ψ = ( u, d, · · · ) ⊤ is a flavor multiplet of quarks, spinor and color indexes are suppressed. D µ = ∂ µ − i gG µ is thecovariant derivative at the zero magnetic field.0The Ward-Takahashi identity [83, 84] for a local transformation ( ¯ ψ ( x ) → ¯ ψ ′ ( x ) = ¯ ψ ( x ) + α ( x ) ¯ ψ ( x ) ,ψ ( x ) → ψ ′ ( x ) = ψ ( x ) + α ( x ) ψ ( x ) . (A.3)is written as, (cid:28) O δ log J δα ( x ) (cid:29) − (cid:28) O δS
QCD δα ( x ) (cid:29) + (cid:28) δOδα ( x ) (cid:29) = 0 . (A.4)Here J is the Jacobian of the transformation, which represents the anomaly. The original relation between renormal-ization factors in QED was derived in [83], and it was generalized as a non-perturbative relation from the symmetryargument [84]. We rely on the Leibniz rule with the covariant derivative, D µ (cid:0) α ( x ) φ ( x ) (cid:1) = (cid:0) D µ α ( x ) (cid:1) φ ( x ) + α ( x ) D µ φ ( x ) (A.5)and associated integrations by parts. If α ( x ) is a neutral scalar function, the covariant derivative reduces into partialderivative, D µ α ( x ) = ∂ µ α ( x ).In the case of nonzero magnetic field we use the Landau gauge for the magnetic field A µ pointing along the z ( x )direction in the infinite volume, A µ ( x ) = ( A , A , A , A ) = (0 , Bx , , D µ → ˜ D µ = ∂ µ − i gG µ − i eA µ Q , (A.6)where Q = 16 σ + 12 σ = 16 + t , (A.7)and σ i is any of three Pauli matrices with i = 1 , , σ = , and Q has the followingcommutation relation with t i ≡ σ i / Q t i = t i Q − i T i , (A.8)where T i = (cid:0) δ i t − δ i t (cid:1) .
1. Chiral rotation for the action
We show the derivation of δS QCD δα i ( x ) in this subsection. Let us firstly introduce a scalar function α i ( x ) which representsan infinitesimal local transformation. The local chiral rotation is then, ( ¯ ψ ( x ) → ¯ ψ ′ ( x ) = ¯ ψ ( x ) + α i ( x ) ¯ ψ ( x ) t i γ ,ψ ( x ) → ψ ′ ( x ) = ψ ( x ) + α i ( x ) t i γ ψ ( x ) . (A.9)The kinetic term in the action is then transformed in the leading order of α i ( x ) as,¯ ψγ µ D µ ψ → ¯ ψ ( x ) γ µ D µ ψ ( x ) + ¯ ψ ( x ) γ µ D µ (cid:0) α i ( x ) t i γ ψ ( x ) (cid:1) + α i ( x ) ¯ ψ ( x ) t i γ γ µ D µ ψ ( x ) . (A.10)The second term in the right hand side of the arrow sign in Eq. A.10 can be decomposed by the Leibniz rule,¯ ψ ( x ) γ µ D µ (cid:0) α i ( x ) t i γ ψ ( x ) (cid:1) = ¯ ψ ( x ) γ µ (cid:0) ∂ µ α i ( x ) (cid:1) t i γ ψ ( x ) + ¯ ψ ( x ) α i ( x ) γ µ D µ (cid:0) t i γ ψ ( x ) (cid:1) . (A.11)In the hereafter derivation we assume a flavor symmetry M = m , the mass term then transforms as,¯ ψ ( x ) ψ ( x ) → ¯ ψ ( x ) ψ ( x ) + 2 α i ( x ) ¯ ψ ( x ) t i γ ψ ( x ) . (A.12)We now derive the variation of the QCD action under the chiral transformation when eB = 0 and eB = 0 as follows1 a. Zero magnetic field For the case of vanishing magnetic field D µ and t i commutes with each other,¯ ψγ µ D µ ψ → ¯ ψ ( x ) γ µ D µ ψ ( x ) − α i ( x ) ∂ µ (cid:0) ¯ ψ ( x ) γ µ t i γ ψ ( x ) (cid:1) . (A.13)The variational part from the action is, δS QCD δα i ( x ) = − ∂ µ (cid:0) ¯ ψ ( x ) γ µ t i γ ψ ( x ) (cid:1) + 2 m ¯ ψ ( x ) t i γ ψ ( x ) . (A.14) b. Non-zero magnetic field We focus on the traceless part, namely we only take i = 1 , ,
3. The kinetic term istransformed as, ¯ ψγ µ ˜ D µ ψ → ¯ ψ ( x ) γ µ ˜ D µ ψ ( x ) + ¯ ψ ( x ) γ µ (cid:0) ∂ µ α i ( x ) (cid:1) t i γ ψ ( x ) − eα i ( x ) A µ ¯ ψ ( x ) γ µ γ T i ψ ( x ) . (A.15)The variational part from the action is, δS QCD δα i ( x ) = − ∂ µ (cid:0) ¯ ψ ( x ) γ µ t i γ ψ ( x ) (cid:1) + 2 m ¯ ψ ( x ) t i γ ψ ( x ) + ∆ i ( x, eB ) , (A.16)where ∆ i ( x, eB ) = − eA µ ( x ) ¯ ψ ( x ) γ µ γ T i ψ ( x ) . (A.17)For the magnetic field in the Landau gauge, i.e. A µ ( x ) = (0 , Bx , , i ( x, eB ) = − eBx ¯ ψ ( x ) γ γ ( δ i t − δ i t ) ψ ( x ) . (A.18)Since the Dirac fields transform as ¯ ψ ( t, ~x ) → ¯ ψ ( t, − ~x ) γ , ψ ( t, ~x ) → γ ψ ( t, − ~x ) under the parity, and the bilinear termtransforms as the odd parity, x ¯ ψγ γ ψ ( t, ~x ) → − x ¯ ψγ γ ψ ( t, − ~x ), ∆ i ( x, eB ) is an odd function of x .
2. Chiral rotation for observables
We show the derivation of δOδα i ( x ) with the chiral rotation for pseudo-scalar and vector operators in this subsectionas follows. a. Chiral rotation for pseudo-scalar operator The pseudo-scalar operator for the SU(2) case is P i ( y ) = ¯ ψ ( y ) γ t i ψ ( y ),where i = 1 , ,
3. It transforms as,¯ ψ ( y ) γ t i ψ ( y ) → ¯ ψ ( y ) γ t i ψ ( y ) + 12 α i ( y ) ¯ ψ ( y ) ψ ( y ) . (A.19)Thus the variation is, δP i ( y ) δα j ( x ) = δ (cid:0) ¯ ψ ( y ) γ t i ψ ( y ) (cid:1) δα j ( x ) = 12 ¯ ψ ( y ) ψ ( y ) δ ( x − y ) δ ij . (A.20) b. Chiral rotation for vector current operator A vector current operator is defined as V µi ( y ) = ¯ ψ ( y ) γ µ t i ψ ( y ). Ittransforms under the chiral transformation as,¯ ψ ( y ) γ µ t i ψ ( y ) → ¯ ψ ( y ) γ µ t i ψ ( y ) − i α j ( y ) ǫ ikj ¯ ψ ( y ) t k γ µ γ ψ ( y ) , (A.21)where ǫ ikj is the anti-symmetric tensor. The corresponding variation term is then expressed as, δV µi ( y ) δα j ( x ) = δ ( ¯ ψ ( y ) γ µ t i ψ ( y )) δα j ( x ) = − i δ ( x − y ) ǫ ikj ¯ ψ ( y ) t k γ µ γ ψ ( y ) . (A.22)
3. Non-anomalous Ward-Takahashi identity
With the variations derived in Appendices A 1 and A 2, we are ready to derive the Ward-Takahashi identities forthe non-anomalous case, (cid:28)
O δS
QCD δα j ( x ) (cid:29) = (cid:28) δOδα j ( x ) (cid:29) , (A.23)for chiral rations with pseudo-scalar and vector operators as follows in the nonzero magnetic field.2 a. Chiral rotation with pseudo-scalar operator We choose an operator as O ( y ) = P j ( y ) = ¯ ψ ( y ) γ t j ψ ( y ). The left hand side of Eq. A.23 is,LHS = − ∂ xµ (cid:10) P j ( y ) (cid:0) ¯ ψ ( x ) γ µ t i γ ψ ( x ) (cid:1)(cid:11) + 2 m (cid:10) P j ( y ) ¯ ψ ( x ) t i γ ψ ( x ) (cid:11) + (cid:10) P j ( y )∆ i ( x, eB ) (cid:11) , (A.24)By integrating over x , and the first term vanishes, Z d x LHS = 2 m Z d x (cid:10) P j ( y ) P i ( x ) (cid:11) + Z d x (cid:10) P j ( y )∆ i ( x, eB ) (cid:11) . (A.25)(A.26)And the right hand side of Eq. A.23 is,RHS = (cid:28) δP j ( y ) δα i ( x ) (cid:29) = 12 (cid:10) ¯ ψ ( y ) ψ ( y ) (cid:11) δ ( x − y ) δ ij . (A.27)Integrating over x we arrive at Z d x RHS = (cid:28) δP j ( y ) δα i ( x ) (cid:29) = 12 (cid:10) ¯ ψ ( y ) ψ ( y ) (cid:11) δ ij . (A.28)Thus the identity becomes,4 m Z d x (cid:10) P j ( y ) P i ( x ) (cid:11) + 2 Z d x (cid:10) P j ( y )∆ i ( x, eB ) (cid:11) = (cid:10) ¯ ψ ( y ) ψ ( y ) (cid:11) δ ij . (A.29)Integrating over y and dividing by four-volume V ,4 m V Z d yd x (cid:10) P i ( y ) P j ( x ) (cid:11) + 2 V Z d y Z d x (cid:10) P j ( y )∆ i ( x, eB ) (cid:11) = 1 V Z d y (cid:10) ¯ ψ ( y ) ψ ( y ) (cid:11) δ ij . (A.30)In our convention the right hand side of the above identity gives a two-flavor chiral condensate, i.e. h ¯ ψψ i u + h ¯ ψψ i d with i = j . At nonzero magnetic field the neutral pion operator is α (¯ uγ u − β ¯ dγ d ) with α + β = 1, thus in the caseof i = j = 3 and α = β = 1 / √ i = 3 for the magnetic field related term ∆ i ( x, eB ) (cf. Eq. A.18) itself vanishes. Theabove identity thus becomes 2 m l χ π = h ¯ ψψ i u + h ¯ ψψ i d . (A.31)The above relation can also be extended to the cases involving the strange quark, m s χ η s = h ¯ ψψ i s , (A.32)( m d + m s ) χ K = h ¯ ψψ i d + h ¯ ψψ i s . (A.33)Here χ π , χ η s and χ K are the space-time sum of the two-point correlation functions for neutral pion, η s ( s ¯ s ) andneutral kaon, and m l = m u = m d .The identity A.31 at zero magnetic field was also obtained using a diagrammatic method [63]. Since the magneticfield represented by the U (1) field can be factored out from the gauge field the extension to nonzero magnetic fieldremains the same, which can be simply observed from the diagrammatic method in [63]. As seen from Fig. 3 and thediscussions in Section III C the identities A.31, A.32 and A.33 hold well in the staggered discretization scheme. b. Chiral rotation with vector current operator The left hand side of Eq. A.23 for the case of the vector current operator is,LHS = − (cid:10) V µj ( y ) ∂ xµ (cid:0) ¯ ψ ( x ) γ µ t i γ ψ ( x ) (cid:1)(cid:11) + 2 m (cid:10) V µj ( y ) ¯ ψ ( x ) t i γ ψ ( x ) (cid:11) + (cid:10) V µj ( y )∆ i ( x, eB ) (cid:11) , (A.34)Integrating over x and using P i ( x ) = ¯ ψ ( x ) γ t i ψ ( x ), Z d x LHS = 2 m Z d x (cid:10) V µj ( y ) P i ( x ) (cid:11) + Z d x (cid:10) V µj ( y )∆ i ( x, eB ) (cid:11) , (A.35)3The right hand side of Eq. A.23 with integration for x can be expressed as Z d x * δV µj ( y ) δα i ( x ) + = − ǫ ikj (cid:10) ¯ ψ ( y ) t k γ µ γ ψ ( y ) (cid:11) , (A.36)By combining them and integrating over y we arrive at2 m Z d xd y (cid:10) V µj ( y ) P i ( x ) (cid:11) + Z d xd y (cid:10) V µj ( y )∆ i ( x, eB ) (cid:11) = − ǫ ikj Z d y (cid:10) ¯ ψ ( y ) t k γ µ γ ψ ( y ) (cid:11) . (A.37) Appendix B: Simulations in the HISQ discretization scheme
In this appendix we describe the implementation of the magnetic field in the lattice QCD simulations using theHISQ action and in particular the procedure to compute the fermion force. The HISQ action is constructed by theKogut-Susskind 1-link action D KS and Naik improvement term D Naik with smeared links, D HISQ [ X ( U ) , W ( U )] ≡ c D KS [ X ( U )] + c D Naik [ W ( U )] , (B.1)where all coefficients are chosen at vanishing magnetic field [85] and D KS [ X ( U )] = X µ η µ ( n ) (cid:2) X µ ( n ) δ n + µ,n ′ − X † µ ( n − ˆ µ ) δ n − µ,n ′ (cid:3) , (B.2)and D Naik [ W ( U )] = X µ η µ ( n ) (cid:2) W µ ( n ) W µ ( n + ˆ µ ) W µ ( n + 2ˆ µ ) δ n +3 µ,n ′ − W † µ ( n − ˆ µ ) W † µ ( n − µ ) W † µ ( n − µ ) δ n − µ,n ′ (cid:3) , (B.3)where η µ ( n ) is the staggered phase, X denotes level-two fat-7 smeared links and W denotes re-unitarized links, whichare constructed by thin links U .The magnetic field on the lattice is represented by U(1) links in the Landau gauge for the electric charge q , u x ( n x , n y , n z , n t ) = ( exp[ − iq ˆ BN x n y ] , ( n x = N x − u y ( n x , n y , n z , n t ) = exp[ iq ˆ Bn x ] ,u z ( n x , n y , n z , n t ) = u t ( n x , n y , n z , n t ) = 1 , where ˆ B = a B with the lattice spacing a . In the HISQ action, the magnetic field can be realized by just replacingall smeared links as D HISQ [ X, W ] → D HISQ [ uX, uW ] . (B.4)while keeping the bare links in the gauge action because gluons do not carry electric charges. Note that we multiply themagnetic field to the smeared links instead of bare links, which could be different from the case in the implementationof imaginary chemical potential as the magnetic field variable u depends on the coordinate. We suppress the Naikterm contribution for simplicity below but the extension is straightforward.We employ the RHMC algorithm to generate gauge configurations. The fermion force is defined as, F µ ( n ) = (cid:20) U µ ( n ) ∂S f [ X ] ∂U µ ( n ) (cid:21) TA (no sum) , (B.5)where S f is rationally approximated pseudofermion action and TA means removing-trace and anti-hermitianizingoperations. In the presence of magnetic fields the force is modified as, F µ ( n ) = (cid:20) U µ ( n ) ∂S f [ uX ] ∂U µ ( n ) (cid:21) TA (no sum) . (B.6)4For the case of zero magnetic field the derivative with respect to bare links can be calculated with the chain rule,which can be symbolically expressed as ∂S f [ X ] ∂U = ∂S f [ X ] ∂X ∂X∂W ∂W∂V ∂V∂U , (B.7)where V is the level one fat-7 link. The first term on the right hand side ∂S f [ X ] /∂X is formally same function form asthe standard staggered force except that smeared links are used instead of thin links. This term can be symbolicallywritten as: ∂S f [ X ] ∂X ∼ X k ψ † k ⊗ Ψ k + · · · (B.8)where Ψ k = ( D † HISQ [ X ] D HISQ [ X ] + β k ) − φ , ψ k = D HISQ [ X ]Ψ k and φ is a pseudo-fermion field and β k represents arational coefficient of order k . For the case of nonzero magnetic field Eq. B.7 becomes, ∂S f [ uX ] ∂U = ∂S f [ uX ] ∂X ∂X∂W ∂W∂V ∂V∂U . (B.9)Since X, W , and V are dummy variables in the chain rule, the above term can be rewritten as, ∂S f [ uX ] ∂U = ∂S f [ ˜ X ] ∂ ˜ X ∂ ˜ X∂W ∂W∂V ∂V∂U , (B.10)where variables with tilde represent U(1) rotated variables ˜ X = uX . Thus the functional form of the first term onthe right hand side of the equal sign is the same as that at zero magnetic field (cf. Eq. B.7). Moreover, in the secondterm, magnetic field is just factored out as ∂ ˜ X/∂W = u∂X/∂W due to the structure of fat-7 links. Finally in thepresence of magnetic fields the fermion force is, F mag-HISQ µ ( n ) = [ u µ ( n ) U µ ( n ) F ] TA (no sum) (B.11)with F = ∂S f ∂X (cid:12)(cid:12)(cid:12) X → uX ∂X∂W ∂W∂V ∂V∂U , (B.12)where ∂S f ∂X (cid:12)(cid:12) X → uX means that the staggered force term with an argument uX instead of X . Thus the force calculationin the molecular dynamics is summarized as follows,1. Prepare smeared links X, W and V from U ,2. Multiply U(1) variable u on the smeared links X ,3. Calculate a part of fermion force, Ψ k = ( D † HISQ [ uX ] D HISQ [ uX ] + β k ) − φ and ψ k = D HISQ [ X ]Ψ k for all k forthe rational approximation,4. Remove U(1) variable u from the smeared links uX ,5. Construct F ,6. Finalize the force F mag-HISQ µ ( n ). [1] D. E. Kharzeev, L. D. McLerran, and H. J. Warringa, The Effects of topological charge change in heavy ion collisions:’Event by event P and CP violation’, Nucl. Phys. A803 , 227 (2008), arXiv:0711.0950 [hep-ph].[2] V. Skokov, A. Yu. Illarionov, and V. Toneev, Estimate of the magnetic field strength in heavy-ion collisions,Int. J. Mod. Phys.
A24 , 5925 (2009), arXiv:0907.1396 [nucl-th].[3] W.-T. Deng and X.-G. Huang, Event-by-event generation of electromagnetic fields in heavy-ion collisions,Phys. Rev.
C85 , 044907 (2012), arXiv:1201.5108 [nucl-th].[4] T. Vachaspati, Magnetic fields from cosmological phase transitions, Phys. Lett.
B265 , 258 (1991). [5] K. Enqvist and P. Olesen, On primordial magnetic fields of electroweak origin, Physics Letters B , 178 (1993).[6] M. D’Elia, S. Mukherjee, and F. Sanfilippo, QCD Phase Transition in a Strong Magnetic Background,Phys. Rev. D82 , 051501 (2010), arXiv:1005.5365 [hep-lat].[7] I. A. Shovkovy, Magnetic Catalysis: A Review, Lect. Notes Phys. , 13 (2013), arXiv:1207.5081 [hep-ph].[8] H.-T. Ding, C. Schmidt, A. Tomiya, and X.-D. Wang, Chiral phase structure of three flavor QCD in a background magneticfield, (2020), arXiv:2006.13422 [hep-lat].[9] G. S. Bali, F. Bruckmann, G. Endrodi, Z. Fodor, S. D. Katz, S. Krieg, A. Schafer, and K. K. Szabo, The QCD phasediagram for external magnetic fields, JHEP , 044, arXiv:1111.4956 [hep-lat].[10] G. S. Bali, F. Bruckmann, G. Endrodi, Z. Fodor, S. D. Katz, and A. Schafer, QCD quark condensate in external magneticfields, Phys. Rev. D86 , 071502 (2012), arXiv:1206.4205 [hep-lat].[11] E. M. Ilgenfritz, M. Muller-Preussker, B. Petersson, and A. Schreiber, Magnetic catalysis (and inverse catalysis) at finitetemperature in two-color lattice QCD, Phys. Rev.
D89 , 054512 (2014), arXiv:1310.7876 [hep-lat].[12] V. G. Bornyakov, P. V. Buividovich, N. Cundy, O. A. Kochetkov, and A. Sch¨afer, Deconfinement transition in two-flavor lattice QCD with dynamical overlap fermions in an external magnetic field, Phys. Rev.
D90 , 034501 (2014),arXiv:1312.5628 [hep-lat].[13] G. S. Bali, F. Bruckmann, G. Endr¨odi, S. D. Katz, and A. Sch¨afer, The QCD equation of state in background magneticfields, JHEP , 177, arXiv:1406.0269 [hep-lat].[14] A. Tomiya, H.-T. Ding, X.-D. Wang, Y. Zhang, S. Mukherjee, and C. Schmidt, Phase structure of three flavor QCD inexternal magnetic fields using HISQ fermions, Proceedings, 36th International Symposium on Lattice Field Theory (Lattice2018): East Lansing, MI, United States, July 22-28, 2018 , PoS
LATTICE2018 , 163 (2019), arXiv:1904.01276 [hep-lat].[15] M. D’Elia and F. Negro, Chiral Properties of Strong Interactions in a Magnetic Background,Phys. Rev. D , 114028 (2011), arXiv:1103.2080 [hep-lat].[16] J. O. Andersen, W. R. Naylor, and A. Tranberg, Phase diagram of QCD in a magnetic field: A review,Rev. Mod. Phys. , 025001 (2016), arXiv:1411.7176 [hep-ph].[17] T. Kojo and N. Su, The quark mass gap in a magnetic field, Phys. Lett. B720 , 192 (2013), arXiv:1211.7318 [hep-ph].[18] F. Bruckmann, G. Endrodi, and T. G. Kovacs, Inverse magnetic catalysis and the Polyakov loop, JHEP , 112,arXiv:1303.3972 [hep-lat].[19] K. Fukushima and Y. Hidaka, Magnetic Catalysis Versus Magnetic Inhibition, Phys. Rev. Lett. , 031601 (2013),arXiv:1209.1319 [hep-ph].[20] M. Ferreira, P. Costa, O. Louren¸co, T. Frederico, and C. Providˆencia, Inverse magnetic catalysis in the (2+1)-flavor Nambu-Jona-Lasinio and Polyakov-Nambu-Jona-Lasinio models, Phys. Rev. D89 , 116011 (2014), arXiv:1404.5577 [hep-ph].[21] L. Yu, H. Liu, and M. Huang, Spontaneous generation of local CP violation and inverse magnetic catalysis,Phys. Rev. D , 074009 (2014), arXiv:1404.6969 [hep-ph].[22] B. Feng, D. Hou, H.-c. Ren, and P.-p. Wu, Bose-Einstein Condensation of Bound Pairs of Relativistic Fermions in aMagnetic Field, Phys. Rev. D , 085019 (2016), arXiv:1512.08894 [hep-ph].[23] X. Li, W.-J. Fu, and Y.-X. Liu, Thermodynamics of 2+1 Flavor Polyakov-Loop Quark-Meson Model under ExternalMagnetic Field, Phys. Rev. D99 , 074029 (2019), arXiv:1902.03866 [hep-ph].[24] S. Mao, From inverse to delayed magnetic catalysis in a strong magnetic field, Phys. Rev.
D94 , 036007 (2016),arXiv:1605.04526 [hep-th].[25] U. G¨ursoy, I. Iatrakis, M. J¨arvinen, and G. Nijs, Inverse Magnetic Catalysis from improved Holographic QCD in theVeneziano limit, JHEP , 053, arXiv:1611.06339 [hep-th].[26] K. Xu, J. Chao, and M. Huang, Spin polarization inducing diamagnetism, inverse magnetic catalyasis and saturationbehavior of charged pion spectra, (2020), arXiv:2007.13122 [hep-ph].[27] M. D’Elia, F. Manigrasso, F. Negro, and F. Sanfilippo, QCD phase diagram in a magnetic background for different valuesof the pion mass, Phys. Rev. D98 , 054509 (2018), arXiv:1808.07008 [hep-lat].[28] G. Endrodi, M. Giordano, S. D. Katz, T. G. Kov´acs, and F. Pittler, Magnetic catalysis and inverse catalysis for heavypions, JHEP , 007, arXiv:1904.10296 [hep-lat].[29] C. Bonati, M. D’Elia, M. Mariti, M. Mesiti, F. Negro, A. Rucci, and F. Sanfilippo, Magnetic field effects on the staticquark potential at zero and finite temperature, Phys. Rev. D94 , 094007 (2016), arXiv:1607.08160 [hep-lat].[30] M. Gell-Mann, R. J. Oakes, and B. Renner, Behavior of current divergences under SU(3) x SU(3),Phys. Rev. , 2195 (1968).[31] P. Boucaud et al. (ETM), Dynamical twisted mass fermions with light quarks, Phys. Lett.
B650 , 304 (2007),arXiv:hep-lat/0701012 [hep-lat].[32] J. Gasser and H. Leutwyler, Chiral Perturbation Theory: Expansions in the Mass of the Strange Quark,Nucl. Phys. B , 465 (1985).[33] M. Jamin, Flavor symmetry breaking of the quark condensate and chiral corrections to the Gell-Mann-Oakes-Rennerrelation, Phys. Lett. B , 71 (2002), arXiv:hep-ph/0201174.[34] J. Bordes, C. Dominguez, P. Moodley, J. Penarrocha, and K. Schilcher, Chiral corrections to the SU (2) × SU (2) Gell-Mann-Oakes-Renner relation, JHEP , 064, arXiv:1003.3358 [hep-ph].[35] J. Bordes, C. Dominguez, P. Moodley, J. Penarrocha, and K. Schilcher, Corrections to the SU ( ) × SU ( ) Gell-Mann-Oakes-Renner relation and chiral couplings L r and H r , JHEP , 102, arXiv:1208.1159 [hep-ph].[36] J. Gasser and H. Leutwyler, Light Quarks at Low Temperatures, Phys. Lett. B184 , 83 (1987).[37] I. A. Shushpanov and A. V. Smilga, Quark condensate in a magnetic field, Phys. Lett.
B402 , 351 (1997),arXiv:hep-ph/9703201 [hep-ph]. [38] N. O. Agasian and I. A. Shushpanov, Gell-Mann-Oakes-Renner relation in a magnetic field at finite temperature,JHEP , 006, arXiv:hep-ph/0107128 [hep-ph].[39] H. T. Ding et al. , Chiral Phase Transition Temperature in ( 2+1 )-Flavor QCD, Phys. Rev. Lett. , 062002 (2019),arXiv:1903.04801 [hep-lat].[40] H.-T. Ding, New developments in lattice QCD on equilibrium physics and phase diagram, in (2020) arXiv:2002.11957 [hep-lat].[41] K. Hattori, T. Kojo, and N. Su, Mesons in strong magnetic fields: (I) General analyses, Nucl. Phys. A951 , 1 (2016),arXiv:1512.07361 [hep-ph].[42] Z. Wang and P. Zhuang, Meson properties in magnetized quark matter, Phys. Rev.
D97 , 034026 (2018),arXiv:1712.00554 [hep-ph].[43] S. Mao, Pions in magnetic field at finite temperature, Phys. Rev.
D99 , 056005 (2019), arXiv:1808.10242 [nucl-th].[44] M. Coppola, D. Gomez Dumm, S. Noguera, and N. N. Scoccola, Neutral and charged pion properties under strong magneticfields in the NJL model, Phys. Rev.
D100 , 054014 (2019), arXiv:1907.05840 [hep-ph].[45] S. S. Avancini, R. L. S. Farias, and W. R. Tavares, Neutral meson properties in hot and magnetized quark matter:a new magnetic field independent regularization scheme applied to NJL-type model, Phys. Rev.
D99 , 056009 (2019),arXiv:1812.00945 [hep-ph].[46] K. Xu, S. Shi, H. Zhang, D. Hou, J. Liao, and M. Huang, Extracting the magnitude of magnetic field at freeze-out inheavy-ion collisions, (2020), arXiv:2004.05362 [hep-ph].[47] J. Chao, Y.-X. Liu, and L. Chang, Light charged pion in ultra-strong magnetic field, (2020), arXiv:2007.14258 [hep-ph].[48] Y. Hidaka and A. Yamamoto, Charged vector mesons in a strong magnetic field, Phys. Rev.
D87 , 094502 (2013),arXiv:1209.0007 [hep-ph].[49] M. N. Chernodub, Superconductivity of QCD vacuum in strong magnetic field, Phys. Rev.
D82 , 085011 (2010),arXiv:1008.1055 [hep-ph].[50] M. N. Chernodub, Spontaneous electromagnetic superconductivity of vacuum in strong magnetic field: evidence from theNambu–Jona-Lasinio model, Phys. Rev. Lett. , 142003 (2011), arXiv:1101.0117 [hep-ph].[51] E. V. Luschevskaya and O. V. Larina, The ρ and A mesons in a strong abelian magnetic field in SU (2) lattice gaugetheory, Nucl. Phys. B884 , 1 (2014), arXiv:1203.5699 [hep-lat].[52] G. S. Bali, B. B. Brandt, G. Endr˝odi, and B. Gl¨aßle, Meson masses in electromagnetic fields with Wilson fermions,Phys. Rev.
D97 , 034505 (2018), arXiv:1707.05600 [hep-lat].[53] E. V. Luschevskaya, O. E. Solovjeva, O. A. Kochetkov, and O. V. Teryaev, Magnetic polarizabilities of light mesons in SU (3) lattice gauge theory, Nucl. Phys. B898 , 627 (2015), arXiv:1411.4284 [hep-lat].[54] E. Luschevskaya, O. Solovjeva, and O. Teryaev, Magnetic polarizability of pion, Phys. Lett. B , 393 (2016),arXiv:1511.09316 [hep-lat].[55] V. D. Orlovsky and Yu. A. Simonov, Nambu-Goldstone mesons in strong magnetic field, JHEP , 136,arXiv:1306.2232 [hep-ph].[56] S. Fayazbakhsh and N. Sadooghi, Weak decay constant of neutral pions in a hot and magnetized quark matter,Phys. Rev. D88 , 065030 (2013), arXiv:1306.2098 [hep-ph].[57] G. S. Bali, B. B. Brandt, G. Endr˝odi, and B. Gl¨aßle, Weak decay of magnetized pions, Phys. Rev. Lett. , 072001 (2018),arXiv:1805.10971 [hep-lat].[58] M. Coppola, D. Gomez Dumm, S. Noguera, and N. N. Scoccola, Pion-to-vacuum vector and axial vector amplitudes andweak decays of pions in a magnetic field, Phys. Rev.
D99 , 054031 (2019), arXiv:1810.08110 [hep-ph].[59] M. Coppola, D. Gomez Dumm, S. Noguera, and N. N. Scoccola, Magnetic field driven enhancement of the weak decaywidth of charged pions, (2019), arXiv:1908.10765 [hep-ph].[60] M. A. Andreichikov and Yu. A. Simonov, Chiral physics in the magnetic field with quark confinement contribution,Eur. Phys. J.
C78 , 902 (2018), arXiv:1805.11896 [hep-ph].[61] H.-T. Ding, S.-T. Li, S. Mukherjee, A. Tomiya, and X.-D. Wang, Meson masses in external magnetic fields with HISQfermions, in (2020) arXiv:2001.05322 [hep-lat].[62] A. Bazavov, S. Dentinger, H.-T. Ding, et al. , Meson screening masses in (2+1)-flavor QCD,Phys. Rev.
D100 , 094510 (2019), arXiv:1908.09552 [hep-lat].[63] G. W. Kilcup and S. R. Sharpe, A Tool Kit for Staggered Fermions, Nucl. Phys.
B283 , 493 (1987).[64] S. Aoki et al. (JLQCD), Pion decay constant for the Kogut-Susskind quark action in quenched lattice QCD,Phys. Rev.
D62 , 094501 (2000), arXiv:hep-lat/9912007 [hep-lat].[65] E. Follana, Q. Mason, C. Davies, K. Hornbostel, G. P. Lepage, J. Shigemitsu, H. Trottier, and K. Wong (HPQCD, UKQCD),Highly improved staggered quarks on the lattice, with applications to charm physics, Phys. Rev.
D75 , 054502 (2007),arXiv:hep-lat/0610092 [hep-lat].[66] H.-T. Ding, F. Karsch, and S. Mukherjee, Thermodynamics of strong-interaction matter from Lattice QCD,Int. J. Mod. Phys.
E24 , 1530007 (2015), arXiv:1504.05274 [hep-lat].[67] M. Al-Hashimi and U.-J. Wiese, Discrete Accidental Symmetry for a Particle in a Constant Magnetic Field on a Torus,Annals Phys. , 343 (2009), arXiv:0807.0630 [quant-ph].[68] A. Bazavov, T. Bhattacharya, M. Cheng, C. DeTar, H.-T. Ding, et al. , The chiral and deconfinement aspects of the QCDtransition, Phys.Rev.
D85 , 054503 (2012), arXiv:1111.1710 [hep-lat].[69] H. Akaike, A new look at the statistical model identification, IEEE Transactions on Automatic Control , 716 (1974).[70] J. E. Cavanaugh, Unifying the derivations for the akaike and corrected akaike information criteria, Statistics & Probability Letters , 201 (1997).[71] C. W. Bernard, M. C. Ogilvie, T. A. DeGrand, C. E. DeTar, S. A. Gottlieb, A. Krasnitz, R. L. Sugar, and D. Toussaint,The Spatial structure of screening propagators in hot QCD, Phys. Rev. Lett. , 2125 (1992).[72] C. W. Bernard, T. Blum, T. A. DeGrand, C. E. Detar, S. A. Gottlieb, A. Krasnitz, R. L. Sugar, and D. Tou-ssaint, Finite size and quark mass effects on the QCD spectrum with two flavors, Phys. Rev. D48 , 4419 (1993),arXiv:hep-lat/9305023 [hep-lat].[73] C. W. Bernard, T. Burch, K. Orginos, D. Toussaint, T. A. DeGrand, C. E. Detar, S. Datta, S. A. Gottlieb, U. M. Heller, andR. Sugar, The QCD spectrum with three quark flavors, Phys. Rev.
D64 , 054506 (2001), arXiv:hep-lat/0104002 [hep-lat].[74] H.-T. Ding, O. Kaczmarek, F. Karsch, S.-T. Li, S. Mukherjee, A. Tomiya, and Y. Zhang, Dirac Eigenvalue spectrum of N f =2+1 QCD towards the chiral limit using HISQ fermions, in (2019) arXiv:2001.05217 [hep-lat].[75] L. Giusti and M. Luscher, Chiral symmetry breaking and the Banks-Casher relation in lattice QCD with Wilson quarks,JHEP , 013, arXiv:0812.3638 [hep-lat].[76] G. Cossu, H. Fukaya, S. Hashimoto, T. Kaneko, and J.-I. Noaki, Stochastic calculation of the Dirac spectrum on the latticeand a determination of chiral condensate in 2+1-flavor QCD, PTEP , 093B06 (2016), arXiv:1607.01099 [hep-lat].[77] Z. Fodor, K. Holland, J. Kuti, S. Mondal, D. Nogradi, and C. H. Wong, New approach to the Dirac spectral densityin lattice gauge theory applications, Proceedings, 33rd International Symposium on Lattice Field Theory (Lattice 2015):Kobe, Japan, July 14-18, 2015 , PoS
LATTICE2015 , 310 (2016), arXiv:1605.08091 [hep-lat].[78] F. X. Lee, L. Zhou, W. Wilcox, and J. C. Christensen, Magnetic polarizability of hadrons from lattice QCD in thebackground field method, Phys. Rev. D , 034503 (2006), arXiv:hep-lat/0509065.[79] C. Adolph et al. (COMPASS), Measurement of the charged-pion polarizability, Phys. Rev. Lett. , 062002 (2015),arXiv:1405.6377 [hep-ex].[80] T. D. Cohen, D. A. McGady, and E. S. Werbos, The Chiral condensate in a constant electromagnetic field,Phys. Rev. C76 , 055201 (2007), arXiv:0706.3208 [hep-ph].[81] E. S. Werbos, The Chiral condensate in a constant electromagnetic field at O(p**6), Phys. Rev.