
The energymomentum tensor plays an important role in QCD thermodynamics. Its
expectation value contains information of the pressure and the energy density
as its diagonal part. Further properties like viscosity and specific heat can
be extracted from its correlation function. Recently a new method based on the
gradient flow was introduced to calculate the energymomentum tensor on the
lattice, and has been successfully applied to quenched QCD.
In this paper, we apply the gradient flow method to calculate the
energymomentum tensor in (2+1)flavor QCD. As the first application of the
method with dynamical quarks, we study at a single but fine lattice spacing
a=0.07 fm with heavy u and d quarks ($m_\pi/m_\rho=0.63$) and approximately
physical s quark. Performing simulations on lattices with Nt=16 to 4, the
temperature range of T=174697 MeV is covered. We find that the results of the
pressure and the energy density by the gradient flow method are consistent with
the previous results using the Tintegration method at T<280 MeV, while the
results show disagreement at T>350 MeV (Nt<8), presumably due to the smallNt
lattice artifact of $O((aT)^2)=O(1/N_t^2)$.
We also apply the gradient flow method to evaluate the chiral condensate
taking advantage of the gradient flow method that renormalized quantities can
be directly computed avoiding the difficulty of explicit chiral violation with
lattice quarks. We compute the renormalized chiral condensate in the MSbar
scheme at renormalization scale $\mu=2$ GeV with a high precision to study the
temperature dependence of the chiral condensate and its disconnected
susceptibility. Even with the Wilsontype quark action, we obtain the chiral
condensate and its disconnected susceptibility showing a clear signal of
pseudocritical temperature at T~190 MeV related to the chiral restoration
crossover.

We propose a general procedure for the detectorresponse correction
(including efficiency correction) of higher order cumulants observed by the
eventbyevent analysis in heavyion collisions. This method makes use of the
moments of the response matrix characterizing the property of a detector, and
is applicable to a wide variety of response matrices such as those having
nonbinomial responses and including the effects of ghost tracks. A procedure
to carry out the detectorresponse correction of realistic detectors is
discussed. In test analyses, we show that this method can successfully re
construct the cumulants of true distribution for various response matrices
including the one having multiplicitydependent efficiency.

The spatial distribution of the stress tensor around the quarkantiquark
($Q\bar{Q}$) system in SU(3) lattice gauge theory is studied. The YangMills
gradient flow plays a crucial role to make the stress tensor welldefined and
practically derivable from the numerical simulations on the lattice. Clear
evidence of the flux tube formation and its transverse structure are presented
from the stresstensor distribution. The linear confining behavior of the
$Q\bar{Q}$ potential at long distances is derived directly from the integral of
the local stress tensor.

We measure correlation functions of the nonperturbatively renormalized
energymomentum tensor in $N_f=2+1$ full QCD at finite temperature by applying
the gradient flow method both to the gauge and quark fields. Our main interest
is to study the conservation law of the energymomentum tensor and to test
whether the linear response relation is properly realized for the entropy
density. By using the linear response relation we calculate the specific heat
from the correlation function. We adopt the nonperturbatively improved Wilson
fermion and Iwasaki gauge action at a fine lattice spacing $=0.07$ fm. In this
paper the temperature is limited to a single value $T=232$ MeV. The $u$, $d$
quark mass is rather heavy with $m_\pi/m_\rho=0.63$ while the $s$ quark mass is
set to approximately its physical value.

Euclidean twopoint correlators of the energymomentum tensor (EMT) in SU(3)
gauge theory on the lattice are studied on the basis of the YangMills gradient
flow. The entropy density and the specific heat obtained from the twopoint
correlators are shown to be in good agreement with those from the onepoint
functions of EMT. These results constitute a first step toward the first
principle simulations of the transport coefficients with the gradient flow.

We derive formulas for the efficiency correction of cumulants with many
efficiency bins. The derivation of the formulas is simpler than the previously
suggested method, but the numerical cost is drastically reduced from the naive
method. From analytical and numerical analyses in simple toy models, we show
that the use of the averaged efficiency in the efficiency correction can lead
to wrong corrected values, which have larger deviation for higher order
cumulants. These analyses show the importance of carrying out the efficiency
correction without taking the average.

