aainstitutetext: Universität Bielefeld,
Universitätsstraße 25, 33615 Bielefeld, Germany

# The Chiral Separation Effect from lattice QCD at the physical point

B. B. Brandt a G. Endrődi a E. Garnacho-Velasco a and G. Markó
<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

In this paper we study the Chiral Separation Effect by means of first-
principles lattice QCD simulations. For the first time in the literature, we
determine the continuum limit of the associated conductivity using 2+1 flavors
of dynamical staggered quarks at physical masses. The results reveal a
suppression of the conductivity in the confined phase and a gradual
enhancement toward the perturbative value for high temperatures. In addition
to our dynamical setup, we also investigate the impact of the quenched
approximation on the conductivity, using both staggered and Wilson quarks.
Finally, we highlight the relevance of employing conserved vector and
anomalous axial currents in the lattice simulations.

##  1 Introduction

The non-trivial topological structure of the QCD vacuum manifests itself in a
variety of fascinating aspects that can be studied theoretically and probed
experimentally. In this context, anomalous transport phenomena have attracted
the most attention recently. These effects appear due to the interplay between
quantum anomalies and electromagnetic fields or vorticity, explaining their
intimate relation to the topological nature of QCD and representing the impact
of event-by-event CP-violation in QCD Kharzeev:2007jp . The prime example
among these phenomena is the Chiral Magnetic Effect (CME): the generation of a
vector current due to magnetic fields and a chiral imbalance Fukushima:2008xe
. In the last decade, major experimental campaigns have been launched in order
to measure the CME in different setups ranging from condensed matter systems
Li:2014bha  to heavy ion collision experiments STAR:2013ksd ; STAR:2014uiw ;
STAR:2021mii .

Another anomalous transport phenomenon is the Chiral Separation Effect (CSE):
the generation of an axial current in a dense and magnetized system Son:2004tq
; Metlitski:2005pr . Although highly analogous to the CME in its formulation,
the equilibrium interpretation of the two effects turns out to be quite
different. Together with the CME, it is expected to form the Chiral Magnetic
Wave (CMW), a collective excitation emerging in a magnetized environment
Kharzeev:2010gd , which is also the subject of intense experimental searches
Huang:2015oca . The key parameter, describing the strength of the effect is
the conductivity $C_{\text{CSE}}$, to be defined in detail below.

The analytical approaches to assess the CSE span a broad range and include for
example chiral kinetic theory Avdoshkin:2017cqp , effective models of QCD
Gorbar:2011ya , holography Jimenez-Alba:2014pea  or the field correlator
method Zubkov:2023vvb . The first calculations were performed for massless
quarks in the absence of gluonic interactions, giving
$C_{\text{CSE}}=1/(2\pi^{2})$ Son:2004tq . It was later recognized that this
value is not fixed by the axial anomaly but rather depends on the quark mass
$m$ as well as the temperature $T$. In fact, as we show in App. A, for non-
interacting quarks, an analytical treatment of this problem gives

$C_{\text{CSE}}^{\text{\,free}}=\frac{1}{2\pi^{2}}\int_{0}^{\infty}\differential
p\,\left[1+\cosh(\sqrt{p^{2}+(m/T)^{2}})\right]^{-1}\,,$ (1)

in agreement with Son:2004tq ; Metlitski:2005pr .

The resulting curve is shown in Fig. 1, revealing a suppression of the
conductivity for heavy quarks (or, equivalently, low temperatures), see also
Refs. Sheng:2017lfu ; Lin:2018aon .

Figure 1: The CSE conductivity as a function of $m/T$ for non-interacting
quarks. The marked points correspond to $m/T=0,1,4$, which will be the values
we compare to in Sec. 3.1.

Given the sensitivity of the CSE to the mass and the temperature, it is
natural to expect that the impact of gluonic interactions on $C_{\text{CSE}}$
will also be substantial (see also Ref. Gorbar:2013upa ). In the
phenomenologically interesting region around the finite temperature QCD
crossover Aoki:2006we , these interactions are non-perturbative. The most
successful non-perturbative tool, based on the first principles of the
underlying quantum field theory, is lattice QCD. However, at nonzero quark
density – which is required for the discussion of the CSE – standard lattice
simulations break down due to the complex action problem, see e.g. the recent
review Nagata:2021ugx . For this reason, many of the existing lattice
investigations Puhr:2016kzp ; Buividovich:2020gnl  were not performed in full
QCD, but either in the two-color version of it, or in the so-called quenched
approximation – both of which are free of the complex action problem. Even
with this caveat in mind, some of the findings in the literature gave
different predictions about whether QCD interactions at intermediate
temperatures affect the value of $C_{\text{CSE}}$ Buividovich:2020gnl  or not
Puhr:2016kzp . The lattice regularization has also been employed to study the
CSE for free fermions Khaidukov:2017exf ; Buividovich:2013hza  and with
classical-statistical real-time simulations Muller:2016jod ; Mace:2016shq .

In this work, we demonstrate for the first time how the CSE conductivity can
be constructed using a Taylor-expansion around zero density, thereby enabling
its treatment in full QCD. We determine $C_{\text{CSE}}$ using dynamical
staggered quarks with physical masses and carry out its continuum
extrapolation for a broad range of temperatures in the transition region as
well as deep in the confined phase. We shed further light on the impact of the
quenched approximation and determine $C_{\text{CSE}}$ in this setup using both
staggered and Wilson quarks. In this setup, we also discuss the effect of
increasing the quark masses away from the physical point. Finally, we also
determine the CSE conductivity in the absence of gluonic interactions for both
the staggered and the Wilson fermion discretizations – this provides a useful
benchmark of the appropriate lattice currents.

The present effort will be useful for studying further anomalous transport
effects like the CME, for which most lattice simulations so far
Buividovich:2009wi ; Buividovich:2009my ; Yamamoto:2011ks ;
Buividovich:2013hza ; Bali:2014vja ; Astrakhantsev:2019zkr  were either based
on indirect approaches or are yet to be performed in full QCD at the physical
point.

##  2 Conductivities and lattice methods

As introduced above, the chiral separation effect amounts to the generation of
an axial current $J_{\nu 5}$ in the presence of a background magnetic field
and nonzero quark density. The latter is parameterized by a chemical potential
$\mu$ that couples to the corresponding conserved number density. The magnetic
field $B$ is considered to be constant and homogeneous and, without loss of
generality, pointing along the third spatial direction. Moreover, the magnetic
field is measured in units of the elementary electric charge $e>0$, so that we
can work with the renormalization group invariant combination $eB$.

###  2.1 Currents and chemical potentials

For a single quark flavor, the definition of the axial current and the
chemical potential is unambiguous. For several quark flavors, like in full
QCD, there are different options which – as we will find below – give
different results for the strength of the effect. One may consider the axial
current

$J_{\nu
5}^{i}=\sum_{f}c^{i}_{f}\,\frac{T}{V}\int\differential^{4}x\,\bar{\psi}_{f}(x)\gamma_{\nu}\gamma_{5}\psi_{f}(x)\,,$
(2)

with different quantum numbers $i=\mathcal{B},\mathcal{Q}$ or $\mathcal{L}$ in
flavor space,

$c_{f}^{\mathcal{B}}=\frac{1}{3}\,,\qquad
c_{f}^{\mathcal{Q}}=\frac{q_{f}}{e}\,,\qquad
c_{f}^{\mathcal{L}}=\begin{cases}1/3\,,&f=u,d\,,\\\
0\,,&\textmd{otherwise}\,,\\\ \end{cases}$ (3)

where $f=u,d,s,\ldots$ labels the quark flavors and $q_{f}$ are the
corresponding electric charges. The currents (2) are rendered intensive by
dividing by the Euclidean four-volume $V/T$. The above choices correspond to
the currents coupling to baryon number $\mathcal{B}$, electric charge
$\mathcal{Q}$ or the so-called light baryon number $\mathcal{L}$
Borsanyi:2012cr . The latter serves to mimic heavy-ion collision setups.

Analogously to the axial currents, the quark density can also be controlled
either by a baryon chemical potential $\mu_{\mathcal{B}}$, an electric charge
chemical potential $\mu_{\mathcal{Q}}$ or a light baryonic chemical potential
$\mu_{\mathcal{L}}$, which couples to the corresponding vector currents,

$J_{\nu}^{i}=\sum_{f}c_{f}^{i}\,\frac{T}{V}\int\differential^{4}x\,\bar{\psi}_{f}(x)\gamma_{\nu}\psi_{f}(x)\,.$
(4)

The chemical potentials for the individual quark flavors are
$\mu_{f}=\mu_{\mathcal{B}}/3$ in the baryonic and
$\mu_{f}=\mu_{\mathcal{Q}}\,q_{f}/e$ in the charge case, while they are set as
$\mu_{u}=\mu_{d}=\mu_{\mathcal{L}}/3$, $\mu_{f}=0$ ($f\neq u,d$) for the light
baryon chemical potential.

The most common choice in the flavor space corresponds to the quantum number
$\mathcal{Q}$ for the current and $\mathcal{B}$ for the chemical potential. In
this case, the CSE conductivity is defined via the leading order-behavior for
the expectation value

$\langle
J^{\mathcal{Q}}_{35}\rangle=C_{\text{CSE}}\,\mu_{\mathcal{B}}\,eB+\mathcal{O}(\mu_{\mathcal{B}}^{3},B^{3})\,.$
(5)

Here we use the short notation $C_{\text{CSE}}\equiv
C_{\textmd{CSE}}^{\mathcal{Q}\mathcal{B}}$ and we analogously define
$C_{\textmd{CSE}}^{ij}$, where the first superscript
$i\in\\{\mathcal{B},\mathcal{Q},\mathcal{L}\\}$ corresponds to the current and
the second one $j\in\\{\mathcal{B},\mathcal{Q},\mathcal{L}\\}$ to the chemical
potential.

Full lattice QCD simulations suffer from the complex action problem at
$\mu_{\mathcal{B}}\neq 0$, $\mu_{\mathcal{Q}}\neq 0$ or $\mu_{\mathcal{L}}\neq
0$. To circumvent this issue, we consider the Taylor-expansion of the current
expectation value in the chemical potential. To leading order, this gives the
first chemical potential derivative of the CSE current, which, using Eq. (5),
yields

$\evaluated{\partialderivative{\langle
J_{35}^{\mathcal{Q}}\rangle}{\mu_{\mathcal{B}}}}_{\mu_{\mathcal{B}}=0}=C_{\text{CSE}}\,eB\,,$
(6)

and analogously for the other flavor combinations. This, leading-order Taylor
expansion of the current only requires simulations at zero density, free of
the sign problem. Finally, $C_{\text{CSE}}$ can be extracted via numerical
differentiation with respect to $eB$.

For each choice of the flavor structures, there is an overall proportionality
constant involving the quark baryon numbers or quark charges and the number of
colors $N_{c}=3$, which can be factored out of the conductivity. Since the
coupling to the magnetic field always occurs via the quark charges, these
overall factors read

