The QCD topological charge and its thermal dependence: the role of the η′

We analyze the contribution of the η′ (958) meson in the first two non-trivial moments of the QCD topological charge distribution, namely, the topological susceptibility and the fourth-order cumulant of the vacuum energy density. We perform our study within U(3) Chiral Perturbation Theory up to next-to-next-to-leading order in the combined chiral and large-Nc expansion. We also describe the temperature dependence of these two quantities and compare them with previous analyses in the literature. In particular, we discuss the validity of the thermal scaling of the topological susceptibility with the quark condensate, which is intimately connected with a Ward Identity relating both quantities. We also consider isospin breaking corrections from the vacuum misalignment at leading order in the U(3) framework.


Introduction and motivation
The rich topological structure of the QCD vacuum is encoded in the θ-angle dependence of the vacuum energy density the QCD generating functional in a θ-vacuum, β = 1/T is the inverse of the temperature, V the volume of the system and the QCD action is written as with L QCD the QCD Lagrangian at θ = 0 and with c 2 = χ top the topological susceptibility and c 4 the fourth-order cumulant, which in Euclidean space-time read with T dx = β 0 dτ d 3 x. The topological susceptibility is meant to be connected with the η mass through the U A (1) anomaly. More specifically, in the quenched approximations where quarks loops are absent, χ quenched top is related to the mass of the η (or rather its singlet part η 0 ) for N f massless quarks as [1,2] where M 0 denotes the anomalous contribution to the η mass and F is the pion decay constant in the chiral m u = m d = m s = 0 limit. The quenched approximation is formally valid in the N c → ∞ limit, where meson-loop contributions are suppressed and the mass of the η 0 becomes of the same order than the other members of the Nambu-Goldstone boson (NGB) meson octet, since the U(1) A anomaly scales with 1/N c . However, even when the fluctuation of the winding number is directly linked to the U(1) A anomaly, meson-loop corrections are indeed very relevant and one of the reasons why the QCD topological charge can be analyzed using low-energy effective field theories. This can be seen by looking at the leading-order (LO) low-energy chiral prediction for the topological susceptibility, a well-known result in Chiral Perturbation Theory (ChPT) that for three light flavors N f = 2 + 1 reads [3][4][5]  where Σ = B 0 F 2 = − qq is the single-flavor quark condensate in the chiral limit and B 0 = M 2 0π ± /(m u + m d ), with M 0π ± the mass of the charged pions. The result in (1.7) can be easily extended to a larger number of flavors [3,5], while the two-flavor result can be recovered by taking the m s m u,d limit. One of the main consequences of (1.7) is that the topological susceptibility is linearly proportional to the quark mass and hence, unlike the pure gluonic quenched result (1.6), it vanishes in the chiral limit. Thus, meson loop corrections generate terms that cancel the contribution in (1.6). This fact is actually crucial to understand why the QCD topological charge can be reliable described just by including the lightest degrees of freedom.
In addition, (1.7) also shows that at LO in the chiral expansion the topological susceptibility is proportional to the quark condensate. This property has been used in several lattice analyses to extract the value of qq in the chiral limit from χ top , both for N f = 2 [6,7] and N f = 2 + 1 [8,9]. In these analyses, χ top is determined at different quark-mass values JHEP11(2019)086 and then fitted using ChPT predictions to extract the quark condensate. This relation has been emphasized in [10] to study the strange mass paramagnetic suppression of the three-flavor condensate. Recent values of the quark condensate obtained with this method can be found in [11], while direct lattice measurements of χ top are provided in [12][13][14] for N f = 2 + 1, in [15] for N f = 2 and in [16] within the framework of two fermion families. Therefore, obtaining explicit analytic expressions for χ top and studying its quark mass dependence within the effective Lagrangian framework under different approximations is of the utmost importance from the point of view of lattice analyses.
Nevertheless, higher order corrections change the simple linear dependence of χ top with the quark condensate given in (1.7). In this sense, an important result that hints towards more complicated dependencies is the existence of the following family of Ward Identities (WI) connecting quark condensates with the topological and pseudoscalar susceptibilities in the isospin limit [17][18][19][20] where qq l = ūu +dd , m q = m u = m d and χ ll P , χ ss P are the pseudoscalar susceptibilities corresponding to the quark bilinear operators η l = i ūγ 5 u +dγ 5 d and η s = iss, respectively.
The identities (1.8) and (1.9) have been verified in U(3) ChPT up to next-to-nextleading order (NNLO) in [18]. Furthermore, the identity (1.8) has been used in lattice works to determine χ top indirectly at finite temperature [21,22] and in [23] to justify that U(1) A restoration approaches the O(4) phase transition when the chiral symmetry is exactly restored. Now, note that for N f = 2 in the isospin limit, the first term in the r.h.s. of (1.8) corresponds to the LO ChPT expression (1.7) when m s m q . Therefore, the term proportional to χ ll P necessarily includes higher order corrections in the chiral series and/or terms suppressed as m q /m s for N f = 2 + 1 flavors. An immediate conclusion is that χ top is not necessarily proportional to the light quark condensate to all orders. In fact, as we will see here in detail, the dependence with temperature of the two terms in the r.h.s. of (1.8) is completely different. Namely, the condensate term drops with T signaling chiral restoration while the χ ll P term shows a much smoother behavior. Higher order corrections to the topological susceptibility within the chiral Lagrangian framework have been obtained in [24] to next-to-leading order (NLO) (one-loop) in SU(N f ) ChPT. This result has been incorporated in the lattice analysis of [7] to extract the quark condensate and in [25,26] to provide a numerical estimate of χ top in terms of the lowenergy constants (LECs) of the O(p 4 ) effective Lagrangian. Higher order isospin-breaking corrections as well as an analysis of χ top within the so called resummed ChPT can be found in [10], whereas in the recent work [26] NNLO and electromagnetic corrections in SU (2) ChPT have been computed.

