# Many flavor approach to study the nature of chiral phase transition of two-flavor QCD

###### Abstract

We perform lattice numerical simulations to study the phase transition of QCD at finite temperature to clarify the nature of the transition of massless two flavor QCD. We investigate QCD with two light and heavy quarks instead of two-flavor QCD, and focus on the light quark mass dependence of the critical heavy mass, below which the transition is of first order. The heavy quarks are incorporated into two flavor configurations in the form of the hopping parameter expansion through the reweighting technique. The nature of the transition is identified by the shape of the constraint effective potential at the critical temperature. Our result indicates that the critical heavy mass remains finite in the chiral limit of the two flavors, suggesting the phase transition of massless two-flavor QCD is of second order.

Many flavor approach to study the nature of chiral phase transition of two-flavor QCD

Ryo Iwami

Graduate School of Science and Technology, Niigata University, Niigata 950-2181, Japan

E-mail: iwami@muse.sc.niigata-u.ac.jp

\abstract@cs

## 1 Introduction

Since the outstanding work by Pisarski and Wilczek [1], clarifying the nature of chiral phase transition of massless two flavor QCD has been one of the central problems in QCD. The nature crucially depends on the presence of the flavor singlet axial symmetry at [1], and the nature itself and the effective restoration of above are under active investigation [2, 3, 4, 5, 6, 7, 8]. Clarifying this point is important not only for understanding the phase structure of QCD but also for the axion dark matter scenario. The axion abundance is essentially determined by the temperature dependence of the topological susceptibility , which vanishes when symmetry is effectively and fully restored. If vanishes very rapidly right above , too much axions would be produced, and the axion dark matter scenario becomes hard or is even excluded, depending on how rapidly it happens [9]. The lattice studies of the axion dark matter scenario recently began in the quenched approximation [9, 10, 11].

(a) 2+1 flavor QCD | (b) 2+ flavor QCD |

In this work, we apply the approach proposed in Ref. [12]. The essential idea is depicted in Fig. 1. Figure 1 (a) shows the so-called Columbia plot for 2+1 flavor QCD. There are two distinct first order regions lying around the quenched limit () and the chiral limit of three flavor QCD (), respectively. In the following, we focus only on the latter. The question is whether the massless two flavor QCD point ( and ) is inside the first order region or not. If we could trace the critical line (solid or dotted curve), the question is resolved. However, it is difficult as the critical line is located in the small quark mass region.

The situation becomes tractable by adding extra flavors. Figure 1 (b) represents the expected Columbia plot for 2+ flavor QCD, where denotes the mass of extra flavors. In this case, the critical line moves upward, and for sufficiently large the critical line becomes reachable by the hopping parameter expansion (HPE) [12]. If the critical heavy mass turns out to remain finite in the limit, it immediately means that massless two flavor QCD corresponding to the point is the outside of the first order region. An important remark is that, in the limit of or , both of 2+1 and 2+ flavor QCD end up with the same two flavor QCD. Thus, the original question is simplified to whether the critical heavy mass in the chiral limit of two flavors stays finite or not. In the following, we briefly describe our approach and the result. For details, pleasee see the full paper [13].

## 2 method

Two flavor configurations at finite temperatures are generated following the standard Hybrid Monte Carlo method at four values of two flavor mass. The effective potential is obtained from the probability distribution function (PDF) for the generalized plaquette [14, 15, 16, 12] as

(2.0) | |||||

(2.0) | |||||

(2.0) |

where and are the lattice actions for the gauge field and two flavors of light quarks, respectively. and denote the averaged plaquette and rectangle, respectively, and and satisfying are the improvement coefficients for lattice gauge action. is the quark matrix for heavy flavors. When measuring the PDF, flavors of extra heavy quarks are introduced in the form of the hopping parameter expansion (HPE) via the reweighting method as shown later.

We separate the effective potential into two parts for convenience as

(2.0) |

Here and hereafter, terms independent of are ignored, because we are interested in the differentiation of the potential with regard to . The first term is defined by