$\displaystyle C_{\rm dof}$ $\displaystyle\equiv C_{\rm
dof}^{\mathcal{Q}\mathcal{B}}=C_{\rm
dof}^{\mathcal{B}\mathcal{Q}}=\frac{N_{c}}{3}\sum_{f}\left(\frac{q_{f}}{e}\right)^{2},\quad\quad
C_{\rm
dof}^{\mathcal{Q}\mathcal{Q}}=N_{c}\sum_{f}\left(\frac{q_{f}}{e}\right)^{3},$
$\displaystyle C_{\rm dof}^{\mathcal{B}\mathcal{B}}$
$\displaystyle=\frac{N_{c}}{9}\sum_{f}\frac{q_{f}}{e}\,,\qquad C_{\rm
dof}^{\mathcal{Q}\mathcal{L}}=\frac{N_{c}}{3}\sum_{f=u,d}\left(\frac{q_{f}}{e}\right)^{2}\,,\qquad
C_{\rm
dof}^{\mathcal{B}\mathcal{L}}=\frac{N_{c}}{9}\sum_{f=u,d}\frac{q_{f}}{e}\,,$
(7)

and similarly for the remaining ones. For massless, non-interacting quarks,
these factors appear as overall normalizations in $C_{\text{CSE}}$. (Note that
$C_{\rm dof}^{\mathcal{B}\mathcal{B}}=0$ for three light quarks.) Note,
however, that for massive, non-degenerate quarks each flavor contributes
differently, depending on the choice of the quantum numbers. In the
interacting case, disconnected diagrams are also different for $\mathcal{B}$,
$\mathcal{Q}$ or $\mathcal{L}$. This implies that the individual
$C_{\textmd{CSE}}^{ij}$ are not equivalent despite scaling out the above
factors. Since these overall factors can always be restored, from this point
we rescale all our results by $C_{\text{dof}}$ unless explicitly stated
otherwise.

Throughout this paper, we are working in Euclidean space-time. Taking into
account the relation between the Minkowski (indicated by the superscript $\rm
M$) and Euclidean Dirac matrices ($\gamma_{0}^{\rm M}=\gamma_{4}$,
$\gamma_{i}^{\rm M}=i\gamma_{i}$), we find that the Minkowski-space observable
reads

$\real\evaluated{\partialderivative{\langle
J_{35}^{\mathcal{Q}}\rangle}{\mu_{\mathcal{B}}}}_{\mu_{\mathcal{B}}=0}^{\rm
M}=-\imaginary\evaluated{\partialderivative{\langle
J_{35}^{\mathcal{Q}}\rangle}{\mu_{\mathcal{B}}}}_{\mu_{\mathcal{B}}=0}\,.$ (8)

Therefore we consider, for the rest of the text, the imaginary part of the
Euclidean two-point function calculated on the lattice. The real part of the
Euclidean observable was checked to be consistent with zero.

###  2.2 Lattice setup – staggered quarks

We begin by describing the setup with rooted staggered quarks, which we used
to study the CSE in full dynamical QCD with up, down and strange quark
flavors. Here, the partition function $\mathcal{Z}$ is written using the
Euclidean path integral over the gluon links $U$ as

$\mathcal{Z}=\int\mathcal{D}U\exp[-\beta S_{g}]\,\prod_{f}\quantity[\det
M_{f}(U,q_{f},m_{f})]^{1/4},$ (9)

with $\beta=6/g^{2}$ denotes the inverse gauge coupling and $m_{f}$ the quark
masses for each flavor $f=u,d,s$. In Eq. (9), $S_{g}$ is the gluonic action,
which in our setup is the tree-level Symanzik action, and $M_{f}$ is the
massive staggered Dirac operator with twice stout-smeared links. The Dirac
operator contains the quark charges $q_{u}/2=-q_{d}=-q_{s}=e/3$. The quark
masses are tuned to the physical point as a function of the lattice spacing
$a$ Borsanyi:2010cj . The simulations are carried out using four-dimensional
lattice geometries with $N_{s}$ spatial and $N_{t}$ temporal points. The
physical spatial volume is given by $V=L^{3}=(aN_{s})^{3}$ and the temperature
by $T=(aN_{t})^{-1}$. The lattice sites are labeled by the four integers
$n=\\{n_{1},n_{2},n_{3},n_{4}\\}$ and $\hat{\nu}$ denotes the unit vector in
the $\nu$ direction.

The magnetic field $B$ is included as a classical background field, so we do
not consider dynamical photons in our setup. The electromagnetic potential
enters the Dirac operator for the quark flavor $f$ in the same way as the
$\mathrm{SU}(3)$ links, i.e. as a parallel transporter between two lattice
points $u_{f\nu}=\exp(iaq_{f}A_{\nu})$. These electromagnetic links are chosen
such that they represent a homogeneous magnetic field in the $x_{3}$
direction. The flux of the field is quantized as $eB=6\pi N_{b}/(aN_{s})^{2}$
with the integer flux quantum $N_{b}\in\mathds{Z}$ Bali:2011qj . Moreover, the
chemical potential is also included in the exponential form, multiplying the
temporal links by $\exp(\pm a\mu_{f})$. Even though our simulations are at
$\mu_{f}=0$, we will need this form to calculate the Taylor-expansion
coefficients.

Within the staggered formalism, the quark fields $\psi_{f}$ are transformed in
coordinate and Dirac space to the so-called staggered fields $\chi_{f}$ in
order to partially diagonalize the Dirac operator. To proceed, we need to
express the vector and axial currents in terms of $\chi_{f}$. The prescription
that maintains a conserved vector current and an anomalous (flavor-singlet)
axial current are the point-split bilinears Sharatchandra:1981si ,

$J_{\nu}^{f}(n)=\bar{\chi}_{f}(n)\Gamma^{f}_{\nu}(n,m)\chi_{f}(m)\,,\qquad
J^{f}_{\nu 5}(n)=\bar{\chi}_{f}(n)\Gamma^{f}_{\nu 5}(n,m)\chi_{f}(m)\,,$ (10)

involving the staggered counterparts of the Dirac matrices Durr:2013gp ,

$\displaystyle\Gamma^{f}_{\nu}(n,m)$
$\displaystyle=\dfrac{\eta_{\nu}(n)}{2}\left[U_{\nu}(n)u_{f\nu}(n)\,e^{a\mu_{f}\delta_{\nu
4}}\delta_{n+\hat{\nu},m}+U^{\dagger}_{\nu}(n-\hat{\nu})u^{*}_{f\nu}(n-\hat{\nu})\,e^{-a\mu_{f}\delta_{\nu
4}}\delta_{n-\hat{\nu},m}\right],$ $\displaystyle\Gamma^{f}_{\nu 5}$
$\displaystyle=\dfrac{1}{3!}\sum_{\rho,\alpha,\beta}\epsilon_{\nu\rho\alpha\beta}\,\Gamma^{f}_{\rho}\Gamma^{f}_{\alpha}\Gamma^{f}_{\beta}\,.$
(11)

Here, $\eta_{\nu}(n)=(-1)^{\sum_{\rho<\nu}n_{\rho}}$ are the staggered phases
and $\epsilon_{\nu\rho\alpha\beta}$ the totally antisymmetric four-index
tensor with the convention $\epsilon_{1234}=+1$. We mention that these Dirac
structures depend explicitly on the links as well as on the magnetic field and
the chemical potential. In particular, note that $\Gamma^{f}_{35}$ involves
hoppings in the temporal direction and thus depends on $\mu_{f}$ (but not on
the chemical potentials for the other quark flavors $\mu_{f^{\prime}}$ with
$f\neq f^{\prime}$).

With these definitions, we can now perform the Grassmann integral over the
staggered fields, giving the expectation value for the (volume-averaged) axial
current,

$\langle
J_{35}^{\mathcal{Q}}\rangle=\frac{T}{V}\frac{1}{4}\sum_{f}\frac{q_{f}}{e}\left\langle\Tr\left(\Gamma^{f}_{35}M_{f}^{-1}\right)\right\rangle\,,$
(12)

where $\Tr$ refers to a trace in color space and a summation over the lattice
coordinates. The factor $1/4$ results from rooting. Next, we can calculate the
derivative (6) required for the CSE conductivity. This derivative generates
the usual disconnected and connected terms, together with an additional
tadpole term arising due to the derivative of $\Gamma^{f}_{35}$ with respect
to $\mu_{\mathcal{B}}$,

$\begin{split}C_{\rm CSE}\,eB=\evaluated{\partialderivative{\langle
J_{35}^{\mathcal{Q}}\rangle}{\mu_{\mathcal{B}}}}_{\mu_{\mathcal{B}}=0}=\dfrac{T}{V}\Bigg{[}\dfrac{1}{16}&\sum_{f,f^{\prime}}c_{f}^{\mathcal{Q}}c_{f^{\prime}}^{\mathcal{B}}\expectationvalue{\text{Tr}\quantity(\Gamma^{f}_{35}M_{f}^{-1})\text{Tr}\quantity(\Gamma^{f^{\prime}}_{4}M_{f^{\prime}}^{-1})}\\\
-\dfrac{1}{4}&\sum_{f}c_{f}^{Q}c_{f}^{\mathcal{B}}\expectationvalue{\text{Tr}\quantity(\Gamma^{f}_{35}M_{f}^{-1}\Gamma^{f}_{4}M_{f}^{-1})}\\\
+\dfrac{1}{4}&\sum_{f}c_{f}^{\mathcal{Q}}c_{f}^{\mathcal{B}}\expectationvalue{\text{Tr}\quantity(\dfrac{\partial\Gamma^{f}_{35}}{\partial\mu_{f}}M_{f}^{-1})}\Bigg{]}\,,\end{split}$
(13)

where the flavor coefficients (3) entered. We emphasize that the tadpole term
– the last line in Eq. (13) – appears due to the exponential form of
introducing the chemical potential. This form is required to maintain gauge
invariance on the lattice and to avoid chemical potential-dependent
ultraviolet divergences Hasenfratz:1983ba . The tadpole contribution is
analogous to the one that appears in the calculation of quark number
susceptibilities for staggered quarks, see e.g. Ref. Endrodi:2011gv . We also
mention that in deriving (13) we exploited the charge conjugation symmetry of
the $\mu_{\mathcal{B}}=0$ system, i.e. that the expectation value of the
baryon density and that of the axial current vanish at $\mu_{\mathcal{B}}=0$
(we refer back to this point in Sec. 3.2).

The expectation values in Eq. (13) are to be evaluated at
$\mu_{\mathcal{B}}=0$ but nonzero magnetic field $eB$. The conductivities
$C_{\rm CSE}^{ij}$ for the other flavor quantum numbers can be calculated
analogously, and differ from Eq. (13) by factors of
$c_{f}^{i}c_{f^{\prime}}^{j}$ in the disconnected and $c_{f}^{i}c_{f}^{j}$ in
the connected and tadpole terms under the flavor sums.

