2.1 Numerical approach

A number of numerical simulations have revealed that the magnitudes of the magnetic fields in heavy-ion collisions can be very large [5-16]. In Fig. 1, the different spatial components of the magnetic fields at r= 0 (i.e., the center of the overlapping region) and t=0 (i.e., the time when the two colliding nuclei completely overlap) are presented as functions of the impact parameter, b. The fields are calculated by using the Lienard-Wiechert potential on the event-by-event basis where the events are generated by the HIJING event generater [10]. Please find the discussions on the quantum correction to the otherwise classical Leinard-Wiechert potential in Refs. [16, 93]. The most important information from Fig. 1 are three folds.

Figure 1:Full-screen View

The magnetic fields at t=0 and r= 0 as functions of the impact parameter b where $\langle \cdots \rangle$ denotes average over events. (Reproduced from Ref. [10].)

First, it shows that heavy-ion collision experiments at RHIC and LHC generate very strong magnetic fields. The event-averaged magnetic fields at Au + Au collisions at $\sqrt{s}=200$ GeV can reach about 10¹⁹ Gauss, and at Pb + Pb collisions at $\sqrt{s}=2.76$ TeV can reach about 10²⁰ Gauss. They are much stronger than the magnetic fields of neutron stars, including the magnetars which may have surface magnetic fields of the order of 10¹⁴-10¹⁵ Gauss [94, 95]. They are also much larger than the masses squared of electron, $m_e^2$, and light quarks, $m_u^2\mathrm{,}{m}_d^2$, and thus are capable of inducing significant quantum phenomena.

Second, although we do see that only the y-component (y axis is set to be perpendicular to the reaction plane) of the magnetic fields remains after the event average (due to the left-right symmetry of the noncentral collisions), we can find that the event-by-event fluctuation of the nuclear distribution results in nonzero x-component of the magnetic fields as well. These are reflected in the averaged absolute values of the fields in Fig. 1 and is most evident for central collisions [8, 10].

Third, the fields in LHC Pb + Pb at 2.76 TeV are roughly 13.8=2.76/0.2 times larger than that in RHIC Au + Au collisions at 200 GeV. In fact, as investigated in Refs. [8, 10], to high precision, the magnitudes of the magnetic fields linearly depend on the collision energy, $\sqrt{s}$. This can be understood from the fact that the fields are proportional to the Lorentz gamma factor of the beam motion which is equal to $\sqrt{s}/(2m_{\text{N}})$, where m_N is the mass of a nucleon.

The simulations depicted in Fig. 1 are only for fields at t=0. As the hot medium generated in heavy-ion collisions is expanding fast, a natural question is: how do the magnetic fields evolve in time in the expanding medium? This is still an unsolved problem. We summarize the recent progress briefly here and in the next subsection.

As we explained, the strong initial magnetic fields are consequences of the smallness of the system and the relativistically fast motion of the colliding nuclei. The main contributors to the initial magnetic fields are the spectator protons (which do not participate into the collision) in the colliding uclei. They leave the collision zone very fast in a time scale characterized by the Lorentz-suppressed radius of the nucleus, i.e., ${\tau }_B=R_A/\gamma \approx 2m_N{R}_A/\sqrt{s}$. Thus in relativistic heavy-ion collisions, one can expect that the magnetic fields generated by the spectators will decay very fast after a short time scale, τB. This fact is indeed verified by numerical simulations. In Fig. 2, one can find the numerical results of the time evolutions of the magnetic fields at r=0 in collisions with b=10 fm for Au + Au collision at $\sqrt{s}=200$ GeV and for Pb + Pb collision at $\sqrt{s}=2.76$ TeV [10]. We can see that the magnetic fields decay very fast after the collision reflecting the fact that the spectators are leaving the collision region very fast. Once the spectators are all far away from the collision region, the remnant charge carriers, which participate in the collisions and thus moves much slower than the spectators along the beam direction, become the dominant sources for the magnetic fields, and they essentially slow down the decay of the magnetic fields.

Figure 2:Full-screen View

The time evolution of the magnetic fields at r=0 with impact parameter b=10 for Au + Au collision at $\sqrt{s}=200$ GeV and Pb + Pb collision at $\sqrt{s}=2.76$ TeV. The effect due to the finite conductivity of quark-gluon matter is not considered. (Reproduced from Ref. [10].)

2.2 Analytical approach

To understand the basic picture of the creation of the magnetic fields, we compared the numerical results shown in the previous section with a simple analytical model proposed in Ref. [28]. With the relativistic beam energy, the momenta of the colliding nuclei are so large that the nuclear stopping power is not really strong enough to stop them: the incident nuclei pass each other and recede from the collision point with almost the same rapidity as the incident beam. Therefore, in this simple model, we assume that a magnetic field is created by protons freely streaming with the beam rapidity, and thus do not distinguish the spectators and participants. From the comparison between the numerical and analytic modeling, we will also discuss effects of the nuclear stopping power.

As we have specified with the classical motion of charged particles, the magnitude of the magnetic field can be estimated by using the Lienard-Wiechert potential

where Y_beam is the absolute value of the beam rapidity. The positive and negative signs correspond to the nuclei moving in the positive and negative z-direction, respectively. The coordinate system is illustrated in Fig. 3 (Left). We have defined the coordinate vector $\tilde{\text{x}}$ from the the center of the nucleus (±b/2,0,v^±t) to the observer at x, and its length in the rest frame $r^{\pm }=\sqrt{{(x\pm b/2)}^2+y^2+{\gamma }^2{(z-v^{\pm }t)}^2}$ with the impact parameter b, the velocity of the nucleus v^±=± tanh Y_beam, and the gamma factor $\gamma =1/\sqrt{1-{(v^{\pm })}^2}$.

Figure 3:Full-screen View

Effective charges contributing to the creation of the magnetic field in the rest and the center-of-mass (COM) frame of the collision (Left). Impact parameter dependence of the effective charge of a nucleus at the origin of the COM coordinate (Right).

The Lienard-Wiechert potential in Eq. (1) is obtained by boosting the electro-static potential in the rest frame of the nuclei (see Fig. 3). As is well known, the net electro-static potential can be calculated from the effective point-like charge sitting at the center of the sphere. Therefore, the effective charge Z_eff appearing in Eq. (1) should be the total amount of the charges inside the sphere, i.e.,

where ρ(r′) is the charge distribution inside the nuclei. Note that the upper boundary of the integral is given by the distance r^± from the center to the observer. When the observer is in the exterior of the nucleus r^±≥RA, one simply has the total charge of the nucleus Z_eff(r^±>RA)=Z. The impact parameter and time dependences of the effective charge is shown in the right panel of Fig. 3. We find significant dependence on both the impact parameter and the time collapsed after the collision, which will be important to reproduce the results from the numerical simulation shown in the previous section.

