Asymptotic behavior of the Fourier coefficients and the analytic structure of the QCD equation of state

In this paper we study the universal properties of the baryon chemical potential Fourier coefficients in Quantum Chromodynamics. We show that by following a well-defined strategy, the Fourier coefficients can be used to locate Yang-Lee edge singularities associated with chiral phase transition (and by extension with the Roberge-Weiss) in the complex chemical potential plane. We comment on the viability of performing this analysis using lattice QCD data.


I. INTRODUCTION
Due to the sign problem, the lattice Monte-Carlo methods cannot directly probe QCD thermodynamics at non-zero real values of the baryon chemical potential.One thus has to utilize indirect approaches 1 .One of them is the evaluation of QCD thermodynamics at imaginary chemical potential with the goal of either analytic continuation [4,5], to the real axis or studying the analytic structure of the QCD equation of state in the complex chemical potential plane [6].In this paper, we follow the latter path.
Consider the QCD equation of state at zero chemical potential above the possible chiral critical point but below the Roberge-Weiss (RW) phase transition [7,8], that is, at temperatures in the range (T c , T RW ).For physical quark masses, the equation of state in this regime shows a smooth non-critical transition from hadrons to quark and gluon degrees of freedom.This transition is colloquially referred to as the "QCD crossover".The apparent smoothness of the QCD crossover obscures non-trivial critical behavior at complex values of the baryon chemical potential, where the remnants of the critical point reside, known by the name of Yang-Lee edge singularities [9,10].
These singularities are continuously connected to the associated critical points [11][12][13][14]: when two YLE singularities merge and pinch a physical axis of the corresponding thermodynamic variable (for the case of the chiral critical point, the baryon chemical potential), the critical point with corresponding critical scaling emerges.Thus, locating and especially tracking the Yang-Lee edge singularities as a function of temperature may reveal the existence and the location of the QCD critical point. 2 These singularities can be treated as standard critical points, with the exception that in contradistinction to standard criticality, there is only one (not two) relevant variable, and thus, there is only one independent critical exponent: the edge critical exponent σ YLE = 1 δ YLE .The bootstrap approach provides the most accurate estimate for the value of the critical exponent σ YLE = 0.085(1) [15].It is independent of the symmetry of the underlying universality class.
In this paper, using the universal properties of the YLE singularities, we predict 1 Despite the recent developments in complex Langevin and contour deformation methods (see e.g. Ref. [1][2][3]), direct calculations at non-zero real baryon chemical potential are still not practical. 2Mere presence of the YLE singularities does not necessitate the existence of a critical point at finite temperature.For instance, there are YLEs in the classic one-dimensional Ising model, but no critical point for T > 0.
the asymptotic behavior of the Fourier coefficients of the baryon chemical potential.
In contrast to the study in Ref. [16], our analysis is focused on the cross-over region T c < T < T RW .Using this asymptotic result as a fitting prior allows us to locate the YLE singularities in a low-energy QCD model.The model calculations are done in both mean-field approximation (in the mean-field, σ YLE = 1  2 ) and by accounting for fluctuations in the local potential approximation of the Functional Renormalization Group approach (in this case, the resulting critical exponent σ YLE = 1  5 approximate better the actual value).
The outline of the paper is as follows.In Section II, we discuss the Fourier coefficients and their asymptotic behavior based on the universal properties.In Section III, we apply our analysis to the low energy model and demonstrate that even with a limited number of the Fourier coefficients, the position of the closest YLE singularity to the imaginary axis can be accurately extracted.We conclude with Section IV.Moreover, in Appendices, we discuss the finite size effects, the required precision of calculating the baryon number at imaginary values of the chemical potential to extract k-th order Fourier coefficient, and, finally, we comment on the numerical quadrature method in computing the Fourier coefficients on the available lattice QCD data.