JHEP11(2019)086
Another important motivation for the study of the topological susceptibility is its relation with the Peccei-Quinn axion [27,28] and hence with various cosmological and astrophysical implications. The axion mass is directly proportional to χ top , which allows for numerical estimates based on chiral Lagrangians [25,26,29]. In addition, the fourthorder self-coupling of the axion field can be obtained from the fourth-order cumulant c 4 of the vac (θ) expansion in (1.4), which has been thoroughly analyzed within SU(N f ) ChPT up to NLO [24,30,31] and computed in the lattice [9,[32][33][34].
In addition, the large-N c behavior of the topological susceptibility and the fourthorder cumulant have been recently analyzed in [35] and compared with the lattice results provided at different large-N c values in [32].
The thermal dependence of the vacuum energy density (θ) is important for several reasons: as stated above χ top (T ) plays an important role in the relation between chiral and U(1) A restoration. In addition, the connection of χ top with the light-and strangequark condensate as a function of temperature is relevant for lattice analyses and for the understanding of the temperature dependence of the different contributions in the WI (1.8)-(1.9). Thermal corrections to the axion potential and its mass are also of importance in the cosmological and astrophysical context [25]. A ChPT analysis for χ top (T ) for N f = 2 has been performed in [25]. In that paper, the fact that for N f = 2 the one-loop corrections to χ top can be encoded in the physical pion mass and decay constant is used to establish the scaling χ top (T )/χ top (0) = qq l (T )/ qq l (0), valid at that order. Actually, that scaling law is nothing but the first term in the r.h.s. of (1.8), which opens the question of how relevant are the additional corrections provided by the second term in that identity. As we have already mentioned, we will discuss that particular aspect in detail in the present work.
As for lattice results at finite temperature, direct measurements of the topological susceptibility have typically large errors. Results can be found e.g. in [12][13][14]36]. As commented above, indirect measurements can be obtained precisely through (1.8). Higher order cumulants at finite temperature in the lattice have been analyzed in [12,37].
With the above motivation in mind, we will carry out here a ChPT-based analysis of the topological susceptibility and the fourth-order cumulant, concentrating in particular in the following aspects: • We will provide a NNLO calculation of the topological susceptibility and the fourthorder cumulant within the formalism of U(3) ChPT, in which the singlet η 0 is incorporated as a ninth pseudo-Goldstone boson within the large-N c framework [1,38,39]. This approach allows us to study η meson effects and to assess its contribution. For instance, the inclusion of η will help to understand the relevance of meson-loop corrections to the quenched result (1.6) and their role on the vanishing of χ top in the chiral limit, as already noticed in [5]. Our calculation will also allow us to estimate numerically the effect of the additional U(3) corrections in terms of the LECs involved. As mentioned above, many lattice analyses of the quark condensate rely on the chiral expansion of χ top . In that respect, studying the influence of an additional heavier degree of freedom is important. The η case is especially significant due to its direct connection with the axial anomaly and its relevant role in the Ward identities JHEP11(2019)086 described above. In addition, our U(3) analysis will provide a natural way to establish contact with recent large-N c analysis of the topological susceptibility and the fourth-order cumulant. Our study will also have the advantage of including explicitly the dependence of those quantities on the η − η mixing angle.
• We will also calculate the leading isospin-breaking m u = m d corrections within the U(3) formalism, thus extending previous SU(3) works. The importance of isospin breaking corrections to the topological susceptibility and the fourth-order cumulant stems from the vacuum misalignment induced by the combined m u = m d and θ = 0 effects [3,4,24,26,31]. It implies corrections proportional to (m u − m d )/(m u + m d ), hence much larger than the typical isospin-breaking correction in other quantities, such as quark condensates or the π 0 η mixing [40], which are proportional to (m u − m d )/m s . Recent estimates within the SU(2) formalism show that these isospin contributions give rise to around a 4% correction to χ 1/4 top . As a natural extension to those analysis, we will include isospin breaking in the LO U(3) correction to vac (θ) to estimate its numerical effect.
• We will extend our analysis to finite temperature. In this way, we will account for corrections both within the SU(3) and U(3) formalisms for the different contributions in the Ward identities (1.8) and (1.9). This will also allow us to test the robustness of the N f = 2 scaling performed in [25] when corrections from m q /m s and η, η loops are properly incorporated. As commented before, the evolution of quark condensates and susceptibilities towards chiral restoration makes it interesting to clarify their relation with the topological susceptibility as the temperature grows, within the context of chiral and (1) A restoration. In addition, we recall that the finite temperature dependence of the η mass has been analyzed within the U(3) ChPT formalism in [41,42], in fermion models [43] and in the lattice [44], confirming the U(1) A restoring behavior, as in the recent analysis of the η − η mixing angle [20].
For that purpose, the paper is organized as follows. In section 2 we will discuss the U(3) ChPT calculation of the topological susceptibility, providing explicit analytic expressions up to NNLO. We will also analyze its various limits of interest. In section 3 we will extend the calculation to the fourth order cumulant. In section 4 we will provide numerical estimates and compare to previous approaches. The isospin breaking corrections to the U(3) susceptibility and the fourth-order cumulant will be calculated in section 5. In section 6 we will study in detail their finite temperature dependence and their connection with chiral and U(1) A restoration. Some of the explicit analytic U(3) expressions will be collected in appendix A.

