††thanks<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
#  New interpretation of chiral phase transition: Violation of trilemma in QCD

Chuan-Xin Cui  Center for Theoretical Physics and College of Physics, Jilin
University, Changchun, 130012, China  Jin-Yang Li  Center for Theoretical
Physics and College of Physics, Jilin University, Changchun, 130012, China
Shinya Matsuzaki  Center for Theoretical Physics and College of Physics, Jilin
University, Changchun, 130012, China  Mamiya Kawaguchi  School of Nuclear
Science and Technology, University of Chinese Academy of Sciences, Beijing
100049, China  Akio Tomiya  RIKEN BNL Research center, Brookhaven National
Laboratory, Upton, NY, 11973, USA  Department of Information Technology,
International Professional University of Technology in Osaka, 3-3-1 Umeda,
Kita-Ku, Osaka, 530-0001, Japan

###### Abstract

We find that the chiral phase transition (chiral crossover) in QCD at physical
point is triggered by big imbalance among three fundamental quantities
essential for the QCD vacuum structure: susceptibility functions for the
chiral symmetry, axial symmetry, and the topological charge. The balance,
dobbed the QCD trilemma, is unavoidably violated when one of the magnitudes
among them is highly dominated, or suppressed. Based on a three-flavor Nambu-
Jona-Lasinio model, we explicitly evaluate the amount of violation of the QCD
trilemma at physical point, and show that the violation takes place not only
at vacuum, but even in a whole temperature regime including the chiral
crossover epoch. This work confirms and extends the suggestion recently
reported from lattice QCD with 2 flavors on dominance of the axial and
topological susceptibilities left in the chiral susceptibility at high
temperatures. It turns out that the imbalance is essentially due to the flavor
symmetry violation of the lightest three flavors, and the flavor breaking
specifically brings enhancement of the axial anomaly contribution in the
chiral order parameter, while the the strength of the axial breaking and the
transition rate of the topological charge are fairly insensitive to the flavor
symmetry. The violation of QCD trilemma and its flavor dependence can be
tested by lattice simulations with 2 + 1 flavors in the future, and would also
give a new guiding principle to explore the flavor dependence of the chiral
phase transition, such as the Columbia plot, including possible extension with
external fields.

##  I Introduction

The chiral phase transition is of importance to comprehend the QCD vacuum, and
is also essential to figure out the origin of mass in a view of thermal
history of the universe. Plenty of studies on the chiral phase transition have
extensively been worked out so far through the nonperturbative analysis in
lattice simulations, and also chiral effective models of QCD. However, as
argued in the literature [1, 2], still, it is not well understood whether the
chiral symmetry breaking is the most dominant source of the origin of mass,
even in presence of contamination with the $U(1)_{A}$ anomaly, and by what
mechanism it is restored faster than the $U(1)_{A}$ symmetry at high
temperature.

The order parameter of the chiral symmetry is given by the quark condensate,
which can alternatively be signaled by difference of meson correlation
functions for the chiral partners: the latter is referred to as an indicator
of the chiral breaking strength. Though being so simple and well-defined, the
chiral order parameter at physical point is actually involved due to finite
quark masses, which explicitly break the chiral symmetry. Indeed, the chiral
symmetry is restored at high temperature only in part, refereed to as the
chiral crossover [3, 4]. Actually it gets more intricate because the chiral
order parameter (the indicator of the chiral breaking strength) couples with
the indicator of the axial breaking strength and topological features of the
QCD vacuum via finite quark masses. The latter tagging is captured by a robust
relation between the indicators for the chiral $SU(2)_{L}\times SU(2)_{R}$
symmetry and $U(1)_{A}$ axial symmetry, which is constructed from a set of
generic anomalous Ward identities for the three-flavor chiral $SU(3)_{L}\times
SU(3)_{R}$ symmetry [5, 6] (for more details, see also the next section):

$\displaystyle\chi_{\eta-\delta}=\chi_{\pi-\delta}+\frac{4}{m_{l}^{2}}\chi_{\rm
top}\,,$ (1)

where $m_{l}=m_{u}=m_{d}$ is the isospin-symmetric mass for the lightest up
and down quarks; $\chi_{\eta-\delta}\equiv\chi_{\eta}-\chi_{\delta}$ and
$\chi_{\pi-\delta}\equiv\chi_{\pi}-\chi_{\delta}$ are differences of meson
susceptibilities related to the partners for the chiral symmetry
($\chi_{\eta}$ and $\chi_{\delta}$) — an indicator for the strength of the
chiral $SU(2)$ symmetry breaking — and axial symmetry ($\chi_{\delta}$ and
$\chi_{\pi}$)

— an indicator for the strength of the $U(1)_{A}$ axial breaking; $\chi_{\rm
top}$ is the topological susceptibility related to the transition rate of the
topological charge carried by the QCD $\theta$ vaccua. By the chiral $SU(2)$
and axial rotations, the meson susceptibilities exchange their partners:
$\chi_{\eta}\leftrightarrow\chi_{\delta}$ (chiral) and
$\chi_{\pi}\leftrightarrow\chi_{\delta}$ (axial), hence $\chi_{\eta-\delta}=0$
and $\chi_{\pi-\delta}=0$ are signals of restorations of the associated
symmetries. ($\chi_{\rm top}<0$ and other susceptibilities are positive in our
sign convention. See also the next section.) Thus Eq.(1) dictates coherence of
the chiral $SU(2)$ symmetry breaking and $U(1)_{A}$ breaking, linked with the
transition rate of the topological charge, where all the breaking is
controlled by nonzero quark masses. This anomalous Ward identity takes the
same form even in the decoupling limit of strange quark, i.e., in the lightest
two-flavor limit.

Equation (1) plays the essential role to comprehend how the effective
restoration of the chiral symmetry is correlated with that of the axial
symmetry and the temperature dependence of topological susceptibility in real-
life QCD. This gives a new guideline in a sense of exploring the chiral phase
transition constrained by Eq.(1), and would provide crucial clues to answer
the questions posed above.

The lattice QCD simulations with 2 + 1 flavors at the physical point have
revealed a faster drop of $\chi_{\eta-\delta}$ , than $\chi_{\pi-\delta}$
around and above the pseudo-critical temperature of the chiral crossover [4].
In the case of 2 flavors at the chiral limit, the effective restoration of the
chiral and axial symmetry has also been discussed through the meson
susceptibilities [7, 8, 9]

However, those are based on independent measurements of two terms,
$\chi_{\eta-\delta}$ and $\chi_{\pi-\delta}$, with the constraint of Eq.(1)
disregarded.

Measurements of $\chi_{\rm top}$ on lattice QCD with 2 + 1 flavors at around
physical point and its temperature dependence have been reported [10, 11, 12],
in light of detecting the effective restoration of the $U(1)_{A}$ symmetry
assuming the much faster restoration of the chiral $SU(2)$ symmetry. However,
those are also individual observations, basically separated from measurements
of the chiral and axial indicators. Therefore, it is yet uncovered how the
temperature dependence of $\chi_{\rm top}$ would correlate with the other two,
with reflecting the constraint of Eq.(1).

In a view of the coherence in Eq.(1), a recent lattice study with two lightest
flavors has for the first time shown significant contributions from the axial
and topological susceptibilities ($\chi_{\pi-\delta}$ and $\chi_{\rm top}$)
left in the chiral susceptibility ($\chi_{\eta-\delta}$) in the chiral
crossover domain [2]. This would imply that the faster chiral crossover is
triggered by a sizable cancellation between axial and topological
susceptibilities, the two terms in the right-hand side of Eq.(1).

To quantify the magnitude of such a cancellation, we may define an ideal case
with no preference among three susceptibilities in magnitude in Eq.(1), so
that the Ward identity acts like a balance equation. We dob this ideal
situation as “QCD trilemma”, and depict a triangle cartoon in Fig. 1. The
degree of formation of QCD trilemma can be evaluated via the following
quantity:

$\displaystyle R$ $\displaystyle\equiv\frac{\frac{4}{m_{l}^{2}}\chi_{\rm
top}+\chi_{\pi-\delta}}{\chi_{\eta-\delta}-\frac{4}{m_{l}^{2}}\chi_{\rm
top}}=1+\frac{\frac{4}{m_{l}^{2}}\chi_{\rm top}}{\chi_{\pi-\delta}}\,.$ (2)

By using this $R$ the Ward identity in Eq.(1) is rewritten as

$\displaystyle\chi_{\eta-\delta}$ $\displaystyle=R\cdot\chi_{\pi-\delta}\,,$
$\displaystyle{\rm or}\qquad-\frac{4}{m_{l}^{2}}\chi_{\rm top}$
$\displaystyle=(1-R)\cdot\chi_{\pi-\delta}\,,$ (3)

so that $R$ measures the size of gap in magnitude between the chiral
($\chi_{\eta-\delta}$) and axial ($\chi_{\pi-\delta}$) susceptibilities , or
the topological $(\chi_{\rm top})$ and axial ($\chi_{\pi-\delta}$)
susceptibilities. Equation (3) tells us that
$\frac{\chi_{\eta-\delta}}{\chi_{\pi-\delta}}=R$ and
$\frac{\left(-\frac{4}{m_{l}^{2}}\chi_{\rm
top}\right)}{\chi_{\pi-\delta}}=1-R$, so one may then quantify the amount of
balance to keep the QCD trilemma, by saying that the three susceptibilities
are balanced when

$\displaystyle{\rm balanced}:\qquad 0.1<R<0.9\,,$ (4)

otherwise imbalanced along with a big gap by more than one order of magnitude
between two of three susceptibilities. An ideal and optimized trilemma is thus
realized when $R=0.5$. We shall dub this $R$ as the trilemma estimator. $R$
becomes $\ll 1$ when axial and topological susceptibilities get close each
other with different sign in Eq.(2), which would signal the faster effective
restoration of the chiral symmetry than that of the axial one through Eq.(3).

Figure 1:  Illustration of QCD trilemma and its violation. The QCD vacuum
structure is built upon the “Chiral SU(2)”, “U(1) Axial”, and “Topological”
features, which are related each other by a balance relation in Eq.(1)), where
the “Chiral SU(2)”, “U(1) Axial” and “Topological” are monitored by
$\chi_{\eta-\delta}$, $\chi_{\pi-\delta}$, and $(-4/m_{l}^{2})\cdot\chi_{\rm
top}$, respectively. Left panel: the QCD vacuum is “balanced” and holds the
trilemma by forming the equilateral triangle with the same order of the weight
amplitudes denoted by blobs. Right panel: the trilemma is violated
(imbalanced) when a big cancellation between “U(1) Axial”, and “Topological”
takes place in Eq.(1), which is represented by the isosceles triangle with one
blob significantly reduced, keeping Eq.(1) and the corresponding two sides
stretched out. As it will turn out in the text, real-life QCD is “imbalanced”.

The aforementioned evidence observed by lattice simulations [4] on the faster
drop of $\chi_{\eta-\delta}$ than $\chi_{\pi-\delta}$ indicates $R\ll 1$ in a
view of Eq.(3). The result from the recent lattice study with two lightest
flavors in [2] can be rephrased as $R\ll 1$ in both Eqs.(2) and (3). Though
not explicitly addressed and restricted only around the crossover regime, this
imbalance could also be read off from the existing lattice QCD data with 2 + 1
flavors in [4] and also [2] with taking into account possible finite volume
effects and statistical errors.