Plugging the effective charge (2) to the Lienard-Wiechert potential (1), we obtain the magnitude of the magnetic field shown in Figs. 4, 5 and 6. Here, we have taken the Woods-Saxon form for the charge distribution ρ(r)=NZ/{1+exp[(r-RA)/a]} with $N_Z=Z{\{4\pi {\displaystyle {\int }_0^{R_A}}\text{}d{r}^{\text{'}}{r}^{\text{'}2}{\rho }_{}(r^{\text{'}})\}}^{-1}$ where Z = 79, RA = 6.38 fm, and a = 0.535 fm for Au, and Z = 82, RA = 6.62 fm, and a = 0.546 fm for Pb. We do not consider the event-by-event fluctuation of the charge distribution in our analytical modeling so that we in Figs. 4, 5 and 6 compare the analytical modeling with the numerical simulations of event-averaged magnetic fields.

Figure 4:Full-screen View

The analytic model (solid lines) reproduces well the time dependence of the magnetic field at x = 0 from the numerical simulation (dots) shown in Fig. 2. The mismatches seen at late time are due to the lack of considering nuclear stopping power in the analytical modeling.

Figure 5:Full-screen View

The analytic model (solid lines) reproduces well the impact parameter dependence of the event-averaged magnetic field at t=0 and x = 0 from the numerical simulation (dots) shown in Fig. 1. The pointlike nucleus model (dashed lines) fails to reproduce the numerical results when the nuclei have an overlap with each other, i.e., $b\lesssim 2R_A$.

Figure 6:Full-screen View

Beam-energy dependence of the magnetic field from the analytic model (left) and the numerical simulation (right) taken from Ref. [7].

Figure 4 shows the time dependence of the magnetic field obtained from the analytic model (solid lines) and the event-averaged numerical simulations (dots). The latter has been already shown in Fig. 2. We find a good agreement between the results from these two models in the time range before and after the collision except for the late time. The origin of the deviation in the late time is, as can be guessed from the asymmetric time dependence in the numerical result with respect to the origin of time, attributed to the effect of the nuclear stopping power. Because of the deceleration by the nuclear stopping power, participant nucleons (also called remnants) will stay in the collision region longer than in the free streaming case so that they will maintain the magnetic field for a longer time.

Figure 5 shows the impact parameter dependence of the magnetic field. Again, the sold lines and dots show the results from the analytic model and the event-averaged numerical simulations, respectively, and we find a good agreement between them. The slight deviation may be originated from the energy loss due to the nuclear stopping power in the numerical simulation. In addition to these results, we show the magnetic field created by point-like nuclei carrying the total charge Z, i.e., Z_eff (t,x)=Z in Eq. (1). We find that this simplest model completely fails to reproduce the numerical results when there is an overlap between the nuclei ($b\lesssim 2R_A$), indicating the importance of taking into account the effective charge (2).

Finally, Fig. 6 shows the beam energy dependence of the magnetic field from the analytic model (left panel) and the numerical simulation (right panel) [7]. While we find a good agreement, the nuclear stopping power will play an important role in the low energy collisions such as those in the beam energy scan program in RHIC, FAIR, and J-PARC.

In the above discussions, we omit one important issue which may substantially modify the time evolution of the magnetic fields, that is, the finite electric conductivity, σ, of the produced quark-gluon matter. If σ is large, the Faraday induction effect may be strong as well and the lifetime of the magnetic fields can be significantly prolonged. As a matter of fact, the equilibrium QGP is a very good conductor according to the theoretical and lattice QCD studies [96-105]; see also Ref. [106] for review.

Now let us understand how the large σ induces the strong Faraday induction and prolongs the lifetime of the magnetic fields; more discussions can be found in review [93]. Let us write down the Maxwell’s equations in the background of a charged flow characterized by velocity field v,

where J is the electric current which according to the Ohm’s law is determined by

Substituting J into Eq. (4) and eliminating E by using Eq. (3), one can obtain an equation that governs the time evolution of the magnetic field

where we have regarded the QGP as locally neutral (this is a good assumption for high-energy collisions) so that ∇⋅E=ρ=0.

The Faraday prolongation of the lifetime of the magnetic fields is most easily seen in two limits. The first limit is the no-flow limit, i.e., the limit when v can be neglected. In this case, if the electric conductivity is so large that $\sigma \gg 1/t_{\text{c}}$ with t_c being the characteristic time scale over which the field strongly varies, one can neglect the second-order time derivative and render Eq. (6) a diffusion equation [79]:

Clearly, in this case, the decay of the magnetic field is due to diffusion, and the diffusion time is given by

where L characterizes the length scale of the QGP system. For sufficiently large σ, the diffusion time, t_D, could be larger than τ_B [79].

The second limit is the ideal magnetohydrodynamic limit, i.e., when σ is very large and at the same time the flow velocity, v, is also large so that we can neglect the term divided by σ in Eq. (6) and obtain

This equation is known to lead to the so-called frozen-in theorem, i.e., the magnetic lines are frozen in the ideally conducting plasma. Thus in this limit, the decay of the mangnetic fields is totally due to the expansion of the QGP. For high-energy heavy-ion collisions, the longitudinal expansion is well described by the Bjorken flow, vz=z/t, which makes the area of cross section of the overlapping region on the reaction plane grow linearly in time. Thus the frozen-in theorem leads us to the conclusion that in the ideal magnetohydrodynamic limit, before the transverse expansion play significant role, the magnetic fields will decay inversely proportional to time, By(t)∝(t₀/t)By(t₀), which is much slower than the 1/t³-type decay in the insulating case. The transverse expansion of the system can happen because of the transverse gradient of the fluid pressure. When this transverse expansion becomes fast, the decay of the fields will become fast as well, see discussions in Refs. [10, 93].

So far, we discussed two special limiting cases from which we can analytically extract the effect of Faraday induction. However, the realistic environment produced by heavy-ion collisions is not the case as described by the two limits. Moreover, the system is essentially non-equilibrium. Thus, the precise time evolution of the magnetic fields is really a hard task to achieve, and some recent progresses can be found in Refs. [10, 15, 79, 82, 107-113].

In the above discussions, we focus only on collisions of the same spherical nuclei, i.e., the Au + Au collisions and Pb + Pb collisions. The Heavy-ion program at RHIC also performed collisions of different nuclei, like Cu + Au collisions, as well as collisions of prolate nuclei like the U + U collisions. Of course, these collisions can also generate strong magnetic fields; see for example, the studies in Refs [12, 114-118] and recent review [93]. It is worth mentioning that, owing to the charge number asymmetry between Cu and Au nuclei, the Cu + Au collisions can generate a strong electric field pointing from Au nucleus to Cu nucleus [114-117, 119]. Such E field may lead to a splitting between the directed flow v₁ of positively and negatively charged hadrons, which has been recently test by STAR Collaboration at RHIC [120].

3.1 The chiral magnetic effect

The CME represents the generation of vector current induced by an external magnetic field in chirality-imbalanced (P-odd) medium [5, 17]. The CME composes a very general class of anomalous transport phenomena spreading in a wide contexts in physics, ranging from astrophysics, condensed matter physics, nuclear physics, to particle physics. We apologize for not being able to mention all the relevant works. We refer the readers to other review articles as mentioned above for more information. The CME can be expressed as

where e is the charge of the fermion, JV is the vector current defined by $J_V^{\mu }=e\langle \overline{\psi }{\gamma }^{\mu }\psi \rangle$, and μA is the chiral chemical potential which parameterizes how large the chirality imbalance of the medium is. Note that, the above expression for CME current is for each species of massless fermions. For QGP at high temperature where the masses of up and down quarks can be neglected, the total CME conductivity is given by ${\sigma }_{\text{C}ME}=N_C{\mu }_A{\displaystyle {\sum }_f}{q}_f^2/(2{\pi }^2)$ with qf the charge of quark of flavor, f, and N_c=3 the number of color.

Although the CME stems from the profound chiral anomaly in gauge theory, its physical origin can be understood in an intuitive way by adopting the Landau quantization picture. We know that the electrically charged fermions in a magnetic field occupy a tower of Landau levels (i.e., the energy spectrum of charged fermions in a magnetic field is discrete) in such a way that the lowest Landau level contains fermions with only one spin polarization (i.e., the spin of positively (negatively) charged fermions is polarized to be along (opposite to) the direction of the magnetic field) while higher Landau levels are occupied by equal numbers of up-polarized and down-polarized fermions. Now, let us consider the situation that the number of right-handed (RH) and the number of left-handed (LH) fermions are unequal, say, for example, NR>NL. Thus there are more RH fermions than LH fermions confined in the lowest Landau level which means that a current consisting of positively charged fermions will flow along the magnetic field as RH fermions lock their momenta with their spins, while a current consisting of negatively charged fermions will flow opposite to the magnetic field. The net effect is the appearance of an electric current along the magnetic field. This is the CME current. Note that the higher Landau levels do not contribute to the CME current as higher Landau levels are occupied symmetrically by spin-up and spin-down fermions and their contribution to the CME current cancel.

Now let us quantify the above intuitive argument to see how we can arrive at Eq. (10). To simplify our lives, we consider a system of massless Dirac fermions with only positive charge, e, in an constant magnetic field, B, pointing along z direction which is supposed to be so strong that all the fermions are confined in the lowest Landau level. The number densities of RH and LH fermions are given by

where the prefactor eB/(2π) is the density of states of the lowest Landau level in the transverse directions and pF denotes the Fermi momentum (which is the momentum that separates the empty states and occupied states). Moreover, chiral fermions on the lowest Landau level move at the speed of light in a way correlated to their spins which tell us that the RH fermions generate a current along z direction, JR=enR, while LH fermions generate a current along -z direction, JL=-enL. Therefore the net vector current reads

which, once we equate $p_F^{R/L}$ with μR/L and μA with (μR-μL)/2, leads us to Eq. (10) and Eq. (11). Here we use the fact that for massless fermions, the chemical potential is equal to the Fermi momentum.

We emphasize that although the above argument is for a non-interacting system and relies on Landau quantization picture, the CME conductivity σ_CME=e²μA/(2π²) is actually fixed by the chiral anomaly equation and universal in the sense that it does not involve concrete microscopic scattering among fermions. This is in parallel to the fact that the chiral anomaly equation itself is universal so that it does not receive perturbative correction from scattering between fermions.

Apparently, the occurrence of CME requires an environmental parity violation as being characterized by μA. A nonzero μA may be induced by different mechanisms in different physical systems. For example, the electro-weak plasma can have a nonzero μA because the weak interaction does not respect parity symmetry [127]. Recently, it is realized in condensed matter experiments the so-called Weyl or Dirac semimetals, whose band structures permit level crossing around the low energy excitations, are Weyl or Dirac fermions. The CME was observed in these semimetals by applying parallel electric and magnetic fields, which induce a difference between the chemical potentials of RH and LH fermions via the chiral anomaly [128-131]; see also discussions in Refs. [92, 93, 126, 132]. In the hot quark-gluon matter produced in heavy-ion collisions, μA may appear due to the topological transition between the degenerate vacua of the gluonic sector of QCD characterized by different winding numbers [5, 133, 134]. The rate of such topological transition is exponentially suppressed at zero temperature as two neighbor vacua are separated by a finite energy barrier of the order of QCD confinement scale Λ_QCD∼200 MeV. However, at high temperature, the topological transition can be induced by another classical and thermal excitation called sphaleron [135, 136] which, instead of tunneling through the barrier, can take the vacuum over the barrier at a large rate. At equilibrium, the probability of having a topological transition that drives the winding number to increase is equal to the probability of having an opposite topological transition that decreases the winding number. Thus, what really matters is the fluctuation of the topological transition above global equilibrium which results in local domains in the quark-gluon matter within which the gluonic configuration has a nonzero winding number. The quarks coupled with such gluons via the chiral anomaly can acquire nonzero chirality and thus a nonzero μA. Furthermore, μA in heavy-ion collisions may also be induced by other mechanisms. For example, the parallel chromoelectric and chromomagnetic fields, Ea·Ba≠0, in the initial glasma can also induce a nonzero μA [137-139].

3.2 Experimental signals of the chiral magnetic effect

If the CME occurs in the quark-gluon matter produced in the heavy-ion collisions, what will be its experimental signals? We have already seen that the magnetic fields in heavy-ion collisions are, on average over events, perpendicular to the reaction plane; thus the CME is expected to induce an electric current which on average is also perpendicular to the reaction plane. Such an electric current will drive a charge dipole in the fireball which, after the hydrodynamic evolution of the fireball, will be converted to a charge dipole in momentum space. However, as μA can flip the sign from event to event, the charge dipole could also flip its direction from event to event, and thus the dipolar distribution in different events cancel each other after event average and cannot be measured directly. Instead, one should look for an observable which measures the dipolar fluctuation over events. Such an observable can be chosen as the following charge-dependent two-particle azimuthal correlation [133, 140]:

where the indices α and β denote the charge of the hadrons, φα and φβ are the azimuthal angles of hadrons of charge α and β, respectively, Ψ_RP is the reaction plane angle, and <…> denotes average over events; see Fig. 7 for illustration. It is easy to understand that a charge separation with respect to the reaction plane, will give a positive γ_+- (and γ_-+) and a negative γ₊₊ or γ_–. We will refer to γ_+-=γ_-+ as γ_OS and γ₊₊≈γ_– as γ_SS where "OS" stands for "opposite-sign" and "SS" stands for "same-sign". In real experiments, the measurements are done with three-particle correlations where the third hadron (of arbitrary charge) is used to reconstruct the reaction plane.

Figure 7:Full-screen View

Illustration of a typical noncentral collision.

The correlation γαβ was measured by STAR Collaboration at RHIC for Au + Au and Cu + Cu collisions at $\sqrt{s}=200$ GeV [141-143] (see Fig. 8), by PHENIX Collaboration at RHIC for Au + Au collisions at $\sqrt{s}=200$ GeV [144], and by ALICE Collaboration at LHC for for Pb + Pb collisions at $\sqrt{s}=2.76$ TeV [145]. At mid-central collisions, these measurements show clear positive opposite-sign correlation and negative same-sign correlation, as we expect from the picture of CME. Furthermore, the observed correlation increases from zero at central collisions to peripheral collisions. This is consistent with what we learned from the simulation in last section that the averaged magnetic field increases with the centrality. More recently, the STAR Collaboration performed the measurement of γαβ at different beam energies [146, 147] and observed that γαβ persists as long as the beam energy is larger than 19.6 GeV; for further lowered beam energies, the difference between γ_OS and γ_SS steeply falls down, which may be understood by noticing that at lower energies the system is probably in a hadronic phase where the chiral symmetry is broken and the CME is strongly suppressed.

Figure 8:Full-screen View

The measured correlation γαβ by STAR Collaboration. This figure is from Ref. [141].

Although the experimental data supports the presence of the CME, the interpretation of the data is not unique owing to possible background effects that are not related to the CME. First, the transverse momentum conservation (TMC) can contribute to γαβ [148-150]. The momentum conservation enforces the sum of the momenta projected to the transverse plane (the plane that is perpendicular to the beam direction) of all the particles to be zero at any time during the evolution of the system. Thus, the detection of a hadron of transverse momentum, p_⊥, requires the presence of other produced hadrons whose transverse momenta must be summed to -p_⊥. This enforces an intrinsic back-to-back two-particle correlation:

where we used that ${\displaystyle {\sum }_{i=1}^N\cos }({\varphi }_i-{\Psi }_{\text{RP}})={\displaystyle {\sum }_{i=1}^N\sin }({\varphi }_i-{\Psi }_{\text{RP}})=0$ and N is the total multiplicity. From the above relation, we can estimate how large the TMC contribution is. For example, let us consider RHIC Au + Au at 200 GeV, v₂∼0.1 while N∼1000 which gives that ${\gamma }_{\text{S}S}^{\text{TMC}}={\gamma }_{\text{O}S}^{\text{TMC}}\sim -{10}^{-4}$ which is roughly at the same order as the experimental result shown in Fig. 8. Obviously, the TMC is charge blind. Thus we can subtract the contribution from TMC by making a difference between opposite-sign and same-sign correlations:

However, there are other background effects contributing to Δγ. The most notable one may be the local charge conservation (LCC) [148, 151, 152] which assumes that the pairs of opposite charges are created locally (i.e., within a small volume) in space at the late stage of the heavy-ion collisions. Then, due to the elliptic flow, opposite-sign pairs are strongly correlated while the same-sign pairs are not correlated; thus LCC gives a finite contribution to Δγ. The detailed analysis reveals that the LCC contribute to Δγ the amount of [93]

where M is the average number of hadrons in a typical local neutral cell. As M is expected to be not very large, so we expect the LCC contribution is at the order of 10^-4 or even larger in, e.g., RHIC Au + Au collisions at 200 GeV. This is comparable to the experimental results in Fig. 8.

The background effects mask the possible CME signals. To disentangle the possible CME signal and the flow-related backgrounds, one can setup experiments to either vary the backgrounds with the signal fixed or vary the signal with the backgrounds fixed (It is also important to setup delicate analysis methods to extract the elliptic flow on the event-by-event basis [153, 154]). The former approach can be carried out by using the prolate shape of the uranium nuclei [155]. In central U + U collisions, there would be sizable v₂ due to the prolate shape of the colliding nuclei while the magnetic field would be very tiny [12]. The STAR Collaboration collected 0-1% most central events from U + U collisions at $\sqrt{s}=197$ GeV in 2012, and indeed found sizable v₂ while Δγ is consistent with zero [146]. However, it was found later that the total multiplicity is far less correlated to the number of binary collisions than expected [156], so that the shape selection in experiments is very hard; see Refs. [118, 157] for recent discussions. The latter approach was proposed to be achieved by using the collisions of the isobaric nuclei, such as ${}_{44}^{96}$Ru and ${}_{40}^{96}$Zr [155, 158, 159]. It is expected that Ru + Ru and Zr + Zr collisions at the same beam energy and same centrality will produce similar elliptic flow but a 10% difference in the magnetic fields. Thus if the observed γ correlation contains contribution from CME, the isobaric collisions could have a chance to see it by comparing the same γ correlation in Ru + Ru and Zr + Zr collisions. The detailed numerical simulation was reported in Ref [158] where the authors demonstrated that the two collision types at ps $\sqrt{s}=200$ GeV have more than 10% difference in the CME signal and less than 2% difference in the elliptic-flow-driven backgrounds for the centrality range of 20-60% by assuming that the CME contribution to the γ correlation is 1/3. Such a difference is feasible in current experimental setups of heavy-ion collisions at RHIC. Therefore, the isobaric collisions at RHIC would be very valuable in isolating the chiral magnetic effect from the background sources.

Recently, the CMS Collaboration at LHC reported the measurement of γαβ in high multiplicity events in p + Pb and Pb + Pb collisions at $\sqrt{s}=5.02$ TeV [160]. They found that the γαβ’s in p + Pb and in Pb + Pb as functions of multiplicity behave similarly at the same multiplicities. These results give useful information about the CME search in heavy-ion collisions. First, although the magnetic field in p + Pb can be very strong, the magnetic-field direction should not be firmly correlated to the v₂ plane, and the lifetime of the magnetic field should be shorter at higher collision energy, thus the magnetic field is not expected to drive strong γαβ in p + Pb. The CMS measurement thus gives at least two very different possibilities. One is that, although it was not expected before, the CME still occurs in p + Pb collisions and the magnetic field, although weakly correlated to the v₂ plane, is still able to drive a finite signal which combined with the background effects gives the observed results. The other one is that there is no CME at $\sqrt{s}=5.02$ TeV (or very weak CME) in both p + Pb and Pb + Pb collisions, the γαβ in both p + Pb and Pb + Pb collisions are due to background effects. This is actually somehow compatible with the recent measurement of STAR collaboration in which the so-called H correlation in measured at various energies through the Beam Energy Scan (BES) program [147]. In the H correlation, the background effects are subtracted under reasonable assumptions; see details in Ref. [147]. If one extrapolates the H correlation of STAR Collaboration to $\sqrt{s}=5.02$ TeV, one can find that H correlation is actually small there. This supports that at higher energy the γαβ is more dominated by background contributions. Nevertheless, more measurements are definitely needed. Especially the CMS collaboration should also measure the so-called δ-correlation from which we can deduce the H correlation. In addition, the charge-dependent v₂ measurements (see next section) in p + Pb collisions will be also very important.

3.3 Chiral separation effect and chiral magnetic wave

By looking at the CME (10), one may wonder whether there is a quantum-anomaly induced axial current in response to the external magnetic field and under what condition such current can appear. The answer is yes, and such a current can be caused by the chiral separation effect (CSE) which is the dual effect to the CME [18, 19, 161]:

where the axial current is defined by $J_A^{\mu }=e\langle \overline{\psi }{\gamma }^{\mu }{\gamma }_5\psi \rangle$, and μV is the vector chemical potential. Note that, unlike the CME, the occurrence of the CSE does not require a parity-violating environment but a charge-conjugation odd environment.

We can understand the CSE from an intuitive picture. In parallel to the discussion of CME in Sec.3.1, let us consider a massless fermionic system under a strong magnetic field, and let us focus on the lowest Landau level on which the spin of positively charged fermions are polarized to be parallel to the magnetic field (we call its direction the "up" direction and its opposite direction the "down" direction) while the spin of negatively charged fermions are polarized to be anti-parallel to the magnetic field. If a positively charged fermion on the lowest Landau level moves upward then it contributes to the RH current while if it moves downward it contributes to the LH current; thus no matter which direction it moves it contributes to an axial current along the magnetic field. Similarly, it is easy to see that the negatively charged fermions contributes an axial current in the direction opposite to the magnetic field. Now, if the system contains more positively charged fermions than negatively charged fermions, namely, if μV>0, the downward-moving axial current cannot compensate the upward-moving axial current and thus we get a net axial current along the direction of the magnetic field. This is the CSE current. The fermions occupying the higher Landau levels do not contribute to the CSE current as their spins are not polarized.

Let us formulate the above argument in a more precise way. For simplicity, we consider a system of massless Dirac fermions with only positive charge, e, in an constant magnetic field, B, in z direction which is supposed to be so strong that all the fermions are confined in the lowest Landau level. The number densities of RH and LH fermions are given by Eq. (12). As the massless fermions on the lowest Landau level move at the speed of light in a way correlated to their spins which means that the RH fermions generate a current along z direction, JR=enR, while LH fermions generate a current along -z direction, JL=-enL. Therefore the net axial current reads

which, once we equate $p_F^{R/L}$ with μR/L and μV with (μR+μL)/2, leads us to Eq. (18) and Eq. (19). More rigorous derivation of CSE and various aspects regarding CSE can be found in Refs. [18, 19, 161-170].

The CME and CSE couple together the the vector and axial densities and currents in the presence of the magnetic field and lead to collective modes called chiral magnetic waves (CMWs). This is most easily seen when we express CME and CSE in the chiral basis:

where the RH and LH quantities are defined as JR/L=JV±JA)/2 and μR/L=μV±μA. Substituting JR and JL in to the continuation equations and considering a small fluctuation in the density $J_{R/L}^0$, we obtain

where $\chi =\partial J^0/\partial \mu$ is the number susceptibility (we assume that the RH and LH chiralities have the same susceptibility). Therefore, we can identify two collective wave modes, one propagating along B and one propagating opposite to the direction of B with the same velocity vχ=e²B/(4π²χ). They are the CMWs [20] (see also Ref. [171] for a derivation of CMW based on the kinetic theory).

3.4 Experimental signals for chiral magnetic wave

The CMWs are able to transport both chirality and electric charge and can induce an electric quadrupole in the QGP [21, 22]; see Fig. 9 for an illustration. Owing to strong in-plane pressure gradient, such an electric quadrupole can cause a splitting between the v₂ of π⁺ and π^- [22]; see also Refs. [139, 172-174] for more advanced numerical implementations. The CMW-induced charge dependent elliptic flows can be parameterized as

Figure 9:Full-screen View

The illustration of the CMW-induced electric quadrupole in QGP and the charge-dependent elliptic flow of pions. Figure is from Ref. [93].

where v₂ is the average elliptic flow and A_ch=(N₊-N_-)/(N₊+N_-) is the net charge asymmetry parameter which is usually small (the events that can generate large A_ch are rare). The slop parameter, r, is predicted to be proportional to <(eB)²cos[2(ψB-Ψ_RP)]> where ψB is the azimuthal direction of the magnetic field [11, 22, 170].

Recently, the CMW has been tested by STAR Collaboration [146, 175, 176] at RHIC and by ALICE Collaboration [177] at LHC. The data shows an elliptic-flow difference, v₂(π^-)-v₂(π⁺), linear in A_ch with a positive slope whose centrality dependence can be fitted by the CMW estimation. Furthermore, the slope parameter, r, displays no obvious trend of the beam energy dependence for 10-60% centrality at $\sqrt{s}=20-200$ GeV. This can be understood by noticing that the strength of magnetic field times is roughly beam energy independent (Fig. 10).

Figure 10:Full-screen View

(Left) The v₂ splitting between π^- and π⁺ as a function of charge asymmetry parameter, A_±, measured by STAR Collaboration. (Right) The measured slop parameter, r, as well as the UrQMD and CMW predictions are functions of the centrality. Figures are from Ref. [176].

However, ambiguities exist in the interpretation of the experimental results as there are background effects which can also lead to v₂ splitting between π⁺ and π^-. A variety of possible background effects were proposed, see Refs. [10, 178-183] and review articles, Refs. [93, 126], for discussions. It should be noted that none of the above backgrounds can successfully explain all the features of the experimental data. It is quite plausible that various effects (including CMW) are coupled to produce the measured pattern. Until we can successfully subtract the background effects, we cannot make a conclusive claim. On the experimental side, it will be important to perform more detailed measurements, such as the measurement of the v₂ splitting versus transverse momentum and rapidity, the measurement for other hadrons like kaons, the measurement of splitting in other harmonic flows like v₃, and the measurement of the cross-correlation of different observables. On the theoretical side, an anomalous hydrodynamic simulation or a simulation based on kinetic theory for the CMW that incorporate both the dynamically evolving magnetic fields and the realistic initial charge distribution will be very useful.

Before we end this section, we note that, in addition to the magnetic-field induced anomalous transports, the electric field and vorticity can also induce certain anomalous transport phenomena and corresponding collective modes. For example, electric field can induce an axial current in the medium with both μA and μV nonzero which is called chiral electric separation effect (CESE) and, once coupled with the normal Ohm’s current (the vector current induced by electric field), CESE leads to chiral electric wave and axial and vecor density waves [184-187]. The CESE may lead to observable effect in, Cu + Au collisions where persistant electric fields from Au to Cu nuclei exist [188]. The fluid vorticity can induce the so-called chiral vortical effects (CVEs), which are the vortical analogues of CME [189-192], and the corresponding collective mode, chiral vortical wave (CVW) [193]. The CVE and CVW may have important consequences in noncentral heavy-ion collisions where nonzero fluid vorticity can be generated as a result of angular momentum conservation [194-197]. Phenomenologically, the CVE can induce baryon number separation with respect to the reaction plane [198], and CVW may cause a v₂ splitting between Λ and $\overline{\Lambda }$ [193]. Recently, the STAR Collaboration has reported the first measurement of CVE in Au + Au collisions at $\sqrt{s}=200$ GeV [199], and the data is consistent with the expectation of CVE. But, similar with the case of CME and CMW detections, the data contains contributions from background effects like the transverse momentum conservation and the local baryon number conservation, and we need more efforts to understand the data.

4.1 Quarkonia in magnetic fields

4.1.1 Spin mixing and level repulsion

We shall begin with the simplest description of effects of magnetic fields on quarkonium spectra. In a weak magnetic field, the dominant effects of the magnetic fields can be described by the hadronic degrees of freedom. This contrasts the case of strong magnetic fields where modifications of nonperturbative vacuum structure and of the internal structure of hadrons are important (see Refs. [58, 52, 62, 68-70, 203-207] and references therein for the ongoing discussions). As long as the magnitude of the magnetic field is much smaller than the heavy quark mass, these effects will be subdominant contributions.

In the hadronic degrees of freedom, effects of magnetic fields appear in the mixing of spin eigenstates [25, 27, 29, 30, 33, 38, 39]. This mixing is caused by the breaking of the spatial rotation symmetry in magnetic fields: there remains only the azimuthal rotational symmetry with respect to the direction of the external magnetic field. This indicates that only the spin state along the magnetic field can persist as a good quantum number of the mesons, and thus that there is a mixing between the spin single and triplet states, S_total, Sz) = (0,0) and (1,0), which are, for example, ηc and the longitudinal component of J/ψ.

We examine general mixing patterns among the pseudoscalar, P, vector, Vμ, scalar, S, and axial-vector, Aμ, quarkonia. Possible interaction vertices among those fields are informed from the Lorentz invariance and the parity and charge-conjugation symmetries. The vertices relevant for interactions among static charmonia are found to be [29, 30]

${\mathcal{L}}_{{\gamma }_{\text{P}V}}=\frac{g_{{}_{\text{P}V}}}{m_0}e{\tilde{F}}_{\mu \nu }({\partial }^{\mu }P)V^{\nu }\mathrm{,}$ (25) ${\mathcal{L}}_{{\gamma }_{\text{V}A}}=i{g}_{{}_{\text{V}A}}e{\tilde{F}}_{\mu \nu }{V}^{\mu }{A}^{\nu }\mathrm{,}$ (26) ${\mathcal{L}}_{{\gamma }_{\text{S}A}}=\frac{g_{{}_{\text{S}A}}}{m_1}e{\tilde{F}}_{\mu \nu }({\partial }^{\mu }S)A^{\nu }\mathrm{,}$ (27)

with m₀=(m_P+m_V)/2, m₁=(mS+m_A)/2, and dimensionless effective coupling constants g_PV, g_VA, and g_SA. The vertices proportional to $F^{\mu \nu }{\partial }_{\mu }$ do not play any role for static charmonia, because they pick up the vanishing spatial derivatives and the temporal component of the vector or axial-vector fields. The nonvanishing vertices are responsible for, radiative decay modes of quarkonia such as J/ψ→ηc + γ [208, 209]. It is also worth mentioning that the three-point vertex (25) induces the conversion of ηc to J/ψ in a magnetic field. The contribution to the J/ψ production from this conversion process was examined in Ref. [24].

From the above vertices, one can identify two pairs of mixing partners: the mixing between P and V and between S and A are possible when they are at rest in external magnetic fields. Since the P and V are not mixed with the S and A, we shall focus on the mixing between S and V. As mentioned in the beginning of this section, a mixing effect arises only in the longitudinal component of the vector field, V_||, and the other two transverse components persist as the energy eigenstates. This can be easily confirmed by identifying the non-vanishing components in Eq. (25). With a vanishing spatial momentum q^μ=(ω, 0, 0, 0) and the vertex (25), the coupled equations of motion are found to be

$\left(\begin{array}{ll}-{\omega }^2+m_{\text{P}}^2\hfill & -i\frac{g_{{}_{\text{P}V}}}{m_0}\omega eB\hfill \\ i\frac{g_{{}_{\text{P}V}}}{m_0}\omega eB\hfill & -{\omega }^2+m_{\text{V}}^2\hfill \end{array}\right)\left(\begin{array}{l}P\hfill \\ V_{\parallel }\hfill \end{array}\right)=0.$ (28)

Following from the equations of motion (28), we obtain the physical mass eigenvalues as

${\overline{m}}_{V\mathrm{,}P}^2=\frac{1}{2}(M_{+}^2+\frac{{\gamma }^2}{m_0^2}\pm \sqrt{M_{-}^4+\frac{2{\gamma }^2{M}_{+}^2}{m_0^2}+\frac{{\gamma }^4}{m_0^4}})\mathrm{,}$ (29)

where $M_{+}^2=m_{\text{P}}^2+m_{\text{V}}^2\mathrm{,}{M}_{-}^2=m_{\text{V}}^2-m_{\text{P}}^2$, and γ=g_PVeB. The upper and lower signs are for the vector and pseudoscalar states, respectively. Expanding the right-hand side of Eq. (29) up to the second order in γ and the leading order in $\frac{1}{2}(m_{\text{V}}-m_{\text{P}})/m_0$, we find

${\overline{m}}_{P\mathrm{,}V}^2(B)=m_{V\mathrm{,}P}^2\pm \frac{{\gamma }^2}{M_{-}^2}\mathrm{,}$ (30)

with eigenvectors given by

$\begin{array}{l}|P{\rangle }_B=(1-\frac{1}{2}\frac{{\gamma }^2}{M_{-}^4})|P\rangle -i\frac{\gamma }{M_{-}^2}|V\rangle \mathrm{,}\\ |V{\rangle }_B=-i\frac{\gamma }{M_{-}^2}|P\rangle +(1-\frac{1}{2}\frac{{\gamma }^2}{M_{-}^4})|V\rangle \mathrm{.}\end{array}$ (31)

The mass eigenvalue (30) clearly shows the level repulsion in the presence of the mixing effect.

In case of mixing between ηc and J/ψ, the coupling constant, g_PV, can be read off from the measured radiative decay width, J/ψ→ηc + γ [29, 30]. From the vertex (25), the invariant amplitude of this decay process is written as

$\begin{array}{l}{\mathcal{M}}_{{}_{\text{P}V}}=\langle \gamma \text{P}\left|{\mathcal{L}}_{\gamma \text{P}V}\right|\text{V}\rangle \\ =-\frac{e{g}_{{}_{\text{P}V}}}{m_0}{ϵ}_{\mu \nu \alpha \beta }{k}_{\gamma }^{\mu }{ϵ}_{\gamma }^{\nu }{p}_{\text{V}}^{\alpha }{ϵ}_{\text{V}}^{\beta }\mathrm{,}\end{array}$ (32)

where $k_{\gamma }^{\mu }$ and ${ϵ}_{\gamma }^{\nu }$ ($p_{\text{V}}^{\alpha }$ and ${ϵ}_{\text{V}}^{\beta }$) denote the polarization and momentum of the photon (J/ψ). Summing the polarizations of the photon and averaging those of the vector meson, we find

$\Gamma [\text{V}\to \gamma \text{P}]=\frac{1}{12}\frac{e^2{g}_{{}_{\text{P}V}}^2{\tilde{p}}^3}{\pi m_0^2}\mathrm{.}$ (33)

We have also integrated over the phase-space volume in the two-body final state. The center-of-mass momentum in the final state is given by $\tilde{p}=(m_{\text{V}}^2-m_{\text{P}}^2)/(2m_{\text{V}})$. Using the measured radiative decay width, Γ_exp[J/ψ→γηc]=1.579 keV, the coupling strength is estimated to be g_PV=2.095. A consistent value of the coupling strength is also obtained from computation in the quark degrees of freedom with help of the heavy quark expansion and diagrammatic methods in an external magnetic field (see Appendix. 6).

Figure 11 shows the mass shifts in the presence of the mixing effect. We find the level repulsion between the ηc and longitudinal J/ψ mass spectra. While the mass of ηc decreases as eB increases, the mass of the longitudinal J/ψ (denoted by J/ψ_∥) increases. This result clearly provides the simple picture of the effect of weak magnetic fields.

Figure 11:Full-screen View

Level repulsion in the mass spectra of the ηc and longitudinal J/ψ. (Reproduced from Refs. [30])

The level repulsion occurs most strongly between a pair of adjacent states. Therefore, the other excited states give rise to only subdominant contributions to the mass shifts of ηc and the longitudinal J/ψ in a weak magnetic field. However, as the mass of the longitudinal J/ψ keeps increasing with an increasing magnetic field, it will approaches the other excited states. In such a strong magnetic field, many states will get involved in the mixing effects, which will result in an intriguing mixing patterns of the spectra as discussed in the next section.

Moreover, the mixing effects induced by the minimal coupling (25) to a magnetic field will not be valid in the strong-field regime. Since the strong magnetic field probes the inner structure of hadrons and modifies the QCD vacuum, we will need to look into the dynamics of the fundamental degrees of freedom. To investigate these effects, one can invoke on the Cornell potential models and QCD sum rules with the help of lattice QCD simulations. We will proceed to exploring these frameworks on top of the simplest framework for the mixing effect discussed in this section.

4.1.2 Cornell potential model in magnetic fields

To study the nonrelativistic bound states of the heavy quark and antiquark on the basis of the fundamental degrees of freedom, one can use the Cornell potential model [210, 211] with the help of the lattice QCD simulations. The basic framework for the 1S states of ηc, J/ψ, ηb, and Υ, was first elaborated by Alford and Strickland [27]. The authors included the interaction with the background magnetic field through the magnetic moments of constituents and the kinetic terms, as well as the mutual interactions by the linear and Coulomb potentials in vacuum. In this framework, the interactions through the magnetic moments give rise to the mixing between the spin singlet and triplet states discussed in the previous section. Based on this model, a systematic spectroscopy for excited states was recently performed by Suzuki and Yoshida owing to an improvement of computational method [38, 39]. Effects of strong magnetic fields on nonperturbative QCD vacuum were investigated by Bonati et al. who computed the modification of the quark-antiquark potential by lattice QCD simulation [32, 34] and plugged it into the bound-state problem [33]. These studies suggest an opportunity to explore the nature of QCD confinement in the presence of a strong magnetic field. In the following, we briefly summarize the framework for the potential model in a magnetic field and the numerical results for the mass spectra.

The Hamiltonian for the two-body bound state in a magnetic field is given by [25, 27, 38, 39, 33]

$\mathcal{H}=2m_Q+{\displaystyle \sum _{i=\mathrm{1,2}}}{\mathcal{H}}_0^{(i)}+{\mathcal{H}}_{\overline{q}q}+{\mathcal{H}}_{\mu }$ (34)

where

${\mathcal{H}}_0^{(i)}=\frac{1}{2m_Q}{\left(p^{(i)}-q^{(i)}A(r^{(i)})\right)}^2$ (35) ${\mathcal{H}}_{\overline{q}q}={\sigma }_{s}r-\frac{a}{r}+b(r)\left({\sigma }^{(1)}\cdot {\sigma }^{(2)}\right)$ (36) ${\mathcal{H}}_{\mu }=-\left({\mu }^{(1)}+{\mu }^{(2)}\right)\cdot {\rm B}$ (37)

with the gauge potential, A(r), for the external magnetic field, B, the Pauli matrix, σ, and the heavy-quark mass, mQ. The indices i=1,2 denote the quark and antiquark, respectively, and thus q⁽¹⁾=-q⁽²⁾(≡q) for quarkonia. The mutual interaction is described by the Cornell potential, ${\mathcal{H}}_{\overline{q}q}$, which includes the linear and Coulomb potentials as well as the spin-spin interaction. For the nonrelativistic quark and antiquark, the magnetic moments read

${\mathcal{H}}_{\mu }=-{\mu }_Q\frac{g}{2}({\sigma }_z^{(1)}-{\sigma }_z^{(2)})B$ (38)

where μQ=q/(2mQ) and $g\simeq 2$. We took the magnetic field in +z direction without loss of generality. The eigenstates of the Pauli matrices are spanned by the direct product $|\pm \pm \rangle \equiv |\pm \rangle \otimes |\pm \rangle$ where the left and right states in the product stand for the eigenstates of ${\sigma }_z^{(1)}$ and ${\sigma }_z^{(2)}$, respectively. In the total-spin basis, the spin singlet state is given by $|\mathrm{0,0}\rangle =(|+-\rangle -|-+\rangle )/\sqrt{2}$, while the triplet states by is |1,±1〉=|±±〉 and $|\mathrm{1,0}\rangle =(|+-\rangle +|-+\rangle )/\sqrt{2}$. Therefore, the operation of ${\mathcal{H}}_{\mu }$ on these states results in

${\mathcal{H}}_{\mu }|\mathrm{1,}+1\rangle ={\mathcal{H}}_{\mu }|\mathrm{1,}-1\rangle =0\mathrm{,}$ (39) ${\mathcal{H}}_{\mu }\left(\text{}\text{}\begin{array}{ll}|\mathrm{0,0}\rangle \hfill & \hfill \\ |\mathrm{1,0}\rangle \hfill & \hfill \end{array}\text{}\text{}\right)=\left(\begin{array}{ll}0\hfill & 2{\mu }_{Q}B\hfill \\ 2{\mu }_{Q}B\hfill & 0\hfill \end{array}\right)\text{}\text{}\left(\text{}\text{}\begin{array}{ll}|\mathrm{0,0}\rangle \hfill & \hfill \\ |\mathrm{1,0}\rangle \hfill & \hfill \end{array}\text{}\text{}\right)\mathrm{.}$ (40)

Consistent with the discussion in the previous section, we find that, while the transverse modes |1,±1〉 are eigenstates of the two-body Hamiltonian, there is a mixing between |0, 0〉 and |1, 0〉 states which need to be diagonalized to the energy eigenstates.

A specific issue of the dynamics in magnetic fields is the non-conservation of the transverse momentum. The gauge potential for an external magnetic field will depend on the transverse coordinate as $A_{\text{e}xt}^{\mu }=(\mathrm{0,}{A}_{\perp }^i({\text{x}}_{\perp })\mathrm{,0})$, and apparently breaks the translational invariance in the transverse plane. Whereas the momentum along the magnetic field is conserved, the transverse momentum is not a conserved quantity and is moreover gauge dependent. For example, one of the transverse components is not conserved in the Landau gauge, while neither of them are conserved in the symmetric gauge. Nevertheless, it is clear that a constant magnetic field has a translational invariance in the transverse plane so that there must be a conserved momentum which works as the generator of the translation. This momentum is called the pseudomomentum, and its generalization to the two-body system is defined as (see Refs. [27, 204] and references therein)

$\text{K}={\displaystyle \sum _{i=\mathrm{1,2}}}\left(p^{(i)}+q^{(i)}{A}^{(i)}\right)\mathrm{,}$ (41)

where the labels denote the quark and antiquark as before. This momentum is related to the center-of-mass kinetic momentum as 〈P_kinetic〉=K-qB×〈r〉 with q≡q⁽¹⁾=-q⁽²⁾ and the relative coordinate r=r₁-r₂. The same issue of the momentum nonconservation is seen in diagrammatic calculations in the presence of a magnetic field. In such calculations, the gauge- and translational-symmetry breaking parts are encoded in the so-called Schwinger phase appearing from the gauge-dependent wave function of fermions in a magnetic field [212]. Effects of the Schwinger phase on the Bethe-Salpeter equation for light mesons were discussed in Ref. [62].

The Schrödinger equation for the two-body system (34) can be solved numerically for a fixed value of the conserved momentum, K. The parameters for the Cornell potential are so tuned to reproduce the vacuum spectrum in Ref. [27, 38, 39], and are imported from the lattice simulation in a finite magnetic field in Ref. [33]. In those papers, the numerical codes were carefully tested to reproduce the analytic solution for a harmonic potential.

Figure 12 shows the charmonium spectra from Ref. [27] with the potential in vacuum. At the vanishing momentum 〈P_kinetic〉=0 GeV, one sees that the masses of the ηc and longitudinal J/ψ decreases and increases, respectively, as expected from the mixing effect. On the other hand, the modification in the transverse J/ψ is weak without the mixing effect. The masses of the moving charmonia, which were obtained by subtracting the center-of-mass kinetic energy, tend to increase with an increasing in its momentum. The bottomoium spectra shown in Ref. [27] have qualitatively the same behavior as in charmonium spectra, with suppression in magnitude due to the smaller magnetic moments of the bottom quark and antiquark.

Figure 12:Full-screen View

Ground-state charmonium spectra from the Cornell potential model. (Reproduced from Ref. [27].)

This study was extended to the excited states of quarkonia, which was made possible by using the cylindrical Gaussian expansion method proposed in Ref. [38]. Since this method does not resort to the imaginary time evolution τ → ∞ to extract the spectra [27], it would be more suitable for investigating the excited states. The result of the systematic spectroscopy is shown in Fig. 13. In a weak magnetic field, the level repulsions between the nearest mixing partners give rise to the dominant modification of the mass spectra. As the magnetic field becomes stronger, there appear the level-crossing points where the wave functions are strongly mixed. The density plots in the small windows show the transition of the wave functions. One can see the squeezing of the wave function along the magnetic field in the strong-field regime as expected from the shrink of cyclotron orbits. Including the contributions of the excited states, the ηc state is, as it is the lowest-lying state, pushed down additively by all of the mixing partners, and its mass monotonically decreases. This systematic computation was very recently extended to the spectra of bottomonia and open heavy flavors [39].

Figure 13:Full-screen View

Charmonium spectroscopy by the Cornell potential model. Small windows exhibit the spatial extensions of the wave functions in the transverse and longitudinal directions with respect to the B which are shown as the horizontal and vertical axes, respectively. (Reproduced from Ref. [38].)

So far, we have not considered the modification of the Cornell potential by the magnetic field. Effects of a magnetic field on the potential was investigated in Ref. [32] by lattice QCD simulation. The authors extracted the strengths of the potential in two cases where the separation between the quark and antiquark is in parallel and transverse to the magnetic field. The left panel in Fig. 14 shows the results from the lattice simulation. Compared to the potential in vacuum, the strength of the potential is enhanced and suppressed in the transverse (XY) and longitudinal (Z) directions, respectively. This behavior may be understood from the fact that the screening due to the vacuum polarization of light-quark pairs is coupled to only the temporal and parallel (z) components of the gluon field. The full azimuthal-angle and temperature dependences were recently studied by the lattice QCD simulation [34], and those effects were also discussed by the strong couping method [35, 36, 213].

Figure 14:Full-screen View

The $q\overline{q}$ potential from the lattice QCD simulation as a function of the $q\overline{q}$ separation (left) and of the azimuthal angle (contour plot on right). The separation is specified by the lattice coordinate ns and spacing a, and the magnetic field is applied in the z-direction. The potential without a magnetic field (b=0) is compared with the anisotropic potentials in the transverse (XY) and the longitudinal (Z) directions, as well as the average of them (XYZ), at a finite magnetic field with b=24 that is the discretized magnetic flux on the periodic lattice. (Reproduced from Refs. [32, 33].)

The right panel in Fig. 14 shows the contour plot of the potential which was obtained by interpolating the values in the transverse (XY) and longitudinal (Z) directions with respect to the azimuthal-angle dependence. In Ref. [33], the quarkonium spectra were investigated by including this anisotropic potential as well as the mixing effects. Figure 15 shows the bottomonium spectra with and without the inclusion of the anisotropic potential. The authors argued that effects of the anisotropic potential can be seen in the bottomonium spectra more strongly than in the charmonium spectra, because the magnetic moment of the bottom quark is suppressed by the larger value of the quark mass and the smaller electric charge. Without the anisotropic potential, the spectra barely change due to the suppression of the mixing effects while the effects of the anisotropic potential act to increase the bottomonium masses in all three of the channels. The effect of the anisotropic potential is also much stronger than the momentum dependence of the spectra.

Figure 15:Full-screen View

Bottomonium spectra with and without including the anisotropic potential from the lattice QCD. (Reproduced from Ref. [33].)

On the basis of the potential model, several phenomenological consequences have also been discussed. Prior to moving to the next section, we briefly summarize these studies.

• In the presence of the mixing effect, the wave function of the mass eigenstate is the mixture of the two components so that the mass eigenstate shares the decay modes with both of the mixing partners. This means, for example, that the physical ηc and ηb are allowed to decay into the dilepton while the dilepton decay channels of the physical longitudinal J/ψ and Υ are suppressed [27]. As mentioned in Ref. [27], a similar modification of the positronium decay mode is known.

• The Lorentz forces acting on a moving quark and antiquark pair are oriented in the opposite direction. Therefore, the force will act to break up the bound state implying a possible dissociation mechanism in a strong magnetic field [23, 27]. This effect may be also understood as an effect of an electric field appearing in the rest frame of a moving quarkonium, which is connected to the original magnetic field by the Lorentz boost. The electric field will stretch the separation between the quark and antiquark pair and induces a tunneling through the binding potential [23].

• Mass shifts of D mesons were also examined with use of the potential model [25]. In turn, these mass shifts induce the shift of the $D\overline{D}$ threshold above the quarkonium states. Based on the color evaporation model, modification of the charmonium production rate was discussed [25].

• A time dependence of the magnetic field was introduced in the potential model [37]. Based on the time-dependent model, the authors investigated the anisotropic production rate of high- pT charmonia. The feed-down effect from the excited states was also examined.

4.1.3 QCD sum rules in magnetic fields

Effects of the magnetic field on the meson spectra were also studied by using QCD sum rule (QCDSR) for open heavy flavor mesons [26, 31] and quarkonia [29, 30]. While we will focus here on the QCD sum rule for charmonia, extension to bottomonia is straightforward and extension to other systems will be discussed in the last of this section.

In the QCD sum rules [214-217], the spectra of bound states are extracted from the time-ordered current correlator defined by

${\Pi }^J(q)=i{\displaystyle \int }\text{}\text{}{d}^{4}x{e}^{iqx}\langle 0|T[J(x)J(0)]|0\rangle \mathrm{,}$ (42)

where the superscript J specifies a channel. The dispersion relation relates the hadronic spectral function in the time-like region (q²=-Q²≥0) to the correlator in the deep Euclidean region (Q² → ∞) as

${\Pi }^J(Q^2)=\frac{1}{\pi }{\displaystyle {\int }_0^{\infty }}\frac{\text{I}m{\Pi }^J(s)}{s+Q^2}ds+(\text{subtractions})\mathrm{.}$ (43)

One may evaluate the correlator on the left-hand side on the basis of the Operator Product Expansion (OPE) in the deep Euclidean region. Here, we shall focus on the charmonia created by the currents:

$J_P=i\overline{c}{\gamma }^{5}c\mathrm{,}{J}_V^{\mu }=\overline{c}{\gamma }^{\mu }c\mathrm{.}$ (44)

Up to the dimension-4 operators, the OPE can be written down as

$\Pi =C_0+C_{G^2}\langle G_{\mu \nu }^a{G}^{a\mu \nu }\rangle +C_{B^2}^{\mu \nu \alpha \beta }\langle F_{\mu \nu }^{\text{e}xt}{F}_{\alpha \beta }^{\text{e}xt}\rangle \mathrm{.}$ (45)

The Wilson coefficients, C₀, C_G², and C_B², are responsible for the short-distance processes, so that one may perturbatively calculate the diagrams for the hard quark loops with insertions of the soft external gluon/magnetic fields (see Refs. [29, 30] for explicit expressions). This can be done by using standard methods in the Fock-Schwinger gauge, which is also as known as the fixed-point gauge [216, 218].

On the other side of the dispersion relation (43), we assume an Ansatz for the functional form of the spectral function ρ(s)=π^-1Im∏J(s) to extract the lowest-lying pole. In conventional analyses, the lowest-lying pole is simply separated from the continuum as

${\rho }_{\text{v}ac}(s)=f_0\delta (s-m^2)+\frac{1}{\pi }\text{I}m{\Pi }_{\text{pert}}^J(s)\theta (s-s_{\text{th}})\mathrm{.}$ (46)

The parameters s_th and f₀ are the threshold of the continuum and the coupling strength between the current and the lowest-lying pole, respectively.

We, however, need to carefully reexamine the Ansatz (46) in the presence of a magnetic field, because effects of a magnetic field should be consistently taken into account on both sides of the dispersion relation (43) [26, 29-31]. The effect in the hadronic degrees of freedom is nothing but the mixing effect discussed in the previous sections. In the language of the spectral function, the mixing effect can be seen as the emergence of a new pole for the mixing partner, which is the longitudinal J/ψ in the pseudoscalar channel and the ηc in the vector channel. A schematic picture in Fig. 16 shows the emergence of the longitudinal J/ψ in the pseudoscalar channel, which induces the repulsive shifts of the poles and the reduction of the residue of ηc. It is important to include these effects in order to properly maintain the normalization of the spectral function.

Figure 16:Full-screen View

Sketch of the spectral function in the pseudoscalar channel. The longitudinal J/ψ_|| pole emerges which induces the repulsive shifts of the poles and the reduction of the residue of the ηc pole. Similar events happen in the longitudinal vector channel.

In the pseudoscalar channel, the second-order term in eB can be obtained by first converting the current to the psedoscalar meson with the strength f₀ and then using the results for the mixing effect shown in Eqs. (30) and (31). Therefore, the correlator can be written in the hadronic degrees of freedom as [29]

${\Pi }_{\text{2}nd}^P(q^2)=f_0\left[\frac{|{\langle P|{\eta }_c\rangle }_B{|}^2}{q^2-m_{{\eta }_c}^2}+\frac{|{\langle P|J/\psi \rangle }_B{|}^2}{q^2-m_{J/\psi }^2}\right]\mathrm{.}$ (47)

Inserting Eqs. (30) and (31) into Eq. (47) and expanding up to the second order in eB, we find

${\Pi }_{\text{2}nd}^P(q^2)=f_0\frac{{\gamma }^2}{M_{-}^4}\left[\frac{1}{q^2-m_{\text{V}}^2}-\frac{1}{q^2-m_{\text{P}}^2}-\frac{M_{-}^2}{{(q^2-m_{\text{P}}^2)}^2}\right]\mathrm{.}$ (48)

Interpretation of the terms in Eq. (48) are as follows. The first term corresponds to production of the on-shell longitudinal J/ψ from the pseudoscalar current via an off-shell ηc. The second term with a negative sign is needed to preserve the normalization because the coupling of ηc to the current must be reduced to balance the emergence of the coupling to J/ψ. This is confirmed in Eq. (47), where these two terms come from overlaps between the properly normalized unperturbated and perturbated states obtained as $|(P|\eta c{)}_B{|}^2\sim 1-{(\gamma /M_{-}^2)}^2$ and $|(P|J/\psi {)}_B{|}^2\sim {(\gamma /M_{-}^2)}^2$. The third term has a double pole on the ηc mass with a factor $M_{-}^2$ corresponding to a virtual transition to the longitudinal J/ψ state between on-shell ηc states, which is nothing but the origin of the mass shift due to the mixing effect. In Eq. (47), this term comes from an expansion with respect to the mass correction shown in Eq. (30). Clearly, if one includes this mixing term in the phenomenological spectral function, its effect is subtracted out from the total mass shift obtained from the QCDSR and thus can be separated from the residual effects of B-fields, not described in the hadronic level. A similar calculation can be applied to the longitudinal vector channel. On the phenomenological side, those terms should be included as well as the vacuum contributions (46) in both the pseudoscalar and longitudinal vector channels.

To make the dispersion integral (43) dominated by the pole contribution and to improve the convergence of the OPE, we use the Borel transform defined by

$\mathcal{B}[f(Q^2)]\equiv \text{}\text{}\underset{\begin{array}{c}{Q}^2\mathrm{,}n\to \infty \\ Q^2/n=M^2\end{array}}{{\displaystyle \lim }}\frac{{(Q^2)}^{n+1}}{n!}{\left(-\frac{d}{d{Q}^2}\right)}^{n}f(Q^2)\mathrm{.}$ (49)

By applying the Borel transform to both sides of Eq. (43), the contribution of the delta function for the lowest-lying pole is picked up to result in a mass formula:

$m^2=-\frac{\partial }{\partial (1/M^2)}\ln [{\mathcal{M}}_{\text{OPE}}-{\mathcal{M}}_{\text{cont}}-{\mathcal{M}}_{\text{mix}}]\mathrm{,}$ (50)

where ${\mathcal{M}}_{\text{OPE}}$, ${\mathcal{M}}_{\text{cont}}$, and ${\mathcal{M}}_{\text{mix}}$ are the Borel transform of the OPE (45), the perturbative continuum, and those induced by the magnetic field (48), respectively. More details and explicit forms are given in Refs. [29, 30].

Figure 17 shows the mass spectra from the QCDSR together with those from the hadronic effective model shown in Fig. 11. First, the results from these two methods are consistent with each other, indicating the dominance of the mixing effect in the shown regime. They are also consistent with the results from the potential model discussed in the previous section [27, 38, 39, 33]. Also, the modification in the transverse components of J/ψ (green line) is much smaller than the others due to the absence of the mixing partner. On top of these observations, one finds deviations between the results from the two methods. Although they are small deviations in the current regime of the magnetic field, they suggest the precursor of other effects, which are not explained by the mixing effect. It would be natural to expect the emergence of such effects in stronger magnetic fields. Possible origin of these deviations were discussed in Ref. [30].

Figure 17:Full-screen View

The ηc and J/ψ mass spectra from the QCD sum rule together with the results from the hadronic effective model shown in Fig. 11. (Reproduced from Ref. [29]).

One caveat is that the modification of the gluon condensate in a magnetic field has not been included in the evaluation of the OPE (45). Since the lattice QCD simulation has provided the values of the gluon condensate in finite magnetic field and temperature [201], it will be interesting to investigate effects of these modifications by using the QCDSR.

Also, applications of the QCDSR to light mesons and open heavy flavor mesons are an interesting topic because these mesons contain a light quark component, which reflects the modification of the quark condensates. In Refs. [26, 31], the OPE was performed by taking into account not only the modification of $\langle \overline{q}q\rangle$ but also the emergence of the tensor-type condensate $\langle \overline{q}{\sigma }^{\mu \nu }q\rangle$, which have been calculated by the lattice QCD simulations [48-50, 200]. Furthermore, the mass spectra of the open heavy flavor mesons would serve as a signature of the magnetically induced QCD Kondo effect [219] which was proposed on the basis of a close analogy to the QCD Kondo effect in the dense QCD [220, 221]. The readers are referred to the literature for the interesting analogy between the (1+1)-dimensional low-energy dynamics at high densisy [222-224] and in strong magnetic fields [43, 53, 225, 226].

4.2 Transport of heavy quarks in a hot medium and magnetic field

As mentioned in the beginning of this section, the heavy quarks are dominantly created by the initial hard scatterings among the partons from the colliding nuclei (see Fig. 4.2). Therefore, the heavy quarks typically have the power law pT spectrum according to the perturbative QCD. This spectrum is however gradually modified by the multiple interactions with the thermal quarks and gluons in the QGP, and the spectrum approaches a thermal distribution, which is realized at a large enough time. However, the QGP has a finite lifetime, and to what extent the thermalization is achieved depends on the transport property of the QGP. Therefore, the heavy quarks will probe the transport property of the QGP.

In the RHIC and LHC experiments, the nuclear modification factor, RAA, and the anisotropic spectrum, v₂, of the open heavy flavors have been measured (see a review [202] and references therein). Both of them are thought to be sensitive to the heavy-quark thermalization and are closely related to each other [227]. However, a consistent theoretical modeling, which simultaneously reproduces RAA and v₂, is still under investigation. Effects of the magnetic field could be new ingredients to resolve the discrepancy between the theoretical and experimental results. In this section, we summarize a recent study on the heavy quark dynamics in a hot medium and magnetic field on the basis of Langevin equation and the perturbative computation of the heavy-quark diffusion constant.

The dynamics of the low-energy heavy particles is modeled as a Brownian motion caused by random kicks by the thermal particles which is described by the Langevin equations [227]:

Since the external magnetic field provides a preferred spatial direction, we have a set of two equations for the heavy quark motions, one parallel and one perpendicular to the magnetic field that is oriented in the z-direction. The random forces are assumed to be white noises;

and these coefficients, κ_∥ and κ_⊥, are related to the drag coefficients, ≤_∥ and ≤_⊥, through the fluctuation-dissipation theorem as

Here, we would like to mention a recent numerical study [42] where the authors investigated the Lorentz force exerting on the heavy quark with the inclusion of the Lorentz force in the Langevin equation (51).

In addition to effects of the Lorentz force exerting on the heavy quark as a direct effect of the magnetic field, the diffusion and drag coefficients will be changed because the light quarks in the hot medium are strongly affected by the magnetic field [40, 41]. At the leading order in gs, the anisotropic momentum diffusion coefficients, κ_∥ and κ_⊥, can be defined by

where q is the amount of the momentum transfer from the thermal particles to the heavy quark, and the static limit (q⁰ → 0) is assumed in the above definitions. The momentum transfer rate $\frac{d\Gamma (q)}{d^{3}q}$ is computed from the gluon-exchange diagrams shown in Fig. 18. Effects of a magnetic field appear in two places (highlighted by red). While the gluons are not directly coupled to the magnetic field, (i) the Debye screening mass and (ii) the dispersion relation of the thermal-quark scatterers will be changed as discussed below. The Debye mass is necessary for cutting off the infrared divergence in the forward scattering.

Fig. 18.Full-screen View

Brownian motion of heavy quarks created by the initial hard processes.

In a strong magnetic field, the fermion wave function is strongly squeezed along the magnetic field, corresponding to the small radius of the cyclotron orbit. Indeed, from the Landau level quantization and the Zeeman effect, the dispersion relation of fermions in the lowest Landau level (LLL) becomes the (1+1) dimensional, i.e., ${ϵ}^2=m^2+p_z^2$ for a massive fermion and ϵ=±pz for massless fermions with right- and left-handed chiralities. To be specific, we focus on the strong-field regime such that the transition from the LLL to the hLL states are suppressed according to a hierarchy T²≪eB.

(1) First, the Debye screening mass is obtained from the gluon self-energy computed by the standard method in finite temperature field theory in the (1+1) dimension. The gluon self-energy is completely factorized into the transverse and longitudinal parts as

where the factor of $q_f\frac{eB}{2\pi }$ comes from the degeneracy factor in the transverse phase space, and the Gaussian is the wave function of the LLL state. The longitudinal part is nothing but the (1+1)-dimensional polarization tensor, which is known as the Schwinger model:

Here, the longitudinal momentum is defined by $q_{\parallel }^{\mu }=(q^0\mathrm{,0,0,}{q}_z)$. There is only one gauge-invariant tensor structure in the (1+1) dimension, so the temperature correction is, if any, contained in a function f(q⁰,qz). However, an important observation is that there is no temperature or density correction in the massless case [40, 228, 229], and f(q⁰,qz) = 1. Accordingly, in the massive case, a temperature correction is proportional to the quark mass and is suppressed by mq/T (see Ref. [40] for explicit expressions). Therefore, in both cases, the Debye mass in strong magnetic fields (eB ≫ T²) is given by the vacuum part of the gluon self-energy [230-232]. Namely, we find $m_D^2\sim q_f\frac{eB}{2\pi }\cdot \frac{g_s^2}{\pi }$ which is much larger than the usual thermal mass squared, ∼(gsT)².

(2) The change in the dispersion relation of fermion scatterers is important for the kinematics of the Coulomb scattering process [40, 233]. We shall consider a massless fermion scatterer in the lowest Landau level (see Fig. 19). For the massless quarks to be on-shell both in the initial and final states, the energy transfer has to be $q^0=\pm (k_z^{\text{'}}-k_z)$ and thus q⁰=±qz. Note that the chirality does not change at any perturbative vertex, so that the signs appear only as the overall ones. When taking the static limit q⁰ → 0 in the definition of the diffusion constant, one can immediately conclude that the longitudinal-momentum transfer is kinematically prohibited in the massless case. On the other hand, this constraint does not apply to the transverse momentum transfer, so the transverse momentum transfer is allowed. This finite contribution to the diffusion constant can be obtained either by the direct computation of the matrix elements in Fig. 19 or by using the cutting rule. At the leading order in αs, the transverse momentum diffusion constant in the massless limit is obtained as [40]

Figure 19:Full-screen View

Coulomb scattering amplitudes contributing to the heavy-quark (HQ) momentum diffusion. A magnetic field acts on the quark loop in the polarization and the thermal-quark scatterers.

Therefore, strong magnetic fields give rise to an anisotropy of the momentum diffusion constant.

Nonvanishing contributions to the longitudinal momentum transfer come from either the finite quark-mass correction or the contribution of the gluon scatterers. In Ref. [40], the contribution of the mass correction was obtained as

while the gluon contribution is obtained by substituting the Debye mass $m_D^2\sim {\alpha }_{s}eB$ for the conventional one (gsT)² in Refs. [227, 234] as

We have assumed a hierarchy αseB≪T²≪eB. Complete expressions, including the prefactors, are found in Ref. [40]. It is instructive to compare the mass correction (58) with the gluon contribution (59). The ratio is written as

While the first two factors are small in our working regime, the last factor can be large. Therefore, the massive contribution ${\kappa }_{\parallel }^{\text{L}O\mathrm{,}massive}$ could be in principle as comparably large as ${\kappa }_{\parallel }^{\text{L}O\mathrm{,}gluon}$, and this happens when $eB\sim {\alpha }_s(T^6/m_q^4)$. However, to be consistent with our assumed regime, αs eB ≪ T², we have a constraint of ${\alpha }_s\ll m_q^2/T^2$, which is not quite likely true in the heavy ion collisions. Hence, the longitudinal momentum diffusion constant is dominated by the gluon contribution ${\kappa }_{\parallel }^{\text{L}O\mathrm{,}gluon}$.

Now that we obtained the leading contributions in Eqs. (57) and (59), the anisotropy in the momentum diffusion constant and the drag force is estimated to be

This large anisotropy is induced by the enhancement of the density of states in the quark contribution ${\kappa }_{\perp }^{\text{L}O}$, which is proportional to eB. On the other hand, a magnetic field does not change the phase space volume of the thermal gluons in the gluon contribution ${\kappa }_{\parallel }^{\text{L}O\mathrm{,}gluon}$, but changes only the Debye screening mass.

Based on the above findings, a toy model for the dynamical modeling was discussed in Ref. [40] by using the Fokker-Planck equation [227]. Since the friction exerted on the heavy quarks is stronger in the direction transverse to the magnetic field (see Fig. 20), the heavy-quark momentum in the final state will acquire an anisotropy. In the small momentum region, the heavy quarks will be dragged by the flow more strongly in the transverse direction while they will feel a larger resistance in the transverse direction as the momentum increases. Therefore, there is a turnover in the anisotropy of the heavy quark spectrum (see the right panel in Fig. 20).

Figure 20:Full-screen View

Schematic pictures of the anisotropic drag force in a magnetic field (left) and the resultant anisotropic heavy-quark spectrum in the final state (right). (Reproduced from Ref. [40])

It will be interesting to implement the anisotropic diffusion and drag coefficient in the numerical simulations. As the magnetic field in the heavy-ion collisions is time dependent, it is also necessary to investigate the weak-field limit as the first step in considering the time dependence. These issues are left as open questions. Finally, the drag force in the presence of the anomalous currents was also investigated in Ref. [233, 235, 236].