The topological susceptibility to NNLO in U(3) chiral perturbation theory
Within U(3) ChPT one follows a similar approach as in standard ChPT [4]. The chiral power counting in terms of momenta and quark masses is used to construct the most JHEP11(2019)086 general effective Lagrangian for SU(N f ) pseudo-Goldstone bosons up to a given order, which ensures renormalizability order by order in the expansion. Nevertheless, in the U(3) formalism the singlet η 0 is also included as the ninth pseudo-Goldstone boson. Given the large η mass value, this can be done consistently only in the large-N c framework, since the winding number charge density ω(θ) defined in (1.3), responsible for the anomalous contribution of the η mass, is suppressed within the large-N c counting [1,3,38,39,45]. Thus, the expansion is performed in terms of a small parameter δ such that {M 2 , p 2 , T 2 , m q , m s , 1/N c } = O(δ), where M and p denote typical meson masses and momenta. In this counting, the tree-level pion decay constant F 2 = O(N c ) = O(1/δ), which hence suppresses loop diagrams. The counting of the different LECs according to their O(N c ) trace structure is given in detail in [38,[46][47][48].
Following the same steps as in [49,50], the topological susceptibility can be calculated by taking functional derivatives with respect to the vacuum angle θ. Thus, taking into account the θ-vacuum coupling in the QCD action defined in (1.2) and (1.3), one can derive expectation values or thermal correlators involving the winding number density. In this way, the topological susceptibility reads where · denotes Euclidean vacuum expectation values for T = 0 or thermal correlators for T = 0 and where we have used that ω(0) = 0. In the effective Lagrangian L eff , θ(x) appears through the operator with U = exp(iΦ/F ) = exp(iλ a φ a /F ) the NGB matrix field including the singlet contribution (i.e., det U = 1), φ i the Goldstone fields, λ a=1,...8 the Gell-Mann matrix and Nevertheless, due to the η − η mixing the flavor eigenstates η 8 and η 0 are not mass eigenstates even at LO in the chiral Lagrangian. Thus, we use the angleθ to describe their mixing at LO where c θ = cosθ and s θ = sinθ. The combination (2.2) ensures that L eff will be invariant under local U(N f ) L ×U(N f ) R transformations. The effective Lagrangians containing the X field entering in χ top and the

JHEP11(2019)086
fourth-order cumulant up to the order needed for our purposes here read where M 2 0 is the contribution to the tree-level mass of the singlet η 0 in the chiral limit, i.e., its anomalous contribution as given in (1.6), and are LECs associated to the η − η mixing. As in the standard SU(3) formalism, the LECs are renormalized to absorb divergences coming from the loops.
The first non-vanishing contribution to the topological susceptibility is O(δ 0 ). On the one hand, the second term in the r.h.s. of (2.1) gives rise at this order to a constant term, which is nothing but the contribution (1.6) with M mq=0 η 0 2 = M 2 0 and N f = 3. On the other hand, the first term in the r.h.s. of (2.1) gives rise to terms of the type , tree level two-point functions at p = 0, where η ( ) denotes generically η or η fields. Gathering the two types of contribution yields in the isospin limit where M 0η and M 0η are respectively the η and η masses at tree level. They depend on the quark (or meson) masses and on M 0 . Their explicit expressions in the isospin limit m u = m d = m q can be found e.g. in [48]. Thus, in terms of quark masses, (2.8) can be recast as: wherem IL denotes the isospin limit value ofm in (1.7), i.e., Note that in this work we use the symbolm, as in (2.11), with a different meaning than the light-quark mass, which we denote here by m q . The result in (2.8) and (2.9) is the extension to N f = 2 + 1 flavors of the result given in [5] for N f degenerated flavors in the large-N c limit. Note that our U(3) LO result can be obtained from the SU(3) one by replacingm withm, which also holds when including isospin breaking corrections at LO, as we will see in section 5. In the M 0 → ∞ limit one recovers the SU(3) expression in (1.7)

JHEP11(2019)086
for m u = m d = m q , sincem IL →m IL in this limit. Furthermore, χ U(3),LO top vanishes in the chiral limit, as it also does in the SU(3) case. Actually, in this limit the second and third terms in the r.h.s. of (2.8)(coming from η, η propagators) cancel the first term, as one can check from the chiral limit behavior of s θ , M 0η , M 0η .
Our present U(3) calculation has also the advantage that one can formally recover the quenched result in (1.6) by taking the limit M 0 m. This limit can also be achieved by re-expanding (2.9) in the 1/N c expansion. Taking into account the 1/N c scaling of the different constants involved, namely which coincides with the result in (1.6) for N f = 3. The last equality in (2.12) reproduces the result given in [2], where M 0i denote tree-level meson masses. The NLO (O(δ)) and NNLO O(δ 2 ) results require including the higher-order effective Lagrangians (2.6)-(2.7) in (2.1). They involve one-loop corrections to the η-and η -meson propagators at zero momentum, which include the LECs L 6 , L 7 , L 25 , C 19 , C 31 , Λ 2 and v (2) 2 [48]. The renormalization of these LECs and the constant B 0 [38,39,48] allows one to absorb all one-loop divergences, rendering the result finite and independent of the renormalization scale µ. In addition, the calculation of χ top , χ ll P , χ ss P and qq l up to NNLO in the U(3) expansion allow one to verify the Ward Identities (1.8)-(1.9). The explicit results for the topological susceptibility at NLO and NNLO in U(3) ChPT are given explicitly in appendix A in the isospin limit. Note that the NLO and NNLO U(3) results do not correspond to simply perform the replacementm →m as it happened at LO.
In order to compare with the SU(3) calculation in ChPT, as given for instance in [10,24,25], we recall that the NLO order ChPT result is included distributed among the NLO and NNLO U(3) outcome. On the one hand, we have checked that the M 0 → ∞ limit of χ  (3) contributions, proportional to L r 6 , L 7 and log M 2 π,K,η /µ 2 . On the other hand, the surviving term proportional to log M 2 η /µ 2 is absorbed in L r 6 . Finally, we remark that our present U(3) formalism allows us to study systematically the large N c corrections in (2.12), corresponding to the Witten-Veneziano result. Performing the 1/N c expansion on the different orders in the U(3) ChPT expansion we obtain the O(1/N c ) correction to (2.12), namely The above result is consistent with the large-N c scaling analysis of the topological susceptibility provided in [35].  Figure 1. LO topologies for each of the five different connected contributions of the fourth-order cumulant in (3.1). Following the notation in [30], θ-induced vertices are denoted by its number of derivatives with respect to θ applied to the effective Lagrangian L eff . Note that we only refer to the LO for each different connected topology, irregardless of their corresponding counting in the δ expansion.

Fourth-order cumulant to NNLO in U(3) ChPT
The fourth-order cumulant is defined in (1.5) and involves the difference between the fourand two-point function square of the winding number density. Similarly to χ top , it can be computed by taking functional derivatives with respect to the θ-vacuum angle. Taking into account again the θ-angle coupling in the QCD action as defined in (1.2) and (1.3), the fourth-order cumulant can be written as where we have used once more that ω(0) = 0. The last terms of each integral are associated to the square of χ top and hence they provide disconnected contributions, i.e., terms proportional to the Euclidean four-dimensional volume that should cancel out exactly with the disconnected contributions coming from the four-point function. Thus, the calculation of the fourth-order cumulant involves five connected contributions associated to all possible combinations of a total even number of derivatives with respect to θ applied to the effective lagrangian. The LO topologies for each of these five connected contributions are depicted in figure 1. 2 Following the notation given in [30], the θ-induced vertices are

JHEP11(2019)086
classified in terms of their number of derivatives with respect to the θ-vacuum angle. Meson lines are always coupled to one θ-induced source coming from a single, second or third derivative with respect to the vacuum angle, 3 which in the isospin limit involve η and η meson propagators at zero momentum. The four different θ-induced vertices contributing to the topologies in figure 1 read where the ellipses denote further terms involving π's and K's fields, which do not contribute to connected contributions. The diagram displayed in figure 1 (a) is generated from the fourth-order derivative induced vertex in (3.5). It implies a contact term proportional to Furthermore, since diagram (f) is the only one contributing at O(1) and O(δ), it can be also dressed at higher orders in two different ways. On the one hand, the fourpoint interaction vertex can be dressed by including the higher-order diagrams depicted in figure 2. Namely, the LO four-point vertex can be replaced by its NLO counterpart in the U(3) expansion, i.e., by including a LEC, figure 2(a). Since the propagators are evaluated at p = 0 this diagram only involves the LECs L 6 , L 7 , L 8 from the L 4 Lagrangian and C 19 , C 31 from L 6 . In addition, one might add a loop to the vertex either by including a π, K or η ( ) tadpole, figure 2(b), or a rescattering diagram figure 2(c). Note that in these two cases off-shell momentum-dependent terms enter. In addition, the one-loop contribution in figure 2(c) is evaluated at p = 0, which leads to the function M 2n 0i log M 2 0i /µ 2 , where M 0i denote a π, K, η or η and n = 0, 1, 2 depending on the fields and derivatives running in the loop. Terms proportional to the LEC L 8 are O(δ). The remaining LECs and loop topologies enter only at O(δ 2 ). On the other hand, one or more of the η ( ) propagators can be dressed at NLO in the U(3) expansion. The η ( ) self energies were computed in [48] in U(3) ChPT. At zero momentum they involve the LECs Λ 2 , L 8 at O(δ) and L 6 , figure 3(a), as well as further π, K or η ( ) tadpole diagrams entering at O(δ 2 ), figure 3(b). Finally, the product of the topology in figure 2(a) proportional to L 8 and in figure 3(a) multiplying Λ 2 and L 8 also contribute at O(δ 2 ).
All together, the LO contribution in the δ expansion to the fourth-order cumulant in the isospin limit reads Like the topological susceptibility, c 4 vanishes in the chiral limit m q → 0 pointing out once more that it is a chiral quantity. The SU(3) result for the fourth-order cumulant at LO in the chiral expansion can be recovered by taking the limit M 0 → ∞. In that casem IL →m IL and one retrieves the results in [24,30]. In addition, one can also study the opposite limit, i.e., the quenched approximation for which M 0 m. In that case one obtains The NLO O(δ) results are always proportional to the O(N c ) LECs Λ 2 and L 8 , while NNLO O(δ 2 ) results involve the remaining pieces. Namely, terms proportional to the LECs 2 , L 2 8 , L 8 Λ 2 , Λ 2 2 and meson logarithms. The renormalization of the LECs [38,39,48] render a finite and scale-independent result. The explicit formulas are too long to be displayed here and they are provided as supplementary material. In addition, in the SU(3) M 0 → ∞ limit we recover the NLO results in the chiral expansion given in [30].
Finally, we can also study the large-N c expansion of the fourth-order cumulant. While the LO N c behavior of χ top was well established time ago [1,2,5], the large-N c behavior of c 4 is still under debate. On the one hand, in [32,34] it was argued from the large-N c structure of the vacuum energy density that the fourth-order cumulant should scale as O(1/N 2 c ). On the other hand, an explicit calculation based on U(N f ) ChPT at NLO in the chiral expansion for degenerate quark masses [35] where, as explained, the first term in the r.h.s. is O 1/N 2 c and the rest of the displayed terms are O 1/N 3 c .

Numerical results
In top , (−c 4 ) 1/4 and b 2 = c 4 /(12χ top ) (the latter is defined following standard lattice analyses, see below) calculated in ChPT for SU(2), SU(3) and U(3) at different orders in the chiral or δ expansion. The numerical values for the parameters involved, i.e., F , M 0 , meson masses and the LECs L r 6 , L r 7 , L r 8 , Λ 1 , Λ 2 , C 19 and C 31 are taken from [48]. Note that the U(3) LECs do not correspond to the usual SU(3) ChPT quantities, but larger differences might be expected between them. More precisely, we consider the values of the NNLOFit-B, i.e., their best fit to lattice results for the η and η masses. However, the constants L r 25 , v (0) 4 and v 2 are not included in [48]. were estimated in [42] in an additional fit to lattice data for the η and η parameters and they enter both in χ 2 to zero. Furthermore, we will neglect the NNLO U(3) corrections coming from Λ 2 2 since they are of the same other than the v  [26] for typical LECs values of SU(2) and SU(3) [52]. Note that the NNLO calculated for SU (2) in [25] includes O(p 6 ) = O(δ 3 , δ 4 ) corrections, which are beyond our present analysis. In particular, we do not recover those corrections when taking the M 0 → ∞ and m s → ∞ limits. For that reason, we do not include those results in table 1 where, as explained, we aim to discuss the different contributions that are obtained as limiting cases of our present approach. In any case, it is pointed out in [25], despite the uncertainties in the O(p 6 ) LECs, the NNLO SU(2) numerical corrections are one order of magnitude smaller than the NLO ones. In the case of the b 2 coefficient, the LO SU(2) and SU(3) expressions are LEC independent and therefore, they are given without theoretical uncertainty in table 1. The cancellation of the F 2 dependence in the U(3) LO for b 2 , e.g. from (2.9) and (3.7), explains also its smaller error compared to higher orders.
The results in table 1 are obtained in the isospin limit. In the case of the topological susceptibility they have to be compared with the lattice result χ latt top 1/4 = 73(9) MeV in that case [12]. As for lattice results on c 4 , customarily given in terms of the b 2 coefficient, they are provided only for pure gauge SU(N ) theories in [32][33][34] and for domain-wall N f = 2 + 1 fermions with large light quark masses m l /m s ≥ 0.25 [9]. Remarkably, the value of b 2 seems to be quite stable under those different approximations and close to the simplest SU(2) ChPT value, as pointed out in [12]. We quote for reference the isospin-limit value for N c = 3 gluodynamics, b 2 = −0.0216 (15) [34]. Nevertheless, more accurate lattice determinations for c 4 and b 2 for the physical N f = 2 + 1 case would be needed to make further claims.
From the results in table 1 we also observe the following features. First, in the three theoretical frameworks, the perturbative corrections remain reasonably under control. Second, although the SU(2) approach already reproduces the main contribution, the η meson and mixing angle corrections that we are including in the present work are actually comparable to the kaon and η ones introduced in the SU(3) approach. For χ top , those corrections lower the central value and get closer to the lattice prediction and so on for |b 2 |. Actually, we see that the full U(3) calculation for both observables remains compatible with the lattice results within the range provided by the LECs uncertainties, which in the case of χ top holds also for all the different approximations collected in table 1. The latter confirms that these are are good chiral quantities in the sense that they can be accurately described within ChPT.

Isospin breaking corrections to χ top and c 4
As mentioned in the introduction, isospin breaking corrections can become important for the topological charge distribution. The main reason is that for θ = 0 and m u = m d , the constant field configuration that minimizes the vacuum energy density is not U 0 = 1, but generally U 0 = diag e φ 1 , . . . , e φ N f [3,4,24,31]. In the SU (3)  det U 0 = 1 leads to j φ j = 0. However, within our present U(3) formalism, such constraint does not hold, since the determinant of U 0 is an additional degree of freedom [3]. The U 0 configuration should be such that the vacuum energy density vac (θ) defined in (1.1) is minimized, which is indeed achieved for a constant value of U 0 . Therefore, to LO in the chiral expansion, we have to consider just the usual NGB mass term plus the M 2 0 term in (2.5) for the Euclidean action. Namely, with M = diag(m u , m d , m s ) the quark mass matrix and φ j (θ) are such that they minimize LO vac (θ), i.e., Note that the solution to (5.2) is equivalent to encode the θ dependence in a complex quark mass matrix M exp(−iθ/N f ) as in [24,31] with the change of variable φ j → φ j + θ/N f . Now, following the same procedure as in [24], in order to solve the minimization problem we expand cos φ j in powers of φ j . The reason for this is that we are only interested in the power expansion of vac (θ) around θ = 0 and hence around the solution φ j = 0 of (5.2) for θ = 0. For N f = 2 and M 0 = 0 the solution of (5.2) can be found in [24,31]. The In order to provide a numerical estimate of the isospin breaking effect, we write (5.5) as The result (5.8) for the fourth-order cumulant corresponds again to the SU(3) one in [24,30,31] with the replacementm →m. As in the previous case, we also can write in terms of measurable meson parameters: ) 1/4 =53 MeV, so that the correction lies also within the 5% level when compared with the isospin-limit values in table 1.