(2.0) |

and represents the constraint effective potential in two flavor QCD system. The second term of eq. (2) is defined by

(2.0) | |||||

(2.0) |

where denotes the ensemble average over two flavor configurations generated with and . Then, assuming the unimproved Wilson Dirac operator for with a sufficiently small and , in eq. (2) can be approximated as

(2.0) | |||||

(2.0) |

where

(2.0) |

In the above, the following quantities have been introduced,

(2.0) |

with the real part of Polyakov loop averaged over the spatial volume. Although eqs. (2) and (2) are algebraically identical, the equality is not necessarily trivial in numerical data because the function is approximated by

(2.0) |

Thus, both expressions are examined to check the consistency. At this order of the HPE, the number of extra heavy flavors () and their mass parameter () appear only in a single parameter . The nature of the phase transition is identified by the shape of the constraint effective potential, i.e. single- or double-well, at . By scanning , we determine the critical value , at which the first order and the crossover regions are separated.

## 3 results

Following Ref. [17], we take the Iwasaki gauge action () and the -improved Wilson fermion action with the perturbatively improved for two flavors of light quarks. Simulations are performed on lattices with 25 to 32 values at each , and 10,000 to 40,000 trajectories have been accumulated at each simulation point. Four light quark masses are ranging from to 0.1505, corresponding to being from .

In the approximation of , we take two values of = 0.0001 and 0.00025 to see the stability of the results, and the discrepancy arising from different choice is taken as a systematic uncertainty. In the following plots, the results with =0.0001 are shown unless otherwise stated. As an example, the dependence of are shown in Fig. 2. The statistical errors are invisible in this scale.

We fit the data of and to polynomial functions of . The fits are made over three fit ranges, three different polynomial orders, two values of . Since not all the fits are successful, we only keep the fit results satisfying in the following analysis. Furthermore, it turns out that five or six parameters are enough to fit the data well, thus we omit the fits with the sixth order polynomial. Once the fit parameters are determined, it is straightforward to calculate the curvature of the potential.

The results for the curvature are plotted in Figs. 3 for (solid curves) and (dashed curves), where the results for and 0.4 are shown as examples. The difference in the curvature between and turns out to be reasonably small at all .

The curvature of the first term of eq. (2) can be calculated using the averaged value and the susceptibility of at each [13], which is shown in Fig. 3 with the dotted curves obtained from the fit to a fifth order of polynomial. It is seen that, independently of , is always positive as expected. The figure shows that and have a peak at slightly below the value at which takes the minimum, which indicates that, in many flavor system, the phase transition or rapid crossover occurs at smaller than the two flavor case.

To determine , we iterated the calculation with varying in steps of 0.02. It turns out that the resulting depends on the details of the fitting procedures, though not by much. We adopt all those results as long as , and the spread is taken as the systematic uncertainty. Figure 4 shows the light quark mass dependence of , where the error bars represent the systematic uncertainties.

Since no clear tendency is observed in the light quark mass dependence, we fit the data to a constant (solid lines), yielding in the chiral limit of the two flavor mass. Nonzero value of in the chiral limit excludes the possibility of the first order transition of massless two flavor QCD.

## 4 Summary and discussion

We have studied the finite temperature phase transition of QCD with two light and many heavy quarks at zero chemical potential to identify the critical line separating the continuous crossover and the first order regions on the - plane. In other words, two flavor QCD with a finite mass is enforced to undergo a first order transition by adding extra quarks. We then try to see whether those extra quarks are necessary to keep the first order transition down to the chiral limit of two flavor mass.

The nature of the transition is identified by the shape of the constraint effective potential for the generalized plaquette. It is that stays constant within the systematic uncertainty in the range of two flavor mass we have studied (), which suggests that the critical heavy mass remains finite in the chiral limit of the two flavors and hence the phase transition of massless two flavor QCD is of second order.