II. FOURIER COEFFICIENTS
The QCD partition function along with thermodynamic quantities derivable from it, is periodic in baryon chemical potential μ = µ/T : This is why analyzing the data obtained for purely imaginary values of baryon chemical potential is natural in terms of the corresponding Fourier coefficients [7,8,[17][18][19][20][21].
Specifically, from a lattice QCD perspective, it is convenient to compute the Fourier transformation of the baryon number density where we explicitly took into account the symmetry property of the baryon density μ In actual lattice QCD, the coefficients b k suffer not only from statistical error originating from the Monte-Carlo nature of the simulation3 but also from the discretization errors of the numerical integration in Eq. ( 2).In our analysis, we assume that an effort is taken to minimize both types of errors and ignore them.We review the effect of finite volume Appendix A and effect of the statistical errors in Appendix B. We propose a specific numerical quadrature for the evaluation of Eq. ( 2) in Appendix C.
The question is, then, what can we learn from analyzing the Fourier coefficients?
The answer is quite a bit.Consider the integral (2) in the complex plane, as illustrated in Fig. 1.In this figure, we took into account the presence of the YLE singularities associated with the chiral and Roberge-Weiss phase transitions.Since the continuous deformation of the contour does not change the value of the integral as long as it does not cross the singularities, one can reduce the integration to that over both sides of the cuts, as shown in Fig. 1.
To proceed further it is convenient to consider just a single singularity in the right half-plane of the complex baryon chemical potential.We will do so in the next subsection.
A. Asymptotic behavior of the Fourier coefficients for a function with a branch point Consider an odd function n B (μ) periodic in an imaginary argument having brunch points in the complex plane located at ±μ br , where μbr = μbr r + i μbr i .Here, we assume with σ > −1 and θ c > 0. In the context of the YLE singularity, θ c is the confluent critical exponent (not to be confused with θ).The regular part of the function is encoded in the coefficients a n .Our goal is to find the asymptotic behavior of the Fourier coefficients for such a function.
From the definition of the Fourier coefficients where instead of Im nB (μ = iθ) we use a more appropriate −i nB (μ = iθ).Next, using the property nB (μ) = −n B (−μ) we further simplify This is the form convenient for our analysis.To compute the integral, we will deform the contour as shown in Fig.
Therefore, we get In the second line, the contribution is due to the analytic part.
The integral over the segment (cd) is identical to the expression above except for the 2π rotation around the branch point and an extra minus sign due to the direction of the segment: Adding both integrals together cancels the analytic part to yield Absorbing k independent factors into constants A and B we finaly have Now we are ready to generalize our result for two branch points corresponding to RW and chiral YLE singularities.

B. Asymptotic behavior of the Fourier coefficients in QCD
Generalizing the result obtained in Sec.II A to the case when both YLE and RW singularities are present and accounting for the analytical properties of the partition function (singularities comes in complex conjugate pairs), we obtain Here the coefficients ÃYLE,RW and BYLE,RW are complex numbers in general.
Taking into account that Im μRW = π, we have The confluent critical exponent θ c = ν c ω = σ+1 3 ω is about 0.64 and thus leads to an appreciable suppression of the corrections.It is thus safe to drop them: When the chiral YLE singularity approaches Im μ = π/2, the sign of b k starts to change very rapidly.Indeed for a special case of μYLE i = π/2 and zero phase, one get cos πk 2 , which is zero for any odd k, positive for k = 0, 4, 8, .. and negative k = 2, 6, 10.In this case it is convenient to multiply b k by (−1) k to get  Note that although we plotted the coefficients at small values of k, Eq. ( 13) is only formally valid for asymptotically large k.
This shifts the frequency of the oscillations to lower values.
To illustrate the behavior of the Fourier coefficients, we consider the scaling assumption for the position of the Yang-Lee edge singularities.For that we will use the input from Refs.[23][24][25] on the universal locations, well known critical exponents for Z (2) and O(4) universality classes and approximation for the scaling variables The non-universal coefficients (T RW , T c , κ B 2 , and z 0 for both transitions) are taken from lattice QCD parametrizations, see Refs.[26].
Using this asymptotic formula and the scaling assumption for the locations of the singularities (the position is displayed in Fig. 3) we plotted the qualitative behavior of the Fourier coefficients in Fig. 4, see also a detailed description of the approach in application to the chiral phase transition in Ref. [27].
Lets examine the figures.Consider Fig. 4a, the period of the oscillations (sign changes in k) is about P k = 4; we thus can easily estimate the value of μYLE .This is to be compared to the input value of 0.74 π.For higher temperatures, the sign change is absent in (−1) k b k demonstrating the fact that the singularities are on the line μi = π.
The real value of the location of the YLE singularity (at smaller T ) and RW singularity (at larger T ) can be read off from the slopes in the figures.We thus see that a simple analysis of the Fourier coefficients provide a good estimate for the location of the singularities.

III. MODEL RESULTS
In this section we perform the calculations of the Fourier coefficient in Quark-meson model in both mean-field approximation and by accounting for fluctuations in the local potential approximation of the Functional Renormalization Group approach (for the later, the critical exponent at YLE singularity is numerically closer to the actual one, see Introduction).We will demonstrate that the analysis of a limited number of Fourier coefficients complemented by the asymptotic expression derived in the previous subsection is sufficient to extract the location of the Yang-Lee edge singularity with a good precision.

A. Mean-field quark-meson model
In this section, we follow the notation of Ref. [28].The Euclidean action for a quark-meson model with N f = 2 degenerate quark flavors and N c = 3 colors is given by γ µ are the Euclidean gamma matrices, τ T = (1, iγ 5 ⃗ τ ) with the Pauli matrices ⃗ τ , µ is the quark chemical potential, and ) is the effective meson potential.An explicit symmetry breaking is introduced through the source h to get massive pions in the low-temperature phase.
Assuming a homogeneous mean field with the non-zero expectation value for the iso-singlet meson (σ) we obtain the following thermodynamic potential where Here q = d 3 q (2π) 3 and we defined the thermal contribution to the quark determinant as The vacuum contribution J 0 is ultraviolet-divergent.Its finite piece depends on the meson field (see Ref. [29]) For the O(4)-symmetric part of the meson potential we use where σ0 is the solution of the equation of motion Further details of the model can be found in Ref. [29].We used the following input values to fix the parameters of the model m σ = 600 MeV, m π = 140 MeV, f π = 93 MeV and the Yukawa coupling is fixed to be g = 3.6.