To enable a direct comparison to existing lattice results in the literature,
we also consider QCD in the quenched approximation. This amounts to dropping
the fermion determinant from the partition function (9), while leaving the
fermionic operators in the observables unchanged. This results in a non-
unitary theory where quarks behave differently in operators as in loop
diagrams. It is a commonly employed approximation that simplifies simulation
algorithms significantly. For completeness, in the quenched case we employ
both the staggered and the Wilson discretization of fermions.

Finally, we also calculate $C_{\text{CSE}}$ in the absence of gluonic
interactions, both for Wilson and for staggered quarks. To this end, we
calculated the eigenvalues and eigenvectors of the staggered Dirac operator
exactly and constructed the necessary traces from these (see App. B), while
for Wilson fermions we relied on stochastic techniques to estimate the traces.

###  2.3 Lattice setup – Wilson quarks

Next we describe our setup involving Wilson quarks. This discretization we
only consider for free fermions and in the quenched approximation –
simulations with non-degenerate light Wilson quarks (differing in their
electric charges) would be a computationally much more challenging setup. In
this case, the currents are constructed from the bispinor fields $\psi_{f}$,

$J_{\nu}^{f}(n)=\bar{\psi}_{f}(n)\Gamma^{f}_{\nu}(n,m)\psi_{f}(m)\,,\qquad
J^{f}_{\nu 5}(n)=\bar{\psi}_{f}(n)\Gamma^{f}_{\nu 5}(n,m)\psi_{f}(m)\,.$ (14)

Here the Dirac structure is the same as in the continuum, (apart from the
Wilson term), but the conserved vector and anomalous axial currents again
involve a point-splitting, just as for staggered quarks Karsten:1980wd

$\displaystyle\Gamma^{f}_{\nu}(n,m)$
$\displaystyle=\frac{1}{2}\Big{[}(\gamma_{\nu}-r)U_{\nu}(n)u_{\nu}(n)\,e^{a\mu_{f}\delta_{\nu
4}}\delta_{m,n+\hat{\nu}}$
$\displaystyle\qquad\quad+(\gamma_{\nu}+r)U^{\dagger}_{\nu}(n-\hat{\nu})u^{*}_{\nu}(n-\hat{\nu})\,e^{-a\mu_{f}\delta_{\nu
4}}\delta_{m,n-\hat{\nu}}\Big{]}\,,$ $\displaystyle\Gamma^{f}_{\nu 5}(n,m)$
$\displaystyle=\frac{1}{2}\Big{[}\gamma_{\nu}\gamma_{5}U_{\nu}(n)u_{\nu}(n)\,e^{a\mu_{f}\delta_{\nu
4}}\delta_{m,n+\hat{\nu}}$
$\displaystyle\qquad\quad+\gamma_{\nu}\gamma_{5}U^{\dagger}_{\nu}(n-\hat{\nu})u^{*}_{\nu}(n-\hat{\nu})\,e^{-a\mu_{f}\delta_{\nu
4}}\delta_{m,n-\hat{\nu}}\Big{]}\,,$

where $r$ is the coefficient of the usual Wilson term in the Dirac operator
(taken to be one in our setup). Notice the presence of the Wilson term in the
vector current and its absence in the axial current111The former follows
because gauge invariance requires the chemical potential to enter as parallel
transporter in all hopping terms in the action, and the vector current can be
written as the $\mu_{f}$-derivative of the action. The latter is quite
nontrivial: the axial anomaly equation is satisfied by this choice because the
terms proportional to $r$ (breaking the $\mathrm{U}_{A}(1)$ symmetry
explicitly) approach exactly the topological term
$\propto\epsilon_{\nu\rho\alpha\beta}G_{\nu\rho}G_{\alpha\beta}$ in the
continuum limit, where $G_{\alpha\beta}$ is the gluon field strength tensor
Karsten:1980wd ..

Using the above definitions, the expectation value of the axial current reads,

$\langle
J_{35}^{\mathcal{Q}}\rangle=\frac{T}{V}\sum_{f}\frac{q_{f}}{e}\left\langle\Tr\left(\Gamma^{f}_{35}M_{f}^{-1}\right)\right\rangle\,.$
(15)

Since $\Gamma^{f}_{35}$ this time only involves hoppings in the $3$ direction,
it is independent of the chemical potential. Thus, in comparison to the
staggered case, our observable takes a simpler form,

$\begin{split}C_{\rm CSE}\,eB=\evaluated{\partialderivative{\langle
J^{\mathcal{Q}}_{35}\rangle}{\mu_{\mathcal{B}}}}_{\mu_{\mathcal{B}}=0}=\dfrac{T}{V}\Bigg{[}&\sum_{f,f^{\prime}}c_{f}^{\mathcal{Q}}c_{f^{\prime}}^{\mathcal{B}}\expectationvalue{\text{Tr}\quantity(\Gamma^{f}_{35}M_{f}^{-1})\text{Tr}\quantity(\Gamma^{f^{\prime}}_{4}M_{f^{\prime}}^{-1})}\\\
-&\sum_{f}c_{f}^{\mathcal{Q}}c_{f}^{\mathcal{B}}\expectationvalue{\text{Tr}\quantity(\Gamma^{f}_{35}M_{f}^{-1}\Gamma^{f}_{4}M_{f}^{-1})}\Bigg{]}\,,\end{split}$
(16)

involving only disconnected and connected terms. The conductivities for the
other flavor quantum numbers follow similarly.

In the literature, it is also customary to consider a local vector current

$J_{\nu}^{f,\rm loc}(n)=\bar{\psi}_{f}(n)\gamma_{\nu}\psi_{f}(n)$ (17)

instead. Even though it is not conserved, it has the same quantum numbers as
$J^{f}_{\nu}$ and is often employed in Wilson fermion simulations, such as for
the study of various aspects of hadron physics at zero and non-zero
temperature, for instance. It was also used to study the CME in quenched and
dynamical QCD Yamamoto:2011ks . It is also possible to introduce a chemical
potential $\mu_{\mathcal{B}}^{\rm loc}$ that couples to these local currents
in the action. The analogue of (16) now reads

$\begin{split}C_{\rm CSE}^{\mathcal{Q}\mathcal{B},\rm
loc}\,eB=\evaluated{\partialderivative{\langle
J^{\mathcal{Q}}_{35}\rangle}{\mu^{\rm loc}_{\mathcal{B}}}}_{\mu^{\rm
loc}_{\mathcal{B}}=0}=\dfrac{T}{V}\Bigg{[}&\sum_{f,f^{\prime}}c_{f}^{\mathcal{Q}}c_{f^{\prime}}^{\mathcal{B}}\expectationvalue{\text{Tr}\quantity(\Gamma^{f}_{35}M_{f}^{-1})\text{Tr}\quantity(\gamma_{4}M_{f^{\prime}}^{-1})}\\\
-&\sum_{f}c_{f}^{\mathcal{Q}}c_{f}^{\mathcal{B}}\expectationvalue{\text{Tr}\quantity(\Gamma^{f}_{35}M_{f}^{-1}\gamma_{4}M_{f}^{-1})}\Bigg{]}\,.\end{split}$
(18)

This introduction of a chemical potential is a naive extension of the
continuum formulation to the discretized theory and it is well known that this
leads to new ultraviolet divergences Hasenfratz:1983ba . Therefore we use this
setup with the non-conserved, local current merely for comparison, and we will
test the impact of these ultraviolet divergences on the observable (18) below.

##  3 Results

To estimate the traces appearing in Eqs. (13), (16) and (18), we use the
standard noisy estimator technique.

We performed our numerical simulations using $100$ Gaussian noise vectors,
which was found to be sufficient to reliably calculate the observables. For
free staggered fermions, we calculated the traces using the exact eigensystem
instead, see App. B for details.

To determine $C_{\text{CSE}}$, we calculated the current derivative for
different magnetic fields, and extracted the coefficient from its slope with
respect to $eB$. We considered a linear fit function with no offset term: a
one parameter fit whose optimal value is found via a usual $\chi^{2}$
minimization method. The error analysis for each simulation is performed using
the jackknife method with 10 bins, which is also used to propagate the error
to the fit. This is what we refer to as the statistical error of
$C_{\text{CSE}}$. Furthermore, we repeat the fit, successively eliminating
data points at the largest value of $eB$, until we are left with only one
point. The largest difference in the slopes between the original fit and the
different repetitions is what we consider the systematical error of
$C_{\text{CSE}}$.

Figure 2: Derivative of the axial current with respect to the chemical
potential as a function of the magnetic field for a $24^{3}\times 6$ lattice
with staggered quarks for free fermions (left) and in full QCD with $2+1$
flavors and physical quark masses at $T=305$ MeV (right). The dashed line
represents the optimal fit, while the band is our estimate for its
uncertainty, obtained by adding the statistical and systematical errors in
quadrature.

In Fig. 2 we show an example of the obtained results and the fits, both for
free fermions and full QCD with staggered quarks. In these plots, we can see
that the expected linear behavior is confirmed. Now we turn to the precise
analysis of $C_{\text{CSE}}$ in different situations.

###  3.1 Free quarks

Since the free case is analytically solvable, we can use this simple situation
as a check of our setup. We consider one color-singlet fermion with charge $q$
and mass $m$ – this implies that the flavor quantum numbers are irrelevant and
it is sufficient to treat a single chemical potential $\mu$ and axial current
$J_{35}$. The overall proportionality constant is thus $C_{\rm
dof}=(q/e)^{2}$, which is used to normalize the results.

In Fig. 3 we show the results for free staggered fermions at different $m/T$
values. To compare to the analytic formula, both the continuum limit ($a\to
0$), as well as the thermodynamic limit ($L\to\infty$) need to be taken. The
different panels of Fig. 3 correspond to the continuum extrapolations for
different aspect ratios $LT$ by increasing $N_{s}$ and $N_{t}$ (with $LT$ and
$m/T$ kept constant). For small values of $m/T$, finite-size effects are
expected to be sizeable, and the results confirm this expectation. In
particular, $LT=4$ is already found to agree with the infinite volume limit
for all masses. A cross-check with Wilson fermions is also shown for the
largest value of $m/T$.

The main message that we get from this calculation is that $C_{\text{CSE}}$
approaches the value given by Eq. (1) when the volume is sufficiently large
and the continuum limit is taken. In the case of Wilson fermions, we also
learn that the correct setup is to consider a conserved vector current and the
anomalous axial current. If instead we use the local vector current (17), the
analytical result is no longer recovered and the continuum limit shows a
possibly divergent behavior. A very similar phenomenon occurs for staggered
fermions if we exclude the tadpole term of (13), again leading to a non-
conserved (in fact ultraviolet divergent) vector current. Having cross-checked
the consistency of our setup with the analytical result, we can proceed toward
the physical setting by turning on the gluonic interactions.

Figure 3: CSE conductivity for free fermions at different values of $m/T$
using the staggered discretization. For $m/T=4$, Wilson fermions as well as
discretizations with a non-conserved vector current are also included. The
dashed lines correspond to the marked values from Fig 2.

###  3.2 Quenched QCD

