††institutetext:  Department of Physics, Indian Institute of Technology
Indore, Simrol, Indore 453552, India

# Electric field induction in quark-gluon plasma due to thermoelectric effects

Kamaljeet Singh  Jayanta Dey, and  Raghunath Sahoo<EMAIL_ADDRESS><EMAIL_ADDRESS>Corresponding Author<EMAIL_ADDRESS>
###### Abstract

Relativistic heavy-ion collisions produce quark-gluon plasma (QGP), which is
locally thermalized. QGP with higher thermal conductivity advances towards
global thermalization. Being made of electrically charged particles (quarks),
QGP exhibits interesting thermoelectric phenomena during its evolution,
resulting in the induction of an electric field in the medium. For the first
time, we estimate the induced electric field in the QGP due to thermoelectric
effect. This can be seen even in the QGP produced in the head-on collisions.
The Seebeck coefficient is essential in determining the induced field.
However, spectator current can produce a magnetic field in peripheral heavy-
ion collisions. This breaks the isotropy of the thermoelectric coefficient
matrix and introduces magneto-Seebeck and Nernst coefficients that contribute
to the induced electric field in peripheral collisions. We have taken care of
the temperature evolution of QGP with different hydrodynamic cooling rates to
calculate the transport coefficients. The induced electric field is estimated
with the cooling rate obtained from Gubser hydrodynamic flow. We estimated the
space-time profile of the induced field and found that it is zero at the
center and increases as we go away from the center. At the early time of
evolution, the electric field can reach a maximum value of $eE\approx
1~{}m_{\pi}^{2}$, and the strength decreases in time. Moreover, we estimated
the transport coefficients in the presence of the external time-varying
magnetic field. The effect of the intensity and decay parameter of the
magnetic field on the induced electric field is also explored.

††arxiv: 0000

##  1 Introduction

The success of heavy-ion collision experiments is not limited to the discovery
of phenomena of strong interaction at very high energy but also contributes to
understanding fundamental physics applicable to many different areas, such as
condensed matter physics. The quark-gluon plasma (QGP), which is a deconfined
and locally thermalized state of quarks and gluons, produced in relativistic
heavy-ion collisions (RHICs), undergoes many phases during its evolution
BREWER2021136189 ; Busza:2018rrf . Some phenomena occurred during the
evolution Elfner:2022iae ; PhysRevC.106.014905  significantly impact final
state observables to be distinguished, and some might not yet be discovered
LIAO201699 . It is found that relativistic dissipative hydrodynamics is the
most effective theory that can successfully explain final state observables
such as spectra and flow PhysRevLett.87.182301 ; PhysRevC.77.054901 ; article
. The thermodynamic PhysRevD.100.076009 ; PhysRevC.99.035210 ; Goswami:2023eol
and transport Gavin:1985ph ; PhysRevD.100.114004 ; PhysRevD.102.016016 ;
PhysRevD.106.034008 ; PhysRevD.108.094007  properties of nearly baryon-free
QGP matter created in the Large Hadron Collider (LHC) and that of QGP with
finite baryon-chemical potential created in Relativistic heavy-ion Collider
(RHIC) are explored very well. Moreover, going from head-on to peripheral
collisions, an intense transient electromagnetic field is created that can be
understood from laws of classical electrodynamics Tuchin:2015oka ;
PhysRevD.90.034018 ; PhysRevLett.125.022301 ; KHARZEEV2008227 ;
PhysRevX.14.011028 . The change in energy and centrality not only modifies the
thermodynamic equation of state (EoS) Kumar2024 ; Brandt2023 ; Borsanyi2012 ;
Bali2014  of the created medium but many different phenomena, such as chiral
magnetic effect KHARZEEV2008227 ; PhysRevD.97.056024 , magnetic and inverse
magnetic catalysis PhysRevLett.73.3499 ; PhysRevD.55.6504 , and plausible
transition through the quantum chromodynamics (QCD) critical point
PhysRevD.85.091901 ; Endrodi2011 .

The electromagnetic (EM) field produced in peripheral RHICs due to the moving
charged spectators magnetize the QGP medium Wang:2021eud ; MCINNES2016173 .
The magnitude of the field generated at the early stages of collisions could
be as high as $10-15$ m${}_{\pi}^{2}$ ($\approx 10^{18}$ Gauss), which
submerged the QCD energy scale $\Lambda_{\rm QCD}\sim 1.5$ mπ Simonov:2021eyt
. Here, mπ is the mass of a pion, which is around $140$ MeV. For a physical
comparison, the Earth magnetic field at its surface is 0.25-0.65 Gauss. The
created EM field decays in time, which can be explained by the medium’s
electrical conductivity Dey:2021fbo  and the flux conservation as it expands.
QGP being made of electrically charged particles (quarks), a non-zero EM field
can be produced in the medium irrespective of the spectator current depending
on the hydrodynamic flow and charge density profile. This is explored only
recently by A. Dash and A. K. Panda in Ref. DASH2024138342 . They have found
that the produced electric and magnetic fields could be as high as $0.15$
m${}_{\pi}^{2}$ for SPS energy range and $1$ m${}_{\pi}^{2}$ for RHIC energies
respectively. In the present work, we have calculated the electric field
produced in QGP via the thermoelectric effect. The creation of electric
potential due to temperature gradient or vice versa is known as thermoelectric
effect goldsmid2010introduction . In baryonic QGP, the heat flow is governed
by the baryon current carried by the constituent quarks. As the quarks are
electrically charged, a net non-zero electric potential will appear, giving
rise to an electric field. This field can be calculated from thermoelectric
transport coefficients - magneto-Seebeck and Nernst coefficients. In RHICs,
these coefficients are studies for QGP Das2021 ; PhysRevD.102.096011  and
hadronic matter PhysRevD.102.014030 ; PhysRevD.99.014015  in finite and zero
external magnetic field cases. However, in all these studies, a static picture
is considered, i.e., the effect of temperature evolution (or cooling) is not
considered. In our recent work PhysRevD.109.014018 , we have estimated the
effect of QGP cooling in the calculation of thermal conductivity and its
effect on the elliptic flow coefficient. Here, we follow the similar approach
to obtain thermoelectric coefficients of evolving QGP and estimate the
produced electric field. It is important to mention here that unlike the EM
field produced in peripheral heavy-ion collisions due to spectators, this
field can be produced in head-on collisions as well in ultra-relativistic
energies where a non-zero baryon chemical potential is expected. We have also
calculated the thermoelectric coefficients in the presence of an external
magnetic field produced in peripheral collision. For a realistic picture, we
considered a time-dependent magnetic field profile. To the best of our
knowledge, this is the first time an electric field produced purely due to the
thermoelectric effect is estimated. The created spatiotemporal electric field
can give rise to a spatiotemporal magnetic field, which is also discussed with
possible limitations.

The paper is organized in the following manner. First, we briefly derive the
cooling rates for Bjorken and Gubser flow in section (2). In section (3), we
derive the Seebeck coefficient for the head-on collision or $eB=0$ case. Then,
we estimate the Seebeck coefficient and induced electric field with QGP EoS.
Then, we extend our study to peripheral collision in section (3) where the
initial magnetic field is non-zero. We calculate the magneto-Seebeck and
Nernst coefficients for the anisotropic system and calculate the induced
electric field for the same in section (4). Finally, in section (5), we have
summarized our study with possible outlook. Detailed calculations for
thermoelectric coefficients are given in the appendix.

##  2 Relativistic hydrodynamics and cooling rates

The relativistic hydrodynamics is a successful macroscopic theory to
understand the dynamics of relativistic fluid created in RHICs, namely QGP
Rischke:1995mt ; Okamoto2016 . Hydrodynamics and the correct EoS of a medium
will provide us with a detailed dynamical evolution of the medium and its
cooling rate or space-time-dependent temperature profile. Here, to estimate
thermoelectric coefficients and induced electric field, an analytical
expression of cooling rate is required, for which we use Bjorken and Gubser
flow PhysRevD.27.140 ; PhysRevD.82.085027 .

The symmetric energy-momentum tensor ($T^{\mu\nu}$) of a fluid can be
expressed in terms of hydrodynamic degrees of freedom as PhysRevC.103.014905 ,

$\displaystyle{T}^{\mu\nu}=\left(\epsilon+p\right){u}^{\mu}{u}^{\nu}-p{g}^{\mu\nu}+\pi^{\mu\nu}+W^{\mu}u^{\nu}+W^{\nu}u^{\mu}~{}.$
(1)

Where $\epsilon$, $p$, and $u^{\mu}$ are the energy density, thermal pressure,
and four-velocity of fluid, respectively. Dissipative quantities are viscous
stress tensor $\pi^{\mu\nu}$ and energy flow $W^{\mu}$. The metric tensor of
flat spacetime is $g^{\mu\nu}$ = $\rm{diag}(1,-1,-1,-1)$. A hydrodynamic
system is governed by the conservation equations - energy-momentum
conservation and number or charge conservation. Now, the projection of the
energy-momentum conservation equation along the fluid velocity

$\displaystyle{u}_{\mu}{\partial}_{\nu}{T}^{\mu\nu}=0,$ (2)

leads to the continuity equation PhysRevLett.99.172301 . One can obtain the
required cooling rate by solving Eq. (2) with the relevant flow profile and
EoS.

###  2.1 Bjorken flow

The Bjorken flow represents a boost invariant longitudinal expansion of the
medium created in heavy-ion collisions PhysRevD.27.140 . Along with
longitudinal boost invariant, Bjorken considered translation and rotational
invariance symmetry in the transverse plane, which leads to the velocity
profile

$u^{\mu}=\gamma~{}\left(1,~{}0,~{}0,~{}\frac{z}{t}\right)~{}.$ (3)

One can solve the hydrodynamic equations with the flow profile to obtain the
system’s dynamics. Here, we will briefly describe the cooling rate for an
ideal hydrodynamic system with and without a magnetic field.

####  2.1.1 Ideal hydrodynamics at $B=0$

At the initial time of QGP formation, it behaves as a conformal fluid where
all the dissipative quantities vanish, such that conservation of the entropy
current is a good assumption. The energy-momentum tensor of such a perfect
fluid in the absence of external magnetic fields can be expressed as

$\displaystyle{T}^{\mu\nu}=\left(\epsilon+p\right){u}^{\mu}{u}^{\nu}-p{g}^{\mu\nu}~{}.$
(4)

Solving the hydrodynamic Eq. (2) by using Bjorken flow Eq. (3) and ideal
energy-momentum tensor Eq. (4) we get PhysRevD.27.140

$\displaystyle\partial_{\tau}\epsilon+\frac{\epsilon+p}{\tau}$
$\displaystyle=0~{},$ (5)

where $\tau$ is the proper time. Now, using ideal equation of state
$p=\frac{\epsilon}{3}\propto T^{4}$, we can obtain cooling rate for ideal
hydrodynamics as,

$\displaystyle T=T_{0}{\left(\frac{\tau_{0}}{\tau}\right)}^{\frac{1}{3}}~{}.$
(6)

Where, initial condition is $T(\tau_{0})=T_{0}$, with initial temperature
$T_{0}$ at the medium formation time $\tau_{0}$.

####  2.1.2 Ideal magnetohydrodynamics

In peripheral RHICs, a huge magnetic field is created, which can magnetize the
QGP. The energy-momentum tensor for an ideal magnetized fluid can be expressed
as ROY201545 ; PhysRevD.109.014018 ; Panda2021

$\displaystyle{T}^{\mu\nu}=\left(\epsilon+p+{B}^{2}\right){u}^{\mu}{u}^{\nu}-\left(p+\frac{{B}^{2}}{2}\right){g}^{\mu\nu}-{B}^{\mu}{B}^{\nu}~{}.$
(7)

Depending on the medium’s electrical conductivity, the magnetic field created
during the collision will decay with time. Here, we consider an exponential
decay profile Satow:2014lia

$\displaystyle B(\tau)=B_{0}\exp{[-\tau/\tau_{B}]}~{}.$ (8)

where $B_{0}$ is the magnitude of initial magnetic field and $\tau_{B}$ is its
decay parameter. The field is directed along the $y$ axis in the transverse
plane, considering the longitudinal z-axis to align with the beam axis. Now,
solving the hydrodynamic equation with ideal EoS, and after some
rearrangement, we get

$\displaystyle\frac{\partial T}{\partial\tau}+p(\tau)T=q(\tau)T^{-3}~{},$ (9)

where p$(\tau)$ = $\frac{1}{3\tau}$, q$(\tau)$ =
($\alpha\beta-\frac{\alpha}{\tau}$)exp(-2$\beta\tau$ ) with
$\alpha=\frac{B_{0}^{2}}{12a}$, $\beta=\frac{1}{\tau_{B}}$, and
$a=\left(16+\frac{21}{2}N_{f}\right)\frac{\pi^{2}}{90}$. The above equation
can be recognized as the Bernoulli differential equation and can be
analytically solved to obtain

$\displaystyle T$
$\displaystyle=\Big{[}T_{0}^{4}{\left(\frac{\tau_{0}}{\tau}\right)}^{\frac{4}{3}}+\frac{4\alpha}{{(2\beta\tau)}^{\frac{4}{3}}}\Big{\\{}\Gamma{(4/3,2\beta\tau)}-\Gamma(4/3,2\beta\tau_{0})\Big{\\}}$
$\displaystyle~{}~{}~{}~{}-\frac{2\alpha}{{(2\beta\tau)}^{\frac{7}{3}}}\Big{\\{}\Gamma(7/3,2\beta\tau)-\Gamma(7/3,2\beta\tau_{0})\Big{\\}}\Big{]}^{\frac{1}{4}}~{}.$
(10)

###  2.2 Gubser flow

Gubser generalizes the Bjorken flow to include radial expansion of the medium
and longitudinal boost-invariant expansion PhysRevD.82.085027 . Working with
rapidity or the Milne coordinate system is convenient in heavy-ion collisions.
Milne coordinates ($\tau,\eta,r,\phi$), can be related to the Cartesian
coordinates $(t,x,y,z)$ as,

$\displaystyle\tau=\sqrt{t^{2}-z^{2}},~{}~{}~{}~{}\eta=arctanh\frac{z}{t},~{}~{}~{}~{}~{}r=\sqrt{x^{2}+y^{2}},~{}~{}~{}~{}\phi=arctan\frac{y}{x}.$
(11)

In this coordinate system, the four-velocity $u^{\mu}$ constructed from
symmetry consideration with boost, rotation, and reflection invariance
($\eta\rightarrow-\eta$) is given as

$\displaystyle
u^{\tau}=\cosh{\kappa}=\gamma_{r},~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}\frac{u^{r}}{u^{\tau}}=\tanh{\kappa}=v_{r},~{}~{}~{}~{}~{}~{}u^{\eta}=u^{\phi}=0.$
(12)

Where Lorentz factor $\gamma_{r}$ = $\frac{1}{\sqrt{1-v_{r}^{2}}}$ and $v_{r}$
= $\sqrt{\vec{v}_{x}^{2}+\vec{v}_{y}^{2}}$, $v_{r}$ is known as transverse
velocity. $\kappa$ is parameterized as
$\kappa(\tau,r)={arctanh}(\frac{2q^{2}\tau r}{1+q^{2}\tau^{2}+q^{2}r^{2}})$.
Here, we take the case of the viscous medium while preserving the conformal
invariance of the theory, that is, $p=\epsilon/3$ and shear viscosity $\eta$=
$H_{0}\epsilon^{3/4}$. $H_{0}$ is a dimensionless quantity. Using the Gubser
flow, we can solve the hydrodynamic equation for viscous medium to obtain the
temperature profile as PhysRevD.82.085027 ,

$\displaystyle T$ $\displaystyle={1\over\tau
f_{*}^{1/4}}\left[{\frac{{T}_{0}}{(1+g^{2})^{1/3}}}+\frac{{\rm
H}_{0}g}{\sqrt{1+g^{2}}}\left\\{1-(1+g^{2})^{1/6}{}_{2}F_{1}\left({1\over
2},{1\over 6};{3\over 2};-g^{2}\right)\right\\}\right]~{}.$ (13)

Where ${T}_{0}$ is an integration constant and ${}_{2}F_{1}$ denotes a
hypergeometric function with $g={1-q^{2}\tau^{2}+q^{2}r^{2}\over 2q\tau}$. For
numerical estimations, we used semi-realistic numbers for a central gold-gold
collision at $\sqrt{s_{\rm NN}}=200\,{\rm GeV}$ are $\hat{T}_{0}=5.55$ and
${\rm H}_{0}=0.33$ if we choose $1/q=4.3$ fm.

Figure 1: Temperature ($T$) as a function of proper time ($\tau$) for ideal
hydrodynamics, ideal magnetohydrodynamics (MHD) with two decay parameters at
$\tau_{B}=3,~{}7$ fm, and Gubser flow at $r=0,~{}3$ fm.

In Fig. (1), we have plotted the cooling rate (temperature as a function of
proper time) for all the three hydrodynamic cases discussed above. Note that
the initial conditions are free to adjust within the Heisenberg uncertainty
principle. The cooling rates are represented as follows. The ideal
hydrodynamics at $B=0$ with a black solid line, ideal magnetohydrodynamics
(MHD) at $\tau_{B}=3,~{}7$ fm with the red dashed-triple dot line and green
dash-dot line, respectively, and for Gubser flow at $r=0,~{}3$ fm by magenta
dashed and blue dashed line, respectively. Interestingly, the temperature in
Gubser flow falls faster due to the presence of transverse flow even at the
finite viscosity of the medium. The cooling rates plotted here will be helpful
in the analysis of results obtained in the subsequent sections.

##  3 Head-on collision

In this section, we calculate the thermoelectric coefficient for the QGP
medium created in a head-on collision. Due to the absence of spectator
current, no magnetic field is expected to be created in the head-on
collisions. Here, we will show that an electric field can be created or
induced in the medium due to the thermoelectric effect, even though there is
no initial magnetic field during the collision.

###  3.1 Seebeck coefficient in evolving QGP

