Axial $U(1)$ symmetry at high temperature in 2-flavor lattice QCD

We investigate the axial $U(1)_A$ symmetry breaking above the critical temperature in two-flavor lattice QCD. The ensembles are generated with dynamical M\"obius domain-wall or reweighted overlap fermions. The $U(1)_A$ susceptibility is extracted from the low-modes spectrum of the overlap Dirac eigenvalues. We show the quark mass and temperature dependences of $U(1)_A$ susceptibility. Our results at $T=220 \, \mathrm{MeV}$ imply that the $U(1)_A$ symmetry is restored in the chiral limit. Its coincidence with vanishing topological susceptibility is observed.


Introduction
In quantum chromodynamics (QCD) at low temperature, the axial U(1) A symmetry is violated by the quantum (chiral) anomaly, which is the origin of much heavier η meson than other pseudoscalar mesons. The U(1) A symmetry breaking is closely related to topological excitations of the background gauge field, such as instantons. As an observables to characterize the U(1) A symmetry breaking, the U(1) A susceptibility ∆ π−δ is defined by the difference between the correlators of isovectorpseudoscalar (π a ≡ iψτ a γ 5 ψ) and isovector-scalar (δ a ≡ψτ a ψ) operators: where a is the isospin index. (We consider the theory with two degenerate quark flavors.) Above the critical temperature, T > T c , while the (spontaneously broken) chiral symmetry is known to be restored, the restoration/violation of the U(1) A symmetry is a long standing problem, which has been studied using in analytic methods [1][2][3] and effective theories [4,5] as well as lattice QCD simulations at N f = 2 [6][7][8][9] and N f = 2 + 1 [10][11][12][13]. In Ref. [1], Cohen claimed that the U(1) A symmetry of massless N f = 2 QCD is restored when the contributions from the zero modes of Dirac eigenvalues (i.e. nontrivial topological sectors) can be ignored. As a result, the mesonic correlators for π, σ, δ, and η channels become degenerate. More recently, the authors of Ref. [2] proved that the violation of the U(1) A symmetry of massless N f = 2 QCD becomes invisible under some plausible assumptions including the analyticity of the Dirac spectral density near the origin (for an alternative proof, see Ref. [3]). They even suggested a possible modification of the phase diagram in the space of up or down quark mass m u,d and strange quark mass m s (the so-called Columbia plot, see Fig. 1), according to the argument based on the effective theory [5]. Because of the restoration of U(1) A symmetry in the chiral limit for N f = 2, the chiral phase transition at m u,d = 0 may be first-order rather than the (usually expected) secondorder 1 . If this is the case, a nonzero "critical mass," m cri u,d , appears, which separates the first-order region m < m cri u,d and the crossover region m > m cri u,d for N f = 2. The existence of such a critical mass can also affect the phase structure of N f = 3 QCD.
In the previous works by the JLQCD Collaboration [6,9], we observed restored U(1) A symmetry above T c in N f = 2 lattice QCD. Because the U(1) A susceptibility is sensitive to the chiral symmetry on the lattice, we used the lattice fermions that maintain the chiral symmetry, i.e. the overlap (OV) or domain-wall (DW) fermion formalism. In Ref. [6], the U(1) A symmetry was investigated from the Dirac spectrum on gauge configurations generated with the dynamical overlap fermions under a fixed global topology, Q = 0. In Ref. [9], the gauge configurations are generated with dynamical Möbius domain-wall fermions [14,15]. Since the Ginsparg-Wilson (GW) relation [16] for the Möbius domain-wall fermion is slightly violated especially for larger lattice spacings [17], we applied the domain-wall/overlap reweighting [9], where an observable on the gauge ensembles generated with dynamical domain-wall fermions is reweighted to that of overlap fermions. The exact chiral symmetry is thus realized without sacrificing the topology sampling in Ref. [9].
In these proceedings, we further study the U(1) A symmetry above T c in N f = 2 lattice QCD. Our numerical setup is updated, compared to that in the previous paper [9]. A finer lattice spacing, 1/a = 2.64 GeV (a ∼ 0.075 fm), is used, which improves the GW relation of the Möbius domain-wall fermion.