Thus, violation of the QCD trilemma has not been yet explicitly explored at
the physical point for 2 + 1 flavors on the same lattice setting, and it is
still unclear how axial and topological susceptibilities, holding single
Eq.(1) with the chiral one, develop in a whole finite temperature regime and
contributes to achieving the chiral crossover. Even in the context of
effective chiral models, no such discussion along with Eq.(1) has so far been
made together with proper incorporation of the flavor-singlet condition for
$\chi_{\rm top}$ [13, 14, 6] (to the latter point, see also the next section).
Real-life QCD having 2 + 1 flavors at physical point might be imbalanced in
realizing the chiral crossover, through undergoing a big cancellation between
axial and topological susceptibilities in a whole temperature regime.

In this paper, we discuss the violation of QCD trilemma in real-life QCD based
on a Nambu-Jona-Lasinio (NJL) model, and give a qualitative interpretation of
the mechanism of the violation, namely the coherence among the chiral, axial,
and topological susceptibilities, constrained by Eq.(1). Prior to the lattice
simulations, we show that real-life QCD indeed yields $R\ll 1$, i.e., exhibits
the violation of QCD trilemma, in a whole temperature regime including the
chiral crossover regime. We find that the violation of QCD trilemma, and the
related dominance of $\chi_{\pi-\delta}$ and $\chi_{\rm top}$ in the chiral
order parameter at the crossover regime are due to the three-flavor symmetry
violation.

Our findings are shortly testable by lattice simulations, and would help
deeper understanding of the flavor dependence of the chiral phase transition,
mapped onto the so-called Columbia plot [15]. Exploring the chiral (crossover)
phase transition along with the violation of QCD trilemma would lead to clues
toward answering the posed questions: the expected dominance of the chiral
symmetry breaking in the origin of mass, and deeper understanding of the
observed faster (effective) restoration of the chiral symmetry in the presence
of contamination with the $U(1)_{A}$ anomaly.

This paper is organized as follows. In Sec. II, we introduce the preliminaries
relevant to the discussion in the later sections, which include definitions
and generic formulas for susceptibilities, as well as a concise derivation of
the anomalous chiral Ward identity. In Sec. III, the NJL model that we work
throughout this paper is introduced, together with showing qualitative
consistency of the model predictions with the lattice data, which includes the
temperature dependence of the quark condensate, meson susceptibilities, and
topological susceptibilities. In Sec. IV we discuss the QCD trilemma estimator
$R$, in a whole temperature region, including the chiral crossover regime, and
show the violation of the trilemma, imbalance of the real-life QCD vacuum. We
then demonstrate that the violation is due to the three-flavor symmetry. Sec.
V devotes to our conclusion, where several possible applications of the notion
of QCD trilemma are also briefly addressed.

##  II Central formulas: topological susceptibility and anomalous chiral Ward-
identities in QCD

In this section we begin by reviewing the generic expression for the
topological susceptibility $\chi_{\rm top}$ [6] with the flavor singlet
condition properly reflected [13, 14], and introduce the related anomalous
chiral Ward identities in QCD involving pseudoscalar susceptibilities
$\chi_{\pi},\chi_{\eta}$ and $\chi_{\delta}$.

###  II.1 Topological susceptibility: flavor singlet nature

The topological susceptibility $\chi_{\rm top}$ is related to the $\theta$
vacuum configuration of QCD. It is defined as the curvature of the
$\theta$-dependent vacuum energy $V(\theta)$ in QCD at $\theta=0$:

$\chi_{\rm
top}=-\int_{T}d^{4}x\frac{\delta^{2}V(\theta)}{\delta\theta(x)\delta\theta(0)}\Bigg{|}_{\theta=0}\,,$
(5)

where the temperature integral $\int_{T}d^{4}x$ is defined as
$\int_{0}^{1/T}d\tau\int d^{3}x$ with the imaginary time $\tau=ix_{0}$, and
$V(\theta)$ denotes the potential of QCD, which is read off from the
generating functional of QCD (in Euclidean space):

$\displaystyle Z_{\rm QCD}$ $\displaystyle=$
$\displaystyle\int[\Pi_{f}dq_{f}d\bar{q}_{f}][dA]\exp\Biggl{[}-\int_{T}d^{4}x\Biggl{\\{}\sum_{f=u,d,s}\Bigg{(}\bar{q}^{f}_{L}i\gamma^{\mu}D_{\mu}q^{f}_{L}+\bar{q}^{f}_{R}i\gamma^{\mu}D_{\mu}q^{f}_{R}$
(6)
$\displaystyle+\bar{q}^{f}_{L}m_{f}q^{f}_{R}+\bar{q}^{f}_{R}m_{f}q^{f}_{L}\Bigg{)}+\frac{1}{4g^{2}}(F_{\mu\nu}^{a})^{2}+\frac{i\theta}{32\pi^{2}}F_{\mu\nu}^{a}\tilde{F}_{\mu\nu}^{a}\Biggl{\\}}\Biggl{]}.$

Here $q^{f}_{L(R)}$ denote the left- (right-) handed quark fields; the
covariant derivative of the quark field is represented as $D_{\mu}$ involving
the gluon fields $A$; $F_{\mu\nu}^{a}$ is the field strength of the gluon
fields with $g$ being the QCD coupling constant; for simplicity, the quark
masses are taken to be real and positive with the electroweak-induced CP
violation disregarded.

The form of the $\theta$ dependence on the QCD generating functional is
ambiguous, because the $\theta$ parameter can always be shifted by the
$U(1)_{A}$ rotation through the induced $U(1)_{A}$ anomaly as well as the
phase shift in the quark mass term. Thus the QCD-$\theta$ vacuum is shifted by
the $U(1)_{A}$ phase as well. Taking into account this shift, the true QCD
vacuum is determined so as to set the net $\theta$ to zero, i.e. the CP
invariant vacuum, as shown in the context of the Vafa-Witten’s theorem [16].
The topological susceptibility $\chi_{\rm top}$ should then be evaluated at
the true vacuum with the net $\theta=0$. This is how to properly compute
$\chi_{\rm top}$ [14].

Under the $U_{A}(1)$ rotation with the rotation angle $\theta_{f}$, the left-
and right-handed quark fields are transformed as

$\displaystyle q^{f}_{L}$ $\displaystyle\to$
$\displaystyle\exp\left(-i\theta_{f}/2\right)q^{f}_{L},$ $\displaystyle
q^{f}_{R}$ $\displaystyle\to$
$\displaystyle\exp\left(i\theta_{f}/2\right)q^{f}_{R}\,.$ (7)

We then find that the extra phase factor shows up in the QCD generating
functional written in terms of the transformed chiral quark fields:

$\displaystyle\int[\Pi_{f}dq_{f}d\bar{q}_{f}][dA]\exp\Biggl{[}-\int_{T}d^{4}x\Biggl{\\{}\sum_{f=u,d,s}\Bigg{(}\bar{q}^{f}_{L}i\gamma^{\mu}D_{\mu}q^{f}_{L}+\bar{q}^{f}_{R}i\gamma^{\mu}D_{\mu}q^{f}_{R}$
$\displaystyle+\bar{q}^{f}_{L}m_{f}e^{i\theta_{f}}q^{f}_{R}+\bar{q}^{f}_{R}m_{f}e^{-i\theta_{f}}q^{f}_{L}\Bigg{)}+\frac{1}{4g^{2}}(F_{\mu\nu}^{a})^{2}+\frac{i(\theta-\bar{\theta})}{32\pi^{2}}F_{\mu\nu}^{a}\tilde{F}_{\mu\nu}^{a}\Biggl{\\}}\Biggl{]},$
(8)

where
$\bar{\theta}=\sum_{f=u,d,s}\theta_{f}=\theta_{u}+\theta_{d}+\theta_{s}$.

We can take a phase convention in such a way that the $\theta$ dependence goes
away from the topological gluonic term $(F\tilde{F})$:

$\displaystyle\theta=\bar{\theta}=\theta_{u}+\theta_{d}+\theta_{s}\,.$ (9)

Instead, the quark mass term fully carries the $\theta$-dependence, which
makes manifest presence of the strong CP violation. Here, the CP violating
phases $\theta_{u,d,s}$ actually involves redundancy, i.e., being not fully
independent, because the QCD interaction and the QCD vacuum characterized by
the three-flavor symmetric quark condensate are flavor blind (singlet).
Therefore, we must recover the flavor independence for those CP violating
phases. Supposing a small enough $\theta_{f}$ (which is to be consistent with
the observation on $\theta<10^{-10}$ from the electric dipole moment of
neutron), and expanding the quark mass terms in powers of $\theta_{f}$s, we
find that the CP violating coupling terms at the nontrivial leading order
arise with the factor of $m_{f}\theta_{f}$. Thus the desired flavor singlet
condition goes like [13],

$\displaystyle m_{u}\theta_{u}=m_{d}\theta_{d}=m_{s}\theta_{s}\equiv x.$ (10)

By using this flavor singlet condition together with Eq. (9), $\theta_{f}$s
are determined to be proportional to $\theta$ as

$\displaystyle\theta_{u}=\frac{\bar{m}}{m_{u}}\theta,\;\;\;\theta_{d}=\frac{\bar{m}}{m_{d}}\theta,\;\;\;\theta_{s}=\frac{\bar{m}}{m_{s}}\theta\,,$
(11)

where

$\displaystyle\bar{m}=\left(\frac{1}{m_{u}}+\frac{1}{m_{d}}+\frac{1}{m_{s}}\right)^{-1}.$
(12)

Thus the $\theta$-dependent vacuum energy of QCD with the flavor singlet
nature properly reflected is:

$\displaystyle V_{\rm QCD}(\theta)$ $\displaystyle=$
$\displaystyle-\ln\left[\int[\Pi_{f}dq_{f}d\bar{q}_{f}][dA]\exp\left(-\int_{T}d^{4}x{\cal
L}_{\rm QCD}^{(\theta)}\right)\right]\,,$ (13)

where

$\displaystyle{\cal L}_{\rm
QCD}^{(\theta)}=\sum_{f}\left(\bar{q}^{f}_{L}i\gamma^{\mu}D_{\mu}q^{f}_{L}+\bar{q}^{f}_{R}i\gamma^{\mu}D_{\mu}q^{f}_{R}\right)+\bar{q}_{L}{\cal
M}_{\theta}q_{R}+\bar{q}_{R}{\cal
M}_{\theta}^{\dagger}q_{L}+\frac{1}{4g^{2}}(F_{\mu\nu}^{a})^{2}\,,$ (14)

with ${\cal M}_{\theta}$ being the $\theta$-dependent quark matrix,

$\displaystyle{\cal M}_{\theta}={\rm
diag}\left[m_{u}\exp\left(i\frac{\bar{m}}{m_{u}}\theta\right),m_{d}\exp\left(i\frac{\bar{m}}{m_{d}}\theta\right),m_{s}\exp\left(i\frac{\bar{m}}{m_{s}}\theta\right)\right].$
(15)

Based on Eq.(13), we evaluate $\chi_{\rm top}$ in Eq.(5), and find [6]