Before moving on to full QCD, there is an intermediate step where one can
analyze the impact of gluons on the CSE: the quenched approximation.
Furthermore, as we will see, the quenched theory reveals an interesting
feature of the CSE due to the presence of an exact center symmetry. To
appreciate this, next we briefly introduce the notion of the center symmetry
and the Polyakov loop.

In the absence of dynamical fermions, the system undergoes a first-order phase
transition from confined to deconfined matter at around $T_{c}^{q}\approx 270$
MeV Boyd:1996bx , for which the Polyakov loop

$P=\frac{1}{V}\sum_{n}\Tr\quantity[\,\prod_{n_{4}=0}^{N_{t}-1}U_{4}(n)]\,,$
(19)

acts as the order parameter.

This behavior arises from the $\mathrm{Z}(3)$ center symmetry of the gauge
action, which is invariant under transformations $U_{4}(n)\rightarrow
zU_{4}(n)$, with $z\in\mathrm{Z}(3)$, while the Polyakov loop is not,
enforcing a vanishing expectation value of $P$. Below the transition
temperature, center symmetry is intact and on typical configurations $P=0$. In
turn, in the deconfined phase center symmetry is spontaneously broken, and the
theory chooses one of the three vacua with $\arg P=0$, $\pm 2\pi/3$ with equal
probability.

Figure 4: Results for the Polyakov loop (left) and the derivative of the
current with staggered fermions (right) in an ensemble of 100 configurations
on a $32^{3}\times 8$ lattice at $T\approx 400$ MeV and at a magnetic field of
$eB=0.74\text{ GeV}^{2}$. The three sectors can be clearly distinguished at
this temperature, and the contribution of configurations with imaginary
Polyakov loops is observed to be very different as compared to those in the
real sector.

The Polyakov loop sectors at high temperature (where these are spatially
approximately homogeneous) can be mimicked by an imaginary baryon chemical
potential $i\mu_{\mathcal{B}}/T=0,\pm 2\pi/3$.

As we show in App. A, nonzero imaginary chemical potentials give a non-trivial
contribution222In fact, at $i\mu_{\mathcal{B}}\neq 0$, the expectation values
of the baryon density and the axial current are nonzero, giving rise to an
additional term
$\propto\expectationvalue{\text{Tr}(\Gamma_{35}^{f}M_{f}^{-1})}\expectationvalue{\text{Tr}(\Gamma_{4}^{f^{\prime}}M_{f^{\prime}}^{-1})}$
in Eq. (13) for example. This term is included for the comparison in Fig. 4.
to the observable (6). In Fig. 4, we present an example of this behavior at a
temperature of $T\approx 400$ MeV. We observe that the Polyakov loop populates
the three different, well separated sectors (left panel). In the right panel
of Fig. 4, we show the contributions of the individual sectors to the real
part of the observable333That is to say, the imaginary part of the Euclidean
observable, see Eq. (8)., indicating the effect of the effective imaginary
baryon chemical potential. The fluctuations of the current derivative are
observed to be enhanced drastically in the imaginary sectors.

In full QCD, the quarks explicitly break the $\mathrm{Z}(3)$ symmetry and,
consequently, the theory is always in the real Polyakov loop sector –
corresponding to a vanishing imaginary baryon chemical potential. This
indicates that the relevant QCD contribution to the CSE originates from
configurations with real Polyakov loops. Therefore, for the quenched analysis
of $C_{\text{CSE}}$ below, we rotate all our gauge configurations to this
sector by the appropriate center transformation. We note that while this
procedure has a clear impact at high temperatures, where the Polyakov loop
background is homogeneous (and thus equivalent to a nonzero $i\mu_{B}$),
around $T_{c}^{q}$ the non-trivial spatial distribution of center clusters
(see e.g. Endrodi:2014yaa ) complicates the interpretation. This is a
shortcoming of the quenched approximation for the CSE.

Having this caveat in mind, we present the results for $C_{\text{CSE}}$ in the
quenched approximation in Fig. 5. We use configurations generated with the
plaquette gauge action at $\beta=5.845,5.9,6.0,6.2,6.26$ and $6.47$, and use
both staggered and Wilson fermions in the valence sector. We employed these
gauge configurations already in Refs. Bali:2017ian ; Bali:2018sey . For the
staggered discretization, the quark masses were tuned to a pion mass of
$M_{\pi}\approx 415$ MeV, while for Wilson fermions we use $M_{\pi}\approx
710$ MeV. The latter choice is motivated by the CME study Yamamoto:2011ks .

Figure 5: Results for $C_{\text{CSE}}$ with Wilson and staggered fermions in
the quenched approximation, measured on configurations generated using the
plaquette gauge action. The configurations were rotated to the real Polyakov
loop sector (for details see the text). The valence pion mass is set to
$M_{\pi}\approx 415$ MeV for staggered quarks and $M_{\pi}\approx 710$ MeV for
Wilson quarks. The quenched critical temperature $T_{c}^{q}$ is indicated by
the dashed vertical line.

Although here we only present results at a few temperatures, we can clearly
observe two distinct regimes: in the confined phase the CSE is severely
suppressed, reaching zero toward $T\approx 0$, while at temperatures well
above $T_{c}^{q}$, the result approaches the massless free fermion value.
Moreover, the transition between these two regimes appears to occur in the
vicinity of $T_{c}^{q}\approx 270\textmd{ MeV}$. All these hints given by the
quenched results will be confirmed by the full QCD simulations, which we
present next.

###  3.3 Full QCD at physical quark masses

Finally, we present the main result of this study: the conductivity
$C_{\text{CSE}}$ in full QCD, in particular for $N_{f}=2+1$ flavors of
staggered fermions at physical quark masses. This is the first fully non-
perturbative result for $C_{\text{CSE}}$ at the physical point. The
measurements were performed on an already existing ensemble of configurations
for different magnetic fields Bali:2011qj ; Bali:2012zg .

In Fig. 6, we present the dependence of the conductivity on the temperature
using several finite-temperature lattice ensembles $24^{3}\times 6$,
$24^{3}\times 8$, $28^{3}\times 10$, $36^{3}\times 12$ as well as two zero-
temperature ensembles $24^{3}\times 32$ and $32^{3}\times 48$. The impact of
the temperature on $C_{\text{CSE}}$ is very pronounced, as we have already
seen in the quenched results. Above the QCD crossover temperature
$T_{c}\approx 150\textmd{ MeV}$, $C_{\text{CSE}}$ is found to approach the
value corresponding to free massless quarks, in accordance with asymptotic
freedom. In turn, at temperatures below $T_{c}$, the conductivity decreases
until reaching zero. This indicates a suppression of the CSE at low
temperatures, which is consistent with a previous study in two-color QCD
Buividovich:2020gnl .

In the low temperature region, we consider a simple model to describe the
system: a non-interacting gas of hadronic degrees of freedom. For this
particular observable, only electrically charged hadrons contribute that
couple to chirality (i.e. the Dirac indices of $\gamma_{5}$)444We do not
consider the impact of Wess-Zumino-Witten type terms Scherer:2002tk  in our
model.. Thus, in our model we only consider a gas of protons and
$\Sigma^{\pm}$ baryons. The contribution of heavier charged spin-$1/2$ baryons
is negligible in the relevant temperature range, and we do not include
spin-$3/2$ baryons either. The so constructed model also features a strong
suppression of $C_{\text{CSE}}$ in the confined regime and is found to agree
with the lattice results for $T\lesssim 120\textmd{ MeV}$. It is also in
qualitative agreement with other approaches to the CSE Avdoshkin:2017cqp ;
Zubkov:2023vvb .

The continuum limit is carried out using the lattice data in the range of
temperatures $90\textmd{ MeV}\lesssim T\lesssim 400$ MeV. To reliably
extrapolate to the continuum555We note that the flavor-singlet axial vector
current entering $C_{\text{CSE}}$ is subject to multiplicative
renormalization. The corresponding (perturbative) renormalization constant
approaches unity in the continuum limit, both for Wilson and for staggered
fermions Constantinou:2016ieh ; Bali:2021qem . It might be used to reduce
discretization errors, but in this work, we perform the continuum limit
without including multiplicative renormalization constants., we consider a
spline fit procedure combined with the continuum limit. For this, we consider
a $T$-dependent spline fit of all lattice ensembles with $a$-dependent
coefficients. The best fitting surface in the $a-T$ plane is found by
minimizing $\chi^{2}/\text{dof}$. The details of the spline fit can be found
in Endrodi:2010ai . The statistical error of this procedure is calculated
using the jackknife samples, while we estimate the systematic error by
repeating the spline fit removing the coarsest lattice. The maximum difference
between the continuum limits obtained with these two data sets is taken as the
systematic error, which is added in quadrature to the statistical one.

Figure 6: $C_{\text{CSE}}$ for a broad range of temperatures in full QCD with
$2+1$ flavors of staggered fermions. The continuum limit extracted via the
spline fit procedure described in the text, together with its uncertainty, is
indicated by the orange band. The result of our low-energy model involving $p$
and $\Sigma^{\pm}$ baryons is shown by the continuous line at low $T$, while
the parameterization for $C_{\text{CSE}}$ in the whole range of temperatures
is displayed as a dashed-dotted line.

We further provide a complete parameterization of $C_{\text{CSE}}$, to be used
for comparisons with effective theories of QCD and in (anomalous) hydrodynamic
descriptions of heavy-ion collisions. The parameterization, which is also
shown in Fig. 6, smoothly connects the continuum extrapolation of the lattice
results with the perturbative limit and the low-energy model. The details of
the parameterization are given in App. C. In summary, we conclude that the CSE
conductivity is an observable very sensitive to the finite temperature QCD
crossover, and it may even be used to define its characteristic temperature.
In particular, the inflection point of the parameterization is found at
$T_{c}\approx 136(1)$ MeV. As an alternative definition, we also consider the
temperature, where the conductivity assumes half its asymptotic value. This
gives $T_{c}\approx 149(4)$ MeV.

The results above correspond to
$C_{\text{CSE}}=C_{\text{CSE}}^{\mathcal{Q}\mathcal{B}}$, i.e. the response of
the electrically charged axial current to the baryon chemical potential. Next,
we compare results for other combinations of the flavor quantum numbers
$\mathcal{Q},\mathcal{B}$ and $\mathcal{L}$. The corresponding proportionality
factor $C_{\rm dof}$ is different for each (see Sec. 2.1), but even after
normalizing the conductivities by these factors, the results are not
equivalent, since the different quark flavors interact with each other via
gluons. Out of the nine possibilities, in Fig. 7 we present the continuum
extrapolated results for the six combinations that are most interesting from a
phenomenological point of view.

Figure 7: Continuum limits of $C_{\text{CSE}}$ for different flavor quantum
numbers (as defined in Eq. (3)) for the axial current and the chemical
potential. We present the phenomenologically relevant cases involving an
electrically charged axial current (left panel) and a baryon axial current
(right panel).