B. Quark-meson model in LPA FRG
The ingredients of the model are similar to those described in the tion.The model is computed using FRG the local potential approximation.This approximation leads to σ YLE = 1/5, which is closer to the actual value.The model is described in detail in Ref. [30].We used the same values for m σ , m π , f π , and g as in the mean-field model.Additionally, we fixed the UV cut-off (a new ingredient required by the FRG approach) to 1200 MeV.Again, a straightforward visual inspection of the plot allows us to establish the period of the sign change and thus estimate the location for the YLE singularity.Consider T = 150 MeV; the period is between 6 and 7. Let's approximate it by 6.5.Thus μYLE ≈ π/6.5 ≈ 0.48; this is amusingly close the the actual value.

IV. CONCLUSIONS
In this manuscript, we derived the asymptotic behavior of the Fourier coefficients in the crossover regime.Our result differs from the one obtained in Ref. [16], where the authors used the Riemann-Lebesgue lemma and where, we believe, the critical exponents and the amplitudes at YLE singularities were misidentified.In this manuscript, we analyzed the behavior of the coefficients and concluded that the position of the YLE (closest to the real axis) singularity can be extracted knowing about 20 Fourier coefficients.Applying our analysis to lattice QCD has two caveats.First, finite size effects alter the asymptotic behavior of the coefficients.This is straightforward to account for, as we demonstrated in Appendix A. Second, lattice QCD calculations come with unavoidable statistical errors.They might significantly alter the higher-order Fourier coefficients.In Appendix B, we estimated the upper bound on the error to extract the required coefficient order.A dedicated effort might make achieving this level of accuracy feasible.In Appendix C we propose a specific numerical method, which seems most appropriate for our purpose.lattice QCD data that is available from simulations at imaginary chemical potentials [6].However, we note that the error is only decaying polynomially in 1/k, while the absolute value of the Fourier coefficients are expected to decay exponentially.I.e., the relative error will nevertheless start to grow above a certain order.
Here we test the plain piece-wise Filon quadrature for N = 10, 20 and S = 1 for the case of the mean-field quark-meson model at T = 150 MeV, as discussed in Sec.III A.
In Fig. 7 (left) we compare the results with the exact Fourier coefficients, and in Fig.Despite the fact that the absolute value of the error decays as expected, we find that the period of the Fourier coefficients is obscured already at O(N ).We note that more advanced extended Filon-type quadrature methods do exist, which we leave for future investigations.