Simulation setup 2.1 U(1) A susceptibility on the lattice
The U(1) A susceptibility (1) can be written in terms of the spectral density ρ(λ) of Dirac eigenvalues λ for fermions with a mass m. In the continuum theory, it reads where the eigenvalue spectrum is defined by On the lattice with overlap fermion formulation, the U(1) A susceptibility is given by [17] ∆ ov where λ (ov,m) i is the i-th eigenvalue of the massive overlap Dirac operator, and we set the lattice spacing a = 1. When the GW relation is violated, there are additional terms in Eq.
The subtraction of chiral zero modes in Eq. (4) is justified as follows [2]. If the GW relation is exact, the second term of Eq. (4) can be written as 2N 0 /Vm 2 , where N 0 is the number of chiral zero modes.
. Therefore, the contribution from exact zero modes should vanish in the thermodynamic limit: The overlap-Dirac eigenvalues λ (ov,m) i measured on the Möbius domain-wall fermion ensembles may include fictitious (zero and nearzero) modes induced by the partially quenched approximation [9]. They come from a mismatch between the eigenvalues and the fermion determinant, which determines the Boltzmann factor for the particular gauge configuration. After the DW/OV reweighting, the mismatch is resolved by giving negligible small reweighting factor for the gauge configuration that suffer from the fictitious modes. In Section 3, we compare λ (ov,m) i and∆ ov π−δ on the DW ensembles with those on the OV ensembles.

Numerical setup
The parameters in our numerical study are summarized in Table 1. We use the lattice with the spatial size L = 32 and the temporal length L t = 12 and 8 corresponding to T = 220 and 330 MeV, respectively, at the lattice spacing, 1/a = 2.64 GeV (a ∼ 0.075 fm). We take quark masses, am = 0.001−0.01 (2.64 − 26.4 MeV).
where "sgn" is the matrix sign function. It is approximated as tanh(L s tanh −1 (H M )) with the length of the fifth dimension, L s . It becomes the exact sign function in the limit L s → ∞. The kernel operator H M is [14,15] where D W is the Wilson-Dirac operator with a large negative mass −1/a, and we set the scale parameter α to 2.
In this study, we focus on the low-lying overlap Dirac eigenmodes to evaluate the U(1) A susceptibility from Eq. (3) or (4), but the polar approximation for the sign function with L s used in our simulations is insufficient for the low modes to satisfy the GW relation [17]. Instead, we exactly calculate the sign function for the lowest eigenmodes of H M below some threshold λ M th . Namely, we construct the overlap Dirac operator as follows [9,18]: where λ M i is the i-th eigenvalue of the kernel operator (6). In this form, the first and second terms correspond to the (separately treated) lower and higher modes, respectively. Since the approximation of the sign function is sufficiently precise for large λ M i 's, the operator constructed in this way satisfies the GW relation nearly exactly.
An observable O measured on the Möbius DW fermion ensembles may be affected by the violation of the GW relation, but they can be transformed to those on the OV fermion ensembles by the DW/OV reweighting [9]: where · · · DW and · · · ov are the ensemble average with the Möbius DW and reweighted OV fermions, respectively. R is the reweighting factor stochastically estimated on the Möbius DW ensembles.