Except for the $\mathcal{B}\mathcal{B}$ case, all combinations show similar
trends, being zero at low $T$ and approaching the massless free fermion result
at high temperature. Still, differences at intermediate temperatures are
visible, in particular in the cases with a baryonic axial current. The
$\mathcal{B}\mathcal{B}$ setup is exceptional, because the three-flavor
perturbative result in this case vanishes (see the discussion in Sec. 2.1).
Our full QCD result for $C_{\text{CSE}}^{\mathcal{B}\mathcal{B}}$ also
approaches zero for high temperatures. While at intermediate $T$, the
perturbatively expected cancellation between the contributions of the three
flavors is only partial, the result is still found to be much smaller than for
the other quantum numbers. Moreover, the data for
$C_{\text{CSE}}^{\mathcal{B}\mathcal{B}}$ is rather noisy due to the
disconnected contribution in Eq. (13), especially around $T_{c}$. While our
current results do not enable us to resolve this temperature dependence, they
do show that the continuum limit of this conductivity lies within the yellow
band shown in the right panel of Fig. 7.

##  4 Summary and outlook

In this paper, we have presented a comprehensive study of the Chiral
Separation Effect using first-principles lattice QCD simulations. In
particular, we have calculated the conductivity $C_{\text{CSE}}$ using
dynamical staggered quarks with physical quark masses as well as quenched
staggered and quenched Wilson quarks. In addition, we calculated
$C_{\text{CSE}}$ in the absence of gluonic interactions both in the continuum
and on the lattice with staggered and Wilson quarks.

In the free case, our analytic formula (1) reproduces the known results in the
literature Son:2004tq ; Metlitski:2005pr . Moreover, we find that both lattice
discretizations reproduce the analytical result once the continuum limit is
taken and the volume is sufficiently large. We highlight that the use of
conserved vector currents and anomalous axial currents on the lattice is
crucial – non-conserved vector currents are demonstrated to lead to incorrect
results.

Having cross-checked our setup, we moved on to full QCD, where we have
determined the continuum limit of $C_{\text{CSE}}$ for a broad range of
temperatures at physical quark masses. At temperatures higher than $T_{c}$,
the conductivity approaches the prediction for massless free fermions, as
expected, while at temperatures below $T_{c}$, the conductivity goes to zero,
revealing a strong suppression of the CSE in the confined phase. The result
depends on the flavor quantum numbers used for the axial current and the
chemical potential. Even after rescaling the conductivity by the
multiplicative factor $C_{\rm dof}$ containing the degrees of freedom, for
example a charge chemical potential induces slightly different baryonic or
charge axial currents.

A widely employed approximation to take into account gluonic interactions is
quenched QCD, which we also explored in this paper. In the quenched theory, we
highlighted the impact of the Polyakov loop sector on the CSE at high
temperatures, an important subtlety that needs to be treated carefully when
working in this approximation. Only the real Polyakov loop sector delivers the
physical result – configurations in imaginary Polyakov loop sectors deviate
from it and also come with drastically enhanced fluctuations. This issue can
be overcome by rotating the configurations to the real sector by a center
transformation. The results are in qualitative agreement with the picture
based on our full QCD results: a strong suppression of the CSE at low
temperatures and a tendency to the massless free fermion result at
temperatures much higher than $T_{c}$.

Our final results in full QCD are interpolated by a simple parameterization,
explained in App. C, for use in hydrodynamic descriptions of heavy-ion
collisions. Finally, we note that our technique can be generalized to study
the CME, which we plan to carry out in an upcoming publication. Performing the
analogous simulations for that case will contribute to a better theoretical
understanding of anomalous transport phenomena, as a counterpart to the large-
scale experimental efforts being made to detect this effect.

###### Acknowledgements.

This research was funded by the DFG (Collaborative Research Center CRC-TR 211
“Strong-interaction matter under extreme conditions” - project number
315477589 - TRR 211) and by the Helmholtz Graduate School for Hadron and Ion
Research (HGS-HIRe for FAIR). GE would like to express special thanks to the
Mainz Institute for Theoretical Physics (MITP) of the Cluster of Excellence
PRISMA+ (Project ID 39083149) for its hospitality and support. The authors are
moreover grateful for inspiring discussions with Pavel Buividovich, Kenji
Fukushima, Dirk Rischke, Sören Schlichting, Igor Shovkovy and Lorenz von
Smekal.

##  Appendix A Analytical treatment of the CSE coefficient for free quarks

We will now shortly summarize the continuum calculation leading to Eq. (1) for
non-interacting fermions. As mentioned in the main text, it is sufficient to
consider one colorless fermion flavor (with mass $m$ and charge $q$), a single
chemical potential $\mu$ and axial current $J_{35}$. To enable a direct
comparison to the lattice results, we will evaluate the $\mu$-derivative of
the axial current in a homogeneous magnetic field background, in which case
free fermion propagators are known exactly Schwinger:1951nm ; Shovkovy:2012zn
. We will carry out the calculation starting out in Minkowski metric and
extending to finite temperature during the calculation using the Matsubara
formalism. Accordingly, the Dirac matrices are the Minkowski ones fulfilling
$\\{\gamma_{\mu},\gamma_{\nu}\\}=2\eta_{\mu\nu}=2\,\textmd{diag}(1,-1,-1,-1)$
contrary to the main text where Euclidean Dirac matrices are used.

Since the chemical potential, $\mu$ couples to the four-volume integral of the
zeroth component of the vector current $J_{\mu}=\bar{\psi}\gamma_{\mu}\psi$,
the $\mu$-derivative of the $z$-component of the chirality current can be
written as

$\displaystyle\left.\frac{\partial\langle
J_{35}\rangle}{\partial\mu}\right|_{\mu=0}=\int\differential^{4}x\int\differential^{4}y\,\langle\bar{\psi}(x)\gamma_{3}\gamma_{5}\psi(x)\bar{\psi}(y)\gamma_{0}\psi(y)\rangle\,.$
(20)

We regularize in the ultraviolet using the Pauli-Villars (PV) scheme, meaning
that all diagrams are replicated by the regulator fields with coefficients
$c_{s}$ and masses $m_{s}$. The details of the regularization is text-book
knowledge found e.g. in Ref. Itzykson:1980rh . We recall here that we need
three extra fields altogether, and reserving $s=0$ for the physical field
(with physical mass $m$) the parameters will be

$\displaystyle c_{0}$ $\displaystyle=c_{1}=1\,,\quad c_{2}=c_{3}=-1\,,$ (21)
$\displaystyle m_{0}^{2}$ $\displaystyle=m^{2}\,,\quad
m_{1}^{2}=m^{2}+2\Lambda^{2}\,,\quad m_{2}^{2}=m_{3}^{2}=m^{2}+\Lambda^{2}\,,$
(22)

with $\Lambda\to\infty$ to be taken at the end of the calculation. From here
on, the dependence on the regulator is made implicit and in fact, the Pauli-
Villars fields turn out to give zero contribution to the CSE conductivity.
However, their role is found to be crucial for the analogous calculation of
the CME conductivity (which we will include in an upcoming publication).

The derivative according to Eq. (6) is proportional to the CSE coefficient,
and using Wick’s theorem we can rewrite the right hand side of Eq. (20) to
obtain

$\displaystyle C_{\rm
CSE}\,qB=\frac{iT}{V}\sum_{s=0}^{3}c_{s}\int\differential^{4}x\int\differential^{4}y\Tr\left[\gamma_{3}\gamma_{5}S_{s}(x,y)\gamma_{0}S_{s}(y,x)\right]\,,$
(23)

where the PV fields are already taken into account, and $S_{s}$ is the fermion
propagator for the field $s$ in a homogeneous magnetic background. The latter
reads Shovkovy:2012zn

$\displaystyle
S_{s}(x,y)=\Phi(x,y)\int\frac{\differential^{4}p}{(2\pi)^{4}}{\rm\,e\,}^{-ip(x-y)}\widetilde{S}_{s}(p)\,.$
(24)

Here

$\displaystyle\Phi(x,y)=\exp\left[iqB(x_{1}+y_{1})(x_{2}-y_{2})\right]\,,$
(25)

is the Schwinger phase, and

$\displaystyle\widetilde{S}_{s}(p)=\int_{0}^{\infty}\differential
z{\rm\,e\,}^{iz(p_{0}^{2}-m_{s}^{2}-p_{3}^{2})-i\frac{p_{1}^{2}+p_{2}^{2}}{|qB|}\tan{(z|qB|)}}$
$\displaystyle\left[\not{p}+m_{s}+(p_{1}\gamma_{2}-p_{2}\gamma_{1})\tan(zqB)\right]$
$\displaystyle\times\left[1-\gamma_{1}\gamma_{2}\tan(zqB)\right]\,.$ (26)

To extend the calculation to finite temperature, we replace the frequency
component of the integral by a sum over fermionic Matsubara frequencies

$\displaystyle\int\frac{\differential p_{0}}{2\pi}\to
T\sum_{n=-\infty}^{\infty}\,,$ (27)

while also replacing $p_{0}\to i\omega_{n}=i2\pi T(n+1/2)$ in
$\widetilde{S}_{s}(p)$. Plugging this into (23) and evaluating the trace over
Dirac indices leads to

$\displaystyle C_{\rm CSE}\,qB=4\sum_{s=0}^{3}c_{s}T\sum_{n}$
$\displaystyle\int\frac{\differential^{3}p}{(2\pi)^{3}}\int_{0}^{\infty}\differential
z_{1}\,\differential
z_{2}{\rm\,e\,}^{i(z_{1}+z_{2})(-\omega_{n}^{2}-m_{s}^{2}-p_{3}^{2})}\left(m_{s}^{2}-\omega_{n}^{2}+p_{3}^{2}\right)$
$\displaystyle{\rm\,e\,}^{-i\frac{p_{1}^{2}+p_{2}^{2}}{qB}[\tan{(z_{1}qB)}+\tan{(z_{2}qB)}]}\left[\tan(z_{1}qB)+\tan(z_{2}qB)\right]\,.$
(28)

The $p_{1}$ and $p_{2}$ integrals are simple Gaussians, and evaluating them
factorizes the magnetic field dependence, uncovering an explicitly linear
behavior in $B$. The rather simple formula emerges

$\displaystyle C_{\rm CSE}=-\frac{1}{2\pi^{2}}\sum_{s=0}^{3}c_{s}T\sum_{n}$
$\displaystyle\int_{-\infty}^{\infty}\differential
p_{3}\int_{0}^{\infty}\differential z_{1}\,\differential
z_{2}{\rm\,e\,}^{i(z_{1}+z_{2})(-\omega_{n}^{2}-m_{s}^{2}-p_{3}^{2})}\left(m_{s}^{2}-\omega_{n}^{2}+p_{3}^{2}\right)\,,$
(29)

allowing the evaluation of the $z_{1}$ and $z_{2}$ integrals. We find