QGP created in RHICs cools down rapidly as it goes through space-time
evolution with rapid expansion. This rapid cooling affects the medium’s
transport coefficients, such as thermal conductivity and thermoelectric
coefficients. In all the earlier studies of thermoelectric coefficients
Das2021 ; Abhishek:2020wjm ; Kurian:2021zyb ; PhysRevD.102.096011 , a static
picture of QGP is considered where the effect of temperature evolution was not
taken into account, which is a more realistic scenario. In Ref.
PhysRevD.109.014018 , for the first time, we have estimated the effect of
temperature evolution on the thermal conductivity of QGP. Here, we will
calculate the thermoelectric coefficient for evolving QGP using a similar
approach of Ref. PhysRevD.109.014018 . We will also show the difference
between the expressions obtained considering a static QGP.

To study the fluid properties of QGP medium, we consider a system of
relativistic fluid consisting of particles with mass $m_{i}$, momentum
$\vec{k_{i}}$, energy $\omega_{i}=\sqrt{\vec{k}_{i}^{2}+m_{i}^{2}}$, and
chemical potential $\mu_{i}={\rm b}_{i}\mu_{B}$ for ith species, with bi as
baryon quantum number and $\mu_{B}$ as baryon chemical potential. The single
particle distribution function at equilibrium for ith species is

$\displaystyle f^{0}_{i}=\frac{1}{e^{\frac{\omega_{i}-{\rm
b}_{i}\mu_{B}}{T}}\pm 1}~{},$ (14)

where $\pm$ stands for fermion and boson, respectively. The total single-
particle distribution function ($f_{i}$) for a system slightly out of
equilibrium ($\delta f_{i}$) can be written as $f_{i}=f^{0}_{i}+\delta f_{i}$.
When a system is out of equilibrium, net non-zero currents arise in the
system. In kinetic theory, electric current density for such a system can be
expressed as PhysRevD.108.094007

$\displaystyle\vec{j}=\sum_{i}q_{i}g_{i}\int\frac{d^{3}|\vec{k}_{i}|}{(2\pi)^{3}}\frac{\vec{k_{i}}}{\omega_{i}}\delta
f_{i}~{}.$ (15)

Here, $q_{i}$ is the electric charge, and $g_{i}$ is the degeneracy of the ith
species particles. Similarly, the microscopic definition of heat flow
PhysRevD.108.094007 ; PhysRevD.109.014018 ,

$\displaystyle{\vec{I}}_{i}$
$\displaystyle=\sum_{i}g_{i}\int\frac{d^{3}|\vec{k}_{i}|}{(2\pi)^{3}}\frac{\vec{k}_{i}}{\omega_{i}}(\omega_{i}-{\rm
b}_{i}h)\delta f_{i}~{}.$ (16)

Where $h$ is the total enthalpy per particle. To find the expression of
$\delta f_{i}$, we solve the Boltzmann transport equation (BTE) with the help
of relaxation time approximation (RTA) PhysRevD.102.016016

$\displaystyle\frac{\partial f_{i}}{\partial
t}+\frac{\vec{k}_{i}}{\omega_{i}}\cdot\frac{\partial
f_{i}}{\partial\vec{x_{i}}}+q_{i}\vec{E}\cdot\frac{\partial
f_{i}}{\partial\vec{k_{i}}}=-\frac{\delta f_{i}}{\tau^{i}_{R}}~{},$ (17)

where $\tau^{i}_{R}$ is the relaxation time of the particle. $\vec{E}$
represents the external electric field. Here, the deviation of the
distribution function from equilibrium is driven by the electric field and
temperature gradients. With leading order contribution, we can consider an
ansatz of $\delta f_{i}$ as PhysRevD.104.094037

$\displaystyle\delta f_{i}=(\vec{k_{i}}\cdot\vec{\Omega})\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}~{}.$ (18)

As QGP is a rapidly cooling system, to account for all leading order
contributions from gradient forces, a general form of unknown vector
${\vec{\Omega}}$ can be assumed as

$\displaystyle\vec{\Omega}=$
$\displaystyle~{}\alpha_{1}\vec{E}+\alpha_{2}\vec{\nabla}T+\alpha_{3}\vec{\nabla}\dot{T}~{}.$
(19)

The unknown coefficients $\alpha_{j}$ ($j=1,2,3$) determine the strength of
the respective gradient force fields driving the system away from equilibrium.
The coefficient $\alpha_{3}$ corresponds to the term that arises due to the
evolving picture of the medium. For the case of a static picture, where the
effects of cooling are not considered, only $\alpha_{1}$ and $\alpha_{2}$
contribute, whereas $\alpha_{3}$ vanishes in Eq. (19) as in Ref. Das2021 ;
Kurian:2021zyb . We can solve Eq. (17) using Eqs. (18) and (19), and get the
expression of $\alpha_{j}$’s. After getting the expressions for
$\alpha_{j}$’s, the simplified form of $\delta f_{i}$ is (see Appendix A),

$\displaystyle\delta f_{i}$
$\displaystyle=-q_{i}\tau_{R}^{i}(\vec{v_{i}}.\vec{E})+\frac{(\omega_{i}-{\rm
b}_{i}h)\tau_{R}^{i}}{T}\Big{[}\Big{(}\vec{v_{i}}\cdot\vec{\nabla}{T}\Big{)}-\tau_{R}^{i}\Big{(}\vec{v_{i}}\cdot\vec{\nabla}{\dot{T}}\Big{)}\Big{]}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}.$ (20)

Using the expressions of $\delta f_{i}$ in Eq. (15), we can express electric
current as

$\displaystyle j^{l}$
$\displaystyle=\sum_{i}\frac{q_{i}g_{i}}{3}\int\frac{d^{3}|\vec{k}_{i}|}{(2\pi)^{3}}v_{i}^{2}\tau_{R}^{i}\Big{[}-q_{i}E^{l}+\frac{(\omega_{i}-{\rm
b}_{i}h)}{T}\Big{\\{}\frac{\partial{T}}{\partial
x^{l}}-\tau_{R}^{i}\frac{\partial{\dot{T}}}{\partial
x^{l}}\Big{\\}}\Big{]}\frac{\partial f^{0}_{i}}{\partial\omega_{i}}.$ (21)

The components of the electric current in the x,y, and z directions are given
as,

$\displaystyle j^{x}$
$\displaystyle=\sum_{i}\frac{q_{i}g_{i}}{3}\int\frac{d^{3}|\vec{k}_{i}|}{(2\pi)^{3}}v_{i}^{2}\tau_{R}^{i}\Big{[}-q_{i}E^{x}+\frac{(\omega_{i}-{\rm
b}_{i}h)}{T}\Big{\\{}\frac{\partial{T}}{\partial
x}-\tau_{R}^{i}\frac{\partial{\dot{T}}}{\partial
x}\Big{\\}}\Big{]}\frac{\partial f^{0}_{i}}{\partial\omega_{i}}.$
$\displaystyle j^{y}$
$\displaystyle=\sum_{i}\frac{q_{i}g_{i}}{3}\int\frac{d^{3}|\vec{k}_{i}|}{(2\pi)^{3}}v_{i}^{2}\tau_{R}^{i}\Big{[}-q_{i}E^{y}+\frac{(\omega_{i}-{\rm
b}_{i}h)}{T}\Big{\\{}\frac{\partial{T}}{\partial
y}-\tau_{R}^{i}\frac{\partial{\dot{T}}}{\partial
y}\Big{\\}}\Big{]}\frac{\partial f^{0}_{i}}{\partial\omega_{i}}.$
$\displaystyle j^{z}$
$\displaystyle=\sum_{i}\frac{q_{i}g_{i}}{3}\int\frac{d^{3}|\vec{k}_{i}|}{(2\pi)^{3}}v_{i}^{2}\tau_{R}^{i}\Big{[}-q_{i}E^{z}+\frac{(\omega_{i}-{\rm
b}_{i}h)}{T}\Big{\\{}\frac{\partial{T}}{\partial
z}-\tau_{R}^{i}\frac{\partial{\dot{T}}}{\partial
z}\Big{\\}}\Big{]}\frac{\partial f^{0}_{i}}{\partial\omega_{i}}.$ (22)

We introduce the following integrals to write Eqs. (3.1) in compact form as

$\displaystyle
L_{1i}=\frac{g_{i}}{3T}\int\frac{d^{3}|\vec{k}_{i}|}{(2\pi)^{3}}\frac{\vec{k}^{2}_{i}}{\omega_{i}^{2}}\tau_{R}^{i}f^{0}_{i}(1-f^{0}_{i})~{},$
$\displaystyle
L_{2i}=\frac{g_{i}}{3T}\int\frac{d^{3}|\vec{k}_{i}|}{(2\pi)^{3}}\frac{\vec{k}^{2}_{i}}{\omega_{i}^{2}}(\omega_{i}-{\rm
b}_{i}h)\tau_{R}^{i}f^{0}_{i}(1-f^{0}_{i})~{},$ $\displaystyle
L_{3i}=\frac{g_{i}}{3T}\int\frac{d^{3}|\vec{k}_{i}|}{(2\pi)^{3}}\frac{\vec{k}^{2}_{i}}{\omega_{i}^{2}}(\omega_{i}-{\rm
b}_{i}h)\tau_{R}^{i2}f^{0}_{i}(1-f^{0}_{i}).$ (23)

Hence, Eqs. (21) can be written as,

$\displaystyle j^{l}$
$\displaystyle=\sum_{i}q_{i}^{2}L_{1i}E_{l}-\frac{1}{T}\sum_{i}q_{i}\Big{(}L_{2i}-L_{3i}\frac{d}{dt}\Big{)}\frac{dT}{dx^{l}}~{}.$
(24)

In above Eq., by setting $j_{x}$ = $j_{y}$ = $j_{z}$=0, we can get $E_{x}$,
$E_{y}$ and $E_{z}$ in terms of temperature gradients $\frac{dT}{dx}$,
$\frac{dT}{dy}$ and $\frac{dT}{dz}$. Here, we rearrange the Eqs. (24) in the
matrix form as

$\displaystyle\mathbb{\bm{\sigma}~{}E=L~{}X},$ (25)

where,

$\displaystyle\bm{\sigma}=\begin{pmatrix}\sigma_{e}&0&0\\\ 0&\sigma_{e}&0\\\
0&0&\sigma_{e}\\\ \end{pmatrix},~{}\hskip
28.45274pt\mathbb{E}=\begin{pmatrix}E_{x}\\\ E_{y}\\\ E_{z}\end{pmatrix},~{}$
$\displaystyle\mathbb{L}=\frac{1}{T}\begin{pmatrix}L_{23}&0&0\\\ 0&L_{23}&0\\\
0&0&L_{23}\\\
\end{pmatrix},~{}\mathbb{X}=\begin{pmatrix}\frac{dT}{dx}\vspace{2mm}\\\
\frac{dT}{dy}\vspace{2mm}\\\ \frac{dT}{dz}\\\ \end{pmatrix}.$

Hence, after finding the inverse of $\bm{\sigma}$ matrix, we can obtain the
components of the electric field as,

$\displaystyle\mathbf{E}$ $\displaystyle=\mathbf{(\bm{\sigma}^{-1}L)X}$
$\displaystyle\Rightarrow\begin{pmatrix}E_{x}\\\ E_{y}\\\ E_{z}\end{pmatrix}$
$\displaystyle=\begin{pmatrix}S&0&0\\\ 0&S&0\\\ 0&0&S\\\
\end{pmatrix}\begin{pmatrix}\frac{dT}{dx}\vspace{2mm}\\\
\frac{dT}{dy}\vspace{2mm}\\\ \frac{dT}{dz}\\\ \end{pmatrix}~{}.$ (26)

Where $\sigma_{el}$ = $\sum_{i}q_{i}^{2}L_{1i}$ is electrical conductivity.
Further, we defined the integrals $L_{23}$ =
$\sum_{i}q_{i}\Big{(}L_{2i}-L_{3i}\frac{d}{dt}\Big{)}$. Note that it is a
generalized study of thermoelectric coefficients, i.e., we did not consider
any preferred direction of the temperature gradient or a particular cooling
rate. Therefore, we get the Seebeck coefficient($S$) in three directions.
Here, we identify the Seebeck coefficient $S$ as,

$\displaystyle S$
$\displaystyle=\frac{(\frac{L_{23}}{T^{2}})}{(\frac{\sigma_{el}}{T})}~{}.$
(27)

Due to the boost invariance, $T$ is invariant along the longitudinal $z$ axis,
therefore $\frac{\partial T}{\partial z}=0$. Hence, only $x$ and $y$
components of the electric field $\vec{E}$ survive in the transverse plane.

####  3.1.1 Analytical solution for massless non-interacting system

Here, we address the analytic solution of Seebeck coefficient Eq. (27) for
massless ($\omega_{i}=k_{i}$) and non-interacting system. In this case, the
Seebeck coefficient can be solved analytically using the Riemann zeta function
$\zeta(s)$ and Dirichlet eta $\eta(s)$ function:

$\displaystyle\zeta(s)=\frac{1}{\Gamma(s)}\int_{0}^{\infty}\frac{x^{s-1}}{e^{x}-1}dx~{},~{}~{}$
$\displaystyle~{}~{}\eta(s)=\frac{1}{\Gamma(s)}\int_{0}^{\infty}\frac{x^{s-1}}{e^{x}+1}dx~{},$
(28)

where, $\Gamma(s)$ is the Gamma function. Now, using the relation
$\eta(s)=1-2^{(1-s)}\zeta(s)$, Eqs. (3.1) can be write in terms of the
$\zeta(s)$ as,

$\displaystyle L_{1i}$
$\displaystyle=\sum_{i}\frac{g_{i}\tau_{R}^{i}T}{6\pi^{2}}\Big{\\{}2\zeta(2)+\Big{(}\frac{b_{i}\mu_{B}}{T}\Big{)}^{2}\Big{\\}},$
$\displaystyle L_{2i}$
$\displaystyle=\sum_{i}\frac{g_{i}\tau_{R}^{i}T^{2}}{6\pi^{2}}\Big{(}\frac{b_{i}^{3}\mu_{B}^{3}}{T^{2}}+6b_{i}\mu_{B}\zeta(2)-2b_{i}h\zeta(2)-\frac{b_{i}^{3}\mu_{B}^{2}h}{T^{2}}\Big{)},$
$\displaystyle L_{3i}$
$\displaystyle=\sum_{i}\frac{g_{i}\tau_{R}^{i2}T^{2}}{6\pi^{2}}\Big{(}\frac{b_{i}^{3}\mu_{B}^{3}}{T^{2}}+6b_{i}\mu_{B}\zeta(2)-2b_{i}h\zeta(2)-\frac{b_{i}^{3}\mu_{B}^{2}h}{T^{2}}\Big{)}.$
(29)

Hence, using the values of $L_{1i,2i,3i}$ we get

$\displaystyle\sigma_{el}$
$\displaystyle=\sum_{i}\frac{g_{i}q_{i}^{2}\tau_{R}^{i}T}{6\pi^{2}}\Big{\\{}2\zeta(2)+\Big{(}\frac{b_{i}\mu_{B}}{T}\Big{)}^{2}\Big{\\}},$
$\displaystyle L_{23}$
$\displaystyle=\sum_{i}q_{i}~{}\frac{g_{i}\tau_{R}^{i}T^{2}}{6\pi^{2}}\Big{(}\frac{b_{i}^{3}\mu_{B}^{3}}{T^{2}}+6b_{i}\mu_{B}\zeta(2)-2b_{i}h\zeta(2)-\frac{b_{i}^{3}\mu_{B}^{2}h}{T^{2}}\Big{)}\left(1-\tau_{R}^{i}\frac{d}{dt}\right).$
(30)

Note that, $L_{23}$ is proportional to electric charge $q_{i}$. Degeneracy is
the same for all three flavors, $g_{i}=6$ (for a spin and color degeneracy).
$\tau_{R}^{i}=\tau_{R}$ (say) is the same for all the flavors for massless and
non-interacting systems of particles. Therefore, the charge of the three
flavors (up, down, and strange) adds up to zero, i.e., $L_{23}=0$. As a
result, the Seebeck coefficient for the massless non-interacting system is

$\displaystyle S=\frac{(\frac{L_{23}}{T^{2}})}{(\frac{\sigma_{el}}{T})}=0~{}.$
(31)

####  3.1.2 Numerical solution for QGP EoS

Figure 2: Seebeck coefficient ($S$) as a function of proper time ($\tau$) at
$\mu_{B}=0.3$ GeV. Left: Upper panel represents $S$ for Bjorken flow, bottom
panel represents its percentage deviation of static from evolving picture
($=\frac{(S^{e}-S^{s})}{S^{e}}\times 100\Big{)}\%$. Right: Similar to the left
plot but for Gubser flow at two different values of $r=0,~{}3$ fm.

We use a quasiparticle model formulated by Gorenstein and Yang
Gorenstein:1995vm  for the numerical estimation. It is a phenomenological
model where the lattice QCD equation of state for QGP is achieved by
considering the thermal masses of the partons. The thermal mass $m(T)$ arises
from the interactions among the partons. The thermodynamic consistency is
achieved by introducing a bag constant arising from vacuum energy
Srivastava:2010xa . The dispersion relation of the particle having energy
$\omega_{i}$ and momentum $k_{i}$ is
$\omega_{i}^{2}(k_{i},T)=k_{i}^{2}+m_{i}^{2}(T)$. Where $m_{i}$ is the total
effective mass of ith quark flavor and can be parameterized as

$\displaystyle m_{i}^{2}=m_{i0}^{2}+\sqrt{2}m_{i0}m_{iT}+m_{iT}^{2}~{}.$ (32)

$m_{i0}$ and $m_{iT}$ are bare mass and thermal mass of the ith flavor with,

$\displaystyle
m_{iT}^{2}(T,\mu_{B})=\frac{N_{c}^{2}-1}{8N_{c}}\Big{(}T^{2}+\frac{\mu_{B}^{2}}{9\pi^{2}}\Big{)}g^{2}(T).$
(33)

The effective mass of gluon $(m_{g})$ in this model can be represented as

$\displaystyle
m_{g}=\frac{N_{c}}{6}g^{2}(T)T^{2}\left(1+\frac{N_{f}+\frac{\mu_{B}^{2}}{\pi^{2}T^{2}}}{6}\right).$
(34)

$N_{c}$ represents the number of color degrees of freedom and
$g^{2}(T)=4\pi\alpha_{s}(T)$, $\alpha_{s}(T)$ is running coupling constant.
For more details of the model, one can follow references Gorenstein:1995vm ;
Srivastava:2010xa ; Bannur:2006hp ; Peshier:1995ty ; Peshier:1999ww . Note
that, for the estimation of thermoelectric coefficients and the induced
electric field, gluons contribute only to the enthalpy of the system. Now, for
the relaxation time of quarks, we use a momentum-independent expression
obtained for QCD matter Hosoya:1983xm