JHEP11(2019)086
Finally, we recall that at the order we are considering the isospin-breaking corrections, temperature dependence is absent since it only enters through loop contributions.

Finite temperature dependence
From our present U(3) ChPT analysis, we can straightforwardly include the temperature dependence coming from meson loops. In the case of χ top , loop effects only arise at NNLO in the δ expansion from tadpole contributions coming from the Euclidean tree-level propagator G i (x = 0). Its finite part reads where i = π, K, η, η , M 0i are the tree level masses and µ is the renormalization scale.
In the case of c 4 , in addition to tadpoles, which enter again from the Euclidean treelevel propagator G i (x = 0) in figure 3(b) but also from the six-point interaction vertex in figure 2(b), one has to take into account the one-loop function depicted in figure 2(c), which finite part can be written in terms of That said, we want to remark that the U(3) thermal expansion is based on a perturbative calculation and hence, it is only expected to converge at very low temperatures. The fact that thermal corrections arise at NNLO in the δ expansion implies that the series breaks down as soon as thermal effects are sizable. In figure 4, we show the temperature dependence of χ top (T ) and b 2 (T ) for the same parameter values and errors used for the T = 0 results in table 1. In addition, we also plot the lattice data results for the topological susceptibility obtained in [12] and [14] and for b 2 in [12]. In the latter case, lattice errors are larger than for the susceptibility and we have actually not considered the data set for a = 0.0824 fm which has even larger errors than those showed in the figure. . Left: temperature dependence of the topological susceptibility calculated within the U(3) formalism compared to lattice data from [12] and [14]. Here, we have taken T c = 155 MeV. Right: temperature dependence of the b 2 coefficient, with lattice data from [12].
We see that the present U(3) ChPT analysis is consistent with the lattice within uncertainties, even beyond its applicability range, which as mentioned before lies at low temperatures, well below the transition. Consequently, its extrapolation close and above T c has to be taken with care. However, the good agreement with the lattice observed in figure 4 reveals once more that accounting properly for the lightest meson degrees of freedom is crucial for the description of the topological susceptibility and the cumulant.
This implies an important difference with other thermodynamic observables like the quark condensate, which accurate description based on effective theories requires the contribution of many hadronic states, like for instance in the Hadron Resonance Gas approach [53][54][55][56]. Actually, in figure 5 we compare different orders of the finite-T ChPT approach for the topological susceptibility and we see that the SU(3) and U(3) calculations represent rather small deviations from the SU(2) one.
In the case of b 2 , our theoretical U(3) result in figure 4 is quite flat with temperature, and lattice data are less accurate. This indicates that both quantities χ(T ) and c 4 (T ) decrease with T in a roughly similar way. Nevertheless, the agreement between theory and lattice results is also quite remarkable.
As discussed in section 1, another important issue regarding the temperature dependence is to what extent it can be approximated by just the scaling of the quark condensate, i.e., whether or not the second term in the Ward Identity in (1.8) can be ignored. This is actually the case if one sticks to SU(2) at NLO in the chiral expansion, what has been used in [25]. In our present work, we are calculating χ top (T ) in U(3) at NNLO, including SU(3) and SU(2) NLO as special cases. Thus, we can provide a much more accurate analysis in that respect. That issue can be also of relevance for lattice analyses. If the quark condensate terms dominate, the combined use of (1.8) and (1.9) may help to relate χ top with quantities much better determined in the lattice. Namely, is the reduced quark condensate used in lattice calculations to eliminate finite-size divergences appearing in individual condensates [21,36]. The relation (6.5) offers a way to measure indirectly χ top , which is alternative to the usual method based on the WI which stems directly from (1.8), with χ π P = − qq l /m q the pion susceptibility [57]. The reduced quark condensate is a particular example of a well determined quantity in the lattice.
The possible dominance of the quark-condensate term in the WI (1.8) is also relevant for a current topic of discussion, which has been actively studied both theoretically and in the lattice [19][20][21][22][23][58][59][60][61][62][63][64][65][66]. Namely, whether the chiral and U(1) A restoration temperatures are close enough. Since χ top and qq l are meant to vanish at exact U(1) A (asymptotically) and chiral O(4) restoration, respectively, their difference, encoded in χ ll P in (1.8), provides a direct measure of the separation between the two transitions. In this sense, it is useful to recall the behavior of these quantities near the light chiral limit m q → 0 + (M π → 0 + ), where the effects of chiral symmetry restoration are meant to be enhanced. In NNLO U(3) ChPT [20], χ top = O(m q ), while qq l (T ) and χ ll P (T ) are both O(1) quantities in the m q → 0 + limit. Thus, in the equation only the quark condensate contribution survives in the right-hand side of (6.6) in the chiral limit, which supports its dominance at low temperatures. However, near the transition qq l (T ) → 0 + in the chiral limit, while χ ll P (T ) changes much more slowly, since it is controlled by a term proportional to T 2 /M 2 0K . Thus, it brings up the question as to whether near the transition the χ ll P term can become important enough for physical masses. In figure 5, we compare the temperature dependence of the full U(3) topological susceptibility with the scaling of the light quark condensate calculated in the same framework, which would correspond to neglect the χ ll term in the WI (1.8) and (6.6). The values of JHEP11(2019)086 the LECs and other parameters involved are the same as those used in figure 4. In the same figure, we also show the simple NLO SU(2) scaling used in [25], corresponding just to the NLO ChPT quark condensate, namely The results in figure 5 show that the contribution from the additional χ ll term in the WI, although not large, may be significant as T approaches the transition point. The simple SU(2) description remains also close to the full U(3) one, which is a test of its robustness despite its simplicity. However, it is very important to point out that, although the topological susceptibility may be well described within a ChPT analysis including only the lightest degrees of freedom, that is definitively not the case for the quark condensate. To show this explicitly, we have plotted in figure 5 the quark condensate resulting from the HRG approximation provided in [56], which includes hadron states with masses up to 2 GeV and provides a very good fit to lattice condensate data. One can see clearly a much larger departure from the scaling of the topological susceptibility than the one observed with the ChPT expressions for the quark condensate, which happen to remain close to the topological susceptibility. In this sense, we remark that the addition of degrees of freedom is expected to reduce drastically the chiral condensate, as expected from approximate chiral restoration at T c . Nevertheless, that may not be the case for the topological susceptibility, whose behavior is not directly related to chiral symmetry restoration but includes U(1) A restoration features, describing lattice data just with the light degrees of freedom. Thus, in a full description of the hadron gas we do expect large deviations from the quark condensate scaling and hence significant contributions from the second term in (1.8), which becomes large close to the transition point. This also indicates that the U(1) A symmetry is still sizeably broken at the chiral transition for physical quark masses. This analysis should prevent from the use of the topological susceptibility to extract the quark condensate at finite temperature.