Figure 1 :
Figure 1: Complex chemical potential plane with the chiral and RW YLE.The integration path along the imaginary chemical potential axis can be deformed to the integration around the branch point singularities and the associated cuts (Stokes lines).

Figure 2 :
Figure 2: Computation of the Fourier coefficient for the function with one singularity in the right half-plane of complex chemical potential.

2 .
In the figure, we assume that the right-most points are extended to the infinity, i.e.Re µ → ∞.The contribution of the segments (ab) and (gh) cancel each other due to the periodicity of the integrand and the opposite direction of the segments.The contributions from (bc) and (f g) is zero due to the exponential decay of exp(−ikθ) = exp(−k μ) for any k > 0. The integral around the branch point (de) is vanishing due to σ > −1.Thus, the only non-trivial contribution is due to the segments on both sides of the cut (cd) and (ef ) as demonstrated in the figure.To evaluate the contribution of these segments, we consider the parametrization for the segment iθ = s + μbr 1 π (ef ) dθ nB (μ = iθ)e −ikθ = 1 iπ e −µ br k ∞ 0 ds nB (μ = s + μbr )e −ks .(5) Now for an arbitrary power p > −1, ∞ 0 ds s p e −ks = Γ(p + 1) k p+1 .

Figure 3 :
Figure 3: Illustration of the positions of the singularities under the scaling assumptions in QCD.The stars (dots) indicate the positions of RW (chiral) YLE.Numerical labels indicate the corresponding temperatures in MeV.The horizontal dashed lines indicate Imμ = π/2 and π.

where ϕ a
and ϕ b are phases due to non trivial phases of ÃYLE and BYLE and trivial real factors were absorbed into | ÂRW | and | BRW |.This is the final expression.Note that coefficients, b k , are exponentially sensitive to the imaginary values of the positions of the YLE singularities.

Figure 4 :
Figure 4: Illustration of different qualitative behaviour of the Fourier coefficients.
Fig.4b, we get an estimate of π/3 to be compared with the input value of μYLE i

Figure 5
Figure5summarizes our results.The data was fitted using Eq.(13).The results of the fit are as follows.The fit yields μfit YLE = 0.441(2) + i 0.325(3) (μ fit YLE = 0.1156(6) + i 0.9952(5)) for T = 150(180) MeV.The actual location of the singularity is μYLE = 0.412884 + i 0.342187 (μ YLE = 0.118657 + i 1.00256) for T = 150(180) MeV.The fit reproduces the location within 7% precision.The fits are performed starting from b 5 and do not require going to asymptotically large Fourier modes.Note that the imaginary part of the chemical potential can be estimated easily without performing the fits and visually extracting the period of the sign change; for example, for T = 180 MeV, the sign changes every three points, this leads to μYLE ≈ π/3 ≈ 1.05, which is very close to the actual value.

Figure 6
Figure 6 summarizes our result.The fit yields μfit YLE = 1.483(7) + i 0.446(6) (μ fit YLE = 0.949(8) + i 0.675(11)) for T = 150(180) MeV.The actual location of the singularity is μYLE = 1.553 + i 0.4794 (μ YLE = 0.9445 + i 0.6618) for T = 150(180) MeV.Although we could not reliably extract the Fourier coefficient for k > 22 (and thus are limited to relatively modest values of k), the fit accuracy is sufficiently high as they reproduce the locations within 5% precision.As for the mean field model, the fits start from b 5 .

3 Figure 7 :
Figure 7: Comparison of the Fourier coefficients from the plain piece-wise Filon quadrature on N sub-intervals with the exact results for the case of the mean-field quark-meson model at T = 150 MeV.

(
right) we show the decadic logarithm of the absolute value of the error log 10 |b exact k − b Filon k |.The expected asymptotic behaviour of the quadrature is indicated by the solid line.