$\displaystyle\chi_{\rm top}$
$\displaystyle=\bar{m}^{2}\left[\frac{\langle\bar{u}u\rangle}{m_{l}}+\frac{\langle\bar{d}d\rangle}{m_{l}}+\frac{\langle\bar{s}s\rangle}{m_{s}}+\chi_{P}^{uu}+\chi_{P}^{dd}+\chi_{P}^{ss}+2\chi_{P}^{ud}+2\chi_{P}^{us}+2\chi_{P}^{ds}\right]$
$\displaystyle=\frac{1}{4}\left[m_{l}\left(\langle\bar{u}u\rangle+\langle\bar{d}d\rangle\right)+m_{l}^{2}\left(\chi_{P}^{uu}+\chi_{P}^{dd}+2\chi_{P}^{ud}\right)\right]=m_{s}\langle\bar{s}s\rangle+m_{s}^{2}\chi_{P}^{ss}\,,$
(16)

where the pseudoscalar susceptibilities $\chi_{P}^{uu,dd,ud}$, $\chi_{P}^{ss}$
and $\chi_{P}^{us,ds}$ are defined as

$\displaystyle\chi_{P}^{f_{1}f_{2}}$
$\displaystyle=\int_{T}d^{4}x\langle(\bar{q}_{f_{1}}(0)i\gamma_{5}q_{f_{1}}(0))(\bar{q}_{f_{2}}(x)i\gamma_{5}q_{f_{2}}(x))\rangle\,,\qquad{\rm
for}\quad{q_{{}_{f_{1,2}}}=u,d,s}\,.$ (17)

In Eq.(16) we have taken the isospin symmetric limit $m_{u}=m_{d}\equiv
m_{l}$. The signs of the quark masses and condensates are chosen to be
positive and negative, respectively, such that $\chi_{\rm top}<0$. Note that
$\chi_{\rm top}\to 0$, when either of quarks becomes massless ($m_{l}$ or
$m_{s}$ $\to 0$), reflecting the flavor-singlet nature of the QCD vacuum [13,
14].

###  II.2 Anomalous chiral Ward identities

The anomalous Ward identities regarding the chiral $SU(3)_{L}\times SU(3)_{R}$
symmetry are directly read off from chiral variations of the QCD potential in
Eq.(13). The central formula then takes the form

$\displaystyle\langle\delta_{a}{\cal
O}_{b}(0)\rangle=-\int_{T}d^{4}x\langle{\cal
O}_{b}(0)\cdot\bar{q}_{f}(x)i\gamma_{5}\\{T_{a},{M}\\}q_{f}(x)\rangle\,,$ (18)

where $T_{a}=\lambda_{a}/2$ ($a=1,\cdots,8$) are generators of $SU(3)$;
$\delta_{a}$ stands for the infinitesimal variation of the chiral $SU(3)$
transformation associated with the generator $T_{a}$, under which $q_{f}$
transforms as $\delta_{a}q_{f}=i\gamma_{5}T_{a}q_{f}$; ${\cal O}_{b}(0)$
($b=0,\cdots,8$) is an arbitrary operator. In particular, for the pseudoscalar
operators ${\cal O}_{b}=\bar{q}_{f}i\gamma_{5}T_{b}q_{f}$, choosing
$a=1,2,3,8$ and $b=0,8$ with $T_{0}=1/\sqrt{6}\cdot{\bf 1}_{3\times 3}$, we
get [5, 6]

$\displaystyle\langle\bar{u}u\rangle+\langle\bar{d}d\rangle$
$\displaystyle=-m_{l}\chi_{\pi}\,,$
$\displaystyle\langle\bar{u}u\rangle+\langle\bar{d}d\rangle+4\langle\bar{s}s\rangle$
$\displaystyle=-\left[m_{l}\left(\chi_{P}^{uu}+\chi_{P}^{dd}+2\chi_{P}^{ud}\right)-2(m_{s}+m_{l})\left(\chi_{P}^{us}+\chi_{P}^{ds}\right)+4m_{s}\chi_{P}^{ss}\right]\,,$
$\displaystyle\langle\bar{u}u\rangle+\langle\bar{d}d\rangle-2\langle\bar{s}s\rangle$
$\displaystyle=-\left[m_{l}\left(\chi_{P}^{uu}+\chi_{P}^{dd}+2\chi_{P}^{ud}\right)+(m_{l}-2m_{s})\left(\chi_{P}^{us}+\chi_{P}^{ds}\right)-2m_{s}\chi_{P}^{ss}\right]\,,$
(19)

where $\chi_{\pi}$ denotes the pion susceptibility defined as

$\displaystyle\chi_{\pi}$
$\displaystyle=\int_{T}d^{4}x\left[\langle(\bar{u}(0)i\gamma_{5}u(0))(\bar{u}(x)i\gamma_{5}u(x))\rangle_{\rm
conn}+\langle(\bar{d}(0)i\gamma_{5}d(0))(\bar{d}(x)i\gamma_{5}d(x))\rangle_{\rm
conn}\right]\,,$ (20)

with $\langle\cdot\cdot\cdot\rangle_{\rm conn}$ being the connected part of
the correlation function. The form of the anomalous chiral Ward-identities in
Eq.(19) will be intact, as long as only the quark masses gives the leading
order of explicit chiral breaking effects, as evident in the chiral variation
of Eq.(18).

Combining Ward identities in Eq.(19), we find

$\displaystyle\chi_{\rm top}$
$\displaystyle=\frac{1}{2}m_{l}m_{s}\left(\chi_{P}^{us}+\chi_{P}^{ds}\right)=\frac{1}{4}m_{l}^{2}(\chi_{\eta}-\chi_{\pi})\,,$
(21)

where $\chi_{\eta}$ is the eta meson susceptibility, defined as

$\displaystyle\chi_{\eta}$
$\displaystyle=\int_{T}d^{4}x\Bigg{[}\langle(\bar{u}(0)i\gamma_{5}u(0))(\bar{u}(x)i\gamma_{5}u(x))\rangle+\langle(\bar{d}(0)i\gamma_{5}d(0))(\bar{d}(x)i\gamma_{5}d(x))\rangle$
$\displaystyle+2\langle(\bar{u}(0)i\gamma_{5}u(0))(\bar{d}(x)i\gamma_{5}d(x))\rangle\Bigg{]}$
$\displaystyle=\chi_{P}^{uu}+\chi_{P}^{dd}+2\chi_{P}^{ud}\,.$ (22)

The last line of Eq.(21) can be written as

$\displaystyle(\chi_{\eta}-\chi_{\delta})=(\chi_{\pi}-\chi_{\delta})+\frac{4}{m_{l}^{2}}\chi_{\rm
top}\,,$ (23)

where $\chi_{\delta}$ is the susceptibility for the delta meson channel
($a_{0}$ meson in terms of the Particle Data Group identification), defined in
the same way as $\chi_{\pi}$ in Eq.(17) with the factors of $(i\gamma_{5})$
replaced with identity $1$.
$\chi_{\eta-\delta}\equiv\chi_{\eta}-\chi_{\delta}$ and
$\chi_{\pi-\delta}\equiv\chi_{\pi}-\chi_{\delta}$ play the roles of the
indicators to detect the strength of the chiral and axial breaking, which
signal the restorations when those (asymptotically) reach zero.

Equation (23) is our central formula, which has been inferred in Eq.(1), and
will be explored in details in the later section.

##  III A chiral effective model: NJL

In this section we introduce an NJL model that we work on, and give a couple
of preliminaries for discussion on the estimate of the QCD trilemma estimator
$R$ in Eq.(2), with showing consistency of the NJL estimates with currently
available lattice data on 2 + 1 flavors at physical point. Since the
methodology to compute observables and thermodynamic quantities in the model
is standard and fully described in a review [17], we will skip all the
details, and just present the final formulas directly used to the numerical
evaluation of the quark condensates and susceptibilities.

The three-flavor NJL model Lagrangian that we work on is constructed as
follows:

$\displaystyle\mathcal{L}$
$\displaystyle=\bar{q}(i\gamma_{\mu}\partial^{\mu}-{\bf
m})q+\mathcal{L}_{4f}+\mathcal{L}_{\rm KMT}\,,$
$\displaystyle\mathcal{L}_{4f}$
$\displaystyle=\frac{g_{s}}{2}\sum^{8}_{a=0}[(\bar{q}\lambda^{a}q)^{2}+(\bar{q}i\gamma_{5}\lambda^{a}q)^{2}]\,,$
$\displaystyle\mathcal{L}_{\rm KMT}$ $\displaystyle=g_{D}[\mathop{\rm
det}\limits_{i,j}\bar{q_{i}}(1+\gamma_{5})q_{j}+{\rm h.c.}]\,,$ (24)

where $q$ is the $SU(3)$ triplet-quark field, $q=(u,d,s)^{T}$. The current
quark masses are embedded in the mass matrix $\bf m$ of the form ${\bf m}={\rm
diag}(m_{u},m_{d},m_{s})$.

The four-fermion interaction term $\mathcal{L}_{4f}$ is invariant under the
chiral $U(3)_{L}\times U(3)_{R}$ transformation: $q\to U\cdot q$ with
$U=\exp[-i\gamma_{5}\sum_{a=0}^{8}(\lambda^{a}/2)\theta^{a}]$ and the chiral
phases $\theta^{a}$. The mass term in ${\cal L}$ explicitly breaks
$U(3)_{L}\times U(3)_{R}$ symmetry. The determinant term $\mathcal{L}_{\rm
KMT}$, called the Kobayashi- Maskawa-‘t Hooft [18, 19, 20, 21] term, induced
from the QCD instanton configuration, preserves $SU(3)_{L}\times SU(3)_{R}$
invariance (associated with the chiral phases labeled as $a=1,\cdots,8$) but
breaks the $U(1)_{A}$ (corresponding to $a=0$) symmetry, measured by the
effective coupling constant $g_{D}$.

The $U(1)_{A}$ symmetry is anomalous due to not only the quark mass terms, but
also the KMT term reflecting the underlying gluonic anomaly. Thus we have the
anomalous conservation law:

$\displaystyle\partial^{\mu}j_{\mu}^{a=0}$ $\displaystyle=$ $\displaystyle
2i\bar{q}{\bm{m}}\gamma_{5}q-12g_{D}{\rm Im}\left[{\rm
det}\bar{q}_{i}(1-\gamma_{5})q_{j}\right]\,.$ (25)

The matching with the underlying QCD leads to an operator relation:
$Q=-4g_{D}{\rm Im}\left[{\rm det}\bar{q}_{i}(1-\gamma_{5})q_{j}\right]$, where
$Q=g^{2}/(32\pi^{2})\,F_{\mu\nu}^{a}\tilde{F}^{a\mu\nu}$ is the topological
charge. Using this operator relation together with Eq.(5), one could evaluate
$\chi_{\rm top}$ as in the literature, e.g., Refs. [22, 23, 24, 25] in the
framework of the NJL model with the mean field approximation. However, this
procedure makes the flavor-singlet nature of the vacuum nontransparent. The
identification of $\chi_{\rm top}$ defined within the NJL model with the mean
field approximation with that derived directly in QCD requires a careful
separate investigation ensuring the flavor singlet nature of the vacuum in
both theories. Instead, in the present paper the NJL-model is considered as a
reduction of the full QCD allowing the evaluation of the right-hand side of
Eq.(16) in which the flavor singlet nature is manifestly built-in.