Spectral density of overlap Dirac eigenvalues
In Fig. 2, we show the spectral density ρ(|λ|) of the overlap Dirac eigenvalues λ calculated on both the Möbius domain-wall and reweighted overlap fermion ensembles at T = 220 MeV, which is above T c . For the small quark mass (the upper panel in Fig. 2), we find that the low modes are suppressed (or a "gap" opens) due to the thermal effect. Then, one can clearly distinguish the zero modes from other higher modes. From Eq. (4), since the finite value of∆ ov π−δ comes from the nonzero low modes, such suppression of the low modes is expected to decrease the value of∆ ov π−δ . In addition, we note that the zero modes in the spectrum calculated on the DW fermion ensemble (shown by blue bins) are artifacts caused by the disagreement between the valence and sea quarks (partially quenched approximation). After the reweighting (shown by magenta bins), such artificial zero modes disappear (see the leftmost bin).
For the large quark mass (the lower panel in Fig. 2), the gap on the spectra is closed, and we cannot clearly separate the zero modes from finite modes. We expect that, the contributions from the nonzero    Figure 4. Quark mass dependences of U(1) A susceptibilities at T ∼ 220 MeV by coarser lattices with a ∼ 0.11 fm (triangle points) in our previous work [9]. The results by the finer lattice (magenta circles) are the same as Fig. 3.
modes lead to a finite value of∆ ov π−δ . Furthermore, in this case, the zero modes on the DW ensemble survive after the DW/OV reweighting, which means that such zero modes are physical, and they can be related to nonzero topological susceptibility.

U(1) A susceptibility
In the left panel of Fig. 3, we show the quark mass dependence of U(1) A susceptibility∆ ov π−δ . Magenta circles and blue squares correspond to the results on the OV and DW ensembles, respectively. Sincē ∆ ov π−δ on the DW suffers from fictitious modes by the partially quenched approximation, we expect that∆ ov π−δ on the OV is closer to the continuum limit. As expected, at a small quark mass m ∼ 3 MeV, ∆ ov π−δ on the OV nearly vanishes. Therefore, in the chiral limit, m → 0, we expect that the U(1) A symmetry is restored. Furthermore, in the region around m ∼ 10 MeV, we find a sudden increase of ∆ ov π−δ . This behavior suggests the existence of a "critical mass" as suggested in Ref. [2]. For the large quark mass region,∆ ov π−δ has a large value, implying that the U(1) A symmetry is clearly broken. Note that the results on the DW ensembles are consistent with those on the OV, for large m, and the partially quenched approximation works relatively well.
In Fig. 4, we show the results at the almost same temperature, T = 217 MeV, shown by triangle points, but on a coarser lattice (a ∼ 0.11 fm) in Ref. [9]. Although the coarser lattice leads to larger violation of the GW relation for Möbius DW fermion [17], we found that∆ ov π−δ is suppressed for the small quark mass region. Therefore, both the finer (a ∼ 0.075 fm) and coarser lattices lead to the similar suppression of∆ ov π−δ . In the right panel of Fig. 3, we show the results at T = 330 MeV which is much larger than T c . At this temperature the U(1) A susceptibility is highly suppressed at all the quark mass parameters we studied. For instance, at m ∼ 27 MeV,∆ ov π−δ after the reweighting is ∼ 10 −4 GeV 2 . This behavior comes from the suppression of Dirac low modes by the appearance of a gap in the Dirac spectra at higher temperature.

Conclusion and outlook
In this study, we investigated the U(1) A symmetry above the critical temperature from numerical simulations in two-flavor lattice QCD. The quark mass dependence of the U(1) A susceptibility at T = 220 MeV implies the restoration of U(1) A symmetry in the chiral limit, which is consistent with the predictions of [2,3]. The zero modes of the overlap Dirac spectra are closely related to the topological susceptibility, which is reported in another talk [19]. As other observables to see the restoration of U(1) A , the spatial isovector meson correlators from the same gauge configuration are investigated in Ref. [20].
In the future, we will investigate U(1) A and topological susceptibilities at other temperature. It is important to investigate the volume dependence and continuum limits in order to justify the restoration/violation of U(1) A symmetry. The study of N f = 2 + 1 is a next step to be explored. Previous studies [10][11][12][13] suggested the violation of U(1) A symmetry, and we are interested in the results of the rigorous setup of this work.