$\displaystyle C_{\rm
CSE}=-\frac{1}{2\pi^{2}}\sum_{s=0}^{3}c_{s}T\sum_{n=-\infty}^{\infty}\int_{-\infty}^{\infty}\differential
p_{3}\,\frac{m_{s}^{2}-\omega_{n}^{2}+p_{3}^{2}}{\left(\omega_{n}^{2}+p_{3}^{2}+m_{s}^{2}\right)^{2}}\,.$
(30)

The Matsubara sum then evaluates to

$\displaystyle C_{\rm
CSE}=-\frac{1}{2\pi^{2}}\sum_{s=0}^{3}c_{s}\int_{-\infty}^{\infty}\differential
p_{3}\,\frac{\differential n_{F}(E_{p}^{(s)})}{\differential E_{p}^{(s)}}\,,$
(31)

where $n_{F}(x)=({\rm e}^{x/T}+1)^{-1}$ is the Fermi-Dirac distribution and
$E_{p}^{(s)}=\sqrt{p_{3}^{2}+m_{s}^{2}}$. The derivative of the Fermi-Dirac
distribution vanishes for infinite masses, therefore only the $s=0$, physical
term remains from the sum over PV fields,

$\displaystyle C_{\rm CSE}=-\frac{1}{\pi^{2}}\int_{0}^{\infty}\differential
p_{3}\,\frac{\differential n_{F}(E_{p})}{\differential E_{p}}\,,$ (32)

which, after carrying out the derivative, gives Eq. (1) of the main text.

We point out here that generalizing this result to finite chemical potentials
is rather simple and only amounts to replacing $n_{F}(E_{p})$ by
$[n_{F}(E_{p}+\mu)+n_{F}(E_{p}-\mu)]/2$. Specifically, e.g. for imaginary
chemical potentials $\mu=i\mu_{I}$ we arrive at

$\displaystyle C_{\rm CSE}=\frac{1}{2\pi^{2}}\int_{0}^{\infty}\differential
p\,\frac{1+\cosh(\sqrt{p^{2}+(m/T)^{2}})\cos(\mu_{I}/T)}{\left[\cosh(\sqrt{p^{2}+(m/T)^{2}})+\cos(\mu_{I}/T)\right]^{2}}\,.$
(33)

This formula reveals the impact of constant Polyakov loop backgrounds on the
CSE conductivity. The three center sectors correspond to $\mu_{I}=\pm 2\pi/3$
and $\mu_{I}=0$, for which (33) gives different responses, as can be seen in
Fig. 8. Interestingly, averaging over the three center sectors gives the same
result as taking merely the real sector, except for a rescaling of the fermion
mass,

$\frac{1}{3}\sum_{k=0,1,2}C_{\rm CSE}(\mu_{I}/T=2\pi k/3,m/T)=C_{\rm
CSE}(\mu_{I}=0,3m/T)\,.$ (34)

The averaged conductivity is also included in Fig. 8 for comparison.

Figure 8: The CSE conductivity as a function of $m/T$ for non-interacting
quarks, at different values of the imaginary chemical potential $\mu_{I}$. The
averaged result corresponds to averaging over the three center sectors, that
is $\mu_{I}/T=0,2\pi/3,-2\pi/3$.

##  Appendix B CSE with free staggered fermions

In this appendix, we illustrate how to calculate $C_{\text{CSE}}$ for free
staggered quarks using explicitly the eigenvalues and eigenvectors of the
Dirac operator. Therefore, here we turn off the gluon fields and consider one
color component only. Just as in App. A, we discuss a single quark flavor with
mass $m$, one chemical potential $\mu$ and axial current $J_{35}$. Below we
work on an $N_{s}^{3}\times N_{t}$ lattice and write everything in lattice
units, setting $a=1$. In the non-interacting theory, the disconnected term of
Eq. (13) vanishes, so we only need to calculate two terms

$\evaluated{\partialderivative{J_{35}}{\mu}}_{\mu=0}=\dfrac{1}{N_{s}^{3}N_{t}}\Bigg{[}-\dfrac{1}{4}\text{Tr}\quantity(\Gamma_{4}M^{-1}\Gamma_{35}M^{-1})+\dfrac{1}{4}\text{Tr}\quantity(\dfrac{\partial\Gamma_{35}}{\partial\mu}M^{-1})\Bigg{]}\,,$
(35)

and the expectation value indicating the fermion path integral was omitted for
brevity. Furthermore, for convenience we redefine the staggered phases as

$\eta_{1}(n)=1,\quad\eta_{2}(n)=(-1)^{n_{1}},\quad\eta_{3}(n)=(-1)^{n_{1}+n_{2}+n_{4}},\quad\eta_{4}(n)=(-1)^{n_{1}+n_{2}}\,.$
(36)

The massless staggered Dirac operator is antihermitian,
$\not{D}^{\dagger}=-\not{D}$, therefore its eigenvalues are purely imaginary.
Moreover, due to staggered chiral symmetry, $\\{\not{D},\eta_{5}\\}=0$ (with
$\eta_{5}=(-1)^{n_{1}+n_{2}+n_{3}+n_{4}}$), the eigenvalues come in complex
conjugate pairs,

$\not{D}\Psi_{i}=\pm i\lambda_{i}\Psi_{i}\,.$ (37)

In addition, we will need the analogous eigensystem for
$MM^{\dagger}=(\not{D}+m)(\not{D}+m)^{\dagger}=\not{D}\not{D}^{\dagger}+m^{2}$,
so

$MM^{\dagger}\Psi_{i}=(\lambda_{i}^{2}+m^{2})\Psi_{i}\,,$ (38)

where each eigenvalue is doubly degenerate due to (37). Using this eigensystem
as basis, the traces in Eq. (35) can be written as

$\begin{split}\evaluated{\partialderivative{J_{35}}{\mu}}_{\mu=0}=\dfrac{1}{N_{s}^{3}N_{t}}\Bigg{[}&-\dfrac{1}{4}\sum_{i,j}\dfrac{1}{(\lambda^{2}_{i}+m^{2})(\lambda^{2}_{j}+m^{2})}\Psi^{\dagger}_{i}\Gamma_{4}M^{\dagger}\Psi_{j}\Psi^{\dagger}_{j}\Gamma_{35}M^{\dagger}\Psi_{i}\\\
&+\dfrac{1}{4}\sum_{i}\dfrac{1}{\lambda^{2}_{i}+m^{2}}\Psi^{\dagger}_{i}\dfrac{\partial\Gamma_{35}}{\partial\mu}M^{\dagger}\Psi_{i}\Bigg{]}\,,\end{split}$
(39)

where we inserted a complete set of eigenstates
$\sum_{j}\Psi_{j}\Psi_{j}^{\dagger}=\mathds{1}$ in the first term and used
$M^{-1}=M^{\dagger}(MM^{\dagger})^{-1}$.

Next, we make use of the separability of the problem and reduce Eq. (38) to
one two-dimensional and two one-dimensional eigenproblems. This will enable us
to determine the complete spectrum on much larger lattices than with a direct,
four-dimensional approach. To this end, we write the free Dirac operator as
$\not{D}=\not{D}_{12}+\not{D}_{3}+\not{D}_{4}$ with
$\not{D}_{12}=\not{D}_{1}+\not{D}_{2}$ and

$\not{D}_{\nu}(n,m)=\dfrac{\eta_{\nu}(n)}{2}\left[u_{\nu}(n)\delta_{n+\hat{\nu},m}-u^{*}_{\nu}(n-\hat{\nu})\delta_{n-\hat{\nu},m}\right]\equiv\eta_{\nu}(n)D_{\nu}(n,m)\,.$
(40)

Similarly, the staggered Dirac matrices (2.2) simplify to

$\displaystyle\Gamma_{\nu}(n,m)=\dfrac{\eta_{\nu}(n)}{2}\left[u_{\nu}(n)\delta_{n+\hat{\nu},m}+u^{*}_{\nu}(n-\hat{\nu})\delta_{n-\hat{\nu},m}\right]\equiv\eta_{\nu}(n)S_{\nu}(n,m)\,.$
(41)

Notice that the hop operators satisfy the property

$[S_{1},S_{3}]=[S_{1},S_{4}]=[S_{2},S_{3}]=[S_{2},S_{4}]=[S_{3},S_{4}]=0\,,$
(42)

because the $\mathrm{U}(1)$ links only enter in $S_{1}$ and $S_{2}$. Moreover,

$\left.\frac{\partial S_{4}}{\partial\mu}\right|_{\mu=0}=D_{4}\,,$ (43)

due to the way that the chemical potential enters the temporal hoppings.

The squared operator
$MM^{\dagger}=MM^{\dagger}_{12}+MM^{\dagger}_{3}+MM^{\dagger}_{4}$ separates
into three terms that act in the respective subspaces and only depend on the
respective coordinates – thus they commute with each other. Therefore, the
eigenvectors in (38) factorize,

$\Psi_{\\{i_{12},i_{3},i_{4}\\}}(n_{1},n_{2},n_{3},n_{4})=\rho_{i_{12}}(n_{1},n_{2})\,\phi_{i_{3}}(n_{3})\,\xi_{i_{4}}(n_{4})\,,$
(44)

with $0\leq i_{12}<N_{s}^{2}$, $0\leq i_{3}<N_{s}$ and $0\leq
i_{4}<N_{t}$.666Note that the same separation does not hold for the
eigensystem (37), since for example $[\not{D}_{12},\not{D}_{3}]\neq 0$ due to
the staggered phases. This is because $[MM^{\dagger},\eta_{5}]=0$ but
$[\not{D},\eta_{5}]\neq 0$. Below we use a shorthand notation and simply write
$\Psi_{i}=\rho_{i}\phi_{i}\xi_{i}$.

We proceed by expanding the operators appearing in the matrix elements in Eq.
(39) into separate contributions that depend only on $n_{1,2}$, $n_{3}$ or
$n_{4}$,

$\displaystyle\Gamma_{4}$ $\displaystyle=(-1)^{n_{1}+n_{2}}S_{4}\,,$ (45)
$\displaystyle M^{\dagger}$
$\displaystyle=-\not{D}_{12}-(-1)^{n_{1}+n_{2}+n_{4}}D_{3}-(-1)^{n_{1}+n_{2}}D_{4}+m\,,$
(46)

and it can be shown with some algebra and Eq. (42) that

$\Gamma_{35}=S_{12}S_{4}(-1)^{n_{2}}\,,$ (47)

with $S_{12}=\\{S_{1},S_{2}\\}/2$. Finally, for the tadpole term we need the
derivative of $\Gamma_{35}$ with respect to the chemical potential. The latter
only appears in $S_{4}$. Therefore, using (43),

$\displaystyle\dfrac{\partial\Gamma_{35}}{\partial\mu}=S_{12}D_{4}(-1)^{n_{2}}\,.$
(48)

The necessary products, appearing in (39) are therefore