Conclusions
The main conclusions achieved in this work are the following: • We have provided a full calculation of the topological susceptibility and the fourthorder cumulant up to NNLO in U(3) Chiral Perturbation Theory. Our result allows one to consider the effect of the η meson consistently, as well as the η − η mixing angle dependence. As limits of interest, we recover the SU(2) and SU (3)  • We have estimated the η corrections to χ top and c 4 at zero temperature, which turn out to be of the same order as the K and η contributions. Furthermore, it provides JHEP11(2019)086 results compatible with lattice analyses, consistently with the idea that the QCD topological charge is an observable well described by the expansion around the chiral limit provided by ChPT.
• Including the dominant isospin breaking effect in the vacuum misalignment, we have provided both χ top and the fourth cumulant c 4 of the vacuum energy density expansion in the θ parameter to LO in the U(3) ChPT expansion. Numerical corrections due to isospin breaking remain below the 5% level for χ 1/4 top and (−c 4 ) 1/4 .
• We have calculated χ top and c 4 at finite temperature up to NNLO in U(3) ChPT. The temperature dependence obtained for the topological susceptibility is consistent with lattice data, supporting again that this quantity is well described by a gas made only of light mesons, unlike for instance the quark condensate. We have also discussed the relation between these two quantities, which are connected through a Ward Identity valid at all temperatures. Although the quark condensate calculated within the same ChPT formalism seems to scale quite similarly to χ top , we argue that this cannot be the case for the full hadron gas. It reveals a sizable gap between the chiral transition and the U(1) A one for physical quark masses, even though recent theoretical analysis show that those transitions tend to coincide near the chiral limit for exact chiral restoration.
A Results for χ top at NLO and NNLO in U(3) ChPT in the isospin limit We provide here the full results for the topological susceptibility in the U(3) ChPT formalism at NLO and NNLO in the δ expansion.