The photon production rate from the deconfined medium is analyzed with the
photon selfenergy constructed from the quark propagator obtained by the
numerical simulation on the quenched lattice for two values of temperature,
$T=1.5T_{\rm c}$ and $3T_{\rm c}$, above the critical temperature $T_{\rm c}$.
The photon selfenergy is calculated by the SchwingerDyson equation with the
lattice quark propagator and a vertex function determined so as to satisfy the
WardTakahashi identity. The obtained photon production rate exhibits a similar
behavior as the perturbative results at the energy of photons larger than
$0.5$~GeV.

We discuss properties and applications of factorial cumulants of various
particle numbers and for their mixed channels measured by the eventbyevent
analysis in relativistic heavyion collisions. After defining the factorial
cumulants for systems with multiparticle species, their properties are
elucidated. The uses of the factorial cumulants in the study of critical
fluctuations are discussed. We point out that factorial cumulants play useful
roles in understanding fluctuation observables when they have underlying
physics approximately described by the binomial distribution. As examples, we
suggest novel utilization methods of the factorial cumulants in the study of
the momentum cut and rapidity window dependences of fluctuation observables.

We study time evolution of critical fluctuations of conserved charges near
the QCD critical point in the context of relativistic heavy ion collisions. A
stochastic diffusion equation is employed in order to describe the diffusion
property of the critical fluctuation arising from the coupling of the order
parameter field to conserved charges. We show that the diffusion property gives
rise to a possibility of probing the early time fluctuations through the
rapidity window dependence of the secondorder cumulant and correlation
function of conserved charges. It is pointed out that their nonmonotonic
behaviors as functions of the rapidity interval are robust experimental signals
for the existence of the critical enhancement around the QCD critical point.

The energy density and the pressure of SU(3) gauge theory at finite
temperature are studied by direct lattice measurements of the renormalized
energymomentum tensor obtained by the gradient flow. Numerical analyses are
carried out with $\beta=6.287$$7.500$ corresponding to the lattice spacing
$a= 0.013$$0.061\,\mathrm{fm}$. The spatial (temporal) sizes are chosen to be
$N_s= 64$, $96$, $128$ ($N_{\tau}=12$, $16$, $20$, $22$, $24$) with the aspect
ratio, $5.33 \le N_s/N_{\tau} \le 8$. Double extrapolation, $a\rightarrow 0$
(the continuum limit) followed by $t\rightarrow 0$ (the zero flowtime limit),
is taken using the numerical data. Above the critical temperature, the
thermodynamic quantities are obtained with a few percent precision including
statistical and systematic errors. The results are in good agreement with
previous highprecision data obtained by using the integral method.

We calculate the energy gap (latent heat) and pressure gap between the hot
and cold phases of the SU(3) gauge theory at the first order deconfining phase
transition point. We perform simulations around the phase transition point with
the lattice size in the temporal direction $N_t=6,$ 8 and 12 and extrapolate
the results to the continuum limit. The energy density and pressure are
evaluated by the derivative method with nonperturabative anisotropy
coefficients. We find that the pressure gap vanishes at all values of $N_t$.
The spatial volume dependence in the latent heat is found to be small on large
lattices. Performing extrapolation to the continuum limit, we obtain $\Delta
\epsilon/T^4 = 0.75 \pm 0.17$ and $\Delta (\epsilon 3 p)/T^4 = 0.623 \pm
0.056.$ We also tested a method using the YangMills gradient flow. The
preliminary results are consistent with those by the derivative method within
the error.

We study temperature dependence of the topological susceptibility with the
$N_{f}=2+1$ flavors Wilson fermion. We have two major interests in this paper.
One is a comparison of gluonic and fermionic definitions of the topological
susceptibility. Two definitions are related by the chiral WardTakahashi
identity but their coincidence is highly nontrivial for the Wilson fermion. By
applying the gradient flow both for the gauge and quark fields we find a good
agreement of these two measurements. The other is a verification of a
prediction of the dilute instanton gas approximation at low temperature region
$T_{pc}< T<1.5T_{pc}$, for which we confirm the prediction that the topological
susceptibility decays with power $\chi_{t}\propto(T/T_{pc})^{8}$ for three
flavors QCD.