$\begin{split}\Gamma_{4}M^{\dagger}=-&[(-1)^{n_{1}+n_{2}}\not{D}_{12}]_{12}\cdot[\mathds{1}]_{3}\cdot[S_{4}]_{4}\\\
-&[\mathds{1}]_{12}\cdot[D_{3}]_{3}\cdot[S_{4}(-1)^{n_{4}}]_{4}\\\
-&[\mathds{1}]_{12}\cdot[\mathds{1}]_{3}\cdot[S_{4}D_{4}]_{4}\\\
+&m[(-1)^{n_{1}+n_{2}}]_{12}\cdot[\mathds{1}]_{3}\cdot[S_{4}]_{4}\,,\end{split}$
(49)

where $[.]_{12}$, $[.]_{3}$ and $[.]_{4}$ indicate operators that only act in
the respective spaces and only depend on the respective coordinates.
Similarly, we obtain

$\begin{split}\Gamma_{35}M^{\dagger}=-&[S_{12}(-1)^{n_{2}}\not{D}_{12}]_{12}\cdot[\mathds{1}]_{3}\cdot[S_{4}]_{4}\\\
-&[S_{12}(-1)^{n_{1}}]_{12}\cdot[D_{3}]_{3}\cdot[S_{4}(-1)^{n_{4}}]_{4}\\\
-&[S_{12}(-1)^{n_{1}}]_{12}\cdot[\mathds{1}]_{3}\cdot[S_{4}D_{4}]_{4}\\\
+&m[S_{12}(-1)^{n_{2}}]_{12}\cdot[\mathds{1}]_{3}\cdot[S_{4}]_{4}\,,\end{split}$
(50)

and

$\begin{split}\dfrac{\partial\Gamma_{35}}{\partial\mu}M^{\dagger}=&-[S_{12}(-1)^{n_{2}}\not{D}_{12}]_{12}\cdot[\mathds{1}]_{3}\cdot[D_{4}]_{4}\\\
&-[S_{12}(-1)^{n_{1}}]_{12}\cdot[D_{3}]_{3}\cdot[D_{4}(-1)^{n_{4}}]_{4}\\\
&-[S_{12}(-1)^{n_{1}}]_{12}\cdot[\mathds{1}]_{3}\cdot[D_{4}D_{4}]_{4}\\\
&+m[S_{12}(-1)^{n_{2}}]_{12}\cdot[\mathds{1}]_{3}\cdot[D_{4}]_{4}\,.\end{split}$
(51)

Combining everything, Eq. (39) becomes

$\begin{split}\dfrac{\partial
J_{35}}{\partial\mu}\Big{|}_{\mu=0}=-&\dfrac{1}{N_{s}^{3}N_{t}}\dfrac{1}{4}\sum_{i,j}\dfrac{1}{(\lambda^{2}_{i}+m^{2})(\lambda^{2}_{j}+m^{2})}\\\
&\Big{\\{}[A_{ij}C_{ji}]_{12}\cdot[\delta_{ij}]_{3}\cdot[G_{ij}G_{ji}]_{4}+[A_{ij}E_{ji}]_{12}\cdot[\delta_{ij}]_{3}\cdot[G_{ij}H_{ji}]_{4}\\\
&+[A_{ij}E_{ji}]_{12}\cdot[\delta_{ij}J_{ji}]_{4}\cdot[G_{ij}I_{ji}]_{3}-m[A_{ij}F_{ji}]_{12}\cdot[\delta_{ij}]_{3}\cdot[G_{ij}G_{ji}]_{4}\\\
&+[\delta_{ij}C_{ji}]_{12}\cdot[\delta_{ij}]_{3}\cdot[H_{ij}G_{ji}]_{4}+[\delta_{ij}E_{ji}]_{12}\cdot[\delta_{ij}]_{3}\cdot[H_{ij}H_{ji}]_{4}\\\
&+[\delta_{ij}E_{ji}]_{12}\cdot[\delta_{ij}J_{ji}]_{3}\cdot[H_{ij}I_{ji}]_{4}-m[\delta_{ij}F_{ji}]_{12}\cdot[\delta_{ij}]_{3}\cdot[H_{ij}G_{ji}]_{4}\\\
&+[\delta_{ij}C_{ji}]_{12}\cdot[J_{ij}\delta_{ji}]_{3}\cdot[I_{ij}G_{ji}]_{4}+[\delta_{ij}E_{ji}]_{12}\cdot[J_{ij}\delta_{ji}]_{3}\cdot[I_{ij}H_{ji}]_{4}\\\
&+[\delta_{ij}E_{ji}]_{12}\cdot[J_{ij}J_{ji}]_{3}\cdot[I_{ij}I_{ji}]_{4}-m[\delta_{ij}F_{ji}]_{12}\cdot[J_{ij}\delta_{ij}]_{3}\cdot[I_{ij}G_{ji}]_{4}\\\
&-m[B_{ij}C_{ji}]_{12}\cdot[\delta_{ij}]_{3}\cdot[G_{ij}G_{ji}]_{4}-m[B_{ij}E_{ji}]_{12}\cdot[\delta_{ij}]_{3}\cdot[G_{ij}H_{ji}]_{4}\\\
&-m[B_{ij}E_{ji}]_{12}\cdot[\delta_{ij}J_{ji}]_{3}\cdot[G_{ij}I_{ji}]_{4}+m^{2}[B_{ij}F_{ji}]_{12}\cdot[\delta_{ij}]_{3}\cdot[G_{ij}G_{ji}]_{4}\Big{\\}}\\\
-&\dfrac{1}{N_{s}^{3}N_{t}}\dfrac{1}{4}\sum_{i}\dfrac{1}{\lambda^{2}_{i}+m^{2}}\\\
&\Big{\\{}[C_{ii}]_{12}\cdot[\delta_{ii}]_{3}\cdot[K_{ii}]_{4}+[E_{ii}]_{12}\cdot[J_{ii}]_{3}\cdot[N_{ii}]_{4}\\\
&+[E_{ii}]_{12}\cdot[\delta_{ii}]_{3}\cdot[L_{ii}]_{4}-m[F_{ii}]_{12}\cdot[\delta_{ii}]_{3}\cdot[K_{ii}]_{4}\Big{\\}}\,,\end{split}$
(52)

with

$\begin{split}&A_{ij}\equiv\rho^{\dagger}_{i}(-1)^{n_{1}+n_{2}}\not{D}_{12}\rho_{j}\,,\\\
&C_{ij}\equiv\rho^{\dagger}_{i}S_{12}(-1)^{n_{2}}\not{D}_{12}\rho_{j}\,,\\\
&F_{ij}\equiv\rho^{\dagger}_{i}S_{12}(-1)^{n_{2}}\rho_{j}\,,\\\
&G_{ij}\equiv\xi^{\dagger}_{i}S_{4}\xi_{j}\,,\\\
&I_{ij}\equiv\xi^{\dagger}_{i}S_{4}(-1)^{n_{4}}\xi_{j}\,,\\\
&L_{ij}\equiv\xi^{\dagger}_{i}D_{4}D_{4}\xi_{j}\,,\end{split}\qquad\begin{split}&B_{ij}\equiv\rho^{\dagger}_{i}(-1)^{n_{1}+n_{2}}\rho_{j}\,,\\\
&E_{ij}\equiv\rho^{\dagger}_{i}S_{12}(-1)^{n_{1}}\rho_{j}\,,\\\
&J_{ij}\equiv\phi^{\dagger}_{i}D_{3}\phi_{j}\,,\\\
&H_{ij}\equiv\xi^{\dagger}_{i}S_{4}D_{4}\xi_{j}\,,\\\
&K_{ij}\equiv\xi^{\dagger}_{i}D_{4}\xi_{j}\,,\\\
&N_{ij}\equiv\xi^{\dagger}_{i}D_{4}(-1)^{n_{4}}\xi_{j}\,.\end{split}$ (53)

We employ the LAPACK library to separately solve the two-dimensional
eigenvalue problem in the $12$ plane in the presence of the magnetic field, as
well as two one-dimensional problems in the $3$ and $4$ directions, both
without electromagnetic phases. The eigensystem of the former problem gives
Hofstadter’s butterfly Hofstadter:1976zz , the spectrum of a well-known solid-
state physics model. Its relevance for lattice QCD has been pointed out in
Ref. Endrodi:2014vza  and it has been generalized to the presence of gluonic
interactions Bruckmann:2017pft ; Bignell:2020dze  as well as for inhomogeneous
magnetic fields Brandt:2023dir . The final value for the current derivative
can be reconstructed using Eq. (52), and $C_{\text{CSE}}$ can be extracted in
the same way as explained in the main text by calculating the observable at
different values of the magnetic field.

##  Appendix C Parameterization of the CSE conductivity

In this appendix we provide a parameterization of the CSE conductivity in QCD
as a function of temperature. We consider the continuum extrapolated data for
$C_{\text{CSE}}^{\mathcal{Q}\mathcal{B}}$ in Fig. 6 and we smoothly
interpolate to the prediction given by our model of a gas of protons and
$\Sigma^{\pm}$ baryons at low $T$, while at high temperatures we impose an
asymptotic approach to the result for a gas of free massless fermions ($m/T=0$
in Eq. (1)). We found the following parameterization to adequately describe
the data,

$C_{\text{CSE}}(T)=\frac{1}{2\pi^{2}}\cdot\exp(-h/t^{4})\cdot\frac{1+a_{1}/t+a_{2}/t^{2}+a_{3}/t^{3}}{1+a_{4}/t+a_{5}/t^{2}+a_{6}/t^{3}}\,,\qquad
t=\frac{T}{1\text{ GeV}}\,,$ (54)

with the parameters included in Table. 1. As visible in Fig. 6, this function
captures all details of our results – more specifically, the confidence
interval of the parameterization is about $95\%$.

$h$ | $a_{1}$ | $a_{2}$ | $a_{3}$ | $a_{4}$ | $a_{5}$ | $a_{6}$
---|---|---|---|---|---|---
0.000601 | -0.359512 | 0.049819 | -0.001234 | -0.346448 | 0.048787 | -0.002011
Table 1: Parameters of the function (54).