$\displaystyle\tau_{R}=\frac{1}{5.1T\alpha_{s}^{2}\log(1/\alpha_{s})[1+0.12(2N_{f}+1)]}.$
(35)

For numerical estimation, the value of the strong coupling constant is taken
to be fixed at $\alpha_{s}=0.5$.

In Fig. (2), we plot the Seebeck coefficient ($S$) as a function of proper
time ($\tau$) at $\mu_{B}$ = 0.3 GeV in head-on collisions (or $eB=0$) for the
static and evolving case of QGP. Plots in the left panels are obtained in
Bjorken flow with the ideal hydrodynamic cooling rate, which in the right
panels are for Gubser flow. In both cases, a decreasing (magnitude) trend of
$S$ with time is observed, and it approaches the massless non-interacting
(black dotted line) limit. Note that the negative sign in the value of $S$
implies that the direction of the induced electric field is opposite to the
temperature gradient. With time, the temperature of the medium falls off, and
$S$ reduces in both the hydrodynamic pictures. $S$ for Gubser flow reduces
comparatively faster because of its cooling rate. In both the hydrodynamic
pictures, $S$ beyond a certain time where $T<0.150$ GeV is not justified in
the used quasiparticle model. From Eq. (27), $S$ is the ratio of two different
quantities; for the static case where integral $L_{23}$ reduces to $L_{2i}$,
the coefficient $S$ becomes independent of relaxation time $\tau_{R}$. This
happens because we consider the same relaxation time for all the quarks
regardless of their flavor. Whereas, for the evolving picture, $S$ shows its
dependence on $\tau_{R}$ through integral $\L_{23}$. In the bottom panels of
Fig. (2), the percentage deviation of static from the evolving picture is
plotted against $\tau$ to quantify the cooling effect. In both the
hydrodynamic cases, $S$ is enhanced due to temperature evolution. The effect
is high in the early times and decreases with time. Moreover, this is higher
in the Gubser flow case, which shows the impact is sensitive to the cooling
rates.

###  3.2 Induced electric field in isotropic system

Figure 3: Electric field induced in the QGP in the head-on collisions.

To estimate the induced electric field ($\vec{E}$) in the head-on collisions,
we considered baryon chemical potential $\mu_{B}=0.3$ GeV, attained at RHIC
energies. Here, we used the Gubser cooling rate for the temperature gradient.
From Gubser flow, $T=T(\tau,~{}r)$, with $r=\sqrt{x^{2}+y^{2}}$. Due to Gubser
symmetry considerations as mentioned in subsection (2.2), $T$ is symmetric in
the transverse ($x-y$) plane and remains invariant along the longitudinal
direction or $z$ axis. Here, we used the Seebeck coefficient calculated from
the Bjorken flow for the case of evolving QGP. In head-on collisions, there is
no spectator current; the magnetic field is not produced during the collision.
In this paper, by the terms “isotropy” and “anisotropy,” we indicate momentum
isotropy or anisotropy due to magnetic field only. We will discuss the induced
electric field created in an isotropic QGP (without an external magnetic
field) due to the thermoelectric effect.

In Fig. (3), we have shown the vector plot of the induced electric field in
the transverse plane for head-on collision. In an isotropic medium, the
thermoelectric transport coefficient matrix is isotropic as obtained in Eq.
(3.1), only diagonal elements, Seebeck coefficient ($S$) are non-zero. Due to
symmetry consideration, $T$ and its gradient are symmetric in the transverse
($x-y$) plane. Moreover, because of boost invariant symmetry along $z$ axis,
$\frac{\partial T}{\partial z}=0$. Therefore, only $x$ and $y$ components of
the electric field survive, and $\vec{E}$ is symmetric in the transverse
plane. Time evolution of $\vec{E}$ is shown by three panels in the Fig. (3),
from left to right $\tau=0.5,~{}1,~{}3$ fm. The direction and magnitude of
$\vec{E}$ are represented in the figures by the arrow tip and its length,
respectively. Magnitude is also colored for clarity; the color bar on the top
of the figures is given for magnitude in the unit of $m_{\pi}^{2}$. $\vec{E}$
is zero at the center (as $\vec{\nabla}T$ vanishes) and directed radially
outward with increasing magnitude. At $\tau=0.5$ fm, the induced field as high
as $eE=0.6~{}m_{\pi}^{2}$, and its strength decreases as the medium evolves in
time. The direction of $\vec{E}$ can be understood as follows. From Eq. (3.1),
$E_{x}=S\frac{dT}{dx}$, $E_{y}=S\frac{dT}{dy}$. From the right panel of Fig.
(2), $S$ is negative at the early time. $\vec{\nabla}T$ is also negative, as
$T$ decreases as we go away from the center ($x=y=0,~{}{\rm or},~{}r=0$). This
explains the direction of $\vec{E}$. Value of $S$, $\vec{\nabla}T$, and hence
$\vec{E}$ is sensitive to the cooling rate. As an analytical expression of the
cooling rate is required for the calculation of thermoelectric coefficients
and induced field, we used Gubser flow. However, in (3+1)D hydrodynamic
simulations, one can incorporate the current study to estimate induced fields
in QGP produced in heavy-ion collisions.

##  4 Peripheral collision

In this section, we calculate the thermoelectric coefficients of an evolving
relativistic fluid in the presence of an external time-varying magnetic field.
For detailed calculations, see appendix A.

###  4.1 Magneto-Seebeck and Nernst coefficients in evolving QGP

The magnetic field breaks the rotational symmetry of the system, due to which
the thermoelectric coefficient matrix becomes anisotropic with two independent
elements - magneto-Seebeck and Nernst coefficients. Moreover, according to the
Landau quantization, the energy level of charged particles gets quantized in
the presence of the magnetic field. However, in the RTA formalism, Landau
quantization is not included. Therefore, deriving the transport coefficients
in the presence of a magnetic field in the upcoming subsection, we introduce
phase-space quantization with the modified dispersion relation for the Landau
quantization in the subsequent subsection.

####  4.1.1 Without Landau quantization

In the presence of a magnetic field, to find the expression of the deviated
part of distribution function $\delta f_{i}$, which is driven by the
temperature gradient, we use the RTA method with a Lorentz force term in BTE.
In the presence of an external electromagnetic field, BTE under RTA can be
expressed as PhysRevD.102.016016

$\displaystyle\frac{\partial f_{i}}{\partial
t}+\frac{\vec{k}_{i}}{\omega_{i}}\cdot\frac{\partial
f_{i}}{\partial\vec{x_{i}}}+q_{i}\left(\vec{E}+\frac{\vec{k}_{i}}{\omega_{i}}\times\vec{B}\right)\cdot\frac{\partial
f_{i}}{\partial\vec{k_{i}}}=-\frac{\delta f_{i}}{\tau^{i}_{R}}~{},$ (36)

where $\tau^{i}_{R}$ is the relaxation time of the particle. We consider a
time-varying electromagnetic field of the form Satow:2014lia ; Hongo:2013cqa

$\displaystyle B=B_{0}\exp{\left(-\frac{\tau}{\tau_{B}}\right)},~{}~{}$
$\displaystyle~{}~{}E=E_{0}\exp{\left(-\frac{\tau}{\tau_{E}}\right)},$ (37)

where $B_{0},~{}E_{0}$ are the magnitudes of the initial fields having decay
parameters of $\tau_{B}$ and $\tau_{E}$, respectively, and $\tau$ is the
proper time. The deviation of the distribution function from equilibrium is
driven by the electromagnetic field, so leading order contribution in $\delta
f_{i}$ can be PhysRevD.104.094037

$\displaystyle\delta f_{i}=(\vec{k_{i}}\cdot\vec{\Omega})\frac{\partial
f^{0}_{i}}{\partial\omega}~{}.$ (38)

The general form of unknown vector ${\vec{\Omega}}$ can be assumed as (up to
first-order in time derivative)

$\displaystyle\vec{\Omega}=$
$\displaystyle~{}\alpha_{1}\vec{E}+\alpha_{2}\dot{\vec{E}}+\alpha_{3}\vec{\nabla}T+\alpha_{4}\vec{\nabla}\dot{T}+\alpha_{5}(\vec{\nabla}T\times\vec{B})+\alpha_{6}(\vec{\nabla}T\times\dot{\vec{B}})$
$\displaystyle+\alpha_{7}(\vec{\nabla}\dot{T}\times\vec{B})+\alpha_{8}(\vec{E}\times{\vec{B}})+\alpha_{9}(\vec{E}\times{\dot{\vec{B}}})+\alpha_{10}(\dot{\vec{E}}\times{\vec{B}}),$
(39)

where the unknown coefficients $\alpha_{j}$ ($j=1,2,..7$) are helpful to
determine the strength of the respective field in driving the system away of
the equilibrium. The contribution of coefficients $\alpha_{4}$ and
$\alpha_{7}$ is due to the considerations of the cooling effects, i.e., for
the case of the evolving picture of the medium. Whereas for the static
picture, both these coefficients will vanish. Similarly, for the constant
electromagnetic field the coefficients $\alpha_{2}$, $\alpha_{6}$,
$\alpha_{9}$, and $\alpha_{10}$ vanishes. The terms $\vec{B}$, $\dot{\vec{B}}$
do not contribute to the current for the case where the chiral chemical
potential is zero Satow:2014lia .

In the current work, we introduce the cooling rate using hydro-dynamical
theories and study the thermoelectric response of evolving QGP in the presence
of a slowly varying magnetic field, where we can consider that the inverse of
cyclotron frequency is approximately equal to the magnetic field decay
parameter, i.e., $\tau_{B}\approx\frac{\omega_{i}}{q_{i}B}$. Therefore, we can
solve Eq. (36) using Eqs. (38) and (4.1.1), and get the expression of
$\alpha_{j}$’s. After getting the expressions for $\alpha_{j}$’s, the
simplified form of $\delta f_{i}$ is(see Appendix A),

$\displaystyle\delta f_{i}$
$\displaystyle=\frac{-q_{i}\tau_{R}^{i}}{(1+\chi_{i})(1+\chi_{i}+\chi_{i}^{2})}\Big{[}\Big{\\{}(1+\chi_{i})+\frac{\chi_{i}(1+\chi_{i}-\chi_{i}^{2})}{(1+\chi_{i}^{2})}\Big{\\}}(\vec{v_{i}}.\vec{E})\pm\Big{\\{}\chi_{i}(1+\chi_{i})$
$\displaystyle+\frac{\chi_{i}^{2}(2+\chi_{i})}{(1+\chi_{i}^{2})}\Big{\\}}\Big{(}\vec{v_{i}}.(\vec{E}\times\hat{b})\Big{)}\Big{]}+\frac{(\omega_{i}-{\rm
b}_{i}h)\tau_{R}^{i}}{T(1+\chi_{i})(1+\chi_{i}+\chi_{i}^{2})}\Big{[}(1+\chi_{i})\Big{(}\vec{v_{i}}\cdot\vec{\nabla}{T}\Big{)}$
$\displaystyle-\tau_{R}^{i}\frac{(1+\chi_{i}-\chi_{i}^{2})}{(1+\chi_{i}^{2})}\Big{(}\vec{v_{i}}\cdot\vec{\nabla}{\dot{T}}\Big{)}\pm\chi_{i}(1+\chi_{i})\Big{(}\vec{v_{i}}\cdot(\vec{\nabla}{T}\times\hat{b})\Big{)}\mp\tau_{R}^{i}\frac{\chi_{i}(2+\chi_{i})}{(1+\chi_{i}^{2})}$
$\displaystyle\Big{(}\vec{v_{i}}\cdot(\vec{\nabla}{\dot{T}}\times\hat{b})\Big{)}\Big{]}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}.$ (40)

For simplicity we have considered $\tau_{E}=\tau_{B}$ and
$\chi_{i}=\frac{\tau_{R}^{i}}{\tau_{B}}=\frac{\tau_{R}^{i}}{\tau_{E}}$. Here
$\pm$ indicates a positive sign for positively charged particles
(antiparticles) and a negative sign for negatively charged particles
(antiparticles). Using the above expressions of $\delta f_{i}$ in Eq. (15) and
(16), we can express electric current as (see Appendix A)

$\displaystyle j^{l}$
$\displaystyle=\sum_{i}\frac{q_{i}g_{i}}{3}\int\frac{d^{3}|\vec{k}_{i}|}{(2\pi)^{3}}v_{i}^{2}\frac{\tau_{R}^{i}}{(1+\chi_{i}+\chi_{i}^{2})(1+\chi_{i})}\Big{[}-q_{i}\Big{\\{}\Big{(}(1+\chi_{i})+\frac{\chi_{i}(1+\chi_{i}-\chi_{i}^{2})}{(1+\chi_{i}^{2})}\Big{)}\delta^{lk}E^{k}$
$\displaystyle\pm\Big{(}\chi_{i}(1+\chi_{i})+\frac{\chi_{i}^{2}(2+\chi_{i})}{(1+\chi_{i}^{2})}\Big{)}\epsilon^{ljk}h^{j}E^{k}\Big{\\}}+\frac{(\omega_{i}-{\rm
b}_{i}h)}{T}\Big{\\{}(1+\chi_{i})\delta^{lk}\frac{\partial{T}}{\partial
x^{k}}$
$\displaystyle-\tau_{R}^{i}\frac{(1+\chi_{i}-\chi_{i}^{2})}{(1+\chi_{i}^{2})}\delta^{lk}\frac{\partial{\dot{T}}}{\partial
x^{k}}\pm\chi_{i}(1+\chi_{i})\epsilon^{ljk}h^{j}\frac{\partial{T}}{\partial
x^{k}}\mp\tau_{R}^{i}\frac{\chi_{i}(2+\chi_{i})}{(1+\chi_{i}^{2})}\epsilon^{ljk}h^{j}\frac{\partial{\dot{T}}}{\partial
x^{k}}\Big{\\}}\Big{]}\frac{\partial f^{0}_{i}}{\partial\omega_{i}}.$ (41)

Hence, components of the electric current in the three spatial directions are
given as,

$\displaystyle j^{x}=$
$\displaystyle\sum_{i}\frac{q_{i}g_{i}}{3}\int\frac{d^{3}|\vec{k}_{i}|}{(2\pi)^{3}}v_{i}^{2}\frac{\tau_{R}^{i}q_{i}}{(1+\chi_{i}+\chi_{i}^{2})(1+\chi_{i})}\Big{\\{}\Big{(}\frac{(1+\chi_{i}^{2})+\chi_{i}(2+\chi_{i})}{(1+\chi_{i}^{2})}\Big{)}E^{x}\pm\chi_{i}$
$\displaystyle\Big{(}\frac{(1+\chi_{i}^{2})(1+\chi_{i})+\chi_{i}(2+\chi_{i})}{(1+\chi_{i}^{2})}\Big{)}(E^{z}-E^{y})\Big{\\}}\Big{(}-\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}\Big{)}+\sum_{i}\frac{q_{i}g_{i}}{3}\int\frac{d^{3}|\vec{k}_{i}|}{(2\pi)^{3}}v_{i}^{2}$
$\displaystyle\frac{\tau_{R}^{i}(\omega_{i}-{\rm
b}_{i}h)}{T(1+\chi_{i}+\chi_{i}^{2})(1+\chi_{i})}\Big{\\{}(1+\chi_{i})\frac{\partial
T}{\partial
x}-\tau_{R}^{i}\frac{(1+\chi_{i}-\chi_{i}^{2})}{(1+\chi_{i}^{2})}\frac{\partial\dot{T}}{\partial
x}\pm\chi_{i}(1+\chi_{i})\Big{(}\frac{\partial T}{\partial z}-\frac{\partial
T}{\partial y}\Big{)}$
$\displaystyle\mp\tau_{R}^{i}\frac{\chi_{i}(2+\chi_{i})}{(1+\chi_{i}^{2})}\Big{(}\frac{\partial\dot{T}}{\partial
z}-\frac{\partial\dot{T}}{\partial y}\Big{)}\Big{\\}}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}~{}.$ (42)

The components $j^{y}$ and $j^{z}$ can be obtained from the above equation by
changing x, y, and z in cyclic order. By defining the following integrals

$\displaystyle
L_{1i}=\frac{g_{i}}{3T}\int\frac{d^{3}|\vec{k}_{i}|}{(2\pi)^{3}}\frac{\vec{k}^{2}_{i}}{\omega_{i}^{2}}f^{0}_{i}(1-f^{0}_{i})\tau_{R}^{i}\frac{(1+\chi_{i}^{2})+\chi_{i}(2+\chi_{i})}{(1+\chi_{i})\left(1+\chi_{i}^{2}\right)\left(1+\chi_{i}+\chi_{i}^{2}\right)}~{},$
$\displaystyle
L_{2i}=\frac{g_{i}}{3T}\int\frac{d^{3}|\vec{k}_{i}|}{(2\pi)^{3}}\frac{\vec{k}^{2}_{i}}{\omega_{i}^{2}}f^{0}_{i}(1-f^{0}_{i})\tau_{R}^{i}\chi_{i}\frac{(1+\chi_{i})(1+\chi_{i}^{2})+\chi_{i}(2+\chi_{i})}{(1+\chi_{i})\left(1+\chi_{i}^{2}\right)\left(1+\chi_{i}+\chi_{i}^{2}\right)}~{},$
$\displaystyle
L_{3i}=\frac{g_{i}}{3T}\int\frac{d^{3}|\vec{k}_{i}|}{(2\pi)^{3}}\frac{\vec{k}^{2}_{i}}{\omega_{i}^{2}}(\omega_{i}-{\rm
b}_{i}h)f^{0}_{i}(1-f^{0}_{i})\tau_{R}^{i}\frac{1}{(1+\chi_{i}+\chi_{i}^{2})}~{},$
$\displaystyle
L_{4i}=\frac{g_{i}}{3T}\int\frac{d^{3}|\vec{k}_{i}|}{(2\pi)^{3}}\frac{\vec{k}^{2}_{i}}{\omega_{i}^{2}}(\omega_{i}-{\rm
b}_{i}h)f^{0}_{i}(1-f^{0}_{i})\tau_{R}^{i2}\frac{(1+\chi_{i}-\chi_{i}^{2})}{(1+\chi_{i})\left(1+\chi_{i}^{2}\right)\left(1+\chi_{i}+\chi_{i}^{2}\right)}~{},$
$\displaystyle
L_{5i}=\frac{g_{i}}{3T}\int\frac{d^{3}|\vec{k}_{i}|}{(2\pi)^{3}}\frac{\vec{k}^{2}_{i}}{\omega_{i}^{2}}(\omega_{i}-{\rm
b}_{i}h)f^{0}_{i}(1-f^{0}_{i})\tau_{R}^{i}\frac{\chi_{i}}{(1+\chi_{i}+\chi_{i}^{2})}~{},$
$\displaystyle
L_{6i}=\frac{g_{i}}{3T}\int\frac{d^{3}|\vec{k}_{i}|}{(2\pi)^{3}}\frac{\vec{k}^{2}_{i}}{\omega_{i}^{2}}(\omega_{i}-{\rm
b}_{i}h)f^{0}_{i}(1-f^{0}_{i})\tau_{R}^{i2}\frac{\chi_{i}(2+\chi_{i})}{(1+\chi_{i})\left(1+\chi_{i}^{2}\right)\left(1+\chi_{i}+\chi_{i}^{2}\right)}~{},$
(43)