The NJL model itself is a (perturbatively) nonrenormalizable field theory
because ${\cal L}_{4f}$ and ${\cal L}_{\rm KMT}$ describe the higher
dimensional interactions with mass dimension greater than four. Therefore, a
momentum cutoff $\Lambda$ must be introduced to make the NJL model
regularized. We adapt a sharp cutoff regularization for three-dimensional
momentum integration, following the literature [17].

###  III.1 Gap equations

We employ the mean-field approximation, corresponding to the large $N_{c}$
limit, and then derive the gap equation and the thermodynamic potential [17].
There the quark condensates (on thermal average) act as the variable of the
potential and are $T$-dependent, which we define as

$\langle\bar{u}u\rangle\equiv\alpha,\quad\langle\bar{d}d\rangle\equiv\beta,\quad\langle\bar{s}s\rangle\equiv\gamma\,.$
(26)

Searching for the minimum point of the thermodynamic potential with respect to
$\alpha$, $\beta$, and $\gamma$ as variational parameters, we find the
stationary conditions, corresponding to the gap equations [26, 17]:

$\langle\bar{q}_{i}q_{i}\rangle=-2N_{c}\int^{\Lambda}\frac{d^{3}p}{(2\pi)^{3}}\frac{M_{i}}{E_{i}}\big{[}1-2(\exp(E_{i}/T)+1)^{-1}\big{]},$
(27)

where $E_{i}=\sqrt{M_{i}^{2}+p^{2}}$, $N_{c}$ denotes the number of colors to
be fixed to three, and $M_{i}$ are full quark masses including the dynamically
generated terms:

$\displaystyle M_{u}$ $\displaystyle=m_{u}-2g_{s}\alpha-2g_{D}\beta\gamma\,$
$\displaystyle M_{d}$ $\displaystyle=m_{d}-2g_{s}\beta-2g_{D}\alpha\gamma\,$
$\displaystyle M_{s}$ $\displaystyle=m_{s}-2g_{s}\gamma-2g_{D}\alpha\beta\,.$
(28)

###  III.2 Chiral and axial susceptibilities

In this subsection, we introduce susceptibilities for pseudoscalar and scalar
meson channels and give their explicit formulas in the present NJL model.

####  III.2.1 Pseudoscalar meson channel

In the $\eta$ \- $\eta\prime$ coupled channel, the pseudoscalar meson
susceptibility is defined on the generator basis as

$\chi_{P}^{ij}=\int_{T}d^{4}x\langle(i\bar{q}(x)\gamma_{5}\lambda^{i}q(x))(i\bar{q}(0)\gamma_{5}\lambda^{j}q(0))\rangle\,,$
(29)

where $i,j=0,8$. This $\chi_{P}^{ij}$ takes a matrix form

$\chi_{P}=\frac{-1}{1+G_{P}\Pi_{P}(0,0)}\cdot\Pi_{P}(0,0)\,,$ (30)

where $G_{P}$ is the coupling strength matrix and $\Pi_{P}$ is the
polarization tensor matrix, which are given respectively as

$G_{P}=\begin{pmatrix}G_{P}^{00}&G_{P}^{08}\\\
G_{P}^{80}&G_{P}^{88}\end{pmatrix}=\begin{pmatrix}g_{s}-\frac{2}{3}(\alpha+\beta+\gamma)g_{D}&-\frac{\sqrt{2}}{6}(2\gamma-\alpha-\beta)g_{D}\\\
-\frac{\sqrt{2}}{6}(2\gamma-\alpha-\beta)g_{D}&g_{s}-\frac{1}{3}(\gamma-2\alpha-2\beta)g_{D}\end{pmatrix}\,,$
(31) $\Pi_{P}=\begin{pmatrix}\Pi_{P}^{00}&\Pi_{P}^{08}\\\
\Pi_{P}^{80}&\Pi_{P}^{88}\end{pmatrix}=\begin{pmatrix}\frac{2}{3}(2I_{P}^{uu}+I_{P}^{ss})&\frac{2\sqrt{2}}{3}(I_{P}^{uu}-I_{P}^{ss})\\\
\frac{2\sqrt{2}}{3}(I_{P}^{uu}-I_{P}^{ss})&\frac{2}{3}(I_{P}^{uu}+2I_{P}^{ss})\end{pmatrix}\,,$
(32)

with $I_{P}^{ii}(\omega,\boldsymbol{p})$ being the pesudoscalar one-loop
polarization functions [28],

$I_{P}^{ii}(0,0)=-\frac{N_{c}}{\pi^{2}}\int^{\Lambda}_{0}dp\,p^{2}\frac{1}{E_{i}}\left[1-2\left(\exp(E_{i}/T)+1\right)^{-1}\right]\,,\qquad{\rm
for}\quad i=u,d,s\,.$ (33)

By performing the basis transformation, the pseudoscalar susceptibilities
defined in Eq.(17) on the flavor basis are thus obtained as

$\begin{pmatrix}\frac{1}{2}\chi_{P}^{uu}+\frac{1}{2}\chi_{P}^{ud}=\frac{1}{4}\chi_{\eta}\\\
\chi_{P}^{us}\\\
\chi_{P}^{ss}\end{pmatrix}=\begin{pmatrix}\frac{1}{6}&\frac{\sqrt{2}}{6}&\frac{1}{12}\\\
\frac{1}{6}&-\frac{\sqrt{2}}{12}&-\frac{1}{6}\\\
\frac{1}{6}&-\frac{\sqrt{2}}{3}&\frac{1}{3}\end{pmatrix}\begin{pmatrix}\chi_{P}^{00}\\\
\chi_{P}^{08}\\\ \chi_{P}^{88}\end{pmatrix}\,,$ (34)

where we have taken the isospin symmetric limit into account, i.e.,
$\chi_{P}^{uu}=\chi_{P}^{dd}$ and $\chi_{P}^{us}=\chi_{P}^{ds}$.

For $\chi_{\pi}$ defined in Eq.(17), the explicit formula in the NJL model
reads [17]

$\chi_{\pi}=\frac{-1}{1+G_{\pi}\Pi_{\pi}(0,0)}\cdot\Pi_{\pi}(0,0)\,,$ (35)

where $G_{\pi}=g_{s}+g_{D}\gamma$, which is the coupling strength in the pion
channel, and $\Pi_{\pi}$ is the quark-loop polarization function for
$\chi_{\pi}$, which is evaluated by using $I_{P}^{ii}$ in Eq.(33) as

$\Pi_{\pi}=I_{P}^{uu}+I_{P}^{dd}=2I_{P}^{uu}\,.$ (36)

####  III.2.2 Scalar meson channel

The definitions of scalar susceptibilities are similar to those for
pseudoscalars’, which are given just by removing $i\gamma_{5}$ in the
definition of pseudoscalar susceptibilities, and supplying the appropriate
one-loop polarization functions and the corresponding coupling constants.

In the $0$ \- $8$ coupled channel, the scalar susceptibility matrix $\chi_{S}$
is evaluated in the present NJL on the generator basis as

$\chi_{S}=\frac{-1}{1+G_{S}\Pi_{S}(0,0)}\cdot\Pi_{S}(0,0)\,,$ (37)

where $G_{S}$ is the coupling strength matrix,

$G_{S}=\begin{pmatrix}G_{S}^{00}&G_{S}^{08}\\\
G_{S}^{80}&G_{S}^{88}\end{pmatrix}=\begin{pmatrix}g_{s}+\frac{2}{3}(\alpha+\beta+\gamma)g_{D}&\frac{\sqrt{2}}{6}(2\gamma-\alpha-\beta)g_{D}\\\
\frac{\sqrt{2}}{6}(2\gamma-\alpha-\beta)g_{D}&g_{s}+\frac{1}{3}(\gamma-2\alpha-2\beta)g_{D}\end{pmatrix}\,.$
(38)

The scalar polarization tensor matrix $\Pi_{S}$ in Eq.(37) is given by

$\Pi_{S}=\begin{pmatrix}\Pi_{S}^{00}&\Pi_{S}^{08}\\\
\Pi_{S}^{80}&\Pi_{S}^{88}\end{pmatrix}=\begin{pmatrix}\frac{2}{3}(2I_{S}^{uu}+I_{S}^{ss})&\frac{2\sqrt{2}}{3}(I_{S}^{uu}-I_{S}^{ss})\\\
\frac{2\sqrt{2}}{3}(I_{S}^{uu}-I_{S}^{ss})&\frac{2}{3}(I_{S}^{uu}+2I_{S}^{ss})\,,\end{pmatrix}$
(39)

with the integral functions,

$I_{S}^{ii}(0,0)=-\frac{N_{c}}{\pi^{2}}\int^{\Lambda}_{0}p^{2}dp\frac{E_{i}^{2}-M_{i}^{2}}{E_{i}^{3}}\\{1-2[\exp(E_{i}/T)+1]^{-1}\\}\qquad
i=u,d,s\,.$ (40)

By moving on to the flavor base via the base transformation, the scalar
susceptibilities are cast into the form:

$\begin{pmatrix}\frac{1}{2}\chi_{S}^{uu}+\frac{1}{2}\chi_{S}^{ud}\\\
\chi_{S}^{us}\\\
\chi_{S}^{ss}\end{pmatrix}=\begin{pmatrix}\frac{1}{6}&\frac{\sqrt{2}}{6}&\frac{1}{12}\\\
\frac{1}{6}&-\frac{\sqrt{2}}{12}&-\frac{1}{6}\\\
\frac{1}{6}&-\frac{\sqrt{2}}{3}&\frac{1}{3}\end{pmatrix}\begin{pmatrix}\chi_{S}^{00}\\\
\chi_{S}^{08}\\\ \chi_{S}^{88}\end{pmatrix}\,,$ (41)

in which we have read $\chi_{S}^{uu}=\chi_{S}^{dd}$ and
$\chi_{S}^{us}=\chi_{S}^{ds}$. From this, the $\sigma$ meson susceptibility
can be also read off as

$\displaystyle\chi_{\sigma}$
$\displaystyle=\int_{T}d^{4}x\left[\langle(\bar{u}(0)u(0))(\bar{u}(x)u(x))\rangle+\langle(\bar{d}(0)d(0))(\bar{d}(x)d(x))\rangle+2\langle(\bar{u}(0)u(0))(\bar{d}(x)d(x))\rangle\right]$
$\displaystyle=2\chi_{S}^{uu}+2\chi_{S}^{ud}\,.$ (42)

We will not directly evaluate this $\chi_{\sigma}$ in the later section, but
it will be inferred when the ambiguity in subtracting the original form of the
Ward-identity in Eq.(21) by the scalar meson susceptibility, to get the
relation between the chiral and axial indicators, and the topological
susceptibility (in the later section, Sec. 4).

For the $\delta$ meson susceptibility, it is defined as

$\chi_{\delta}=\int_{T}d^{4}x\left[\langle(\bar{u}(0)u(0))(\bar{u}(x)u(x))\rangle_{\rm
conn}+\langle(\bar{d}(0)d(0))(\bar{d}(x)d(x))\rangle_{\rm conn}\right]\,.$
(43)

Similar to $\chi_{\pi}$ in Eq.(35), the explicit formula for $\chi_{\delta}$
reads [17]

