Quark Mass Correction to Chiral Separation Effect and Pseudoscalar Condensate

We derived an analytic structure of the quark mass correction to chiral separation effect (CSE) in small mass regime. We confirmed this structure by a D3/D7 holographic model study in a finite density, finite magnetic field background. The quark mass correction to CSE can be related to correlators of pseudo-scalar condensate, quark number density and quark condensate in static limit. We found scaling relations of these correlators with spatial momentum in the small momentum regime. They characterize medium responses to electric field, inhomogeneous quark mass and chiral shift. Beyond the small momentum regime, we found existence of normalizable mode, which possibly leads to formation of spiral phase. The normalizable mode exists beyond a critical magnetic field, whose magnitude decreases with quark chemical potential.


Introduction and Summary
The chiral magnetic effect (CME) [1,2,3] and chiral separation effect (CSE) [4,5] characterize the response of vector/axial current to the axial/vector chemical potential in external magnetic field. Both effects are manifestation of axial anomaly and are of phenomenological interest in heavy ion collision experiment. In particular, CME leads to charge separation and the interplay of CME and CSE gives rise to chiral magnetic wave (CMW) [6], which leads to charge dependent flow [7]. There have been significant experimental efforts in search of CME [8,9,10] and CMW [11,12], see [13,14,15] and references therein.
While CME and CSE share many similarities, they are known to differ in certain aspects. The chiral magnetic current is known to be independent from quark mass, temperature etc [2]. Correction may arise in dynamical cases, where axial chemical potential is not well defined and the dynamics of axial charge becomes important [16,17,18,19]. The chiral separation current does not suffer from the issue of axial chemical potential, but it does receive correction from quark mass [4,20,21,22,23]. In the static case, the correction to CSE can be derived in an ad-hoc way: with C = − Nce 2 Q 2 2π 2 . In the massless limit, we can writeẼ = −∇µ q , with µ q being the quark chemical potential. Since ∇ ·B = 0, we easily arrive at the celebrated CSE ∇ · j 5 = −∇ · Cµ qB ⇒ j 5 = −Cµ qB . ( To obtain the mass correction to j 5 , we first write the pseudoscalar operator σ 5 ≡ iM qψ γ 5 ψ as the response to quark chemical potential: σ 5 (x) = d 4 yG σ5n (x − y)µ q (y). The structure of the Green's function G σ5n can be deduced from discrete symmetry: σ 5 is odd in both parity and time reversal, thus it should contain magnetic field B, which is odd in time reversal and spatial gradient ∇, which is odd in parity. Therefore, to the lowest order in gradient, we have Following the same steps as (2), we find correction to j 5 , The function g is related to the Green's function as (in momentum space) This relation will be confirmed analytically in model study. Note that g vanishes when the quark mass M q vanishes. We can expand it in small M q regime: assuming µ q ≪ T,B ≪ T 2 . The dimensionless prefactor # is to be determined by dynamics. In fact, the analytic form of g also constraints the response of σ 5 to quark mass M q . Note that we have assumed a spatially inhomogeneous µ q and constant M q . Instead, we can assume an inhomogeneous M q and constant µ q . This should induce vev of σ 5 (x) = d 4 yG σ5σ (x − y)M q (y). Consistency with (3) and (6) indicates the following correction which implies G σ5σ = 2i#M q kBµ q + o(M q ). We will provide clear numerical evidence for this correlator in model study.
The correction to j 5 is more interesting in regime of large µ q andB. WhenB = 0 and µ q large, different instabilities have been discussed in large N c field theory with spontaneous generation of chiral density wave [24,25], current density [26,27,28] and quarkyonic spiral [29,30,30,31] etc. At strong magnetic field, formation of chiral magnetic spiral [32,33,34] is possible. Here we discuss a different type of instability characterized by pseudoscalar condensate. This instability is already identified in [35] see also [36] in low temperature confined phase. We extended the discussion and found it only exists within a window of magnetic field. Formation of this instability leads to spontaneous generation of chiral shift, first introduced in [37], which induces further correction to j 5 .
The paper is organized as follows: In Sec II, we give a brief review of the holographic model and the finite density and magnetic field background; Sec III contains a study of correlators among pseudoscalar condensate, quark condensate and quark number density in small momentum regime; Sec IV extends the study of correlators in arbitrary momentum regime and discussed the instability towards formation of spiral phase. We close the paper in Sec V with some outlooks.
2 A brief review of the model

The finite density background
We use the D3/D7 model to study the effect of finite quark mass. The background consists of N c D3 branes and N f D7 branes. In the probe limit N f ≪ N c , the background is simply given by black hole background sourced by D3 branes, with suppressed backreaction from D7 brane. The D3/D7 model is dual to N = 4 Super Yang-Mills (SYM) fields and N = 2 hypermultiplet fields, which transform in adjoint and fundamental representations of SU (N c ) gauge group respectively. By analogy with QCD, we loosely refer to the N = 4 and N = 2 fields as gluons and quarks respectively. The black hole background of D3 branes is given by [38]: where The temperature of the gluon plasma is given by T = r 0 π . Note that we set AdS radius L = 1. It can be reinstated by dimension. We also explicitly factorize S 5 into S 3 and two additional angular coordinates θ and φ. There is also a nontrivial Ramond-Ramond form The D7 branes share the worldvolume coordinates with D3 branes. In addition, they occupy the coordinates x 4 -x 7 parametrized by the S 3 coordinates. Their position in x 8 -x 9 plane can be parametrized by radius ρ sin θ and polar angle φ. The rotational symmetry in the x 8 -x 9 plane corresponds to U (1) R symmetry in the dual field theory. The D7 branes has an additional U (1) B symmetry carried by its worldvolume gauge field. We will use the U (1) R and U (1) B symmetries as axial and vector symmetries respectively.
We are interested in the field theory state at finite temperature and finite quark chemical potential µ q with background magnetic fieldB. To this end, we introduce worldvolume gauge field A t (ρ) andF xy =B. The embedding function θ(ρ) of D7 branes in D3 background is determined by minimizing the action including a DBI term and WZ term Here T D7 is the D7 brane tension. g ab andF ab are the induced metric and worldvolume field strength respectively. Defining the action simplifies to The asymptotic behavior of θ is given by The coefficients are related to bare quark mass M q and quark condensate ψ ψ [38]: Similarly, the asymptotic behavior of A t determines dimensionless quark chemical potential µ and density n: with physical quark chemical potential and number density given by µ q = r 0 µ 2πα ′ and n q = 4πα ′ N r 3 0 n. The phase diagram of the system has been obtained in [39,40,41,42,43]. There are two possible embeddings with D7 branes crossing/not crossing the black hole horizons, corresponding to meson melting/mesonic phase respectively [44]. We will focus on meson melting phase for studying CSE in quark gluon plasma (QGP).

CSE at finite quark mass
We consider the fluctuation of embedding function φ in the above background. The part of quadratic action containing φ can be written in the following form with M = t, x 1 , x 2 , x 3 , ρ. For the evaluation of CSE, we need Ω = cos 4 θ, κ = 1 8 .
We do not need explicit form of G M N for now. The axial current is defined by [45] J µ R = dρ δL δ∂ µ φ .
Using EOM of φ, we obtain the following non-conservation equation of axial current We will identify J R as axial current. The non-conservation of J R follows from two boundary terms in the integration. The boundary term at ρ = ∞ is related to axial anomaly: Note that the factor r 4 0 follows from dimension of L. In doing this, we have chosen r 0 to set unit and work with dimensionless coordinates t, corresponding to the anomaly term. Therefore the term O η corresponds to the mass term iM qψ γ 5 ψ 1 . The other boundary term at horizon ρ = ρ h is an artifact of the model. Its presence is tied to our modeling of axial symmetry: since we make use of U (1) R symmetry for axial symmetry, the gluon plasma is also charged under axial symmetry. The horizon term represents axial charge exchange between the quarks (fundamental matter) and gluons (adjoint matter). The term is indeed non-vanishing in known examples [45,46]. However, we will study CSE and correlation functions in static limit. We claim the above artifact is absent in those quantities because charge exchange is not possible in static case. This can be checked explicitly. Now we proceed to evaluate CSE, which is the axial current J 3 R . Note that φ = 0 in the background, we obtain We stress that had we assumed φ = 0 at the beginning, we would have obtained a vanishing CSE current. The case m = 0 is trivial. In this case, the embedding function is given by R as a function of m, µ and B. The temperature of QGP is set to T = 300MeV. To guide eyes, we mark phenomenological relevant parameters (strange quark mass M q = 100MeV,B = m 2 π , µ q = 0.5T ) with red dots in corresponding panels. In the upper right and lower panels we use M q = 300MeV to magnify the dependence on µ(µ q ) and B(B). We use α s = 0.3 in determination of λ. For strange quark mass M q = 100MeV, the µ and B dependence is barely visible. This is the standard CSE fixed entirely by anomaly upon restoring units in the last step.
Correction to standard CSE exists for m = 0. In this case, the embedding function θ and gauge potential A t are only known numerically. The corresponding J 3 R can be obtained by numerical integration. We obtain its dependence on m, µ and B in Figure 1. To convert to physical unit, we use with phenomenologically relevant coupling and temperature. We observe that quark mass tends to suppress CSE as expected. Chemical potential and magnetic field both tends to enhance CSE. The qualitative dependence can be understood from (4). The leading term −Cµ qB gives the baseline 1, while the correction is The m(M q ) dependence is quadratic, while the µ q andB dependence is absent at this order.
Assuming higher order terms in M q can be ignored, the dependence in Fig.1 implies the magnitude of the prefactor # drops with growing µ q andB. Note that the prefactor is negative, a suppressed magnitude leads to enhancement of CSE.

Correlators
In this section, we wish to study the correlator among scalar condensate σ ≡ψψ, pseudoscalar condensate σ 5 ≡ iM qψ γ 5 ψ and quark number density n q ≡ψγ 0 ψ. We study the Euclidean correlators at vanishing frequency (in static case when axial charge exchange is absent) and finite momentum.
with µ, ν = σ, n, σ 5 . To this end, we introduce the following fluctuations to the background: The open string metric up to quadratic order in fluctuation is given by with Here we use the following coordinates as D7 brane worldvolume coordinates: t, z, ρ, x, y and Ω, with Ω denoting collectively three angular coordinates on S 3 . They are ordered as they appear in the open string metric. It is straight forward but tedious task to work out the quadratic action of the DBI and WZ terms Here g SS = cos 2 θ, δg We use dot and prime for derivatives with respect to z and ρ respectively. Note that we work with dimensionless z, i.e. ∂ z → r 0 ∂ z . This amounts to setting the scale of spatial momentum by temperature. The rescaling makes the r 0 dependence of S DBI and S WZ appears as an overall r 4 0 factor, thus r 0 drops out completely from the EOM. The EOM following from (30) are given by By observation, we find the ansatz φ(z, ρ) = sin(kz)φ k (ρ), a t (z, ρ) = cos(kz)a t (ρ), δθ(z, ρ) = cos(kz)δθ(ρ), solves the z-dependence of (32). To proceed, we note that a generic set of solution has the following asymptotic expansion where f h = −k 2 f 0 , a h = −k 2 a 0 and t h = −k 2 t 1 . The leading coefficients are the sources to operators σ 5 , δn and δσ respectively. The subleading coefficients are related to their vevs. A holographic renormalization procedure is needed to determine the vevs. We will elaborate this procedure in appendix A. Here we only show results of correlator G ab where we have defined individual responses S ab We proceed to solve (32). Since we have three coupled differential equations, we expect to have three independent solutions. We solve (32) by numerical integration from the horizon to the boundary. The initial condition we impose at the horizon is regularity condition. In practice, we start off the horizon with the following three independent solutions: These solutions give rise to the following asymptotics at the boundary with i = 1, 2, 3 labeling different solutions. In order to calculate individual responses S ab , we need to construct proper solution for which the other two sources vanish. This can be done efficiently in the following way On general ground, we expect the Euclidean correlator to be real and symmetric G ab = G * ab = G ba . Our numerical results confirm that this is indeed the case. We also find the following scaling of all individual responses at small k.
Similarly, the diagonal responses S σσ and S nn are related to ∂c ∂m and ∂n ∂µ 2 . We have compared the results of individual responses at k = 0 and those of thermodynamics, finding expected agreement.
The second line of (40) is of more interest to us. It follows from (35) that G σσ5 , G nσ5 ∼ O(kB). This is consistent with parity (P) and time-reversal (T) symmetry of the corresponding operators: σ 5 is odd under both P and T, while σ and n are even under P and T . The external B and momentum k are odd under T and P respectively. These Euclidean correlators characterize response of the system to external parameter φ: Here we use φ to denote the field theory source coupled to σ 5 . The scaling O(k) can be understood as follows: the source φ enters field theory Lagrangian as M qψ e iφγ 5 ψ. If we perform a chiral rotation, ψ → e −iγ 5 φ/2 ψ, the relevant terms in the Lagrangian is modified as [45] M This implies that only ∂ µ φ appears as physical parameter. In our case, φ only depends on z in field theory coordinates, therefore the physical parameter isφ. Interestingly,φ can be identified as the chiral shift parameter proposed in [37,20]. Assuming the response of σ and n toφ is O(1) at small k, we naturally explain the O(k) scaling of correlators G σσ5 , G nσ5 .
The third line of (40) indicates G σ5σ5 ∼ O(k 2 ). G σ5σ5 is by definition the susceptibility of σ 5 . The susceptibility G σ5σ5 is parity even, thus it scales as even power of k. As we argued above, any response to φ has to start from O(k), the most probable scaling is Let us take a closer look at correlators involving σ 5 . We will present results primarily on these correlators at small k. In fact, we can confirm the linear scaling relation described above by perturbative calculation in k. Note that φ ∼ O(k), a t , δθ ∼ O(1). To the lowest non-trivial order in k, we only need to solve the following equations From the first two equations of (44), we can solve for δθ and a t . Plugging the solution into the third equation and integrating from the horizon to the boundary, we obtain On the left hand side (LHS), the boundary term at the horizon vanishes, the boundary term at infinity is just − m 2 f 2 2 . On the right hand side (RHS), it is related to the sources t 1 and a 0 . There are two independent solutions. We denote their asymptotics as with i = 1, 2. Using (45), each solution give rise to vev of σ 5 ∼ m 2 f (i) 2 . Similar to (39), we obtain the correlators as Note that at the boundaryȧ t ∼ E, δθ ∼ δṁ. The correlator S σ5n (G σ5n ) measures the response of σ 5 to parallel E and B fields. To study the response in more detail, we define a dimensionless ratio Note that we have included a minus sign in the definition of r such that r is always positive.
In terms of correlators, r = S σ5n 2ik . Since S σ5n ∼ O(k), r approaches a constant in the limit k → 0. We plot the m and µ-dependence of r(k → 0, µ = 0) in Figure 2. We find r increases with m, but decreases with µ. The dependence is in qualitative agreement with our discussion before: r ∼ g ∼ # M 2 q T 2 , which grows with m, and drops with µ q from the µ q dependence of the prefactor #. Now we are in a position to confirm the claim (4) in Sec I.
We begin by working in the background with µ = 0 (A t = 0). The corresponding correlator G σ5n becomes particularly simple then Noting that r 0 µ 2πα ′ = µ q , we obtain One the other hand, the correction to CSE can be obtained by performing an integration by part on (22) The first term of (51) corresponds to the standard CSE, while the second term comes from mass correction, which is precisely (50). In this case, r(k → 0) takes the form It is instructive to analyze the coupling dependence of r: we first note that in the case µ = 0, the coupling enters only through m = 2 √ λ Mq T . On the other hand, we have argued in the introduction that r ∼ g ∼ # M 2 q T 2 , which suggests the following dependence r ∼ 1 λ M 2 q T 2 . This seems to imply that stronger interaction leads to weaker response of σ 5 to external electromagnetic fields. This interpretation is misleading for the following reason: in D3/D7 model, the electromagnetic coupling to quark is the same as strong coupling, thus N EB ∼ O(λ), so the actual response of σ 5 is O(λ 0 ). It is also interesting to compare r with the same quantity studied in [47], which is defined in the regime ω → 0, k = 0. In fact, we can show analytically that they do not agree. For monotonic a t , we have This reveals noncommutativity of the limits ω → 0, k → 0 and k → 0, ω → 0 in the response of σ 5 . The correlator G σ5n tells us more than the response of σ 5 . Note that we have G σ5n = G nσ5 by symmetry. G nσ5 characterizes the response of n to chiral shiftφ. The result of r indicates that chiral shift can also induce correction to n in the presence of B, with the correction increases with m, but decreases with µ.
Now we turn to S σ5σ . This correlator measures the response of σ 5 to spatially varying quark mass δm. We plot the m-dependence and µ-dependence of S σ5σ in Figure 3. Indeed, we can see in Fig. 3 that G σ5σ vanishes approximately linearly in µ and m, which is clear evidence for (7). By symmetry S σ5σ = S σσ5 , we also obtain that chiral shiftφ can induce Finally we plot the m and µ dependence of S σ5σ5 in Figure 4. The scaling S σ5σ5 ∼ O(k 2 ) allows for the following parametrization of σ 5 From (54), we easily obtain an induced j 5 in the presence of chiral shiftφ: A similar current from spatial gradient of axion is also discussed in [48]. As we discussed before, ∇φ is just the chiral shift parameter, which couples to the axial current in the Lagrangian. The function h can be viewed as an effective susceptibility. Fig.4   vanishes: We show the momentum k of the mode as a function of B in Figure 5. We find that normalizable modes exist for medium with general nonvanishing µ beyond certain critical magnetic field B c . For each B > B c , there are two normalizable modes with different momenta k. The low momentum branch appears monotonic decreasing function of B, while the high momentum branch is non-monotonic. The normalizable modes we find are numerically consistent with the quasi-normal mode reported in [35]. It is interesting to note that the critical magnetic field corresponds to the point where the two momenta merge. Furthermore, the modes extend to the region of large B, where the state possibly becomes metastable [43]. We do not keep the corresponding mode in Fig. 5. Turning to the µ dependence, we see that as µ is lowered, B c grows. This is qualitatively in agreement with the chiral soliton solution found in [36] in confined phase. We also show k as a function of m in Fig. 5, which clearly shows the low/high momentum branches. As a function of m, the high momentum branch appears monotonic increasing function, while the low momentum branch is non-monotonic.
To have an idea on the magnitude of magnetic field, we convert B c for the case µ = 3 to physical unit. For gluon plasma at temperature T = 300MeV and coupling α s = 0.3.

Outlook
We close this paper by discussing several open questions that we may address based on the results of this paper. Firstly, how does quark mass affect the dynamics of axial and vector charge. As we have seen the pseudoscalar condensate responses to gradient of quark chemical potential etc. This brings in an additional coupling between axial and vector charges. Consequently, it should also modify the dispersion of CMW. For strange quark mass, we expect from the Fig. 1 that the modification in phenomenology is modest. Since our study focuses on Euclidean correlators, the same quantities can be reliably studied on the lattice, which will provide quantitative answers for quark mass effect in real world QCD.
Secondly, the normalizable mode we found at sufficient large B and µ suggests possible formation of spiral phase. To find the true ground state, we need to go beyond the linear analysis. We expect that the true ground state is characterized by the spontaneous generation of chiral shift, which induces further correction to axial current. The correlator G σσ5 and G nσ5 indicates the correction to σ and n as well. It would be interesting to find out detail about this state. We leave this for future work.
Last but not the least, the normalizable mode provides an explicit example of spiral phase in meson melting phase. It would be interesting to extend this work to mesonic phase. Such an example has been found at zero temperature in [36] based on effective field theory models. It would be interesting to study the stability of such state against finite temperature fluctuation.

Note added
When this work was near complete, we learned that Qun Wang et al was about to finish a closely related work [22]. We thank Qun Wang for sharing with us notes of their work before publication.
Plugging (58) into (57), we find the second line always vanishes in the limit Λ → ∞. The first line gives the following contribution We need the following counter terms to remove quadratic divergence in (60) The coefficient of the logarithmic terms (59) is a special case of a more general relation: They do not encode dynamics of the theory. The corresponding logarithmic and finite terms can be removed by the following counter terms with appropriate normalizations N √ −γδθ δθ ln Λ, N √ −γa t a t ln Λ, N √ −γφ φ ln Λ, Adding all counter terms to (57) and dropping contributions like m 2 t 2 1 , which do not encode dynamics of the theory, we obtain the following renormalized on-shell action S ren = N r 4 0 − (64) Now we can do variation of (64) with respect to sources to obtain vev of the corresponding operators. Note that partial derivatives hit twice in each terms in the bracket. We obtain δσ = δS ren δt 1 = N r 3 0 (2πα ′ ) − t 3 2 , δn = δS ren δa 0 = N r 3 0 (2πα ′ )a 2 , We use the δ symbol to indicate that the vev is on top of a nonvanishing background.
Taking the derivatives once more, we obtain correlators shown in (35).