we can write Eq. (4.1.1) as

$\displaystyle j^{x}$
$\displaystyle=\sum_{i}q_{i}^{2}L_{1i}E_{x}-\sum_{i}q_{i}^{2}L_{2i}E_{y}+\sum_{i}q_{i}^{2}L_{2i}E_{z}-\frac{1}{T}\sum_{i}q_{i}\Big{(}L_{3i}-L_{4i}\frac{d}{dt}\Big{)}\frac{dT}{dx}+\frac{1}{T}\sum_{i}q_{i}$
$\displaystyle\Big{(}L_{5i}-L_{6i}\frac{d}{dt}\Big{)}\frac{dT}{dy}-\frac{1}{T}\sum_{i}q_{i}\Big{(}L_{5i}-L_{6i}\frac{d}{dt}\Big{)}\frac{dT}{dz}~{}.$
(44)

Now, by setting $j_{x}$ = $j_{y}$ = $j_{z}$ = 0, i.e., when net current due to
the external field is zero, we get

$\displaystyle j^{x}=0$
$\displaystyle=\sigma_{e}E_{x}-\sigma_{H}E_{y}+\sigma_{H}E_{z}-\frac{L_{34}}{T}\frac{dT}{dx}+\frac{L_{56}}{T}\frac{dT}{dy}-\frac{L_{56}}{T}\frac{dT}{dz}~{},$
$\displaystyle j^{y}=0$
$\displaystyle=\sigma_{H}E_{x}+\sigma_{e}E_{y}-\sigma_{H}E_{z}-\frac{L_{56}}{T}\frac{dT}{dx}-\frac{L_{34}}{T}\frac{dT}{dy}+\frac{L_{56}}{T}\frac{dT}{dz}~{},$
$\displaystyle j^{z}=0$
$\displaystyle=-\sigma_{H}E_{x}+\sigma_{H}E_{y}+\sigma_{e}E_{z}+\frac{L_{56}}{T}\frac{dT}{dx}-\frac{L_{56}}{T}\frac{dT}{dy}-\frac{L_{34}}{T}\frac{dT}{dz}~{}.$
(45)

Where, $\sigma_{el}$ = $\sum_{i}q_{i}^{2}L_{1i}$ and $\sigma_{H}$ =
$\sum_{i}q_{i}^{2}L_{2i}$ are the Ohmic and the Hall components of electrical
conductivity, respectively. Further, we defined the intergrals $L_{34}$ =
$\sum_{i}q_{i}\Big{(}L_{3i}-L_{4i}\frac{d}{dt}\Big{)}$ and $L_{56}$ =
$\sum_{i}q_{i}\Big{(}L_{5i}-L_{6i}\frac{d}{dt}\Big{)}$. To solve the Eqs.
(45), we rearrange them in the matrix form of

$\displaystyle\mathbb{\bm{\sigma}~{}E=L~{}X},$ (46)

where,

$\displaystyle\bm{\sigma}=\begin{pmatrix}\sigma_{e}&-\sigma_{H}&\sigma_{H}\\\
\sigma_{H}&\sigma_{e}&-\sigma_{H}\\\ -\sigma_{H}&\sigma_{H}&\sigma_{e}\\\
\end{pmatrix},~{}\hskip 17.07164pt\mathbb{E}=\begin{pmatrix}E_{x}\\\ E_{y}\\\
E_{z}\end{pmatrix},~{}$
$\displaystyle\mathbb{L}=\frac{1}{T}\begin{pmatrix}-L_{34}&L_{56}&-L_{56}\\\
-L_{56}&-L_{34}&L_{56}\\\ L_{56}&-L_{56}&-L_{34}\\\
\end{pmatrix},~{}\mathbb{X}=\begin{pmatrix}\frac{dT}{dx}\vspace{2mm}\\\
\frac{dT}{dy}\vspace{2mm}\\\ \frac{dT}{dz}\\\ \end{pmatrix}.$

Hence, finding the inverse of $\bm{\sigma}$ matrix, we can obtain the
components of the electric field as,

$\displaystyle\mathbf{E}$ $\displaystyle=\mathbf{(\bm{\sigma}^{-1}L)X}$
$\displaystyle\Rightarrow\begin{pmatrix}E_{x}\\\ E_{y}\\\ E_{z}\end{pmatrix}$
$\displaystyle=\begin{pmatrix}S_{B}&\overline{NB}&NB\\\
NB&S_{B}&\overline{NB}\\\ \overline{NB}&NB&S_{B}\\\
\end{pmatrix}\begin{pmatrix}\frac{dT}{dx}\vspace{2mm}\\\
\frac{dT}{dy}\vspace{2mm}\\\ \frac{dT}{dz}\\\ \end{pmatrix}~{},$ (47)

where the transport coefficients are

$\displaystyle S_{B}$
$\displaystyle=\frac{(\sigma_{e}^{2}+\sigma_{H}^{2})L_{34}+(2\sigma_{e}\sigma_{H})L_{56}}{T(\sigma_{e}^{3}+3\sigma_{e}\sigma_{H}^{2})},$
$\displaystyle\overline{NB}$
$\displaystyle=\frac{\sigma_{H}(\sigma_{e}+\sigma_{H})L_{34}-\sigma_{e}(\sigma_{e}+\sigma_{H})L_{56}}{T(\sigma_{e}^{3}+3\sigma_{e}\sigma_{H}^{2})},$
$\displaystyle NB$
$\displaystyle=\frac{-\sigma_{H}(\sigma_{e}-\sigma_{H})L_{34}+\sigma_{e}(\sigma_{e}-\sigma_{H})L_{56}}{T(\sigma_{e}^{3}+3\sigma_{e}\sigma_{H}^{2})}.$
(48)

Note that this is the generalized case of the magneto-thermoelectric effect,
i.e., we did not consider any preferred direction of the magnetic field.
Therefore, we obtained three independent normalized transport coefficients -
one magneto-Seebeck coefficient $S_{B}$ and two Nernst coefficients $NB$ and
$\overline{NB}$.

Now, in heavy ion collisions, the net magnetic field could possibly have a
preferred direction Hongo:2013cqa . Therefore, for a special case when the
magnetic field is directed along the y-axis (say), i.e., the magnetic field
unit vector $h^{j}\equiv(0,1,0)$, Eq. (4.1.1) gives,

$\displaystyle j^{x}$
$\displaystyle=\sigma_{e}E_{x}+\sigma_{H}E_{z}-\frac{L_{34}}{T}\frac{dT}{dx}-\frac{L_{56}}{T}\frac{dT}{dz}~{},$
$\displaystyle j^{y}$
$\displaystyle=\sigma_{e}E_{y}-\frac{L_{34}}{T}\frac{dT}{dy}~{},$
$\displaystyle j^{z}$
$\displaystyle=-\sigma_{H}E_{x}+\sigma_{e}E_{z}+\frac{L_{56}}{T}\frac{dT}{dx}-\frac{L_{34}}{T}\frac{dT}{dz}~{}.$
(49)

Hence, by setting $j_{x}$ = $j_{y}$ = $j_{z}$ = 0, we can get the components
of the induced electric field as,

$\displaystyle E_{x}$
$\displaystyle=\frac{\Big{\\{}\sigma_{e}L_{34}+\sigma_{H}L_{56}\Big{\\}}\frac{dT}{dx}}{T\Big{\\{}(\sigma_{e})^{2}+(\sigma_{H})^{2}\Big{\\}}}+\frac{\Big{\\{}\sigma_{e}L_{56}-\sigma_{H}L_{34}\Big{\\}}\frac{dT}{dz}}{T\Big{\\{}(\sigma_{e})^{2}+(\sigma_{H})^{2}\Big{\\}}}~{},$
$\displaystyle E_{y}$
$\displaystyle=\frac{L_{34}}{T\sigma_{e}}\frac{dT}{dy}~{},$ $\displaystyle
E_{z}$
$\displaystyle=\frac{\Big{\\{}\sigma_{H}L_{34}-\sigma_{e}L_{56}\Big{\\}}\frac{dT}{dx}}{T\Big{\\{}(\sigma_{e})^{2}+(\sigma_{H})^{2}\Big{\\}}}+\frac{\Big{\\{}\sigma_{e}L_{34}+\sigma_{H}L_{56}\Big{\\}}\frac{dT}{dz}}{T\Big{\\{}(\sigma_{e})^{2}+(\sigma_{H})^{2}\Big{\\}}}~{}.$
(50)

The above equations can be written in matrix form as,

$\displaystyle\begin{pmatrix}E_{x}\\\ E_{y}\\\ E_{z}\end{pmatrix}$
$\displaystyle=\begin{pmatrix}S_{B}&0&NB\\\ 0&\overline{S}_{B}&0\\\
-NB&0&S_{B}\\\ \end{pmatrix}\begin{pmatrix}\frac{dT}{dx}\vspace{2mm}\\\
\frac{dT}{dy}\vspace{2mm}\\\ \frac{dT}{dz}\\\ \end{pmatrix}.$ (51)

Here, we identify the dimensionless transport coefficients - two magneto-
Seebeck coefficients ($S_{B}$, $\overline{S}_{B}$) and one Nernst coefficient
($NB$) as,

$\displaystyle S_{B}$
$\displaystyle=\frac{(\frac{\sigma_{el}}{T})(\frac{L_{34}}{T^{2}})+(\frac{\sigma_{H}}{T})(\frac{L_{56}}{T^{2}})}{(\frac{\sigma_{el}}{T})^{2}+(\frac{\sigma_{H}}{T})^{2}},~{}~{}~{}\overline{S}_{B}=\frac{(\frac{L_{34}}{T^{2}})}{(\frac{\sigma_{el}}{T})},~{}~{}~{}NB=\frac{(\frac{\sigma_{el}}{T})(\frac{L_{56}}{T^{2}})-(\frac{\sigma_{H}}{T})(\frac{L_{34}}{T^{2}})}{(\frac{\sigma_{el}}{T})^{2}+(\frac{\sigma_{H}}{T})^{2}}~{}.$
(52)

Note that, at vanishing Hall conductivity ($\sigma_{H}=0$) the magneto-Seebeck
coefficient ${S_{B}}$ becomes the same as $\overline{S}_{B}$. This can be
realised as follows. The coefficient $\overline{S}_{B}$ is responsible for the
electric field in the y-axis or along the magnetic field, where the Lorentz
force has null contribution. Therefore, the expression of $\overline{S}_{B}$
is similar to that of $eB_{0}=0$ case as in Eq. (27), except the fact that the
integral in $L_{34}$ depends on magnetic decay parameter $\tau_{B}$ at finite
magnetic field. Furthermore, for a static picture, where cooling effects are
absent $\left(\frac{dT}{dt}=0\right)$, the integrals $L_{34}$ and $L_{56}$
reduce to $L_{34}$ = $\sum_{i}q_{i}L_{3i}$ and $L_{56}$ =
$\sum_{i}q_{i}L_{5i}$, then the expression of Seebeck and Nernst coefficients
becomes similar to that obtained in earlier studies PhysRevD.102.096011 ;
Das2021 .

####  4.1.2 With Landau quantization

Now, we study the effect of Landau quantization on the thermoelectric
coefficient of the QGP medium and see how it deviates from classical results.
Basic modifications will occur in the dispersion relation for medium
constituents and the phase space integration. As we considered the magnetic
field along the y direction therefore, momentum quantization will occur in its
perpendicular plane, i.e., the x-z plane. Hence,

$\displaystyle\omega=(\vec{k}^{2}+m^{2})^{1/2}$
$\displaystyle~{}~{}~{}\rightarrow~{}~{}~{}$
$\displaystyle\omega_{l}=(k_{y}^{2}+m^{2}+2l|{q_{i}}|B)^{1/2},$ (53)
$\displaystyle 2\int\frac{d^{3}\vec{k}}{(2\pi)^{3}}$
$\displaystyle~{}~{}~{}\rightarrow~{}~{}~{}$
$\displaystyle\sum_{l=0}^{\infty}\alpha_{l}\frac{|{q_{i}}|B}{2\pi}\int\limits^{+\infty}_{-\infty}\frac{dk_{y}}{2\pi},$
(54)

The factor 2 of spin degeneracy in left hand side of last line will be
converted to $\alpha_{l}$, which will be 2 for all Landau levels $l$, except
lowest Landau level (LLL) $l=0$, where $\alpha_{l}=1$. In general, one can
write $\alpha_{l}=2-\delta_{l,0}$. Under the conditions of an extremely high
magnetic field for which all medium constituents occupy the lowest Landau
energy level $l=0$. It means that the perpendicular motion of medium
constituents completely vanishes as $k_{x}\approx k_{z}\approx 0$ at $l=0$.
However, below that strong magnetic field limt, $l>0$ energy levels might have
some non-negligible contributions, and the lowest Landau level (LLL)
approximation is not sufficient enough PhysRevC.106.044914 ; Dey:2021fbo .
This is the case in our study, where the magnetic field is not constant but
decays with time. Here, we also consider, $k_{x}^{2}\approx
k_{z}^{2}\approx(\frac{k_{x}^{2}+k_{z}^{2}}{2})=\frac{2l|{q_{i}}|B}{2}$, the
integrals in Eq. (4.1.1) can be expressed as,

$\displaystyle
L_{1i}=\frac{g_{i}}{T}\sum_{l=0}^{\infty}\alpha_{l}\frac{|{q_{i}}|B}{2\pi}\int\limits^{+\infty}_{-\infty}\frac{dk_{y}}{2\pi}\frac{{l|{q_{i}}|B}}{\omega^{2}_{l}}f^{0}_{i}(1-f^{0}_{i})\tau_{R}^{i}~{}\frac{(1+\chi_{i}^{2})+\chi_{i}(2+\chi_{i})}{(1+\chi_{i})\left(1+\chi_{i}^{2}\right)\left(1+\chi_{i}+\chi_{i}^{2}\right)}~{},$
$\displaystyle
L_{2i}=\frac{g_{i}}{T}\sum_{l=0}^{\infty}\alpha_{l}\frac{|{q_{i}}|B}{2\pi}\int\limits^{+\infty}_{-\infty}\frac{dk_{y}}{2\pi}\frac{{l|{q_{i}}|B}}{\omega^{2}_{l}}f^{0}_{i}(1-f^{0}_{i})\tau_{R}^{i}~{}\chi_{i}\frac{(1+\chi_{i})(1+\chi_{i}^{2})+\chi_{i}(2+\chi_{i})}{(1+\chi_{i})\left(1+\chi_{i}^{2}\right)\left(1+\chi_{i}+\chi_{i}^{2}\right)}~{},$
$\displaystyle
L_{3i}=\frac{g_{i}}{T}\sum_{l=0}^{\infty}\alpha_{l}\frac{|{q_{i}}|B}{2\pi}\int\limits^{+\infty}_{-\infty}\frac{dk_{y}}{2\pi}\frac{{l|{q_{i}}|B}}{\omega^{2}_{l}}(\omega_{i}-{\rm
b}_{i}h)\tau_{R}^{i}~{}f^{0}_{i}(1-f^{0}_{i})\frac{1}{(1+\chi_{i}+\chi_{i}^{2})}~{},~{}$
$\displaystyle
L_{4i}=\frac{g_{i}}{T}\sum_{l=0}^{\infty}\alpha_{l}\frac{|{q_{i}}|B}{2\pi}\int\limits^{+\infty}_{-\infty}\frac{dk_{y}}{2\pi}\frac{{l|{q_{i}}|B}}{\omega^{2}_{l}}(\omega_{i}-{\rm
b}_{i}h)\tau_{R}^{i2}~{}f^{0}_{i}(1-f^{0}_{i})\frac{(1+\chi_{i}-\chi_{i}^{2})}{(1+\chi_{i})\left(1+\chi_{i}^{2}\right)\left(1+\chi_{i}+\chi_{i}^{2}\right)}~{},~{}$
$\displaystyle
L_{5i}=\frac{g_{i}}{T}\sum_{l=0}^{\infty}\alpha_{l}\frac{|{q_{i}}|B}{2\pi}\int\limits^{+\infty}_{-\infty}\frac{dk_{y}}{2\pi}\frac{{l|{q_{i}}|B}}{\omega^{2}_{l}}(\omega_{i}-{\rm
b}_{i}h)\tau_{R}^{i}~{}f^{0}_{i}(1-f^{0}_{i})\frac{\chi_{i}}{(1+\chi_{i}+\chi_{i}^{2})}~{},~{}$
$\displaystyle
L_{6i}=\frac{g_{i}}{T}\sum_{l=0}^{\infty}\alpha_{l}\frac{|{q_{i}}|B}{2\pi}\int\limits^{+\infty}_{-\infty}\frac{dk_{y}}{2\pi}\frac{{l|{q_{i}}|B}}{\omega^{2}_{l}}(\omega_{i}-{\rm
b}_{i}h)\tau_{R}^{i2}~{}f^{0}_{i}(1-f^{0}_{i})\frac{\chi_{i}(2+\chi_{i})}{(1+\chi_{i})\left(1+\chi_{i}^{2}\right)\left(1+\chi_{i}+\chi_{i}^{2}\right)}.~{}$
(55)