$\chi_{\delta}=\frac{-\Pi_{\delta}(0,0)}{1+G_{\delta}\Pi_{\delta}(0,0)}\,,$
(44)

where $G_{\delta}=g_{s}-g_{D}\gamma$, which is the coupling strength in the
$\delta$ channel, and $\Pi_{\delta}=I_{S}^{uu}+I_{S}^{dd}=2I_{S}^{uu}$ is the
corresponding quark-loop polarization function.

###  III.3 Topological susceptibility in NJL model

Combining Eq. 21, 30, 31, and 34, we get the formula of $\chi_{\rm top}$
written in terms of the present NJL-model parameters:

$\displaystyle\chi_{\rm top}$
$\displaystyle=\frac{1}{2}m_{l}m_{s}(\chi_{p}^{us}+\chi_{p}^{ds})$
$\displaystyle=m_{l}m_{s}\left(\frac{1}{6}\chi_{P}^{00}-\frac{\sqrt{2}}{12}\chi_{P}^{08}-\frac{1}{6}\chi_{P}^{88}\right)$
$\displaystyle=\frac{-m_{l}m_{s}}{\rm{det}(1+G_{P}\Pi_{P})}\bigg{[}\frac{1}{6}(\Pi^{00}_{P}\Pi^{88}_{P}G^{88}_{P}-{\Pi^{08}_{P}}^{2}G^{88}_{P})-\frac{\sqrt{2}}{12}({\Pi^{08}_{P}}^{2}G^{08}_{P}-\Pi^{00}_{P}\Pi^{88}_{P}G^{08}_{P})-\frac{1}{6}(\Pi^{88}_{P}\Pi^{00}_{P}G^{00}_{P}-{\Pi^{80}_{P}}^{2}G^{00}_{P})\bigg{]}$
$\displaystyle=\frac{-m_{l}m_{s}g_{D}\alpha}{6\,\rm{det}(1+G_{P}\Pi_{P})}((\Pi^{80}_{P})^{2}-\Pi^{88}_{P}\Pi^{00}_{P})\,.$
(45)

One can clearly see that $\chi_{\rm top}\propto m_{l}m_{s}$, which is
reflected by the flavor-singlet nature of the QCD vacuum. Note also that
$\chi_{\rm top}\propto g_{D}$. This should be so because the KMT determinant
coupling $g_{D}$ is directly linked with the gluonic anomaly term in the axial
anomaly equation (25).

###  III.4 NJL estimates

In this subsection, we evaluate the temperature dependence of the (subtracted)
quark condensate, scalar and pseudoscalar susceptibilities, and the
topological susceptibility. We also check the consistency with the recent
lattice QCD data on 2 + 1 flavors at physical point, and also with other
effective models of QCD.

####  III.4.1 Parameter setting

In the present NJL model of Eq.(24), we have five parameters that need to be
fixed: the light quark mass $m_{l}$, the strange quark mass $m_{s}$, the
coupling constants $g_{s}$ and $g_{D}$, and the three-momentum cutoff
$\Lambda$. To fix the parameters, we take the following conservative and
empirical input values at $T=0$ in the isospin symmetric limit [17]:

$m_{\pi}=136\,{\rm MeV},\quad f_{\pi}=93\,{\rm MeV},\quad m_{K}=495.7\,{\rm
MeV},\quad m_{\eta\prime}=957.5\,{\rm MeV}\,.$ (46)

For readers convenience, concise NJL formulae are presented in Appendix A. To
fix the remaining one degree of freedom, we follow the literature [29, 17] to
take light quark mass $m_{l}=5.5$ MeV (at the renormalization scale of 1 GeV).
Thus all the model parameters are fixed, which are presented in Table 1.

model parameter description | input value
---|---
isospin symmetric light quark mass $m_{l}$ | 5.5 MeV
strange quark mass $m_{s}$ | 138 MeV
four-fermion coupling constant $g_{s}$ | 0.358 ${\rm fm}^{2}$
six-fermion coupling constant $g_{D}$ |  $-$ 0.0275 ${\rm fm}^{5}$
cutoff $\Lambda$ | 631.4 MeV
Table 1: The model parameter setting, followed from [17].

With the above parameter set #1#1#1 The constituent quark masses are also
estimated as $M_{u}(T=0)=M_{d}(T=0)\simeq 334.2{\rm MeV}$, $M_{s}(T=0)\simeq
530.1{\rm MeV}$, where use has been made of Eqs.(27), (28), together with the
model parameters listed in Table 1. , we estimate the topological
susceptibility

$\displaystyle\chi_{\rm top}\simeq 0.025/{\rm fm}^{4}\,.$ (47)

For this $\chi_{\rm top}$, comparison with the results from the lattice QCD
simulations is available, which are $\chi_{\rm top}=0.019(9)/{\rm fm}^{4}$
[30], and $\chi_{\rm top}=0.0245(24)_{\rm stat}(03)_{\rm flow}(12)_{\rm
cont}/{\rm fm}^{4}$ [12]. Here, for the latter the first error is statistical,
the second one comes from the systematic error, and the third one arises due
to changing the upper limit of the lattice spacing range in the fit. Although
their central values do not agree each other, we may conservatively say that
the difference between them is interpreted as a systematic error from the
individual lattice QCD calculation.

We will not consider intrinsic-temperature dependent couplings, instead, all
the $T$ dependence should be induced only from the thermal quark loop
corrections to the couplings defined and introduced at vacuum. As it will turn
out below, the present NJL shows quite good agreement with lattice QCD results
on the temperature scaling (normalized to the pseudo-critical temperature) for
the chiral, axial, and topological susceptibilities. In this sense, we do not
need to introduce such an intrinsic $T$ dependence for the model parameters in
the regime up to temperatures around the chiral crossover.

####  III.4.2 Subtracted quark condensate

The quark condensate in the NJL model involves a ultraviolet divergence (which
is dominated by a quadratic divergence) due to its vacuum part
($\langle-\bar{q}q\rangle\sim N_{c}m_{q}\Lambda^{2}/(4\pi^{2})$), and is
needed to be renormalized when compared with lattice data. Since the quadratic
divergences in the quark condensate come along with current quark masses (as
above), we use a subtracted quark condensate as the chiral order parameter,
which has been adopted in the lattice simulations:
$\Delta_{l,s}(T)\equiv\langle\bar{l}l\rangle-\frac{m_{l}}{m_{s}}\langle\bar{s}s\rangle$
, where
$\langle\bar{l}l\rangle=\langle\bar{u}u\rangle=\langle\bar{d}d\rangle$.

Figure 2 shows the subtracted quark condensate as a function of temperature
predicted from the present NJL model, in comparison with the 2+1 flavor data
from the lattice QCD at the physical point [31]. The pseudo-critical
temperature $T_{\rm pc}$ is (for the NJL prediction) defined as
$d^{2}\langle\bar{l}l\rangle(T)/dT^{2}|_{T=T_{\rm pc}}=0$. We have found
$T_{\rm pc}|_{\rm NJL}$ $\simeq 188$ MeV, which is compared with the lattice
result $T_{\rm pc}|_{\rm lat.}$ $\simeq 155$ MeV [31, 32, 33, 34, 35].

In the figure, we have normalized $T$ by their $T_{\rm pc}$. The reason to
take the rescaled dimensionless ratios is to make possible systematic errors
for the model selection and calculation reduced. For instance, the NJL model
predicts somewhat larger $T_{pc}$, and $\langle\bar{q}q\rangle$ at any $T$,
and basically all dimensionful quantities tend to be overestimated by about
30%: that is thought of as a systematic trend which could be associated with
the validity of the large $N_{c}$ approximation. Then, the dimensionless
quantities like $T/T_{pc}$ and
$\langle\bar{q}q\rangle_{T}/\langle\bar{q}q\rangle_{T=0}$ can have reduced
systematic errors. So, it would be better to take the rescaled ratios, in
order to comprehend how the current effective model can be compatible with the
lattice result.

Figure 2:  $T/T_{\rm pc}$ dependence of the subtracted quark condensate, in
comparison with data from the lattice QCD with 2 + 1 flavors [31]. The
normalization factor, the pseudo-critical temperature for the chiral crossover
($T_{\rm pc}$) has been set to individual values estimated from the present
NJL model ($T_{\rm pc}|_{\rm NJL}\simeq 188$ MeV) and the lattice simulation
($T_{\rm pc}|_{\rm lat.}\simeq 155$ MeV).

From Fig. 2, we see that the present NJL prediction is consistent with the
lattice data, confirming that the present model describes the chiral crossover
phenomenon quite well.

####  III.4.3 Chiral and axial susceptibility partners

The scalar and pseudoscalar susceptibilities ($\chi_{\eta}$, $\chi_{\pi}$,
$\chi_{\sigma}$, $\chi_{\delta}$) presented in Eqs. (34), (35), (42), and (44)
are correlated with each other by the chiral $SU(2)_{L}\times SU(2)_{R}$ and
$U(1)_{A}$ transformations [36], which can be summarized as the following
cartoon:

${\chi_{\pi}}$${\chi_{\sigma}}$${\chi_{\delta}}$${\chi_{\eta}}$$U(1)_{A}$SU(2)SU(2)$U(1)_{A}$

The chiral and axial partners will be degenerate each other in the symmetric
limits:

$\displaystyle\chi_{\pi}$
$\displaystyle=\chi_{\sigma},\quad\chi_{\delta}=\chi_{\eta}\qquad\text{(chiral
$SU(2)$ symmetric limit)}$ $\displaystyle\chi_{\pi}$
$\displaystyle=\chi_{\delta},\quad\chi_{\sigma}=\chi_{\eta}\qquad\text{($U(1)_{A}$
axial symmetric limit)}\,.$ (48)

Then, observation of null difference between the above partners can
effectively monitor the restoration of the related symmetry.  Figure 3 shows
the plots of $\chi_{\pi-\delta}/T^{2}$ (left panel) and
$\chi_{\eta-\delta}/T^{2}$ (right panel) normalized their values measured at
$T=T_{\rm pc}$, in comparison to the lattice QCD result [4]. Both of the
$T/T_{\rm pc}$ dependence on the normalized susceptibilities show qualitative
agreement of the present NJL estimates with the current lattice data.

Figure 3: The normalized susceptibility differences versus temperature
normalized by the pseudo-critical temperature $T_{\rm pc}$ for the chiral
crossover, in comparison with lattice QCD data for 2 + 1 flavors [4]. The left
panel corresponds to $\chi_{\pi-\delta}/T^{2}$ divided by its value evaluated
at $T=T_{\rm pc}$, and the right panel is the same plots for
$\chi_{\eta-\delta}/T^{2}$. The present NJL yields $T_{\rm pc}|_{\rm NJL}=188$
MeV, and the quoted lattice result predicts $T_{\rm pc}|_{\rm lat}=155$ MeV
(at the central value).

####  III.4.4 Topological susceptibility