Our approach is applicable for any kinds of light quark action. can be considered to be arbitrary small by assuming arbitrary large . Therefore, the above statement is valid independently of the convergence of the hopping parameter expansion.

## 5 Acknowledgments

We would like to thank members of WHOT-QCD Collaboration for discussions. We also thank Ken-Ichi Ishikawa for providing us his simulation codes. This work is in part supported by JSPS KAKENHI Grant-inAid for Scientific Research (B) (No. 15H03669 [NY], 23540295 [SE] ) and (C) (No. 26400244 [SE]), and by the Large Scale Simulation Program of High Energy Accelerator Research Organization (KEK) No. 14/15-23.

## References

- [1] R. D. Pisarski and F. Wilczek, Phys. Rev. D 29, 338 (1984).
- [2] C. Bonati, G. Cossu, M. D’Elia, A. Di Giacomo and C. Pica, PoS LATTICE 2008, 204 (2008) [arXiv:0901.3231 [hep-lat]].
- [3] S. Aoki, H. Fukaya and Y. Taniguchi, Phys. Rev. D 86 (2012) 114512 [arXiv:1209.2061 [hep-lat]].
- [4] G. Cossu, S. Aoki, H. Fukaya, S. Hashimoto, T. Kaneko, H. Matsufuru and J. I. Noaki, Phys. Rev. D 87, no. 11, 114514 (2013) [Phys. Rev. D 88, no. 1, 019901 (2013)] [arXiv:1304.6145 [hep-lat]].
- [5] A. Pelissetto and E. Vicari, Phys. Rev. D 88, no. 10, 105018 (2013) [arXiv:1309.5446 [hep-lat]].
- [6] Y. Nakayama and T. Ohtsuki, Phys. Rev. D 91, no. 2, 021901 (2015) [arXiv:1407.6195 [hep-th]].
- [7] C. Bonati, P. de Forcrand, M. D’Elia, O. Philipsen and F. Sanfilippo, Phys. Rev. D 90, no. 7, 074030 (2014) [arXiv:1408.5086 [hep-lat]];
- [8] T. Sato and N. Yamada, Phys. Rev. D 91, no. 3, 034025 (2015) [arXiv:1412.8026 [hep-lat]].
- [9] R. Kitano and N. Yamada, JHEP 1510, 136 (2015) [arXiv:1506.00370 [hep-ph]].
- [10] E. Berkowitz, M. I. Buchoff and E. Rinaldi, Phys. Rev. D 92, no. 3, 034507 (2015) [arXiv:1505.07455 [hep-ph]].
- [11] S. Borsanyi et al., Phys. Lett. B 752, 175 (2016) [arXiv:1508.06917 [hep-lat]].
- [12] S. Ejiri and N. Yamada, Phys. Rev. Lett. 110, no. 17, 172001 (2013) [arXiv:1212.5899 [hep-lat]].
- [13] S. Ejiri, R. Iwami and N. Yamada, arXiv:1511.06126 [hep-lat].
- [14] S. Ejiri, Phys. Rev. D 77, 014508 (2008).
- [15] S. Ejiri, Phys. Rev. D 78, 074507 (2008) [arXiv:0804.3227 [hep-lat]].
- [16] H. Saito, et al. Y. Maezawa, H. Ohno, and T. Umeda (WHOT-QCD Collaboration), Phys. Rev. D 84, 054502 (2011).
- [17] S. Ejiri, Y. Maezawa, N. Ukita, S. Aoki, T. Hatsuda, N. Ishii, K. Kanaya, and T. Umeda (WHOT-QCD Collaboration), Phys. Rev. D 82, 014508 (2010).
- [18] A. Ukawa, UTHEP-302, C93-06-21.1.
- [19] S. Ejiri, PoS (LATTICE 2008) 002 (2008).
- [20] S. Ejiri, Eur. Phys. J. A 49, 86 (2013) [arXiv:1306.0295 [hep-lat]].
- [21] R. Iwami, S. Ejiri and N. Yamada, these proceedings.