## References

  * (1) D.E. Kharzeev, L.D. McLerran and H.J. Warringa, _The Effects of topological charge change in heavy ion collisions: ’Event by event P and CP violation’_ , _Nucl. Phys. A_ 803 (2008) 227 [0711.0950].
  * (2) K. Fukushima, D.E. Kharzeev and H.J. Warringa, _The Chiral Magnetic Effect_ , _Phys. Rev. D_ 78 (2008) 074033 [0808.3382].
  * (3) Q. Li, D.E. Kharzeev, C. Zhang, Y. Huang, I. Pletikosic, A.V. Fedorov et al., _Observation of the chiral magnetic effect in ZrTe5_ , _Nature Phys._ 12 (2016) 550 [1412.6543].
  * (4) STAR collaboration, _Fluctuations of charge separation perpendicular to the event plane and local parity violation in $\sqrt{s_{NN}}=200$ GeV Au+Au collisions at the BNL Relativistic Heavy Ion Collider_, _Phys. Rev. C_ 88 (2013) 064911 [1302.3802].
  * (5) STAR collaboration, _Beam-energy dependence of charge separation along the magnetic field in Au+Au collisions at RHIC_ , _Phys. Rev. Lett._ 113 (2014) 052302 [1404.1433].
  * (6) STAR collaboration, _Search for the chiral magnetic effect with isobar collisions at $\sqrt{s_{NN}}$=200 GeV by the STAR Collaboration at the BNL Relativistic Heavy Ion Collider_, _Phys. Rev. C_ 105 (2022) 014901 [2109.00131].
  * (7) D.T. Son and A.R. Zhitnitsky, _Quantum anomalies in dense matter_ , _Phys. Rev. D_ 70 (2004) 074018 [hep-ph/0405216].
  * (8) M.A. Metlitski and A.R. Zhitnitsky, _Anomalous axion interactions and topological currents in dense matter_ , _Phys. Rev. D_ 72 (2005) 045011 [hep-ph/0505072].
  * (9) D.E. Kharzeev and H.-U. Yee, _Chiral Magnetic Wave_ , _Phys. Rev. D_ 83 (2011) 085007 [1012.6026].
  * (10) X.-G. Huang, _Electromagnetic fields and anomalous transports in heavy-ion collisions — A pedagogical review_ , _Rept. Prog. Phys._ 79 (2016) 076302 [1509.04073].
  * (11) A. Avdoshkin, A.V. Sadofyev and V.I. Zakharov, _IR properties of chiral effects in pionic matter_ , _Phys. Rev. D_ 97 (2018) 085020 [1712.01256].
  * (12) E.V. Gorbar, V.A. Miransky and I.A. Shovkovy, _Normal ground state of dense relativistic matter in a magnetic field_ , _Phys. Rev. D_ 83 (2011) 085003 [1101.4954].
  * (13) A. Jimenez-Alba and L. Melgar, _Anomalous Transport in Holographic Chiral Superfluids via Kubo Formulae_ , _JHEP_ 10 (2014) 120 [1404.2434].
  * (14) M.A. Zubkov and R.A. Abramchuk, _Effect of interactions on the topological expression for the chiral separation effect_ , _Phys. Rev. D_ 107 (2023) 094021 [2301.12261].
  * (15) X.-l. Sheng, D.H. Rischke, D. Vasak and Q. Wang, _Wigner functions for fermions in strong magnetic fields_ , _Eur. Phys. J. A_ 54 (2018) 21 [1707.01388].
  * (16) S. Lin and L. Yang, _Mass correction to chiral vortical effect and chiral separation effect_ , _Phys. Rev. D_ 98 (2018) 114022 [1810.02979].
  * (17) E.V. Gorbar, V.A. Miransky, I.A. Shovkovy and X. Wang, _Radiative corrections to chiral separation effect in QED_ , _Phys. Rev. D_ 88 (2013) 025025 [1304.4606].
  * (18) Y. Aoki, G. Endrődi, Z. Fodor, S.D. Katz and K.K. Szabó, _The Order of the quantum chromodynamics transition predicted by the standard model of particle physics_ , _Nature_ 443 (2006) 675 [hep-lat/0611014].
  * (19) K. Nagata, _Finite-density lattice QCD and sign problem: Current status and open problems_ , _Prog. Part. Nucl. Phys._ 127 (2022) 103991 [2108.12423].
  * (20) M. Puhr and P.V. Buividovich, _Numerical Study of Nonperturbative Corrections to the Chiral Separation Effect in Quenched Finite-Density QCD_ , _Phys. Rev. Lett._ 118 (2017) 192003 [1611.07263].
  * (21) P.V. Buividovich, D. Smith and L. von Smekal, _Numerical study of the chiral separation effect in two-color QCD at finite density_ , _Phys. Rev. D_ 104 (2021) 014511 [2012.05184].
  * (22) Z.V. Khaidukov and M.A. Zubkov, _Chiral Separation Effect in lattice regularization_ , _Phys. Rev. D_ 95 (2017) 074502 [1701.03368].
  * (23) P.V. Buividovich, _Anomalous transport with overlap fermions_ , _Nucl. Phys. A_ 925 (2014) 218 [1312.1843].
  * (24) N. Müller, S. Schlichting and S. Sharma, _Chiral magnetic effect and anomalous transport from real-time lattice simulations_ , _Phys. Rev. Lett._ 117 (2016) 142301 [1606.00342].
  * (25) M. Mace, N. Mueller, S. Schlichting and S. Sharma, _Non-equilibrium study of the Chiral Magnetic Effect from real-time simulations with dynamical fermions_ , _Phys. Rev. D_ 95 (2017) 036023 [1612.02477].
  * (26) P.V. Buividovich, M.N. Chernodub, E.V. Luschevskaya and M.I. Polikarpov, _Numerical evidence of chiral magnetic effect in lattice gauge theory_ , _Phys. Rev. D_ 80 (2009) 054503 [0907.0494].
  * (27) P.V. Buividovich, M.N. Chernodub, E.V. Luschevskaya and M.I. Polikarpov, _Quark electric dipole moment induced by magnetic field_ , _Phys. Rev. D_ 81 (2010) 036007 [0909.2350].
  * (28) A. Yamamoto, _Lattice study of the chiral magnetic effect in a chirally imbalanced matter_ , _Phys. Rev. D_ 84 (2011) 114504 [1111.4681].
  * (29) G.S. Bali, F. Bruckmann, G. Endrődi, Z. Fodor, S.D. Katz and A. Schäfer, _Local CP-violation and electric charge separation by magnetic fields from lattice QCD_ , _JHEP_ 04 (2014) 129 [1401.4141].
  * (30) N. Astrakhantsev, V.V. Braguta, M. D’Elia, A.Y. Kotov, A.A. Nikolaev and F. Sanfilippo, _Lattice study of the electromagnetic conductivity of the quark-gluon plasma in an external magnetic field_ , _Phys. Rev. D_ 102 (2020) 054516 [1910.08516].
  * (31) S. Borsányi, G. Endrődi, Z. Fodor, S.D. Katz, S. Krieg, C. Ratti et al., _QCD equation of state at nonzero chemical potential: continuum results with physical quark masses at order $\mu^{2}$_, _JHEP_ 08 (2012) 053 [1204.6710].
  * (32) S. Borsányi, G. Endrődi, Z. Fodor, A. Jakovác, S.D. Katz, S. Krieg et al., _The QCD equation of state with dynamical quarks_ , _JHEP_ 11 (2010) 077 [1007.2580].
  * (33) G.S. Bali, F. Bruckmann, G. Endrődi, Z. Fodor, S.D. Katz, S. Krieg et al., _The QCD phase diagram for external magnetic fields_ , _JHEP_ 02 (2012) 044 [1111.4956].
  * (34) H.S. Sharatchandra, H.J. Thun and P. Weisz, _Susskind Fermions on a Euclidean Lattice_ , _Nucl. Phys. B_ 192 (1981) 205.
  * (35) S. Dürr, _Taste-split staggered actions: eigenvalues, chiralities and Symanzik improvement_ , _Phys. Rev. D_ 87 (2013) 114501 [1302.0773].
  * (36) P. Hasenfratz and F. Karsch, _Chemical Potential on the Lattice_ , _Phys. Lett. B_ 125 (1983) 308.
  * (37) G. Endrődi, Z. Fodor, S.D. Katz and K.K. Szabó, _The QCD phase diagram at nonzero quark density_ , _JHEP_ 04 (2011) 001 [1102.1356].
  * (38) L.H. Karsten and J. Smit, _Lattice Fermions: Species Doubling, Chiral Invariance, and the Triangle Anomaly_ , _Nucl. Phys. B_ 183 (1981) 103.
  * (39) G. Boyd, J. Engels, F. Karsch, E. Laermann, C. Legeland, M. Lutgemeier et al., _Thermodynamics of SU(3) lattice gauge theory_ , _Nucl. Phys. B_ 469 (1996) 419 [hep-lat/9602007].
  * (40) G. Endrődi, C. Gattringer and H.-P. Schadler, _Fractality and other properties of center domains at finite temperature: SU(3) lattice gauge theory_ , _Phys. Rev. D_ 89 (2014) 054509 [1401.7228].
  * (41) G.S. Bali, B.B. Brandt, G. Endrődi and B. Gläßle, _Meson masses in electromagnetic fields with Wilson fermions_ , _Phys. Rev. D_ 97 (2018) 034505 [1707.05600].
  * (42) G.S. Bali, B.B. Brandt, G. Endrődi and B. Gläßle, _Weak decay of magnetized pions_ , _Phys. Rev. Lett._ 121 (2018) 072001 [1805.10971].
  * (43) G.S. Bali, F. Bruckmann, G. Endrődi, Z. Fodor, S.D. Katz and A. Schäfer, _QCD quark condensate in external magnetic fields_ , _Phys. Rev. D_ 86 (2012) 071502 [1206.4205].
  * (44) S. Scherer, _Introduction to chiral perturbation theory_ , _Adv. Nucl. Phys._ 27 (2003) 277 [hep-ph/0210398].
  * (45) M. Constantinou, M. Hadjiantonis, H. Panagopoulos and G. Spanoudes, _Singlet versus nonsinglet perturbative renormalization of fermion bilinears_ , _Phys. Rev. D_ 94 (2016) 114513 [1610.06744].
  * (46) RQCD collaboration, _Masses and decay constants of the $\eta$ and $\eta$’ mesons from lattice QCD_, _JHEP_ 08 (2021) 137 [2106.05398].
  * (47) G. Endrődi, _Multidimensional spline integration of scattered data_ , _Comput. Phys. Commun._ 182 (2011) 1307 [1010.2952].
  * (48) J.S. Schwinger, _On gauge invariance and vacuum polarization_ , _Phys. Rev._ 82 (1951) 664.
  * (49) I.A. Shovkovy, _Magnetic Catalysis: A Review_ , _Lect. Notes Phys._ 871 (2013) 13 [1207.5081].
  * (50) C. Itzykson and J.B. Zuber, _Quantum Field Theory_ , International Series In Pure and Applied Physics, McGraw-Hill, New York (1980).
  * (51) D.R. Hofstadter, _Energy levels and wave functions of Bloch electrons in rational and irrational magnetic fields_ , _Phys. Rev. B_ 14 (1976) 2239.
  * (52) G. Endrődi, _QCD in magnetic fields: from Hofstadter’s butterfly to the phase diagram_ , _PoS_ LATTICE2014 (2014) 018 [1410.8028].
  * (53) F. Bruckmann, G. Endrődi, M. Giordano, S.D. Katz, T.G. Kovács, F. Pittler et al., _Landau levels in QCD_ , _Phys. Rev. D_ 96 (2017) 074506 [1705.10210].
  * (54) R. Bignell, W. Kamleh and D. Leinweber, _Pion magnetic polarisability using the background field method_ , _Phys. Lett. B_ 811 (2020) 135853 [2005.10453].
  * (55) B.B. Brandt, F. Cuteri, G. Endrődi, G. Markó, L. Sandbote and A.D.M. Valois, _Thermal QCD in a non-uniform magnetic background_ , 2305.19029.