Hence, the thermoelectric coefficients obtained in Eq. (52) will also be
modified under the conditions of Landau quantization.

####  4.1.3 Analytical solution for massless non-interacting system

To find the analytical expressions for the massless non-interacting system, we
follow the same procedure mentioned in Sec. (3.1.1). Hence, for the case of a
massless non-interacting system, the integral in Eq. (4.1.1) turns into the
following forms:

$\displaystyle L_{1i}$
$\displaystyle=\sum_{i}\frac{g_{i}\tau_{R}^{i}}{6\pi^{2}}\Big{(}2\zeta(2)T^{2}+b_{i}^{2}\mu_{B}^{2}\Big{)}\frac{(1+\chi_{i}^{2})+\chi_{i}(2+\chi_{i})}{(1+\chi_{i})\left(1+\chi_{i}^{2}\right)\left(1+\chi_{i}+\chi_{i}^{2}\right)}~{},$
$\displaystyle L_{2i}$
$\displaystyle=\sum_{i}\frac{2g_{i}\tau_{R}^{i}}{3\pi^{2}}~{}\eta(1)b_{i}\mu_{B}T\frac{(1+\chi_{i})(1+\chi_{i}^{2})+\chi_{i}(2+\chi_{i})}{(1+\chi_{i})\left(1+\chi_{i}^{2}\right)\left(1+\chi_{i}+\chi_{i}^{2}\right)}~{},$
$\displaystyle L_{3i}$
$\displaystyle=\sum_{i}\frac{g_{i}\tau_{R}^{i}}{6\pi^{2}}\Big{(}b_{i}^{3}\mu_{B}^{3}+6b_{i}\mu_{B}\zeta(2)T^{2}-2b_{i}h\zeta(2)T^{2}-b_{i}^{3}\mu_{B}^{2}h\Big{)}\frac{1}{(1+\chi_{i}+\chi_{i}^{2})}~{},$
$\displaystyle L_{4i}$
$\displaystyle=\sum_{i}\frac{g_{i}\tau_{R}^{i}}{6\pi^{2}}\Big{(}b_{i}^{3}\mu_{B}^{3}+6b_{i}\mu_{B}\zeta(2)T^{2}-2b_{i}h\zeta(2)T^{2}-b_{i}^{3}\mu_{B}^{2}h\Big{)}\frac{(1+\chi_{i}-\chi_{i}^{2})}{(1+\chi_{i})\left(1+\chi_{i}^{2}\right)\left(1+\chi_{i}+\chi_{i}^{2}\right)}~{},$
$\displaystyle L_{5i}$
$\displaystyle=\sum_{i}\frac{g_{i}\tau_{R}^{i}}{2\pi^{2}}\Big{(}3T^{3}\zeta(3)+2b_{i}^{2}\mu_{B}^{2}T\eta(1)-\frac{4}{3}hb_{i}^{2}\mu_{B}T\eta(1)\Big{)}\frac{\chi_{i}}{(1+\chi_{i}+\chi_{i}^{2})}~{},$
$\displaystyle L_{6i}$
$\displaystyle=\sum_{i}\frac{g_{i}\tau_{R}^{i2}}{2\pi^{2}}\Big{(}3T^{3}\zeta(3)+2b_{i}^{2}\mu_{B}^{2}T\eta(1)-\frac{4}{3}hb_{i}^{2}\mu_{B}T\eta(1)\Big{)}\frac{\chi_{i}(2+\chi_{i})}{(1+\chi_{i})\left(1+\chi_{i}^{2}\right)\left(1+\chi_{i}+\chi_{i}^{2}\right)}.$
(56)

Here, the integrals $L_{34}$ and $L_{56}$ in Eq. (52) are proportional to
$q_{i}$, hence for the massless non-interacting system $L_{34}$ and $L_{56}$
both individually summed upto zero. Therefore, dimensionless transport
coefficients - two magneto-Seebeck coefficients ($S_{B}$, $\overline{S}_{B}$)
and one Nernst coefficient ($NB$) in Eq. (52) vanishes for the massless non-
interacting system, hence

$\displaystyle S_{B}$ $\displaystyle=\overline{S}_{B}=NB=0.$ (57)

####  4.1.4 Numerical solution for QGP EoS

(a)

(b)

(c)

(d)

Figure 4: Magneto-Seebeck coefficient ($S_{B}$) as a function of proper time
($\tau$) at $\mu_{B}=0.3$ GeV. (a) $S_{B}$ for static, evolving, with and
without (w/o) Landau quantization (LQ), four cases are plotted here in the
upper panel for $\tau_{B}=3,~{}7$ fm for Bjorken flow at $1~{}m_{\pi}^{2}$.
The lower panel represents their percentage deviation at $\tau_{B}=7~{}$fm.
(b) $S_{B}$ in Gubser flow for static and evolving case at $r=0$ and $3$ fm
with LQ is plotted at $\tau_{B}=3~{}$fm (upper panel) and $\tau_{B}=7~{}$fm
(lower panel) at $eB_{0}=1~{}m_{\pi}^{2}$. (c) and (d) are the same as (a) and
(b), respectively, but at $eB_{0}=5~{}m_{\pi}^{2}$.

Fig. (4) represents the magneto-Seebeck coefficient ($S_{B}$) as a function of
proper time ($\tau$) at $\mu_{B}=0.3$ GeV in the presence of time-varying
magnetic field with initial value $eB_{0}=1~{}m_{\pi}^{2}$ (a) and (b) and
$eB_{0}$ = 5 $m_{\pi}^{2}$ (c) and (d) with its decay parameter $\tau_{B}$ = 3
fm, 7 fm. In all the plots, we observe similar behavior of $S_{B}$. For the
case of Bjorken flow, at the early times, the coefficient $S_{B}$ is negative,
which means that the direction of the produced electric field is opposite to
the temperature gradient of the medium. Later, for a particular value of
$\tau$, the $S_{B}$ becomes zero and then increases with positive value as
$\tau$ increases. Here, we take the static and evolving pictures for both the
cases along with Landau quantization (LQ) and without Landau quantization (w/o
LQ). For the Bjorken flow case, the ideal MHD cooling rate is used as
mentioned in Eq. (2.1.2). The black solid line represents the static picture
without Landau quantization, and the green dash-dot-dash line represents the
static picture with Landau quantization. On the other hand, the blue dash-dot
line represents the evolving picture w/o LQ, and the red dashed line
represents the evolving picture with LQ. $S_{B}$ is highly sensitive to
magnetic decay parameter $\tau_{B}$. In the lower panels of (a) and (c), the
percentage deviation without LQ from LQ is plotted with a magenta dash-dotted
line. As expected, the effect of LQ improves with the strength of the magnetic
field ($eB_{0}$). The cyan dashed line represents the percentage deviation
from static to evolving case. In both plots, this deviation is nearly 40%-45%
in the early evolution and decreases as the medium evolves. As the value of
$S_{B}$ becomes positive, this deviation shows a peak by shifting from a
negative value to a positive value. Once this percentage deviation achieves
its maximum value, later on, it starts to decrease with proper time. However,
the validity of the peak-like structure in the deviation plots is
questionable. Because as we mentioned before, the temperature falls beyond
$0.150$ GeV at the corresponding times, where the quasiparticle model is not
justified. In Fig. (4) (b) and (d), the upper and bottom panels correspond to
decay parameter $\tau_{B}$ = 3 fm and $\tau_{B}$ = 7 fm, respectively. Here,
we use the cooling rate obtained from Gubser flow as mentioned in Eq. (13).
Here, we also study the effect of radial expansion of medium on coefficient
$S_{B}$. The red solid line corresponds to $r=0$ fm, and the green dashed line
corresponds to $r=3$ fm in the static picture. On the other hand, The blue
dash-triple dot line corresponds to $r=0$ fm, and the black dashed line
corresponds to $r=3$ fm in the evolving picture. $S_{B}$ follows the similar
trend in both the plots. The numerical value of $S_{B}$ increases as a
function of $\tau$. The effect of $\tau_{B}$ is identical as observed in the
Bjorken flow case.

(a)

(b)

(c)

(d)

Figure 5: Magneto-Seebeck coefficient ($\overline{S}_{B}$) as a function of
proper time ($\tau$) at $\mu_{B}=0.3$ GeV. (a) $\overline{S}_{B}$ for static,
evolving, with and without (w/o) Landau quantization (LQ), four cases are
plotted here in the upper panel for $\tau_{B}=3,~{}7$ fm for Bjorken flow at
$1~{}m_{\pi}^{2}$. The lower panel represents their percentage deviation at
$\tau_{B}=7~{}$fm. (b) $\overline{S}_{B}$ in Gubser flow for static and
evolving case at $r=0$ and $3$ fm with LQ is plotted at $\tau_{B}=3~{}$fm
(upper panel) and $\tau_{B}=7~{}$fm (lower panel) at $eB_{0}=1~{}m_{\pi}^{2}$.
(c) and (d) are the same as (a) and (b), respectively, but at
$eB_{0}=5~{}m_{\pi}^{2}$.

Fig. (5) is the same as Fig. (4), but for $\overline{S}_{B}$. From Eq. (52),
at vanishing Hall like components $\sigma_{H}\rightarrow 0$ and
$I_{56}\rightarrow 0$, $\overline{S}_{B}=S_{B}$. Expression of
$\overline{S}_{B}$ is similar to that of Seebeck $S$ as in $eB=0$ case, except
the fact that in magnetic field effective relaxation time modifies which is a
function of $\tau_{R}$ and $\tau_{B}$. The Hall-like components
($\sigma_{H},I_{56}$) are approximately one order smaller than their leading
order counterparts ($\sigma_{e},I_{34}$). Therefore, at the early evolution,
$\overline{S}_{B}$ is close to $S_{B}$. In the later time, as the magnetic
field significantly reduces, $\overline{S}_{B}$ approaches the Seebeck
coefficient ($S$) as for the $eB=0$ case.

(a)

(b)

(c)

(d)

Figure 6: Dimensionless Nernst coefficient ($NB$) as a function of proper time
($\tau$) at $\mu_{B}=0.3$ GeV. (a) $NB$ for static, evolving, with and without
(w/o) Landau quantization (LQ), four cases are plotted here in the upper panel
for $\tau_{B}=3,~{}7$ fm for Bjorken flow at $1~{}m_{\pi}^{2}$. The lower
panel represents their percentage deviation at $\tau_{B}=7~{}$fm. (b) $NB$ in
Gubser flow for static and evolving case at $r=0$ and $3$ fm with LQ is
plotted at $\tau_{B}=3~{}$fm (upper panel) and $\tau_{B}=7~{}$fm (lower panel)
at $eB_{0}=1~{}m_{\pi}^{2}$. (c) and (d) are the same as (a) and (b),
respectively, but at $eB_{0}=5~{}m_{\pi}^{2}$.

In Fig. (6), dimensionless Nernst coefficient ($NB$) as a function of proper
time ($\tau$) is plotted. Representation is same as in Fig. (4) and (5).
Unlike magneto-Seebeck coefficients, $NB$ is positive and increases with
evolution. The effect of the magnetic field is similar as in $S_{B}$ and
$\overline{S}_{B}$. In the Gubser flow case, $NB$ increases as we go away from
the center, as we see for $r=0,~{}3$ fm. $NB$ enhances in the evolution
picture. Therefore, all thermoelectric coefficients improved or enhanced due
to the temperature evolution, and the effect is prominent in the early
evolution.

###  4.2 Induced electric field in anisotropic system

Figure 7: Electric field induced in the QGP in peripheral collisions with
$eB_{0}=5~{}m_{\pi}^{2}$, $\tau_{B}=3$ fm.

Figure 8: Magnitude of electric field in the transverse plane with
$eB_{0}=5~{}m_{\pi}^{2}$.

Figure 9: Induced electric field vector (left figure) and magnitude (right
figure) for $eB_{0}=1$ m${}_{\pi}^{2}$.

Here, we estimate the induced electric field in the peripheral collisions for
the different magnetic fields at $\mu_{B}=0.3$ GeV. To have a more realistic
picture of the induced field, we used evolving and Landau quantized values of
thermoelectric coefficients. In peripheral collisions, the initial magnetic
field created by moving spectators can break the momentum isotropy of the
medium. As a result, in Eq. (51), we obtained an anisotropic coefficient
matrix, which further produces an anisotropic electric field due to the
thermoelectric effect.

In Fig. (7), we plotted the space-time profile of electric field $\vec{E}$,
induced in QGP created in peripheral collision with the initial magnetic field
of $eB_{0}=5~{}m_{\pi}^{2}$. Upper panels are for magnetic field decay
parameter $\tau_{B}=3$ fm at two different times $\tau=0.5,~{}3$ fm, and lower
panels are for $\tau_{B}=7$ fm. Unlike the head-on collision case, $E_{z}\neq
0$ in the peripheral collision. Due to boost invariant symmetry, $\vec{E}$ is
symmetric along the $z$ axis. However, symmetry in the transverse plane is now
broken due to the magnetic field along the $y$ axis. This can be easily seen
from the magnitude plot in Fig. (8). In the early time during the QGP
evolution, the electric field could go as high as $eE\approx 1~{}m_{\pi}^{2}$,
close to the QCD scale. We see that increasing $\tau_{B}$ from $3$ to $7$ fm,
the strength of $\vec{E}$ reduces at an early time; however, at a later time,
it is quite the same for $\tau_{B}=3$ and $7$ fm. This signifies that in the
current formalism, the electric field induced by the thermoelectric effect
increases for faster decay of the initially produced magnetic field. This can
be understood as follows. For QGP with lower electrical conductivity, the
decay of the magnetic field is comparatively fast. Now, from Eq. (52)
thermoelectric coefficients varies inversely with electrical conductivity
($S_{B}/~{}\bar{S}_{B}/~{}NB\sim\frac{1}{\sigma_{e}}$). Therefore, smaller
$\sigma_{e}$ means fast decay of magnetic field or smaller value of
$\tau_{B}$, which enhances Seebeck and Nernet coefficients. As a result,
$\vec{E}$ increases. In Fig. (9), we plotted $\vec{E}$ for
$eB_{0}=1~{}m_{\pi}^{2}$ to see the effect of different magnetic field
strengths. There is a minor change in the induced field upon changing the
initial magnetic field strength from $eB=5$ to $1~{}m_{\pi}^{2}$. A similar
thing was already observed in the thermoelectric coefficients in Figs. (4, 5),
and 6.

It is important to note that we did not account for the effect of the magnetic
field on the Gubser cooling rate. A more realistic picture of the cooling rate
should be calculated with a finite and time-varying magnetic field. Moreover,
according to Maxwell’s electromagnetic relation, the induced space-time-
dependent electric field will give rise to a magnetic field. However, this
would require the calculation of permittivity and permeability of the evolving
QGP, which is subjected to further research and beyond the scope of the
current work. In the present study, for the first time, we have estimated the
induced electric field in QGP due to the thermoelectric effect. We used the
hydrodynamic theory, where analytical solutions are achievable. This study can
be quickly extended to incorporate in (3+1)D hydrodynamics or
magnetohydrodynamics such as ECHO QGP by dividing temperature into space-time
grids.

##  5 Summary and Conclusion

In summary, we have calculated thermoelectric coefficients for evolving QGP
medium created in RHICs. Thermoelectric transport phenomena can induce an
electric field in the medium. We have estimated the same for the first time
using Gubser hydrodynamics. We used the kinetic theory-based RTA approach to
calculate thermoelectric coefficients. Unlike the static picture of QGP
considered in all the earlier calculations Das2021 ; Abhishek:2020wjm ;
Kurian:2021zyb ; PhysRevD.102.096011 , we solved the Boltzmann equation for
evolving QGP by taking care of temperature evolution. Deriving the general
expressions of the transport coefficients, we used different cooling rates
(ideal hydrodynamic, ideal MHD, Bjorken flow) for their estimation. We used a
quasiparticle model to reproduce the lattice QCD EoS for QGP. Moreover,
calculations and results are obtained for head-on and peripheral collisions
where we expect a magnetic field due to spectator currents. All the results
are at finite baryon chemical potential $\mu_{B}=0.3$ GeV.

In a fluid of charged particles with finite thermal conductivity,
thermoelectric phenomena can create or induce an electromagnetic field. QGP
consists of quarks, which are electrically charged and have finite baryon
quantum numbers. Therefore, they are capable of electrical (charge) and
thermal (baryon) conduction. As a result, thermoelectric transport in QGP at
finite baryon chemical potential leads to the induction of an electric field.
QGP medium is governed by the QCD interaction. Studies show that the EM field
of the order of $\Lambda_{\rm QCD}$ scale can interfere with the interactions,
contributing to many interesting phenomena PhysRevD.78.074033 ;
PhysRevLett.73.3499 ; PhysRevD.56.5251 . It is widely accepted that such a
huge EM field was created in peripheral collisions Tuchin2013ie . This work
showed that an EM field can also be created during the QGP evolution due to
the thermoelectric effect, even in head-on collisions. The maximum electric
field produced in head-on and peripheral collisions is $\sim
0.6~{}m_{\pi}^{2}$ and $0.9~{}m_{\pi}^{2}$, respectively, which are two to one
times smaller than the QCD scale ($\sim 2~{}m_{\pi}^{2}$). However, this can
contribute to many interesting phenomena, as mentioned before. Nevertheless,
the impact might be small. It is essential to mention here that we have taken
care of temperature evolution in calculating thermoelectric coefficients.
i.e., the effect of cooling is considered, which was missing in all the
earlier studies. Results for static and evolving pictures are compared, and we
found that the impact of evolution is significant though non-monotonic. In
peripheral collisions, we considered an external time-varying magnetic field.
The effect of the initial magnetic field and its decay parameter are also
studied. The induced electric field falls off with time. At the early time of
evolution, the magnetic decay parameter has a significant impact; however, the
strength of the initial magnetic field has a minor effect. Considering RHIC
and LHC matter, we considered the external magnetic field along the $y$
direction for the numerical estimation. However, we have also generalized the
calculations, where the direction of the magnetic field was not preferred, and
we obtained the thermoelectric coefficient matrix for the same.