We numerically evaluate $\chi_{\rm top}$ in Eq.(16), with the present NJL
estimates on the quark condensates and pseudoscalar susceptibilities, as a
function of temperature. In Fig. 4, we plot the temperature dependence of the
unnormalized topological susceptibility $\chi_{\rm top}^{1/4}$, where we have
taken the absolute value of $\chi_{\rm top}$. Comparison with the dilute
instanton gas approximation (DIGA) [37, 38], the linear sigma model result
(denoted as CJT in the figure) [6] and the result from lattice simulation in
the continuum limit [10, 11, 12] have also been displayed. The DIGA prediction
has been quoted from the literature [10]. For the way of error bars associated
with the DIGA, see the cited reference. The temperature is normalized by the
pseudo-critical temperature in the figure, where we have taken $T_{\rm
pc}|_{\rm NJL}=188$ MeV for the NJL case, $T_{\rm pc}|_{\rm CJT}=215$ MeV for
the linear sigma model case, and $T_{\rm pc}|_{\rm lat}=155$ MeV for the
lattice.

Figure 4: $T/T_{\rm pc}$ dependence of topological susceptibilities, compared
with lattice data [10, 11, 12] and other models as described in the text. The
vertical line at $T/T_{\rm pc}=1.6$ indicates the theoretical validity of the
present NJL model (denoted as the seemingly-cutoff), above which regime the
present model description may not be reliable because of lack of the gluonic
degrees of freedom (see also footnote #2). This limit needs to be taken into
account in referring to the range of $T/T_{\rm pc}$, when compared to the
lattice data displayed in the figure.

Figure 4 shows good qualitative agreement between the NJL analysis and lattice
result. We see that even though the NJL description without gluonic
contribution may not be reliable for $T/T_{\rm pc}>$ 1.6, the good agreement
keeps in the whole range of the available lattice data, $T/T_{\rm pc}\sim
0-4$.

At $T<T_{\rm pc}$, all the results fit perfectly with each other, including
the linear sigma model estimate. In contrast, when $T>T_{\rm pc}$, we see
substantial deviation for the linear sigma model prediction from the NJL’s and
lattice results #2#2#2 Even within the linear sigma model description, it has
been also shown that improved treatment of the thermal characteristics of the
linear sigma model based on the application of the functional renormalization
group technique to the effective potential reproduces quantitatively correctly
lattice results for the pseduocritical temperature [39]. . In the literature
[6], the pseudoscalar susceptibility terms were not able to evaluate, because
the authors did not include the higher order terms in the current quark
masses, and therefore, performing the second order derivative on the mass
parameter to obtain pseudoscalar susceptibility would not be worked out. Thus,
their $\chi_{\rm top}$ only includes the quark condensate terms. The present
NJL model is able to give the pseudoscalar susceptibility contribution to
$\chi_{\rm top}$, to achieve an improved estimate on the quark condensate. The
better qualitative agreement of the NJL with the lattice result may thus imply
the importance of contributions from the pseudoscalar susceptibilities at
higher temperatures, though the model estimate may not rigorously be valid
beyond the seemingly-cutoff temperature ($T>1.6T_{\rm pc}$) as noted above.

The topological susceptibility $\chi_{\rm top}$ has been discussed based on
the NJL model descriptions similar to ours [22, 23, 24, 25]. However, the
anomalous Ward-identity in Eq.(1) and the flavor-singlet condition necessary
in deriving the proper $\chi_{\rm top}$ in Eq.(16) have been ignored there.
Absence of the former led to miss-identification of the restoration of the
axial symmetry, which was played by $\chi_{\rm top}$ in the literature, while
the latter missing factor made improper temperature dependence of $\chi_{\rm
top}$.

##  IV Evaluation of QCD trilemma estimator

###  IV.1 Violation of QCD trilemma at physical point in a whole temperature
regime

Figure 5 shows values of the trilemma estimator $R$ evolved with $T$, allowing
$m_{s}$ off the physical point with $m_{l}$ kept physical. See the middle-
solid curve with $m_{s}=138$ MeV, which corresponds to real-life QCD.
Comparison with the available $2+1$ flavor-lattice QCD data (with
$m_{\pi}=135$ MeV) on $R$ [4] — reconstructed from the data on
$\chi_{\pi-\delta}$ and $\chi_{\eta-\delta}$ through the relation Eq.(3) — has
also been displayed (in the zoomed-in window), which shows good agreement
including the error bars, for $140\,{\rm MeV}\lesssim T\lesssim 200\,{\rm
MeV}$. The reconstructed data of $R$ include large errors, which is mainly due
to the large uncertainty of the lattice measurement on $\chi_{\eta-\delta}$
(See Fig. 3).

Figure 5:  Plots showing that real-life QCD is imbalanced, which is monitored
by the trilemma estimator $R$ defined in Eq.(2). Estimates have been done
based on the NJL model described as in the text. Comparison with the $2+1$
flavor-lattice QCD data (with $m_{\pi}=135$ MeV) in the available $T$ range
has also been displayed with the error bars [4] (in the zoomed-in window). The
curve with $m_{s}=138$ MeV points to real-life QCD with three flavors, while
the quenched-strange quark limit has been achieved by taking $m_{s}=50$ GeV,
corresponding to the two-flavor limit. Another curve with $m_{s}=5.5$ MeV
denotes a conjectured prediction in the three-flavor symmetric limit. The
“balanced” and “imbalanced” regimes are defined in Eq.(4).

Remarkably, in a whole temperature regime including the chiral crossover
regime, real-life QCD stays outside the “balanced” region defined as in
Eq.(4). We have observed $R\simeq 0.05$ at around $T$ covering the crossover
point ($T_{\rm pc}|_{\rm NJL}\simeq 188$ MeV: $140\,{\rm MeV}\lesssim
T\lesssim 200\,{\rm MeV}$), consistently with the lattice data, and $R\lesssim
0.01$ at $T\gtrsim 300$ MeV. Namely, the amount of imbalance is slightly
amplified by thermal loop effects as $T$ develops from zero #3#3#3 Above
$T\sim 300$ MeV corresponding to the typical scale of the constituent quark
mass, the NJL description as the effective theory of QCD may be somewhat
unreliable because the deconfining color degrees of freedom and thermal
gluonic contributions would be significant. .

One might note that subtraction by $\chi_{\delta}$ in Eq.(1) is ambiguous, and
can be replaced by another chiral susceptibility in the sigma meson channel
($\chi_{\sigma}$). We have checked that this replacement does not alter our
main conclusion that real-life QCD involves big imbalance. We have also found
that $\chi_{\pi-\sigma}\gg\chi_{\eta-\sigma}$ at $T=0$,
$\chi_{\pi-\sigma}\ll\chi_{\eta-\sigma}$ at around the chiral crossover, then
$\chi_{\pi-\sigma}$ will get close to $\chi_{\pi-\delta}$, and finally go to
zero. The latter trend is consistent with the currently available lattice data
[4].

Although the present model parameters are fixed at the physical point, we may
deduce some conjectures on the violation of QCD trilemma in a view of the
quark mass difference. Extrapolating off real-life QCD, one can then observe
that the “imbalanced” domain still covers the two-flavor limit case with
$m_{s}=50$ GeV (bottom-dot-dashed curve), where strange quark is decoupled,
and the amount of imbalance is greater than that in the real-life QCD case.
Taking the three-flavor symmetric limit $m_{s}=m_{l}$ with $m_{l}$ fixed to
the physical value, we find “balanced” QCD (top-dashed curve), which keeps
almost constant $R$ at any finite $T$ within the “balanced” interval in
Eq.(4). This implies that the three-flavor symmetry would be related to the
relaxation of the QCD trilemma.

###  IV.2 QCD trilemma and flavor symmetry

Since the order of magnitude for $R$ tends to be almost fixed at $T=0$, we may
focus only on $R$ at $T=0$, and look into the flavor-symmetry dependence on
$R$, by varying $m_{s}$ in a wide range, with fixed $m_{l}$ to the physical
value. Figure 6 shows plots on $R$ as a function of $m_{s}/m_{l}$, together
with the “balanced” interval in Eq.(4). As $m_{s}$ goes off the flavor
symmetric limit in the “balanced” domain to be smaller, $R$ tends to get
larger, to flow into the “imbalanced” domain with gigantically suppressed
$\chi_{\rm top}$. The figure clearly shows that “balanced” QCD should have had
some approximate three-flavor symmetry for up, down and strange quarks with
$0.06\lesssim m_{s}/m_{l}\lesssim 6$.

Figure 6:  Plots on the QCD trilemma estimator $R$ at $T=0$ as a function of
$m_{s}/m_{l}$, along with the “balanced” interval defined in Eq.(4). The
shaded domain surrounded by the real-life QCD point implies confidence level
intervals for the model prediction, where the model parameters, except for
$m_{s}$, have been fixed at the physical point, as noted in the text. The
thinner-shaded regions should be understood as indefinitely extrapolated
results with somewhat poor reliability.

We shall investigate the dependence of the flavor-symmetry violation on the
imbalanced QCD trillemma in more details.  First of all, we may simply suppose
that the scalar and pseudoscalar susceptibilities are scaled with the
associated meson masses #4#4#4 The susceptibilities correspond to meson-
correlation functions at zero momentum transfer. This is in contrast to the
conventional meson correlators depending on the transfer momentum, from which
meson masses are read off. Furthermore, the susceptibilities involve contact
term contributions independent of momenta, which could be sensitive to a high-
energy scale physics, while the conventional meson correlators are dominated
by the low-lying meson mass scale. Nevertheless, the degeneracy of the chiral
or axial partners at high temperatures, similar to those detected in the
susceptibility, can also be seen in the mass difference or equivalently the
degeneracy of the conventional meson correlators for the partners, which is
simply because the mass difference plays an alternative indicator of the
chiral or axial breaking as observed in the lattice simulations [40, 41]. ,
like $\chi_{\delta}\propto 1/m_{\delta}^{2}$, $\chi_{\eta}\propto
1/m_{\eta}^{2}$ and $\chi_{\pi}\propto 1/m_{\pi}^{2}$, and consider the light
quark mass $m_{l}$ to generically differ from the strange quark mass $m_{s}$,
including the real-life QCD case with the three-flavor symmetry broken. Among
the susceptibilities, $\chi_{\pi}$ is most sensitive to the current mass of
the light quarks ($m_{l}$), because the pion is the pseudo Nambu-Goldstone
boson of spontaneous breaking of $SU(2)_{L}\times SU(2)_{R}$ symmetry carried
by the light quarks. The $\chi_{\pi}$ thus monotonically gets smaller
(larger), as $m_{l}$ gets larger (smaller), by following $\chi_{\pi}\propto
1/m_{\pi}^{2}\sim 1/m_{l}$. On the other hand, the other pseudoscalar
susceptibility $\chi_{\eta}$ significantly involves the U(1)-axial anomaly
contribution in $m_{\eta}$, so it almost keeps constant in $m_{l}$. The scalar
susceptibility $\chi_{\delta}$, free from the Nambu-Goldstone boson nature,
also keeps constant with the change of $m_{l}$. Besides, the topological
susceptibility $\chi_{\rm top}$ also simply scales with $m_{l}$, respecting
the flavor-singlet condition: $\chi_{\rm top}\to 0$ as $m_{l}\to 0$, and will
be completely constant in $m_{l}$ for $m_{l}>m_{s}$ due to decoupling of the
“light” quarks. Thus the difference in magnitude of susceptibilities are
simply originated from the scaling properties with respect to the current mass
of the light quarks.  We plot those $m_{l}$ scaling behaviors (at $T=0$) in
Fig. 7. The light quark mass is allowed to vary from $10^{-2}$eV to the cutoff
scale of the presently employed NJL model (631.4 MeV), since the result from
$m_{l}$ above the cutoff scale would be of poor reliability. From the figure,
the $m_{l}$ dependence is read off and the susceptibilities are found to take
simple power laws when $m_{l}\lesssim m_{s}$:

$\displaystyle\chi_{\pi}\sim m_{l}^{-1}\,,$ $\displaystyle\chi_{\eta}\sim{\rm
constant}\,,\quad{\rm for}\qquad m_{l}<m_{s}\,,$
$\displaystyle\chi_{\delta}\sim{\rm constant}\,,\quad{\rm for}\qquad
m_{l}<m_{s}\,,$ $\displaystyle\chi_{\rm top}\sim m_{l}\quad{\rm for}\,,\qquad
m_{l}<m_{s}\,.$ (49)

Figure 7:  (Without the three-flavor symmetry): the $m_{l}$ dependence on
$\chi_{\pi},\chi_{\eta}$, $\chi_{\delta}$ [MeV2] (left panel), and $\chi_{\rm
top}$ [MeV4] (right panel) at $T=0$. The value of strange quark mass $m_{s}$
is fixed to be equal to $m_{l}$ in the plots.

Next, consider the three-flavor symmetric limit, where $m_{l}$= $m_{s}$. In
this case QCD is balanced as noted above. It also turns out that the scaling
law of $\chi_{\eta}$ in Eq. (49) is broken: the Ward identity Eq. (21) tells
us that the difference between $\chi_{\pi}$ and $\chi_{\eta}$ is controlled by
the $4\frac{m_{s}}{m_{l}}\chi_{P}^{ls}$ term (where
$\chi_{P}^{ls}=\chi_{P}^{us}=\chi_{P}^{ds}$). Since no preference among quark
flavors is present in the flavor symmetric case, $\chi_{P}^{ls}$ should be on
the same order of magnitude as that of $\chi_{\pi}$, which we have indeed
numerically confirmed. See Fig, 8. Since $m_{l}=m_{s}$, there is no extra
power scaling of $1/m_{l}$ which is present in the flavor asymmetric case and
leads to big enhancement of the ($4\frac{m_{s}}{m_{l}}\chi_{P}^{ls}$) part to
destructively interfere with $\chi_{\pi}$, yielding a much suppressed
$\chi_{\eta}$ compared to $\chi_{\pi}$ (See Fig. 7). Thus the scaling law of
$\chi_{\eta}$ is the same as that of $\chi_{\pi}$, i.e., $\chi_{\eta}$ $\sim$
$m_{l}^{-1}$, while others take the same scaling laws as in Eq.(49), namely,

$\displaystyle\chi_{\pi}\sim\chi_{\eta}\sim m_{l}^{-1}\quad{\rm for}\qquad
m_{l}=m_{s}\,,$ $\displaystyle\chi_{\delta}\sim{\rm constant}\quad{\rm
for}\qquad m_{l}=m_{s}\,,$ $\displaystyle\chi_{\rm top}\sim m_{l}\quad{\rm
for}\qquad m_{l}=m_{s}\,,$ (50)

as depicted in Fig.8.

This scaling violation in the flavor symmetric case can also be understood as
a big suppression of the $U(1)_{A}$ anomaly contribution, coupled to the
flavor violation, to $m_{\eta}^{2}$, which dominates in $\chi_{\eta}$ in the
flavor asymmetric case: in the flavor symmetric case we have
$\chi_{\pi}=\chi_{P}^{88}$, and $\chi_{\eta}=\chi_{\pi}+4\chi_{P}^{ls}$ with
$\chi_{P}^{ls}=1/6(\chi_{P}^{00}-\chi_{P}^{88})$. Straightforward numerical
evaluation reveals that $\chi_{P}^{88}\gg\chi_{P}^{00}$ for small $m_{l}$.
Then, we find $\chi_{\eta}\approx\chi_{\pi}/3\sim 1/m_{l}$ for small $m_{l}$.
In particular, note that $\chi_{P}^{88}=\chi_{\pi}$ does not include the
$U(1)_{A}$ anomaly effect, and is now much larger than the $U(1)_{A}$-anomaly
affected $\chi_{P}^{00}$ part, which implies the $U(1)_{A}$ anomaly
contribution is much suppressed in $\chi_{\eta}$, hence in $m_{\eta}$ as well.

Figure 8:  (With the three-flavor symmetry): The same plots as those in Fig.
7, but with the three-flavor symmetry, where $m_{s}=m_{l}$.

Thus the three-flavor symmetry brings the $m_{l}$ scaling universal among
susceptibilities for the pion and eta meson, by suppressing the $U(1)_{A}$
anomaly contribution only in the eta meson channel. This, on the other side of
the same coin, indicates that the amplification of the $U(1)_{A}$ anomaly
effect in $\chi_{\eta-\delta}$ so that the chiral symmetry is made faster
restored than the axial symmetry with holding the anomalous chiral Ward
identity Eq.(1).

##  V Conclusion

In conclusion, real-life QCD is required to relax the trilemma $(R\ll 1)$,
meaning that the much smaller strength of the chiral breaking is given by a
big cancellation between the strength of the axial breaking and the transition
rate of the topological charge. This is schematically depicted in Fig. 1. This
is “imbalance” of the QCD vacuum, present in a whole temperature regime of
thermal QCD including the vacuum at $T=0$. This imbalance or violation of QCD
trilemma is triggered due to the three-flavor symmetry violation for up, down
and strange quarks, in particular the lightness of up and down quarks. The
violation of the three-flavor symmetry specifically brings enhancement of the
$U(1)_{A}$ anomaly contribution in the chiral $SU(2)$ breaking to be
insensitive to the smallness of the light quark mass, while the axial
indicator and the transition rate of the topological charge are fairly
insensitive to the flavor symmetry — in other words, the flavor symmetry
dramatically suppresses the $U(1)_{A}$ anomaly effect in the chiral $SU(2)$
breaking. This implies that in the QCD vacuum with the trilemma realized, all
the three-flavor octet pseudoscalar mesons act like identical pseudo Nambu-
Goldstone bosons, so that the chiral and axial breaking together with the
topological charge contribute to the QCD vacuum on the same order of
magnitude.

The present work confirms and extends the suggestion recently reported from
lattice QCD with 2 flavors on dominance of the axial and topological
susceptibilities left in the chiral susceptibility at high temperatures over
the chiral crossover [2].

The violation of QCD trilemma would be crucial to deeply pursue the expected
dominance of the chiral symmetry breaking in the origin of mass, and provides
the understanding of mechanism how the faster (effective) restoration of the
chiral symmetry in the presence of contamination with the $U(1)_{A}$ anomaly
is achieved: it is understood as a big cancellation between the axial and
topological susceptibilities, due to the three-flavor violation among up, down
and strange quarks. Our findings can directly be tested on lattice QCD with 2
+ 1 flavors at physical point and also off physical point, in the future.

It would be nice if the violation of QCD trilemma could be evaluated directly
using lattice QCD, but since the quark masses need to be varied, the lattice
calculation is so costly. As the first step, it would be appropriate to
perform the calculation using a conventional effective model like the present
NJL. In fact, it has been shown that the present NJL model matches the lattice
results on a couple of observables relevant to the chiral crossover regime,
within a deviation of about 30% (consistent with the expected theoretical
uncertainty of the large $N_{c}$/mean field approximation on that the present
NJL is built). This shows that the present NJL model, though being based on
the rough mean field approximation, indeed is even quantitatively valid as an
effective model of QCD, as good as the lattice QCD. Thus the present work
gives motivation to the research in other QCD-compatible theories, such as
lattice and/or functional renormalization group approach.

In closing, we give several comments related to possible applications of
concept of QCD trilemma:

  * •

The notion of QCD trilemma and its violation would also provide us with a new
guiding principle to explore the flavor dependence of the chiral phase
transition, such as the Columbia plot [15].

  * •

It would be also intriguing to study the violation of QCD trilemma, by means
of the nonperturbative renormalization group.

  * •

Correlations between the violation of QCD trilemma and deconfinement-
confinement phase transition can be addressed, when the present NJL model is
extended by including the Polyakov loop terms.

  * •

Since the form of the anomalous chiral Ward identity in Eq.(1) will be intact
as long as the quark mass terms are only the leading source to explicitly
break the three-flavor chiral symmetry, the presently described argument would
possibly be applicable also to dense QCD, and thermomagnetic QCD, where in the
latter case the electromagnetic interactions contribute as the subleading
(loop) order to the chiral breaking, as in Eq.(18).

Those will deserve to another publication.

## Acknowledgement

We thank Tetsuo Hatsuda for giving us a crucial comment on the
susceptibilities, and are also grateful to Hidenori Fukaya for fruitful
discussion. This work was supported in part by the National Science Foundation
of China (NSFC) under Grant No.11747308, 11975108, 12047569, and the Seeds
Funding of Jilin University (S.M.). The work of A.T. was supported by the
RIKEN Special Postdoctoral Researcher program and partially by JSPS KAKENHI
Grant Number JP20K14479.

##  Appendix A NJL formulae

In this Appendix we list the NJL formulas used for the inputs observables in
Eq.(46), which can be found in Ref. [17]:

  * •

The pion decay constant $f_{\pi}$: it is computed by directly evaluating quark
loop contributions to the spontaneously broken $SU(2)$ axial current
$J_{\mu}^{A,a}=\bar{l}\gamma_{\mu}\gamma_{5}(\sigma^{a}/2)l$ ($a=1,2,3$),
through the definition of $f_{\pi}$, $\langle
0|J_{\mu}^{A\,a}(0)|\pi^{b}(p)\rangle=-ip_{\mu}f_{\pi}\delta^{ab}$. In the NJL
model with the large $N_{c}$ limit taken (summing up the ring diagrams), we
thus have

$f_{\pi}=G_{\pi
q}(m_{\pi})M_{u}\frac{2N_{c}}{\pi^{2}}\int^{\Lambda}_{0}dp\frac{p^{2}}{\sqrt{M_{u}^{2}+p^{2}}[4(M_{u}^{2}+p^{2})-m_{\pi}^{2}]}\,,$
(51)

where $G_{\pi q}(p=m_{\pi})$ is the pion wavefunction renormalization
amplitude evaluated at the onshell,

$G^{2}_{\pi
q}(m_{\pi})=\left(\frac{N_{c}}{2\pi^{2}}\int^{\Lambda}_{0}p^{2}dp\frac{\sqrt{M_{u}^{2}+p^{2}}}{(p^{2}+M_{u}^{2}-\frac{m_{\pi}^{2}}{4})^{2}}\right)^{-1}\,.$
(52)

  * •

The pion mass $m_{\pi}$ is computed by extracting the pole of the pion
propagator dynamically generated by the quark loop contribution in the NJL
with a resummation technique (Random Phase Approximation) applied [17]. The
pole position is thus detected as

$1+G_{\pi}\Pi_{\pi}(m_{\pi}^{2})=0\,.$ (53)

This pion mass is actually related to the light quark condensate, via the the
low energy theorem (the so-called Gell-Mann-Oakes-Renner relation):

$m_{\pi}^{2}=-2m_{l}\langle\bar{l}l\rangle/f_{\pi}^{2}\,.$ (54)

  * •

The kaon mass $m_{K}$ is calculable in the same way as in the case of
$m_{\pi}$ above:

$1+G_{K}\Pi_{K}(m_{K}^{2})=0\,,$ (55)

where

$G_{K}=g_{s}+g_{D}\langle\bar{u}u\rangle\,,$ (56)
$\Pi_{K}(w)=2F(w;u,s)+2F(-w;s,u)\,,$ (57)
$F(w;i,j)=-\frac{N_{c}}{4\pi^{2}}\int^{\Lambda}_{0}p^{2}dp\left[\frac{1}{E_{i}}f_{ij}(w)+\frac{1}{E_{j}}f_{ji}(w)\right]\,,$
(58) $f_{ij}(w)=2\frac{M_{i}(M_{j}-M_{i})-E_{i}w}{E_{j}^{2}-(E_{j}+w)^{2}}\,.$
(59)

  * •

The $\eta^{\prime}$ mass $m_{\eta^{\prime}}$ is identified as the highest mass
eigenvalue arising from the mass mixing in the $0-8$ channel. Similarly to the
pion and kaon cases, $m_{\eta^{\prime}}$ is then extracted by the highest pole
of the mixed propagator in the $0-8$ channel, $D(q^{2}),$ as

$D_{\eta^{\prime}}(m_{\eta^{\prime}}^{2})=0\,,$ (60)

where

$D(q^{2})=-G_{P}^{-1}\left(\frac{1}{1+G_{P}\Pi^{P}(q^{2})}\right)\equiv\begin{pmatrix}A(q^{2})&B(q^{2})\\\
B(q^{2})&A(q^{2})\end{pmatrix}\,,$ (61)

through the diagonalization process like

$T_{\theta}D^{-1}(q^{2})T_{\theta}^{-1}=\begin{pmatrix}D_{\eta^{\prime}}^{-1}(q^{2})&0\\\
0&D_{\eta}^{-1}(q^{2})\,\end{pmatrix}\,.$ (62)

Here $\Pi^{P}(q^{2})$ is a function given as the generalization of Eq. (32)
with the replacement of the loop function $I_{ii}^{P}(w^{2})$:

$I_{ii}^{P}(w^{2})=-\frac{4N_{c}}{\pi^{2}}\int^{\Lambda}_{0}p^{2}dp\frac{E_{i}}{4E_{i}^{2}-w^{2}}\,,$
(63)

and $T_{\theta}$ is the diagonalization matrix,

$T_{\theta}=\begin{pmatrix}\cos(\theta)&-\sin(\theta)\\\
\sin(\theta)&\cos(\theta)\end{pmatrix},\quad\tan(2\theta)=\frac{2B(q^{2})}{C(q^{2})-A(q^{2})}\,.$
(64)

## References

  * [1] E. V. Shuryak, Comments Nucl. Part. Phys. 21, no.4, 235-248 (1994) [arXiv:hep-ph/9310253 [hep-ph]].
  * [2] S. Aoki et al. [JLQCD], PTEP 2022, no.2, 023B05 (2022) doi:10.1093/ptep/ptac001 [arXiv:2103.05954 [hep-lat]].
  * [3] Y. Aoki, G. Endrodi, Z. Fodor, S. D. Katz and K. K. Szabo, Nature 443, 675-678 (2006) doi:10.1038/nature05120 [arXiv:hep-lat/0611014 [hep-lat]].
  * [4] T. Bhattacharya, M. I. Buchoff, N. H. Christ, H. T. Ding, R. Gupta, C. Jung, F. Karsch, Z. Lin, R. D. Mawhinney and G. McGlynn, et al. Phys. Rev. Lett. 113, no.8, 082001 (2014) doi:10.1103/PhysRevLett.113.082001 [arXiv:1402.5175 [hep-lat]].
  * [5] A. Gómez Nicola and J. Ruiz de Elvira, JHEP 03, 186 (2016) doi:10.1007/JHEP03(2016)186 [arXiv:1602.01476 [hep-ph]].
  * [6] M. Kawaguchi, S. Matsuzaki and A. Tomiya, Phys. Rev. D 103, no.5, 054034 (2021) doi:10.1103/PhysRevD.103.054034 [arXiv:2005.07003 [hep-ph]].
  * [7] T. D. Cohen, Phys. Rev. D 54, R1867-R1870 (1996) doi:10.1103/PhysRevD.54.R1867 [arXiv:hep-ph/9601216 [hep-ph]].
  * [8] T. D. Cohen, [arXiv:nucl-th/9801061 [nucl-th]].
  * [9] S. Aoki, H. Fukaya and Y. Taniguchi, Phys. Rev. D 86, 114512 (2012) doi:10.1103/PhysRevD.86.114512 [arXiv:1209.2061 [hep-lat]].
  * [10] P. Petreczky, H. P. Schadler and S. Sharma, Phys. Lett. B 762, 498-505 (2016) doi:10.1016/j.physletb.2016.09.063 [arXiv:1606.03145 [hep-lat]].
  * [11] C. Bonati, M. D’Elia, G. Martinelli, F. Negro, F. Sanfilippo and A. Todaro, JHEP 11, 170 (2018) doi:10.1007/JHEP11(2018)170 [arXiv:1807.07954 [hep-lat]].
  * [12] S. Borsanyi, Z. Fodor, J. Guenther, K. H. Kampert, S. D. Katz, T. Kawanai, T. G. Kovacs, S. W. Mages, A. Pasztor and F. Pittler, et al. Nature 539, no.7627, 69-71 (2016) doi:10.1038/nature20115 [arXiv:1606.07494 [hep-lat]].
  * [13] V. Baluni, Phys. Rev. D 19, 2227-2230 (1979) doi:10.1103/PhysRevD.19.2227
  * [14] J. E. Kim, Phys. Rept. 150, 1-177 (1987) doi:10.1016/0370-1573(87)90017-2
  * [15] F. R. Brown, F. P. Butler, H. Chen, N. H. Christ, Z. h. Dong, W. Schaffer, L. I. Unger and A. Vaccarino, Phys. Rev. Lett. 65, 2491-2494 (1990) doi:10.1103/PhysRevLett.65.2491
  * [16] C. Vafa and E. Witten, Nucl. Phys. B 234, 173-188 (1984) doi:10.1016/0550-3213(84)90230-X
  * [17] T. Hatsuda and T. Kunihiro, Phys. Rept. 247, 221-367 (1994) doi:10.1016/0370-1573(94)90022-1 [arXiv:hep-ph/9401310 [hep-ph]].
  * [18] M. Kobayashi and T. Maskawa, Prog. Theor. Phys. 44, 1422-1424 (1970) doi:10.1143/PTP.44.1422
  * [19] M. Kobayashi, H. Kondo and T. Maskawa, Prog. Theor. Phys. 45, 1955-1959 (1971) doi:10.1143/PTP.45.1955
  * [20] G. ’t Hooft, Phys. Rev. Lett. 37, 8-11 (1976) doi:10.1103/PhysRevLett.37.8
  * [21] G. ’t Hooft, Phys. Rev. D 14, 3432-3450 (1976) [erratum: Phys. Rev. D 18, 2199 (1978)] doi:10.1103/PhysRevD.14.3432
  * [22] K. Fukushima, K. Ohnishi and K. Ohta, Phys. Rev. C 63, 045203 (2001) doi:10.1103/PhysRevC.63.045203 [arXiv:nucl-th/0101062 [nucl-th]].
  * [23] P. Costa, M. C. Ruivo, C. A. de Sousa, H. Hansen and W. M. Alberico, Phys. Rev. D 79, 116003 (2009) doi:10.1103/PhysRevD.79.116003 [arXiv:0807.2134 [hep-ph]].
  * [24] Y. Jiang and P. Zhuang, Phys. Rev. D 86, 105016 (2012) doi:10.1103/PhysRevD.86.105016 [arXiv:1209.0507 [hep-ph]].
  * [25] Y. Jiang, T. Xia and P. Zhuang, Phys. Rev. D 93, no.7, 074006 (2016) doi:10.1103/PhysRevD.93.074006 [arXiv:1511.06466 [hep-ph]].
  * [26] T. Kunihiro and T. Hatsuda, Prog. Theor. Phys. 71, 1332 (1984) doi:10.1143/PTP.71.1332
  * [27] M. L. Bellac, doi:10.1017/CBO9780511721700
  * [28] T. Kunihiro, Nucl. Phys. B 351, 593-622 (1991) doi:10.1016/S0550-3213(05)80035-5
  * [29] G. A. Miller, B. M. K. Nefkens and I. Slaus, Phys. Rept. 194, 1-116 (1990) doi:10.1016/0370-1573(90)90102-8
  * [30] C. Bonati, M. D’Elia, M. Mariti, G. Martinelli, M. Mesiti, F. Negro, F. Sanfilippo and G. Villadoro, JHEP 03, 155 (2016) doi:10.1007/JHEP03(2016)155 [arXiv:1512.06746 [hep-lat]].
  * [31] Y. Aoki, S. Borsanyi, S. Durr, Z. Fodor, S. D. Katz, S. Krieg and K. K. Szabo, JHEP 06, 088 (2009) doi:10.1088/1126-6708/2009/06/088 [arXiv:0903.4155 [hep-lat]].
  * [32] S. Borsanyi et al. [Wuppertal-Budapest], J. Phys. Conf. Ser. 316, 012020 (2011) doi:10.1088/1742-6596/316/1/012020 [arXiv:1109.5032 [hep-lat]].
  * [33] H. T. Ding, F. Karsch and S. Mukherjee, Int. J. Mod. Phys. E 24, no.10, 1530007 (2015) doi:10.1142/S0218301315300076 [arXiv:1504.05274 [hep-lat]].
  * [34] A. Bazavov et al. [HotQCD], Phys. Lett. B 795, 15-21 (2019) doi:10.1016/j.physletb.2019.05.013 [arXiv:1812.08235 [hep-lat]].
  * [35] H. T. Ding, Nucl. Phys. A 1005, 121940 (2021) doi:10.1016/j.nuclphysa.2020.121940 [arXiv:2002.11957 [hep-lat]].
  * [36] A. Bazavov et al. [HotQCD], Phys. Rev. D 86, 094503 (2012) doi:10.1103/PhysRevD.86.094503 [arXiv:1205.3535 [hep-lat]].
  * [37] R. D. Pisarski and L. G. Yaffe, Phys. Lett. B 97, 110-112 (1980) doi:10.1016/0370-2693(80)90559-6
  * [38] D. J. Gross, R. D. Pisarski and L. G. Yaffe, Rev. Mod. Phys. 53, 43 (1981) doi:10.1103/RevModPhys.53.43
  * [39] G. Fejös and A. Patkos, Phys. Rev. D 105, no.9, 096007 (2022) doi:10.1103/PhysRevD.105.096007 [arXiv:2112.14903 [hep-ph]].
  * [40] B. B. Brandt, A. Francis, H. B. Meyer, O. Philipsen, D. Robaina and H. Wittig, JHEP 12, 158 (2016) doi:10.1007/JHEP12(2016)158 [arXiv:1608.06882 [hep-lat]].
  * [41] B. B. Brandt, O. Philipsen, M. Cè, A. Francis, T. Harris, H. B. Meyer and H. Wittig, PoS CD2018, 055 (2019) doi:10.22323/1.317.0055 [arXiv:1904.02384 [hep-lat]].