The energymomentum tensor and equation of state are studied in
finitetemperature (2+1)flavor QCD with improved Wilson quarks using the
method proposed by Makino and Suzuki based on the gradient flow. We find that
the results of the gradient flow are consistent with the previous results using
the $T$integration method at $T \mathrel{\rlap{\raise 0.511ex
\hbox{$<$}}{\lower 0.511ex \hbox{$\sim$}}} 280$ MeV ($N_t\mathrel{\rlap{\raise
0.511ex \hbox{$>$}}{\lower 0.511ex \hbox{$\sim$}}}10$), while a disagreement is
found at $T \mathrel{\rlap{\raise 0.511ex \hbox{$>$}}{\lower 0.511ex
\hbox{$\sim$}}} 350$ MeV ($N_t \mathrel{\rlap{\raise 0.511ex \hbox{$<$}}{\lower
0.511ex \hbox{$\sim$}}} 8$) presumably due to the small$N_t$ lattice artifact.
We also report on the results on the renormalized chiral condensate and its
disconnected susceptibility using the method of Hieda and Suzuki. The results
show a clear signal of the expected chiral restoration crossover even with
Wilsontype quarks which violate the chiral symmetry explicitly.

We study inmedium spectral properties of charmonia in the vector and
pseudoscalar channels at nonzero momenta on quenched lattices, especially
focusing on their dispersion relation and weight of the peak. We measure the
lattice Euclidean correlation functions with nonzero momenta on the anisotropic
quenched lattices and study the spectral functions with the maximum entropy
method. The dispersion relations of charmonia and the momentum dependence of
the weight of the peak are analyzed with the maximum entropy method together
with the errors estimated probabilistically in this method. We find significant
increase of the masses of charmonia in medium. It is also found that the
functional form of the charmonium dispersion relations is not changed from that
in the vacuum within the error even at $T\simeq1.6T_c$ for all the channels we
analyzed.

We analyze the production rate of photons from the thermal medium above the
deconfinement temperature with a quark propagator obtained from a lattice QCD
numerical simulation. The photonquark vertex is determined gaugeinvariantly,
so as to satisfy the WardTakahashi identity. The obtained photon production
rate shows a suppression compared to perturbative results.

We discuss two topics on the experimental measurements of fluctuation
observables in relativistic heavyion collisions. First, we discuss the effects
of the thermal blurring, i.e. the blurring effect arising from the experimental
measurement of fluctuations in momentum space as a proxy of the thermal
fluctuations defined in coordinate space, on higher order cumulants. Second, we
discuss the effect of imperfect efficiency of detectors on the measurement of
higher order cumulants. We derive effective formulas which can carry out the
correction of this effect for higher order cumulants based on the binomial
model.

We study the effect of thermal blurring caused by the use of (momentumspace)
rapidity as a proxy of coordinatespace rapidity in experimental measurements
of conserved charge fluctuations in relativistic heavy ion collisions. In
theoretical studies assuming statistical mechanics, calculated fluctuations are
those in a spatial volume. Experiments, on the other hand, can measure
fluctuations only in a momentumspace in the final state. In a standard
argument to compare experimental results for a momentum space with theoretical
studies for a coordinate space, rapidities of particles are implicitly regarded
as equivalent to their coordinatespace rapidity. We show that the relation of
two fluctuations is significantly altered by the existence of the thermal
motion, i.e. thermal blurring. We discuss that the thermal blurring can be
regarded as a part of the diffusion process, and the effect can be understood
by studying the rapidity window dependences of fluctuations. Centrality
dependence of the thermal blurring effect is also discussed.

Bulk fluctuations of conserved charges measured by eventbyevent analysis in
relativistic heavy ion collisions are observables which are believed to carry
significant amount of information on the hot medium created by the collisions.
Active studies have been done recently experimentally, theoretically, and on
the lattice. In particular, nonGaussianity of the fluctuations has acquired
much attention recently. In this review, we give a pedagogical introduction to
these issues, and survey recent developments in this field of research.
Starting from the definition of cumulants, basic concepts in fluctuation
physics, such as thermal fluctuations in statistical mechanics and time
evolution of fluctuations in diffusive systems, are described. Phenomena which
are expected to occur in finite temperature and/or density QCD matter and their
measurement by eventbyevent analyses are also elucidated.