Despite the interesting findings, our approach has a few limitations. We need
to have an analytical solution for the cooling rate to calculate the
thermoelectric coefficients for an evolving system. However, the hydrodynamic
equations in ($3+1$) D result in coupled differential equations, which can
only be solved numerically. Therefore, we considered simplistic hydrodynamic
pictures (Bjorken and Gubser flow), where an analytical solution for the
cooling rate is achievable. Nevertheless, the results of thermoelectric
coefficients obtained here can be used in ($3+1$) D hydrodynamic simulation
for a better understanding of the space-time profile of the induced EM field.

##  Appendix A Appendix

###  A Thermoelectric coefficients

####  A.1 Case-I (B = 0)

To find $\delta f_{i}$ in Eq. (15), we solve the Boltzmann transport equation
(BTE) with the help of relaxation time approximation (RTA). In the presence of
temperature gradients, the BTE under RTA can be expressed as
PhysRevD.102.016016

$\displaystyle\frac{\partial f_{i}}{\partial
t}+\frac{\vec{k}_{i}}{\omega_{i}}\cdot\frac{\partial
f_{i}}{\partial\vec{x_{i}}}+q_{i}\vec{E}\cdot\frac{\partial
f_{i}}{\partial\vec{k_{i}}}=-\frac{\delta f_{i}}{\tau^{i}_{R}}~{},$ (A.1)

where $\tau^{i}_{R}$ is the relaxation time of the particle. $\vec{E}$ is an
external electric field. We can solve Eq. (A.1) using Eqs. (18) and (19), and
get the expression of $\alpha_{j}$’s. Here, the first term in left-hand side
(lhs) of the Eq. (A.1) becomes,

$\displaystyle-(\omega_{i}-{\rm b}_{i}\mu)\frac{\dot{T}}{T}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}+\omega_{i}{\vec{v}_{i}}.\Big{\\{}\dot{\alpha_{1}}{\vec{E}}+\alpha_{1}{\dot{\vec{E}}}+\dot{\alpha_{2}}{\vec{\nabla}T}+\alpha_{2}{\vec{\nabla}\dot{T}}\Big{\\}}.$
(A.2)

The second term in lhs of the Boltzmann equation leads to,

$\displaystyle\frac{\partial f_{i}}{\partial x_{i}}$
$\displaystyle=\frac{\partial}{\partial x_{i}}(f^{0}_{i}+\delta
f_{i})=\frac{\partial}{\partial x_{i}}(f^{0}_{i})+\frac{\partial
f_{i}}{\partial x_{i}}(\delta f_{i})=\frac{\partial}{\partial
x_{i}}\Big{(}\frac{1}{1\pm\exp{\big{(}\beta(\omega_{i}-{\rm
b}_{i}\mu_{i})\big{)}}}\Big{)}+0$ $\displaystyle=-(\omega_{i}-{\rm
b}_{i}h)\frac{\vec{\nabla}T}{T}\frac{\partial f^{0}_{i}}{\partial\omega_{i}}.$
(A.3)

The third term in lhs leads to,

$\displaystyle\frac{\partial f_{i}}{\partial k_{i}}$
$\displaystyle=\frac{\partial f^{0}_{i}}{\partial k_{i}}+\frac{\partial\delta
f_{i}}{\partial k_{i}}=\vec{v_{i}}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}+\vec{\Omega}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}.$ (A.4)

Hence,

$\displaystyle q_{i}\vec{E}\cdot\frac{\partial
f_{i}}{\partial\vec{k_{i}}}=q_{i}\vec{v_{i}}\cdot\vec{E}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}+q_{i}E\big{(}\alpha_{1}\vec{E}+\alpha_{2}\vec{\nabla}T+\alpha_{3}\vec{\nabla}\dot{T}\big{)}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}.$ (A.5)

Finally, after the substitution of the above results in both sides of Eq.
(A.1), we get

$\displaystyle-(\omega_{i}-{\rm b}_{i}\mu)\frac{\dot{T}}{T}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}+\omega_{i}{\vec{v}_{i}}.\Big{(}\dot{\alpha_{1}}{\vec{E}}+\alpha_{1}{\dot{\vec{E}}}+\dot{\alpha_{2}}{\vec{\nabla}T}+\alpha_{2}{\vec{\nabla}\dot{T}}\Big{)}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}-(\omega_{i}-{\rm
b}_{i}h)\frac{\vec{\nabla}T}{T}\frac{\partial f^{0}_{i}}{\partial\omega_{i}}$
$\displaystyle+q_{i}\vec{v_{i}}\cdot\vec{E}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}=-\frac{\omega_{i}}{\tau_{R}^{i}}\Big{(}\alpha_{1}\vec{E}+\alpha_{2}\vec{\nabla}T+\alpha_{3}\vec{\nabla}\dot{T}\Big{)}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}.$ (A.6)

In the current analysis, we consider only the terms with first-order
derivatives of the fields and neglect higher-order derivative terms. The
comparison of the coefficients of ${\vec{v}}\cdot{\vec{E}}$,
${\vec{v}}\cdot{\vec{\nabla}T}$ and ${\vec{v}}\cdot{\vec{\nabla}\dot{T}}$ on
both sides of the above equation leads to,

$\displaystyle\dot{\alpha_{1}}$
$\displaystyle=-\Bigg{(}\frac{1}{\tau_{R}^{i}}\alpha_{1}+\frac{q_{i}}{\omega_{i}}\Bigg{)},~{}~{}~{}\dot{\alpha_{2}}=-\bigg{(}\frac{1}{\tau_{R}^{i}}\alpha_{2}-\frac{\omega_{i}-{\rm
b}_{i}h}{\omega_{i}T}\bigg{)},~{}~{}~{}{\alpha_{3}}=-{\tau_{R}^{i}}\alpha_{2},$
(A.7)

Here, $\dot{\alpha_{1}}$ and $\dot{\alpha_{3}}$ from Eq. (A.7) can be
expressed in terms of matrix equation as

$\frac{dX}{dt}=AX+G,$ (A.8)

where the matrices take the following forms

$\displaystyle X=\begin{pmatrix}\alpha_{1}\\\ \alpha_{2}\\\
\end{pmatrix},~{}~{}~{}A=\begin{pmatrix}-\frac{1}{\tau_{R}^{i}}&0\\\
0&-\frac{1}{\tau_{R}^{i}}\\\
\end{pmatrix},~{}~{}~{}G=\begin{pmatrix}-\frac{q_{i}}{\omega_{i}}\\\
\frac{\omega_{i}-{\rm b}_{i}h}{\omega_{i}T}\end{pmatrix}.$ (A.9)

Eq. (A.8) can be solved by using the method of the variation of constants. The
eigenvalues corresponding to matrix $A$ are,
$\lambda=\lambda_{1}=\lambda_{2}$ = $-\frac{1}{\tau_{R}^{i}}.$ Hence, one can
write the linearly independent solutions corresponding to the homogeneous part
of differential Eq. (A.8) in terms of its eigenvectors as

$\displaystyle y_{1}=\begin{pmatrix}e^{\eta}\\\
0\end{pmatrix},~{}~{}~{}y_{2}=\begin{pmatrix}0\\\
e^{\eta}\end{pmatrix};~{}~{}{\rm where}~{}~{}\eta$
$\displaystyle=-\frac{\tau}{\tau_{R}^{i}}.$ (A.10)

Therefore, the fundamental matrix for Eq. (A.8) is

$\displaystyle Y=\begin{pmatrix}e^{\eta}&0\\\ 0&e^{\eta}\\\ \end{pmatrix},$
(A.11)

We seek a particular solution of the equation with a given form of Eq. (A.8)
is

$\displaystyle Y_{p}$ $\displaystyle=YU,$ (A.12)

where $U$ is a column matrix of unknowns.

$\displaystyle U=\begin{pmatrix}u1\\\ u2\\\ .\\\ .\\\ .\\\ u_{n}\\\
\end{pmatrix}.$ (A.13)

From Eq. (A.12) and Eq. (A.8) we can see that $Y_{p}$ is a column matrix with
the same coefficients as that of matrix X. Further, the differentiation of Eq.
(A.12) with respect to time gives us $Y_{p}^{{}^{\prime}}$ =
$Y^{{}^{\prime}}U+YU^{{}^{\prime}}$, where $Y^{{}^{\prime}}=AY$. Hence,
$Y_{p}^{{}^{\prime}}$ = $AY_{p}+YU^{{}^{\prime}}$. Comparison of above
equation with Eq. (A.8) tells us

$\displaystyle G=YU^{{}^{\prime}}.$ (A.14)

The determinant of matrix $Y$ is $e^{2\eta}$. Then,

$\displaystyle u^{{}^{\prime}}_{1}$
$\displaystyle=\frac{1}{e^{2\eta}}det\begin{pmatrix}-\frac{q_{i}}{\omega_{i}}&0\\\
\frac{\omega_{i}-{\rm b}_{i}h}{\omega_{i}T}&e^{\eta}\\\
\end{pmatrix},~{}~{}~{}u^{{}^{\prime}}_{2}=\frac{1}{e^{2\eta}}det\begin{pmatrix}{e^{\eta}}&-\frac{q_{i}}{\omega_{i}}\\\
0&\frac{\omega_{i}-{\rm b}_{i}h}{\omega_{i}T}\\\ \end{pmatrix}.$

After integrating $u^{{}^{\prime}}_{1}$ and $u^{{}^{\prime}}_{2}$ with respect
to time , we get the matrix $U$ as

$\displaystyle
U=\begin{pmatrix}\frac{-q_{i}}{\omega_{i}}\xi_{1}\vspace{2mm}\\\
\frac{\omega_{i}-{\rm b}_{i}h}{\omega_{i}T}\xi_{2}\\\ \end{pmatrix}.$ (A.15)

Where $\xi=\int e^{-\eta d\tau}$.

After substituting the above value of $U$ in Eq. (A.30) one obtain,

$\displaystyle Y_{p}$ $\displaystyle=\begin{pmatrix}{e^{\eta}}&0&\\\
0&e^{\eta}\\\ \end{pmatrix}\begin{pmatrix}\frac{-q_{i}}{\omega_{i}}\xi_{1}\\\
\frac{\omega_{i}-{\rm b}_{i}h}{\omega_{i}T}\xi_{2}\\\ \end{pmatrix},$
$\displaystyle\begin{pmatrix}\alpha_{1}\\\ \alpha_{2}\\\ \end{pmatrix}$
$\displaystyle=\begin{pmatrix}e^{\eta}&0&\\\ 0&e^{\eta}\\\
\end{pmatrix}\begin{pmatrix}u_{1}\\\ u_{2}\\\ \end{pmatrix}.$

Hence,

$\displaystyle\alpha_{1}=u_{1}e^{\eta},~{}~{}~{}\alpha_{2}=u_{2}e^{\eta}.$
(A.16)

The functions $u_{1}(\tau)$ and $u_{2}(\tau)$ can be defined as
$u_{1}=-\frac{q_{i}}{\omega_{i}}\xi_{1}$, $u_{2}=i\frac{\omega_{i}-{\rm
b}_{i}h}{\omega_{i}T}\xi_{2}$. Hence, the alphas in Eq. (A.7) get the forms,

$\displaystyle\alpha_{1}=-\frac{q_{i}\tau_{R}^{i}}{\omega_{i}},~{}~{}~{}\alpha_{2}=\frac{\tau_{R}^{i}(\omega_{i}-{\rm
b}_{i}h)}{T\omega_{i}},~{}~{}~{}\alpha_{3}=-\tau_{R}^{i}\alpha_{3}.$

After getting the expressions for $\alpha_{j}$’s, the simplified form of
$\delta f_{i}$ is,

$\displaystyle\delta f_{i}$
$\displaystyle=-q_{i}\tau_{R}^{i}(\vec{v_{i}}.\vec{E})+\frac{(\omega_{i}-{\rm
b}_{i}h)\tau_{R}^{i}}{T}\Big{[}\Big{(}\vec{v_{i}}\cdot\vec{\nabla}{T}\Big{)}-\tau_{R}^{i}\Big{(}\vec{v_{i}}\cdot\vec{\nabla}{\dot{T}}\Big{)}\Big{]}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}.$ (A.18)

Here $\pm$ indicates a positive sign for positively charged
particles(antiparticles) and a negative sign for negatively charged particles
(antiparticles).

####  A.2 Case-II (B $\not=$ 0)

To find the expression of $\delta f_{i}$ in Eq. (15), we solve the Boltzmann
transport equation (BTE) with the help of relaxation time approximation (RTA).
In the presence of an external electromagnetic field, BTE under RTA can be
expressed as PhysRevD.102.016016

$\displaystyle\frac{\partial f_{i}}{\partial
t}+\frac{\vec{k}_{i}}{\omega_{i}}\cdot\frac{\partial
f_{i}}{\partial\vec{x_{i}}}+q_{i}\left(\vec{E}+\frac{\vec{k}_{i}}{\omega_{i}}\times\vec{B}\right)\cdot\frac{\partial
f_{i}}{\partial\vec{k_{i}}}=-\frac{\delta f_{i}}{\tau^{i}_{R}}~{},$ (A.19)

where $\tau^{i}_{R}$ is the relaxation time of the particle. We can solve Eq.
(A.19) using Eqs. (38) and (4.1.1), and get the expression of $\alpha_{j}$’s.
Here, the first term in left-hand side (lhs) of the Eq. (A.19) becomes,

$\displaystyle\omega_{i}{\vec{v}_{i}}.\Big{\\{}\dot{\alpha_{1}}{\vec{E}}+\alpha_{1}\dot{{\vec{E}}}+\dot{\alpha_{3}}{\vec{\nabla}T}+{\alpha_{3}}{\vec{\nabla}\dot{T}}+\alpha_{5}({\vec{\nabla}T}\times\dot{{\vec{B}}})+\dot{\alpha_{5}}({\vec{\nabla}T}\times{\vec{B}})+{\alpha_{5}}({\vec{\nabla}\dot{T}}\times{\vec{B}})$
$\displaystyle+\alpha_{5}(\vec{\nabla}\dot{T}\times{\vec{B}})+\alpha_{8}({\vec{E}}\times\dot{{\vec{B}}})+\dot{\alpha_{8}}({\vec{E}}\times{\vec{B}})+\alpha_{8}(\dot{\vec{E}}\times{\vec{B}})-(\omega_{i}-{\rm
b}_{i}\mu)\frac{\dot{T}}{T}\Big{\\}}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}.$ (A.20)

The second term in lhs of the Boltzmann equation leads to,

$\displaystyle\frac{\partial f_{i}}{\partial x_{i}}$
$\displaystyle=\frac{\partial}{\partial x_{i}}(f^{0}_{i}+\delta
f_{i})=\frac{\partial}{\partial x_{i}}(f^{0}_{i})+\frac{\partial
f_{i}}{\partial x_{i}}(\delta f_{i})=\frac{\partial}{\partial
x_{i}}\Big{\\{}\frac{1}{1\pm\exp{\big{(}\beta(\omega_{i}-{\rm
b}_{i}\mu_{i})\big{)}}}\Big{\\}}+0$ $\displaystyle=-(\omega_{i}-{\rm
b}_{i}h)\frac{\vec{\nabla}T}{T}\frac{\partial f^{0}_{i}}{\partial\omega_{i}}.$
(A.21)

The third term in lhs leads to,

$\displaystyle\frac{\partial f_{i}}{\partial k_{i}}$
$\displaystyle=\frac{\partial f^{0}_{i}}{\partial k_{i}}+\frac{\partial\delta
f_{i}}{\partial k_{i}}=\vec{v_{i}}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}+\vec{\Omega}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}.$ (A.22)

Hence,

$\displaystyle
q_{i}\left(\vec{E}+\frac{\vec{k}_{i}}{\omega_{i}}\times\vec{B}\right)\cdot\frac{\partial
f_{i}}{\partial\vec{k_{i}}}=q_{i}\vec{v_{i}}\cdot\vec{E}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}+\big{(}\vec{v}\times\vec{B}\big{)}\big{\\{}\alpha_{1}\vec{E}+\alpha_{2}\dot{\vec{E}}+\alpha_{3}\vec{\nabla}T+\alpha_{4}\vec{\nabla}\dot{T}+\alpha_{5}(\vec{\nabla}T\times\vec{B})$
$\displaystyle+\alpha_{6}(\vec{\nabla}T\times\dot{\vec{B}})+\alpha_{7}(\vec{\nabla}\dot{T}\times\vec{B})+\alpha_{8}(\vec{E}\times{\vec{B}})+\alpha_{9}(\vec{E}\times{\dot{\vec{B}}})+\alpha_{10}(\dot{\vec{E}}\times{\vec{B}})\big{\\}}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}.$ (A.23)

Finally, after the substitution of the above results on both the sides of Eq.
(A.19), the left side becomes

$\displaystyle\omega_{i}{\vec{v}_{i}}.\Big{\\{}\dot{\alpha_{1}}{\vec{E}}+\alpha_{1}\dot{{\vec{E}}}+\dot{\alpha_{3}}{\vec{\nabla}T}+{\alpha_{3}}{\vec{\nabla}\dot{T}}+\alpha_{5}({\vec{\nabla}T}\times\dot{{\vec{B}}})+\dot{\alpha_{5}}({\vec{\nabla}T}\times{\vec{B}})+{\alpha_{5}}({\vec{\nabla}\dot{T}}\times{\vec{B}})$
$\displaystyle+\alpha_{5}(\vec{\nabla}\dot{T}\times{\vec{B}})+\alpha_{8}({\vec{E}}\times\dot{{\vec{B}}})+\dot{\alpha_{8}}({\vec{E}}\times{\vec{B}})+\alpha_{8}(\dot{\vec{E}}\times{\vec{B}})-(\omega_{i}-{\rm
b}_{i}\mu)\frac{\dot{T}}{T}\Big{\\}}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}.$ $\displaystyle-(\omega_{i}-{\rm
b}_{i}h)\vec{v_{i}}\cdot\frac{\vec{\nabla}T}{T}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}+q_{i}\vec{v_{i}}\cdot\vec{E}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}+\big{\\{}-\alpha_{1}q_{i}\vec{v_{i}}\cdot(\vec{E}\times\vec{B})-\alpha_{2}q_{i}\vec{v_{i}}\cdot(\dot{\vec{E}}\times\vec{B})$
$\displaystyle-\alpha_{3}q_{i}\vec{v_{i}}\cdot(\vec{\nabla}T\times\vec{B})-\alpha_{4}q_{i}\vec{v_{i}}\cdot(\vec{\nabla}\dot{T}\times\vec{B})+\alpha_{5}q_{i}\\{-(\vec{v_{i}}\cdot\vec{B})(\vec{\nabla
T}\cdot\vec{B})+(\vec{v_{i}}\cdot\vec{\nabla}T)({\vec{B}}\cdot\vec{B})\\}$
$\displaystyle+\alpha_{6}q_{i}\\{-(\vec{v_{i}}\cdot\dot{\vec{B}})(\vec{\nabla
T}\cdot\vec{B})+(\vec{v_{i}}\cdot\vec{\nabla}T)(\dot{\vec{B}}\cdot\vec{B})\\}+\alpha_{7}q_{i}\\{-(\vec{v_{i}}\cdot\vec{B})(\vec{\nabla}\dot{T}\cdot\vec{B})+(\vec{v_{i}}\cdot\vec{\nabla}\dot{T})({\vec{B}}\cdot\vec{B})\\}$
$\displaystyle+\alpha_{8}q_{i}\\{-(\vec{v_{i}}\cdot\vec{B})(\vec{E}\cdot\vec{B})+(\vec{v_{i}}\cdot\vec{E})({\vec{B}}\cdot\vec{B})\\}+\alpha_{9}q_{i}\\{-(\vec{v_{i}}\cdot\dot{\vec{B}})(\vec{E}\cdot\vec{B})+(\vec{v_{i}}\cdot\vec{\dot{E}})({\vec{B}}\cdot\dot{\vec{B}})\\}$
$\displaystyle+\alpha_{9}q_{i}\\{-(\vec{v_{i}}\cdot\dot{\vec{B}})(\vec{E}\cdot\vec{B})+(\vec{v_{i}}\cdot\dot{\vec{E}})({\vec{B}}\cdot\dot{\vec{B}})\\}+\alpha_{10}q_{i}\\{-(\vec{v_{i}}\cdot\vec{B})(\dot{\vec{E}}\cdot\vec{B})+(\vec{v_{i}}\cdot\dot{\vec{E}})({\vec{B}}\cdot\vec{B})\\}\Big{\\}}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}},$ (A.24)

and the right side becomes

$\displaystyle=-\frac{\omega_{i}}{\tau_{R}^{i}}\alpha_{1}\vec{E}+\alpha_{2}\dot{\vec{E}}+\alpha_{3}\vec{\nabla}T+\alpha_{4}\vec{\nabla}\dot{T}+\alpha_{5}(\vec{\nabla}T\times\vec{B})+\alpha_{6}(\vec{\nabla}T\times\dot{\vec{B}})$
$\displaystyle+\alpha_{7}(\vec{\nabla}\dot{T}\times\vec{B})+\alpha_{8}(\vec{E}\times{\vec{B}})+\alpha_{9}(\vec{E}\times{\dot{\vec{B}}})+\alpha_{10}(\dot{\vec{E}}\times{\vec{B}})\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}.$ (A.25)

In the current analysis, we consider only the terms with first-order
derivatives of the fields and neglect higher-order derivative terms. The
comparison of the coefficients of ${\vec{v}}\cdot{\vec{E}}$,
${\vec{v}}\cdot{\dot{\vec{E}}}$, ${\vec{v_{i}}}\cdot{\vec{\nabla}T}$,
${\vec{v_{i}}}\cdot{\vec{\nabla}\dot{T}}$,
${\vec{v_{i}}}\cdot({\vec{\nabla}T\times{\vec{B}}})$,
${\vec{v_{i}}}\cdot({\vec{\nabla}T\times\dot{\vec{B}}})$,
${\vec{v_{i}}}\cdot({\vec{\nabla}\dot{T}\times\vec{B}})$,
${\vec{v_{i}}}\cdot({\vec{E}\times{\vec{B}}})$,
${\vec{v_{i}}}\cdot({\vec{E}}\times\dot{\vec{B}})$ and
${\vec{v_{i}}}\cdot(\dot{\vec{E}}\times{\vec{B}})$ on both sides, gives us
$\dot{\alpha_{1}}$, ${\alpha_{2}}$, $\dot{\alpha_{3}}$, $\alpha_{4}$,
$\dot{\alpha_{5}}$, $\dot{\alpha_{3}}$, $\alpha_{6}$, $\alpha_{7}$,
$\dot{\alpha}_{8}$, $\alpha_{9}$ and $\alpha_{10}$ respectively as

$\displaystyle\dot{\alpha_{1}}$
$\displaystyle=-\Bigg{[}\frac{1}{\tau_{R}^{i}}\alpha_{1}+\frac{q_{i}\big{\\{}(\vec{B}\cdot\vec{B}-\tau_{R}^{i}\vec{B}\cdot\dot{\vec{B}})\alpha_{7}\big{\\}}}{\omega_{i}}+\frac{q_{i}}{\omega_{i}}\Bigg{]},$
$\displaystyle\alpha_{2}$
$\displaystyle=-\tau_{R}^{i}\alpha_{1}-\frac{\tau_{R}^{i}q_{i}\vec{B}^{2}}{\omega_{i}}\alpha_{10},$
$\displaystyle\dot{\alpha_{3}}$
$\displaystyle=-\bigg{\\{}\frac{1}{\tau_{R}^{i}}\alpha_{3}+(\frac{q_{i}(\vec{B}\cdot\vec{B}-\tau_{R}^{i}\vec{B}\cdot\dot{\vec{B}})}{\omega_{i}})\alpha_{5}-\frac{\omega_{i}-{\rm
b}_{i}h}{T\omega_{i}}\bigg{\\}},$ $\displaystyle{\alpha_{4}}$
$\displaystyle=-{\tau_{R}^{i}}\Big{\\{}\alpha_{3}+\frac{q_{i}(\vec{B}\cdot\vec{B})}{\omega_{i}}\alpha_{7}\Big{\\}},$
$\displaystyle\dot{\alpha_{5}}$
$\displaystyle=-\frac{1}{\tau_{R}^{i}}\alpha_{5}+\frac{q_{i}}{\omega_{i}}\alpha_{3},$
$\displaystyle\alpha_{6}$ $\displaystyle=-\tau_{R}^{i}\alpha_{5},$
$\displaystyle\alpha_{7}$
$\displaystyle=-\tau_{R}^{i}\alpha_{5}+\frac{\tau_{R}^{i}q_{i}}{\omega_{i}}\alpha_{4},$
$\displaystyle\dot{\alpha_{8}}$
$\displaystyle=-\frac{1}{\tau_{R}^{i}}\alpha_{8}+\frac{q_{i}}{\omega_{i}}\alpha_{1},$
$\displaystyle\alpha_{9}$ $\displaystyle=-\tau_{R}^{i}\alpha_{8},$
$\displaystyle\alpha_{10}$
$\displaystyle=-\tau_{R}^{i}\alpha_{8}+\frac{\tau_{R}^{i}q_{i}}{\omega_{i}}\alpha_{2}.$
(A.26)

Here, $\dot{\alpha_{1}}$, $\dot{\alpha_{3}}$, $\dot{\alpha_{5}}$ and
$\dot{\alpha_{5}}$ from Eq. (A.2) can be expressed in terms of matrix equation
as

$\frac{dX}{dt}=AX+G,$ (A.27)

where the matrices take the following forms

$\displaystyle X=\begin{pmatrix}\alpha_{1}\\\ \alpha_{3}\\\ \alpha_{5}\\\
\alpha_{8}\end{pmatrix},~{}~{}~{}A=\begin{pmatrix}-\frac{1}{\tau_{R}^{i}}&0&0&-\frac{q_{i}F^{2}}{\omega_{i}}\\\
0&-\frac{1}{\tau_{R}^{i}}&-\frac{q_{i}F^{2}}{\omega_{i}}&0\\\
0&\frac{q_{i}}{\omega_{i}}&-\frac{1}{\tau_{R}^{i}}&0\\\
\frac{q_{i}}{\omega_{i}}&0&0&-\frac{1}{\tau_{R}^{i}},\end{pmatrix},~{}~{}~{}G=\begin{pmatrix}-\frac{q_{i}}{\omega_{i}}\\\
\frac{\omega_{i}-{\rm b}_{i}h}{T\omega_{i}}\\\ 0\\\ 0\end{pmatrix},$

with $F=\sqrt{B(B-\tau_{R}^{i}\dot{B})}$. Eq. (A.27) can be solved by
diagonalizing the matrix $A$ and using the method of the variation of
constants. The eigenvalues corresponding to matrix $A$ are,
$\lambda_{j}$ = $-\frac{1}{\tau_{R}^{i}}+a_{j}i\frac{q_{i}F}{\omega_{i}}$,
with $a_{1}=$ $a_{2}=-1$, $a_{3}=$ $a_{4}=1$.
Hence, one can write the linear independent solutions corresponding to the
homogeneous part of differential Eq. (A.27) in terms of its eigenvectors as

$\displaystyle y_{1}=\begin{pmatrix}-iFe^{\eta_{1}}\\\ 0\\\ 0\\\
e^{\eta_{1}}\end{pmatrix},~{}~{}y_{2}=\begin{pmatrix}0\\\ -iFe^{\eta_{2}}\\\
e^{\eta_{2}}\\\ 0\end{pmatrix},~{}~{}y_{3}=\begin{pmatrix}iFe^{\eta_{3}}\\\
0\\\ 0\\\ e^{\eta_{3}}\end{pmatrix},~{}~{}y_{4}=\begin{pmatrix}0\\\
iFe^{\eta_{4}}\\\ e^{\eta_{4}}\\\ 0\end{pmatrix}.$

Where,

$\displaystyle\eta_{j}$
$\displaystyle=-\frac{\tau}{\tau_{R}^{i}}+a_{j}\frac{q_{i}i}{\omega}\int
Fd\tau.$ (A.28)

Therefore, the fundamental matrix for Eq. (A.27) is

$\displaystyle Y=\begin{pmatrix}-iFe^{\eta_{1}}&0&iFe^{\eta_{3}}&0\\\
0&-iFe^{\eta_{2}}&0&iFe^{\eta_{4}}\\\ 0&e^{\eta_{2}}&0&e^{\eta_{4}}\\\
e^{\eta_{1}}&0&e^{\eta_{3}}&0\\\ \end{pmatrix},$ (A.29)

We seek a particular solution of the equation with a given form of Eq. (A.27)
is

$\displaystyle Y_{p}$ $\displaystyle=YU,$ (A.30)

where $U$ is a column matrix of unknowns.

$\displaystyle U=\begin{pmatrix}u1\\\ u2\\\ .\\\ .\\\ .\\\ u_{n}\\\
\end{pmatrix}.$ (A.31)

From Eq. (A.30) and Eq. (A.27) we can see that $Y_{p}$ is a column matrix with
the same coefficients as that of matrix X. Further, the differentiation of Eq.
(A.30) with respect to time gives us $Y_{p}^{{}^{\prime}}$ =
$Y^{{}^{\prime}}U+YU^{{}^{\prime}}$, where $Y^{{}^{\prime}}=AY$. Hence,
$Y_{p}^{{}^{\prime}}$ = $AY_{p}+YU^{{}^{\prime}}$. Comparison of above
equation with Eq. (A.27) tells us

$\displaystyle G=YU^{{}^{\prime}}.$ (A.32)

The determinant of matrix $Y$ is $4F^{2}e^{\eta}$. Then,

$\displaystyle u^{{}^{\prime}}_{1}$
$\displaystyle=\frac{1}{4F^{2}e^{\eta}}det\begin{pmatrix}-\frac{q_{i}}{\omega_{i}}&0&iFe^{\eta_{3}}&0\\\
\frac{\omega_{i}-{\rm
b}_{i}h}{T\omega_{i}}&-iFe^{\eta_{2}}&0&iFe^{\eta_{4}}\\\
0&e^{\eta_{2}}&0&e^{\eta_{4}}\\\ 0&0&e^{\eta_{3}}&0\\\ \end{pmatrix},$
$\displaystyle u^{{}^{\prime}}_{2}$
$\displaystyle=\frac{1}{4F^{2}e^{\eta}}det\begin{pmatrix}-iFe^{\eta_{1}}&-\frac{q_{i}}{\omega_{i}}&iFe^{\eta_{3}}&0\\\
0&\frac{\omega_{i}-{\rm b}_{i}h}{T\omega_{i}}&0&iFe^{\eta_{4}}\\\
0&0&0&e^{\eta_{4}}\\\ e^{\eta_{1}}&0&e^{\eta_{3}}&0\\\ \end{pmatrix},$
$\displaystyle u^{{}^{\prime}}_{3}$
$\displaystyle=\frac{1}{4F^{2}e^{\eta}}det\begin{pmatrix}-iFe^{\eta_{1}}&0&-\frac{q_{i}}{\omega_{i}}&0\\\
0&-iFe^{\eta_{2}}&\frac{\omega_{i}-{\rm
b}_{i}h}{T\omega_{i}}&iFe^{\eta_{4}}\\\ 0&e^{\eta_{2}}&0&e^{\eta_{4}}\\\
e^{\eta_{1}}&0&0&0\\\ \end{pmatrix},$ $\displaystyle u^{{}^{\prime}}_{4}$
$\displaystyle=\frac{1}{4F^{2}e^{\eta}}det\begin{pmatrix}-iFe^{\eta_{1}}&0&iFe^{\eta_{3}}&-\frac{q_{i}}{\omega_{i}}\\\
0&-iFe^{\eta_{2}}&0&\frac{\omega_{i}-{\rm b}_{i}h}{T\omega_{i}}\\\
0&e^{\eta_{2}}&0&0\\\ e^{\eta_{1}}&0&e^{\eta_{3}}&0\\\ \end{pmatrix}.$ (A.33)

After integrating $u^{{}^{\prime}}_{1}$, $u^{{}^{\prime}}_{2}$,
$u^{{}^{\prime}}_{3}$ and $u^{{}^{\prime}}_{4}$ with respect to time , we get
the matrix $U$ as

$\displaystyle U$
$\displaystyle=\begin{pmatrix}\frac{q_{i}}{2i\omega_{i}}\xi_{1}\\\
i\frac{\omega_{i}-{\rm b}_{i}h}{2\omega_{i}T}\xi_{2}\\\
-\frac{q_{i}}{2i\omega_{i}}\xi_{3}\\\ -i\frac{\omega_{i}-{\rm
b}_{i}h}{2\omega_{i}}\xi_{4}\\\
\end{pmatrix},~{}~{}~{}\xi_{j}=\int\frac{e^{-\eta_{j}}}{F}d\tau.$ (A.34)

After substituting the above value of $U$ in Eq. (A.30) one obtain,

$\displaystyle Y_{p}$
$\displaystyle=\begin{pmatrix}-iFe^{\eta_{1}}&0&iFe^{\eta_{3}}&0\\\
0&-iFe^{\eta_{2}}&0&iFe^{\eta_{4}}\\\ 0&e^{\eta_{2}}&0&e^{\eta_{4}}\\\
e^{\eta_{1}}&0&e^{\eta_{3}}&0\\\
\end{pmatrix}\begin{pmatrix}\frac{q_{i}}{2i\omega_{i}}\xi_{1}\\\
i\frac{\omega_{i}-{\rm b}_{i}h}{2\omega_{i}T}\xi_{2}\\\
-\frac{q_{i}}{2i\omega_{i}}\xi_{3}\\\ -i\frac{\omega_{i}-{\rm
b}_{i}h}{2\omega_{i}}\xi_{4}\\\ \end{pmatrix},$
$\displaystyle\begin{pmatrix}\alpha_{1}\\\ \alpha_{3}\\\ \alpha_{5}\\\
\alpha_{8}\end{pmatrix}$
$\displaystyle=\begin{pmatrix}-iFe^{\eta_{1}}&0&iFe^{\eta_{3}}&0\\\
0&-iFe^{\eta_{2}}&0&iFe^{\eta_{4}}\\\ 0&e^{\eta_{2}}&0&e^{\eta_{4}}\\\
e^{\eta_{1}}&0&e^{\eta_{3}}&0\\\ \end{pmatrix}\begin{pmatrix}u_{1}\\\ u_{2}\\\
u_{3}\\\ u_{4}\\\ \end{pmatrix}.$

Hence,

$\displaystyle\alpha_{1}=-u_{1}iFe^{\eta_{1}}+u_{3}iFe^{\eta_{3}},~{}~{}\alpha_{3}$
$\displaystyle=-u_{2}iFe^{\eta_{2}}+u_{4}iFe^{\eta_{4}},$
$\displaystyle\alpha_{5}=u_{2}e^{\eta_{2}}+u_{4}e^{\eta_{4}},~{}~{}\alpha_{8}$
$\displaystyle=u_{1}e^{\eta_{1}}+u_{3}e^{\eta_{3}}.$ (A.35)

The functions $u_{1}(\tau)$, $u_{2}(\tau)$, $u_{3}(\tau)$ and $u_{4}(\tau)$
can be defined as $u_{1}=\frac{q_{i}}{2i\omega_{i}}\xi_{1}$,
$u_{2}=i\frac{\omega_{i}-{\rm b}_{i}h}{2\omega_{i}T}\xi_{2}$,
$u_{3}=\frac{q_{i}}{2i\omega_{i}}\xi_{3}$, and $u_{4}=-i\frac{(\omega_{i}-{\rm
b}_{i}h)}{2\omega_{i}}\xi_{4}$. For time-varying field, we get
$F=B\sqrt{1+\frac{\tau_{R}^{i}}{\tau_{B}}}$ , $\int Fd\tau$ =
$B\tau\sqrt{1+\frac{\tau_{R}^{i}}{\tau_{B}}}$.
Hence, it leads to the given form,

$\displaystyle\eta_{j}$
$\displaystyle=-\frac{\tau}{\tau_{R}^{i}}+a_{j}i\frac{\sqrt{1+\frac{\tau_{R}^{i}}{\tau_{B}}}}{\tau_{B}}~{}\tau,$
$\displaystyle\xi_{j}$
$\displaystyle=\frac{1}{\sqrt{1+\frac{\tau_{R}^{i}}{\tau_{B}}}B_{0}}\frac{e^{\Big{(}\frac{1}{\tau_{R}^{i}}+\frac{1}{\tau_{B}}-a_{j}i\frac{\sqrt{1+\frac{\tau_{R}^{i}}{\tau_{B}}}}{\tau_{B}}\Big{)}\tau}}{\Big{(}\frac{1}{\tau_{R}^{i}}+\frac{1}{\tau_{B}}-a_{j}i\frac{\sqrt{1+\frac{\tau_{R}^{i}}{\tau_{B}}}}{\tau_{B}}\Big{)}}.$
(A.36)