We calculate the energy gap (latent heat) and pressure gap between the hot
and cold phases of the SU(3) gauge theory at the first order deconfining phase
transition point. We perform simulations around the phase transition point with
the lattice size in the temporal direction Nt=6, 8 and 12 and extrapolate the
results to the continuum limit. We also investigate the spatial volume
dependence. The energy density and pressure are evaluated by the derivative
method with nonperturabative anisotropy coefficients. We adopt a multipoint
reweighting method to determine the anisotropy coefficients. We confirm that
the anisotropy coefficients approach the perturbative values as Nt increases.
We find that the pressure gap vanishes at all values of Nt when the
nonperturbative anisotropy coefficients are used. The spatial volume
dependence in the latent heat is found to be small on large lattices.
Performing extrapolation to the continuum limit, we obtain $ \Delta
\epsilon/T^4 = 0.75 \pm 0.17 $ and $ \Delta (\epsilon 3 p)/T^4 = 0.623 \pm
0.056.$

We derive formulas which connect cumulants of particle numbers observed with
efficiency losses with the original ones based on the binomial model. These
formulas can describe the case with multiple efficiencies in a compact form.
Compared with the presently suggested ones based on factorial moments, these
formulas would drastically reduce the numerical cost for efficiency corrections
when the order of the cumulant and the number of different efficiencies are
large. The efficiency correction with realistic $p_T$dependent efficiency
would be carried out with the aid of these formulas.

We study the parametrization of lattice spacing and thermodynamics of SU(3)
gauge theory on the basis of the YangMills gradient flow on fine lattices. The
lattice spacing of the Wilson gauge action is determined over a wide range
$6.3\le\beta\le7.5$ with high accuracy. The measurements of the flow time and
lattice spacing dependences of the expectation values of the energymomentum
tensor are performed on fine lattices.

A parametrization of the lattice spacing ($a$) in terms of the bare coupling
($\beta$) for the SU(3) YangMills theory with the Wilson gauge action is
given in a wide range of~$\beta$. The YangMills gradient flow with respect to
the flow time~$t$ for the dimensionless observable, $t\frac{d}{dt}t^2\langle
E(t)\rangle$, is utilized to determine the parametrization. With fine lattice
spacings ($6.3\le\beta\le7.5$) and large lattice volumes ($N_{\rm
s}=64$$128$), the discretization and finitevolume errors are significantly
reduced to the same level as the statistical error.

A novel method to study the bulk thermodynamics in lattice gauge theory is
proposed on the basis of the YangMills gradient flow with a fictitious time t.
The energy density (epsilon) and the pressure (P) of SU(3) gauge theory at
fixed temperature are calculated directly on 32^3 x (6,8,10) lattices from the
thermal average of the welldefined energymomentum tensor (T_{mu nu}^R(x))
obtained by the gradient flow. It is demonstrated that the continuum limit can
be taken in a controlled manner from the tdependence of the flowed data.

The dilepton production rate from the deconfined medium is analyzed with the
photon selfenergies constructed from quark propagators obtained by lattice
numerical simulation for two values of temperature $T=1.5T_{\rm c}$ and
$3T_{\rm c}$ above the critical temperature $T_{\rm c}$. The photon selfenergy
is calculated by the SchwingerDyson equation with the lattice quark propagtor
and a vertex function determined so as to satisfy the WardTakahashi identity.
The obtained dilepton production rate at zero momentum exhibits divergences
reflecting van Hove singularity, and is significantly enhanced around
$\omega\simeq T$ compared with the rate obtained by the perturbative analysis.

We study the rapidity window dependences of higher order cumulants of
conserved charges observed in relativistic heavy ion collisions. The time
evolution and the rapidity window dependence of the nonGaussian fluctuations
are described by the diffusion master equation. Analytic formulas for the time
evolution of cumulants in a rapidity window are obtained for arbitrary initial
conditions. We discuss that the rapidity window dependences of the nonGaussian
cumulants have characteristic structures reflecting the nonequilibrium
property of fluctuations, which can be observed in relativistic heavy ion
collisions with the present detectors. It is argued that various information on
the thermal and transport properties of the hot medium can be revealed
experimentally by the study of the rapidity window dependences, especially by
the combined use, of the higher order cumulants. Formulas of higher order
cumulants for a probability distribution composed of subprobabilities, which
are useful for various studies of nonGaussian cumulants, are also presented.