Finally, the alphas in Eq. (A.2) are,

$\displaystyle\alpha_{1}$
$\displaystyle=-\frac{q_{i}}{\omega_{i}}~{}\frac{\tau_{R}^{i}}{(1+\chi_{i}+\chi_{i}^{2})},$
$\displaystyle\alpha_{2}$
$\displaystyle=\frac{-\tau_{R}^{i}(1+\chi_{i}-\chi_{i}^{2})}{(1+\chi_{i})(1+\chi_{i}^{2})}\alpha_{1}$
$\displaystyle\alpha_{3}$ $\displaystyle=\frac{(\omega_{i}-{\rm
b}_{i}h)}{\omega_{i}}~{}\frac{\tau_{R}^{i}}{(1+\chi_{i}+\chi_{i}^{2})},$
$\displaystyle\alpha_{4}$
$\displaystyle=\frac{-\tau_{R}^{i}(1+\chi_{i}-\chi_{i}^{2})}{(1+\chi_{i})(1+\chi_{i}^{2})}\alpha_{3}$
$\displaystyle\alpha_{5}$
$\displaystyle=\frac{\chi_{i}}{B(1+\chi_{i})}\alpha_{3}$
$\displaystyle\alpha_{6}$
$\displaystyle=-\tau_{R}^{i}\frac{\chi_{i}}{B(1+\chi_{i})}\alpha_{3}$
$\displaystyle\alpha_{7}$
$\displaystyle=-\tau_{R}^{i}\frac{\chi_{i}(2+\chi_{i})}{B(1+\chi_{i})(1+\chi_{i}^{2})}\alpha_{3}$
$\displaystyle\alpha_{8}$
$\displaystyle=\frac{\chi_{i}}{B(1+\chi_{i})}\alpha_{1}$
$\displaystyle\alpha_{9}$
$\displaystyle=-\tau_{R}^{i}\frac{\chi_{i}}{B(1+\chi_{i})}\alpha_{1}$
$\displaystyle\alpha_{10}$
$\displaystyle=-\tau_{R}^{i}\frac{\chi_{i}(2+\chi_{i})}{B(1+\chi_{i})(1+\chi_{i}^{2})}\alpha_{1}.$
(A.37)

After getting the expressions for $\alpha_{j}$’s, the simplified form of
$\delta f_{i}$ in Eq. (38) is,

$\displaystyle\delta f_{i}$
$\displaystyle=\frac{-q_{i}\tau_{R}^{i}}{(1+\chi_{i})(1+\chi_{i}+\chi_{i}^{2})}\Big{[}\Big{\\{}(1+\chi_{i})+\frac{\chi_{i}(1+\chi_{i}-\chi_{i}^{2})}{(1+\chi_{i}^{2})}\Big{\\}}(\vec{v_{i}}.\vec{E})\pm\Big{\\{}\chi_{i}(1+\chi_{i})$
$\displaystyle+\frac{\chi_{i}^{2}(2+\chi_{i})}{(1+\chi_{i}^{2})}\Big{\\}}\Big{(}\vec{v_{i}}.(\vec{E}\times\hat{b})\Big{)}\Big{]}+\frac{(\omega_{i}-{\rm
b}_{i}h)\tau_{R}^{i}}{T(1+\chi_{i})(1+\chi_{i}+\chi_{i}^{2})}\Big{[}(1+\chi_{i})\Big{(}\vec{v_{i}}\cdot\vec{\nabla}{T}\Big{)}$
$\displaystyle-\tau_{R}^{i}\frac{(1+\chi_{i}-\chi_{i}^{2})}{(1+\chi_{i}^{2})}\Big{(}\vec{v_{i}}\cdot\vec{\nabla}{\dot{T}}\Big{)}\pm\chi_{i}(1+\chi_{i})\Big{(}\vec{v_{i}}\cdot(\vec{\nabla}{T}\times\hat{b})\Big{)}\mp\tau_{R}^{i}\frac{\chi_{i}(2+\chi_{i})}{(1+\chi_{i}^{2})}$
$\displaystyle\Big{(}\vec{v_{i}}\cdot(\vec{\nabla}{\dot{T}}\times\hat{b})\Big{)}\Big{]}\frac{\partial
f^{0}_{i}}{\partial\omega_{i}}.$ (A.38)

For simplicity we have considered $\tau_{E}=\tau_{B}$ and
$\chi_{i}=\frac{\tau_{R}^{i}}{\tau_{B}}=\frac{\tau_{R}^{i}}{\tau_{E}}$. Here
$\pm$ indicates positive sign for positively charged particles(antiparticles)
and negative sign for negatively charged particles (antiparticles).

###### Acknowledgements.

KS acknowledge the doctoral fellowship from UGC, Government of India. JD and
RS gratefully acknowledge the DAE-DST, Govt. of India funding under the mega-
science project – “Indian participation in the ALICE experiment at CERN”
bearing Project No. SR/MF/PS-02/2021- IITI (E-37123).

## References

  * (1) J. Brewer, L. Yan and Y. Yin, _Adiabatic hydrodynamization in rapidly-expanding quark–gluon plasma_ , _Physics Letters B_ 816 (2021) 136189.
  * (2) W. Busza, K. Rajagopal and W. van der Schee, _Heavy Ion Collisions: The Big Picture, and the Big Questions_ , _Ann. Rev. Nucl. Part. Sci._ 68 (2018) 339 [1802.04801].
  * (3) H. Elfner and B. Müller, _The exploration of hot and dense nuclear matter: introduction to relativistic heavy-ion physics_ , _J. Phys. G_ 50 (2023) 103001 [2210.12056].
  * (4) S. Alzhrani, S. Ryu and C. Shen, _$\mathrm{\Lambda}$ spin polarization in event-by-event relativistic heavy-ion collisions_, _Phys. Rev. C_ 106 (2022) 014905.
  * (5) J. Liao, _Chiral magnetic effect in heavy ion collisions_ , _Nuclear Physics A_ 956 (2016) 99.
  * (6) STAR Collaboration collaboration, _Identified particle elliptic flow in $au+au$ collisions at $\sqrt{{s}_{\mathrm{NN}}}\phantom{\rule{0.0pt}{0.0pt}}=\phantom{\rule{0.0pt}{0.0pt}}130\mathrm{GeV}$_, _Phys. Rev. Lett._ 87 (2001) 182301.
  * (7) STAR Collaboration collaboration, _Centrality dependence of charged hadron and strange hadron elliptic flow from $\sqrt{{s}_{\mathit{NN}}}=200$ gev au+au collisions_, _Phys. Rev. C_ 77 (2008) 054901.
  * (8) L. Del Zanna, V.C. Joshi, G. Inghirami, V. Rolando, A. Beraudo, A. Pace et al., _Relativistic viscous hydrodynamics for heavy-ion collisions with echo-qgp_ , _European Physical Journal C_ 73 (2013) .
  * (9) Z.V. Khaidukov and Y.A. Simonov, _Thermodynamics of a quark-gluon plasma at finite baryon density_ , _Phys. Rev. D_ 100 (2019) 076009.
  * (10) S. Koothottil and V.M. Bannur, _Thermodynamic behavior of magnetized quark-gluon plasma within the self-consistent quasiparticle model_ , _Phys. Rev. C_ 99 (2019) 035210.
  * (11) K. Goswami, D. Sahu, J. Dey, R. Sahoo and R. Stock, _Anisotropy of magnetized quark matter_ , _Phys. Rev. D_ 109 (2024) 074012 [2310.02711].
  * (12) S. Gavin, _TRANSPORT COEFFICIENTS IN ULTRARELATIVISTIC HEAVY ION COLLISIONS_ , _Nucl. Phys. A_ 435 (1985) 826.
  * (13) A. Das, H. Mishra and R.K. Mohapatra, _Transport coefficients of hot and dense hadron gas in a magnetic field: A relaxation time approach_ , _Phys. Rev. D_ 100 (2019) 114004.
  * (14) A. Dash, S. Samanta, J. Dey, U. Gangopadhyaya, S. Ghosh and V. Roy, _Anisotropic transport properties of a hadron resonance gas in a magnetic field_ , _Phys. Rev. D_ 102 (2020) 016016.
  * (15) K.K. Gowthama, M. Kurian and V. Chandra, _Thermal and thermoelectric responses of hot qcd medium in time-varying magnetic fields_ , _Phys. Rev. D_ 106 (2022) 034008.
  * (16) K. Singh, J. Dey, R. Sahoo and S. Ghosh, _Effect of time-varying electromagnetic field on wiedemann-franz law in a hot hadronic matter_ , _Phys. Rev. D_ 108 (2023) 094007.
  * (17) K. Tuchin, _Initial value problem for magnetic fields in heavy ion collisions_ , _Phys. Rev. C_ 93 (2016) 014905.
  * (18) D. Satow, _Nonlinear electromagnetic response in quark-gluon plasma_ , _Phys. Rev. D_ 90 (2014) 034018.
  * (19) ALICE Collaboration collaboration, _Probing the effects of strong electromagnetic fields with charge-dependent directed flow in pb-pb collisions at the lhc_ , _Phys. Rev. Lett._ 125 (2020) 022301.
  * (20) 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”_ , _Nuclear Physics A_ 803 (2008) 227.
  * (21) STAR Collaboration collaboration, _Observation of the electromagnetic field effect via charge-dependent directed flow in heavy-ion collisions at the relativistic heavy ion collider_ , _Phys. Rev. X_ 14 (2024) 011028.
  * (22) Y. Kumar, P. Jain, P. Bangotra, V. Kumar, D.V. Singh and S.K. Rajouria, _Role of time-varying magnetic field on qgp equation of state_ , _Advances in High Energy Physics_ 2024 (2024) 1870528.
  * (23) B.B. Brandt, F. Cuteri and G. Endrődi, _Equation of state and speed of sound of isospin-asymmetric qcd on the lattice_ , _Journal of High Energy Physics_ 2023 (2023) 55.
  * (24) 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_, _Journal of High Energy Physics_ 2012 (2012) 53.
  * (25) G.S. Bali, F. Bruckmann, G. Endrődi, S.D. Katz and A. Schäfer, _The qcd equation of state in background magnetic fields_ , _Journal of High Energy Physics_ 2014 (2014) 177.
  * (26) S. Li and H.-U. Yee, _Shear viscosity of the quark-gluon plasma in a weak magnetic field in perturbative qcd: Leading log_ , _Phys. Rev. D_ 97 (2018) 056024.
  * (27) V.P. Gusynin, V.A. Miransky and I.A. Shovkovy, _Catalysis of dynamical flavor symmetry breaking by a magnetic field in 2 + 1 dimensions_ , _Phys. Rev. Lett._ 73 (1994) 3499.
  * (28) D.-S. Lee, C.N. Leung and Y.J. Ng, _Chiral symmetry breaking in a uniform external magnetic field_ , _Phys. Rev. D_ 55 (1997) 6504.
  * (29) S.S. Avancini, D.P. Menezes, M.B. Pinto and C.m.c. Providência, _Qcd critical end point under strong magnetic fields_ , _Phys. Rev. D_ 85 (2012) 091901.
  * (30) G. Endrődi, Z. Fodor, S.D. Katz, Szabó and K. K., _The qcd phase diagram at nonzero quark density_ , _Journal of High Energy Physics_ 2011 (2011) 1.
  * (31) X. Wang and I. Shovkovy, _Polarization tensor of magnetized quark-gluon plasma at nonzero baryon density_ , _Eur. Phys. J. C_ 81 (2021) 901.
  * (32) B. McInnes, _A rotation/magnetism analogy for the quark–gluon plasma_ , _Nuclear Physics B_ 911 (2016) 173.
  * (33) Y.A. Simonov, _The Fundamental Scale of QCD_ , _Phys. Atom. Nucl._ 84 (2021) 1195 [2103.08223].
  * (34) J. Dey, A. Bandyopadhyay, A. Gupta, N. Pujari and S. Ghosh, _Electrical conductivity of strongly magnetized dense quark matter - possibility of quantum Hall effect_ , _Nucl. Phys. A_ 1034 (2023) 122654.
  * (35) A. Dash and A.K. Panda, _Charged participants and their electromagnetic fields in an expanding fluid_ , _Physics Letters B_ 848 (2024) 138342.
  * (36) H.J. Goldsmid et al., _Introduction to thermoelectricity_ , vol. 121, Springer https://doi.org/10.1007/978-3-662-49256-7.
  * (37) A. Das and H. Mishra, _Thermoelectric transport coefficients of hot and dense qcd matter_ , _The European Physical Journal Special Topics_ 230 (2021) 607.
  * (38) D. Dey and B.K. Patra, _Seebeck effect in a thermal qcd medium in the presence of strong magnetic field_ , _Phys. Rev. D_ 102 (2020) 096011.
  * (39) A. Das, H. Mishra and R.K. Mohapatra, _Magneto-seebeck coefficient and nernst coefficient of a hot and dense hadron gas_ , _Phys. Rev. D_ 102 (2020) 014030.
  * (40) J.R. Bhatt, A. Das and H. Mishra, _Thermoelectric effect and seebeck coefficient for hot and dense hadronic matter_ , _Phys. Rev. D_ 99 (2019) 014015.
  * (41) K. Singh, J. Dey and R. Sahoo, _Thermal conductivity of evolving quark-gluon plasma in the presence of a time-varying magnetic field_ , _Phys. Rev. D_ 109 (2024) 014018.
  * (42) D.H. Rischke, Y. Pursun and J.A. Maruhn, _Relativistic hydrodynamics for heavy ion collisions. 2. Compression of nuclear matter and the phase transition to the quark - gluon plasma_ , _Nucl. Phys. A_ 595 (1995) 383 [nucl-th/9504021].
  * (43) K. Okamoto, Y. Akamatsu and C. Nonaka, _A new relativistic hydrodynamics code for high-energy heavy-ion collisions_ , _The European Physical Journal C_ 76 (2016) 579.
  * (44) J.D. Bjorken, _Highly relativistic nucleus-nucleus collisions: The central rapidity region_ , _Phys. Rev. D_ 27 (1983) 140.
  * (45) S.S. Gubser, _Symmetry constraints on generalizations of bjorken flow_ , _Phys. Rev. D_ 82 (2010) 085027.
  * (46) S. Mitra, _Relativistic hydrodynamics with momentum-dependent relaxation time_ , _Phys. Rev. C_ 103 (2021) 014905.
  * (47) P. Romatschke and U. Romatschke, _Viscosity information from relativistic nuclear collisions: How perfect is the fluid observed at rhic?_ , _Phys. Rev. Lett._ 99 (2007) 172301.
  * (48) V. Roy, S. Pu, L. Rezzolla and D. Rischke, _Analytic bjorken flow in one-dimensional relativistic magnetohydrodynamics_ , _Physics Letters B_ 750 (2015) 45.
  * (49) A.K. Panda, A. Dash, R. Biswas and V. Roy, _Relativistic non-resistive viscous magnetohydrodynamics from the kinetic theory: a relaxation time approach_ , _Journal of High Energy Physics_ 2021 (2021) 216.
  * (50) D. Satow, _Nonlinear electromagnetic response in quark-gluon plasma_ , _Phys. Rev. D_ 90 (2014) 034018 [1406.7032].
  * (51) A. Abhishek, A. Das, D. Kumar and H. Mishra, _Thermoelectric transport coefficients of quark matter_ , _Eur. Phys. J. C_ 82 (2022) 71 [2007.14757].
  * (52) M. Kurian, _Thermoelectric behavior of hot collisional and magnetized QCD medium from an effective kinetic theory_ , _Phys. Rev. D_ 103 (2021) 054024 [2102.00435].
  * (53) K.K. Gowthama, M. Kurian and V. Chandra, _Electromagnetic response of hot qcd medium in the presence of background time-varying fields_ , _Phys. Rev. D_ 104 (2021) 094037.
  * (54) M.I. Gorenstein and S.-N. Yang, _Gluon plasma with a medium dependent dispersion relation_ , _Phys. Rev. D_ 52 (1995) 5206.
  * (55) P.K. Srivastava, S.K. Tiwari and C.P. Singh, _QCD Critical Point in a Quasiparticle Model_ , _Phys. Rev. D_ 82 (2010) 014023.
  * (56) V.M. Bannur, _Comments on quasiparticle models of quark-gluon plasma_ , _Phys. Lett. B_ 647 (2007) 271.
  * (57) A. Peshier, B. Kampfer, O.P. Pavlenko and G. Soff, _A Massive quasiparticle model of the SU(3) gluon plasma_ , _Phys. Rev. D_ 54 (1996) 2399.
  * (58) A. Peshier, B. Kampfer and G. Soff, _The Equation of state of deconfined matter at finite chemical potential in a quasiparticle description_ , _Phys. Rev. C_ 61 (2000) 045203.
  * (59) A. Hosoya and K. Kajantie, _Transport Coefficients of QCD Matter_ , _Nucl. Phys. B_ 250 (1985) 666.
  * (60) M. Hongo, Y. Hirono and T. Hirano, _Anomalous-hydrodynamic analysis of charge-dependent elliptic flow in heavy-ion collisions_ , _Phys. Lett. B_ 775 (2017) 266.
  * (61) J. Dey, S. Samanta, S. Ghosh and S. Satapathy, _Quantum expression for the electrical conductivity of massless quark matter and of the hadron resonance gas in the presence of a magnetic field_ , _Phys. Rev. C_ 106 (2022) 044914.
  * (62) K. Fukushima, D.E. Kharzeev and H.J. Warringa, _Chiral magnetic effect_ , _Phys. Rev. D_ 78 (2008) 074033.
  * (63) V.P. Gusynin and I.A. Shovkovy, _Chiral symmetry breaking in qed in a magnetic field at finite temperature_ , _Phys. Rev. D_ 56 (1997) 5251.
  * (64) K. Tuchin, _Particle production in strong electromagnetic fields in relativistic heavy-ion collisions_ , _Advances in High Energy Physics_ 2013 (2013) 490495.