Generation of strong magnetic fields in old neutron stars driven by the chiral magnetic effect

# Generation of strong magnetic fields in old neutron stars driven by the chiral magnetic effect

Maxim Dvornikov  and Victor B. Semikoz
Pushkov Institute of Terrestrial Magnetism, Ionosphere
and Radiowave Propagation (IZMIRAN),
108840 Troitsk, Moscow, Russia;
Physics Faculty, National Research Tomsk State University,
36 Lenin Avenue, 634050 Tomsk, Russia
maxdvo@izmiran.rusemikoz@yandex.ru
###### Abstract

We suggest the generalization of the Anomalous Magneto-Hydro-Dynamics (AMHD) in the chiral plasma of a neutron star (NS) accounting for the mean spin in the ultrarelativistic degenerate electron gas within the magnetized NS core as a continuing source of the chiral magnetic effect. Using the Parker dynamo model generalized in AMHD, one can obtain the growth of a seed magnetic field up to for an old non-superfluid NS at its neutrino cooling era , while neglecting any matter turbulence within its core and assuming the rigid NS rotation. The application of the suggested approach to the evolution of magnetic fields observed in magnetars, , should be self-consistent with all approximations used in the suggested laminar dynamo, at least, up to the jumps of growing fields.

## 1 Introduction

The axial Ward anomaly known as Adler-Bell-Jackiw (ABJ) anomaly for the axial-vector current in QED was used to formulate the chiral magnetic effect (CME) [1, 2] in different media such as QCD plasma [3, 4] or the hot plasma of the early Universe [5].

There are multiple attempts to apply the chiral magnetic phenomena in astrophysics. First, we mention recent refs. [6, 7], where the Anomalous Magneto-Hydro-Dynamics (AMHD) turbulence is used to generate strong magnetic fields in a protoneutron star (PNS) and explain observed magnetic fields in some compact stars called magnetars [8]. The AMHD is applied in ref. [9] to account for high linear velocities of pulsars. Let us also mark the recent development on laminar and turbulent dynamos in the AMHD, called the chiral magnetohydrodynamics in ref. [10], applied both for the hot plasma in early universe and in PNS in ref. [11].

A seed magnetic field in PNS can be amplified if the CME is accounted for in the system of chiral fermions. This amplification model is based on the appearance of the chiral imbalance at the initial moment , where are the initial densities of right and left electrons, which arises owing to the left electron capture in a corresponding PNS. Nevertheless the application of the AMHD for the generation of a magnetic field in a neutron star (NS), driven by the CME, fails. The failure of this approach is because of the fast spin-flip due to Coulomb collisions of electrons in NS plasma resulting in during  [12].

In the present work, we look for alternative mechanisms to support the CME during a time much longer than accounting both for the chiral vortical effect (CVE) [13], given by the vorticity , and the chiral separation effect (CSE) [14], given by the mean spin in magnetized plasma. In particular, we study the saturation regime, , resulting in . It arises due to the CSE when a magnetic field feeds the chiral imbalance . It happens in the presence of the non-uniform mean spin , where is the electron chemical potential in the NS core.

Then we study the generation of strong magnetic fields driven by the CME in a non-superfluid NS at late stages of its evolution without any turbulence within the NS core. Note that the turbulence (convection) is rather specific for a nascent NS because of the huge neutrino emission rate leading to the generation of magnetic fields in PNS at the first seconds after a supernova (SN) explosion [15].

Thus, instead of the standard MHD dynamo mechanism based on the parameter given by the kinetic helicity, , where is the random (fluctuation) velocity, we apply the CME originated by the presence of the non-uniform chiral imbalance . In this situation, the helicity parameter leads to the instability of a seed magnetic field in NS. Note that such a pseudo-scalar field, resulting from the mean spin, has no relation to an effective axion field discussed in ref. [16].

The work is organized as follows. In section 2, we present the full set of AMHD equations completed by the evolution equation for the chiral anomaly density . In section 3, we derive the complete system of the evolution equations for the chiral imbalance and the magnetic field assuming the rigid NS rotation, . Then, in section 4, we compare the action of the CME versus the CVE and then the CSE versus the axial current given by the vorticity .

In section 5, we consider the saturation regime that allows to reduce the system of the evolution equations in AMHD to the single nonlinear Faraday equation for . In section 5.1, we present the approximate solution of this Faraday equation using -dynamo approach when assuming the one dimensional field as the Chern-Simons wave. Then, in section 5.2, we derive the closed system of the evolution equations for the azimuthal components of the axisymmetric 3D magnetic field: and the potential . For this purpose we use the Parker’s dynamo model [17] modified here in AMHD. In section 5.3, using the low mode approximation in a thin layer just under the NS crust, , we derive finally the closed system of the ordinary differential equations for the four amplitudes: for the azimuthal potential and for the toroidal field . The numerical solution of these equations is illustrated in section 6. Finally, in section 7, we discuss both useful issues and shortcomings of our approach for the generation of strong magnetic fields in such NSs as magnetars.

## 2 AMHD in a neutron star with magnetic field

In this section, we write down the full system of AMHD equations and derive the equation for the magnetic helicity evolution.

The full set of AMHD equations for a plasma in NS, accounting for both the CME and the CVE, is given by [18, 19]

 ρ[∂v∂t+(v∇)v−ν∇2v]= −∇p+J×B, (2.1) ∂ρ∂t+∇(ρv)= 0, (2.2) ∂B∂t= −(∇×E), (2.3) (∇×B)= σ[E+v×B+e22π2σμ5B+eμe2π2σμ5ω], (2.4) ∂n5(x,t)∂t= −∇⋅S(x,t)+2αemπ(E⋅B)−Γfn5(x,t), (2.5)

where is the plasma velocity, is the matter density, and are the electric and magnetic fields, is the plasma pressure, is the electric conductivity, is the viscosity coefficient, is the spin flip rate, and is the fine structure constant. Here, in contrast, e.g., to refs. [18, 19], we presented eq. (2.5) in its local form and modified eq. (2.5) adding a new pseudoscalar term given by the divergence of the mean spin .

Equation (2.5) stems from the statistical averaging of the ABJ anomaly,

 (2.6)

completed by losses for due to the spin-flip through collisions in NS plasma (). In eq. (2.5), is the chiral anomaly density. The densities of right and left electrons are related to the corresponding chemical potential by . Thus, for one obtains that in the NS degenerate electron gas where is the non-uniform chemical potential. Using the results of ref. [20], we get that the conductivity of NS matter at the temperature has the value, .

Finally the new term in kinetic eq. (2.5) contains the mean spin in magnetized plasma,

 S(x,t)=⟨ψ+eΣψe⟩0=−(eμe(r)2π2)B(x,t), (2.7)

which arises due to the plasma polarization through the motion of massless electrons and positrons at the main Landau level along the magnetic field; c.f. eq. (2.2) in ref. [22]. Note that the mean spin in eq. (2.7) does not depend on the plasma temperature. It is known as the CSE (see, e.g., refs. [14, 24]).

Then, integrating local eq. (2.5) over space, , one obtains the master eq. (3.14) in ref. [22] for the ultrarelativistic plasma when 111In ref. [22], the evolution eq. (2.8) stems from the same ABJ anomaly (2.6) without the collision term . The densities , as well as the magnetic helicity density , depend on time only. The sum reduces to the mean spin in eq. (2.7), , since , originated by the mean pseudoscalar , vanishes in the massless limit , cf. ref. [25].,

 ddt(nR−nL+αemπh)=−αemπV∮S([E×A+A0B]⋅n)d2S−∮S(S⋅n)d2SV. (2.8)

The magnetic helicity density in eq. (2.8) evolves as

 dhdt=−2∫d3xV(E⋅B)−∮(n⋅[BA0+E×A])d2SV, (2.9)

which is known result in the standard MHD [23].

## 3 Evolution of the chiral anomaly in NS

In this section, we derive the modified Faraday equation and the equation for the evolution of the chiral imbalance accounting for both the CME and the CVE.

Substituting the electric field from the Maxwell eq. (2.4) into eq. (2.5) and using the chiral imbalance , we can change the term there as

 2αemπ(E⋅B)=2αemπ⎡⎣B⋅(∇×B)σ−(μ2e2σ)(eBμ2e)2n5−(eμe2σ)(ω⋅Bμ2e)n5⎤⎦. (3.1)

Such a modification allows to recast local eq. (2.5) as

 ∂n5(x,t)∂t= −∇⋅S(x,t)+2αemπσ[B⋅(∇×B)] (3.2) −n5(x,t)[Γf+(2αem)2σμ2eB2+αemeμeπσμ2e(ω⋅B)].

Returning to the chiral imbalance , one obtains

 ∂μ5(x,t)∂t=π2μ2e{−∇⋅S+2αemπσ[B⋅(∇×B)]}−Γfμ5[1+B2B20+μe(ω⋅B)eB20], (3.3)

where . The strength of this magnetic field is , which corresponds to . The electric conductivity is  [20] at the temperature and in the NS core222We put stationary and uniform, , except of a dependence of the quantum (spin) contribution on , given by the NS density profile , see in eq. (2.7)..

Substituting the electric field from the Maxwell eq. (2.4) into eq. (2.3) one obtains the induction (Faraday) equation in AMHD,

 ∂B∂t=∇×(v×B)+1σ∇2B+e2π2σ∇×μ5(x,t)[eB+μeω], (3.4)

which completes the system of AMHD eqs. (3.3) and (3.4) for the two functions, and .

Note that we assume the rigid rotation of NS when the dynamo term in the Faraday eq. (3.4) vanishes. Indeed, since one neglects a random velocity , the fluid velocity coincides with the NS rotation as a whole. For rigid rotation , the differential rotation is absent, and . Thus, the first dynamo term in the Faraday eq. (3.4) vanishes. In addition, we neglect below a small vorticity term both in eq. (3.4) and in eq. (3.3). The validity of this approximation is discussed in section 4. For such a case we do not need to involve the Navier-Stokes equation, while the instability of the magnetic field is presented through the last term in eq. (3.4). Here we assume also a non-uniform chiral imbalance for the stationary non-uniform electron density in NS,  [26], where is the central (neutron) density and is the electron abundance.

The formal solution of eq. (3.3) without the small vorticity term takes the form,

 μ5(x,t)=μ5(t0)e−A(t)+(π2μ2e)e−A(t)∫tt0dt′e+A(t′)[−∇⋅S+2αemπσB⋅(∇×B)]. (3.5)

In eq. (3.5), the index in the exponent reads

 A(x,t)=Γf∫tt0dt′(1+B2(x,t′)B20), (3.6)

which is huge333 is the rate of spin-flip [12]. for times .

It means that the initial imbalance vanishes rapidly in eq. (3.5). However, the main question in this paper whether survives in SN or vanishes at all requires to solve, in a self-consistent way, the Faraday eq. (3.4), which depends on the mean spin through the second integral term for in eq. (3.5).

## 4 Comparison of the contributions of the CVE and the CME

In this section, we compare the contributions of the CVE and the CME to the system of AMHD equations.

The total vector current in the Maxwell eq. (2.4), where , includes also the sum of the CME and the CVE [19], . Here is the anomalous current that drives the CME and is given by the fluid vorticity  [13]. For we substituted and , which stems from the sum , where .

Substituting for the NS rotation as a whole with the velocity , where is the NS rotation frequency444We use units for which ., one can compare the values of and . We find the CME prevails over the CVE, since . Substituting , one gets

 eB≫13.2×10−17MeV2=66×10−4G,orB≫10−2G, (4.1)

where the known relations and are used.

Let us stress that both and are the pure vector currents since is the pseudoscalar, , with respect to the space inversion. There is the additional axial current given by the vorticity  [27, 24],

 jA=[T26+12π2(μ2e+μ25)]ω≈μ2e2π2ω, (4.2)

which can not compete with the axial mean spin (the CSE) that contributes to the density imbalance evolution in eq. (2.5). Indeed, if we substitute a fast NS rotation frequency into eq. (4.2), one obtains that , comparing eqs. (2.7) and (4.2), due to the same inequality in eq. (4.1), , or .

Thus, the CVE is negligible for NS plasma, as well as the axial current in eq. (4.2), in comparison with the mean spin (the CSE) contribution to eq. (2.5).

## 5 Saturation regime

In this section, we explore the saturation regime of the chiral imbalance in the magnetic field evolution. It should be noted that the saturation of the -parameter in the solar dynamo was discussed previously. It leads to the algebraic quenching [28].

From eq. (3.3) one gets that, in the case , the stationary reads

 μ(sat)5(x)= (π2μ2eΓf)[(2αem/πσ)B⋅(∇×B)−∇⋅S]1+B2/B20+μe(ω⋅B)/eB20 (5.1) ≈(π2αemB20)[B⋅(∇×B)−(πσ2αem)∇⋅S],

where we use the inequalities in the denominator of the first line in eq. (5.1) neglecting there the vorticity term . Of course, we assume here that the magnetic field reaches the saturation as well, , that should be arisen earlier than settles.

For the negative derivative , the term in eq. (5.1) equals to

 −∇⋅S=+eBr2π2dμe(r)dr=eBr2μ2edne(r)dr≃−10−17(rRNS)(BrMeV2)MeV4, (5.2)

where is the radial component of the magnetic field. In eq. (5.2), we substitute the mean spin from eq. (2.7) and use the electron density profile555We substitute in eq.(5.2), where is the central (neutron) density and is the electron abundance, as well as put for . , proposed in ref. [26].

Substituting the spin term from eq. (5.2) into eq. (5.1), one can estimate the value just below the NS crust, ,

 μ(sat)5(RNS,θ)eV=π2σeV (2αem)2B20[eBr(RNS,θ)2μ2ednedr+2αemπσB⋅(∇×B)] =−3×10−4(Br(RNS,θ)MeV2)+1.4×10−15RNS(B⋅(∇×B)MeV4). (5.3)

Note that the non-uniform remains parity odd with respect to the space inversion analogously to in a uniform matter: , since components , , like the diffusion term in eq. (5), are pseudoscalars for the total axial vector , which is parity even accounting for the unit vectors , . In what follows, we neglect the small diffusion term in eq. (5) since it does not change issues in Parker’s dynamo model [17] applied below in section 5.2.

### 5.1 Is the anomaly saturation in eq. (5) is sufficient to drive the growth of the magnetic field?

First, we mention that it is very difficult to solve the dynamical system of eqs. (3.3) and (3.4). Assuming instead the saturation solution of eq. (3.3) given by eq. (5), where we omit the diffusion term, and substituting such into Faraday eq. (3.4) for the rigid NS rotation, one gets the nonlinear Faraday equation,

 ∂B∂t=1σ∇2B+2αemπσ∇×μ(sat)5(Br)B, (5.4)

where should vary somehow together with .

In order to feel how efficient is the instability driven by the modified , we could assume changing (probe) radial components , or independently of the latitude. This could be implemented, e.g., for the one-dimensional field [29], , where , that is similar to the Chern-Simons wave. The simplified solution of eq. (5.4) rewritten in the standard MHD form,

 ∂B∂t=η∇2B+α(t)∇×B, (5.5)

where , for changing , takes the following form:

 B(t,k,Br)=B0exp[∫tt0(±α(Br)k−ηk2)dt′]⟹B0exp[α2(Br)4η(t−t0)]. (5.6)

Such a solution corresponds to the -dynamo amplification of the initial amplitude valid for the extreme wave number in the interval where the field in eq. (5.6) is growing.

One can easily find from eq. (5) the corresponding wave number

 kext=(|Br|MeV2)114cm. (5.7)

The scale diminishes from at down to for growing . The index in the exponent, entering eq. (5.6), becomes bigger than unity at times for the magnetic field close to the Schwinger value, and for times if the field rises up to observed in magnetars [8].

For extreme , the -dynamo mechanism is efficient starting even from the early NS age . Of course, for we should include the saturation factor presented in the first line in eq. (5.1).

### 5.2 Parker dynamo for the 3D axially symmetric field

To proceed in the analysis of eq. (5.4), we use the spherical coordinates system, that is natural if one deals with a magnetic field in NS. Then, following ref. [30, pg. 373], we decompose the magnetic field into the toroidal and the poloidal components: . Moreover, we introduce the vector potential for .

Using in eq. (5.4) the helicity parameter given by the chiral imbalance in eq. (5), one can get for the toroidal component the following nonlinear equation:

 ∂B∂t=η(∇2B−Br2sin2θ)+1r∂∂r(rαsatBθ)−1r∂∂θ(αsatBr) =η(1r∂2(rB)∂r2+1r2∂∂θ[1sinθ∂(sinθB)∂θ]) −1r∂∂r[αsat∂∂r(rA)]−1r∂∂θ[αsatrsinθ∂∂θ(sinθA)], (5.8)

where is the magnetic diffusion coefficient. Here we do not take into account the turbulence contribution to ; cf. ref. [30, pg. 370]. In the last line of eq. (5.2), we substitute the poloidal components for the axisymmetric field given by the azimuthal potential ,

 Bθ=−1r∂∂r(rA),Br=1rsinθ∂∂θ(sinθA), (5.9)

where .

Finally the dynamo in AMHD is completed by the equation for the azimuthal potential ,

 ∂A∂t=η(1r∂2(rA)∂r2+1r2∂∂θ[1sinθ∂(sinθA)∂θ])+αsatB. (5.10)

Equations (5.2) and (5.10) should be further analyzed and solved numerically.

### 5.3 Low mode approximation in the thin layer r≃RNS

To factorize the -dependence in and , we decompose them using the orthogonal functions [31],

 A(t,θ)=a1(t)sinθ+a2(t)sin3θ+…,B(t,θ)=b1(t)sin2θ+b2(t)sin4θ+…, (5.11)

obeying

 ∫π0dθsinθsin3θ=∫π0dθsin2θsin4θ=0, (5.12)

and normalized as

 2π∫π0dθsin2(lθ)=1,l=1,2,…. (5.13)

Neglecting the radial dependence in a thin layer close to the NS crust, , analogously to ref. [17], we simulate the derivative with respect to as . Thus the parameter simulates a layer width where the magnetic field changes significantly.

Then, we define the dimensionless time , where the diffusion time exceeds the Universe age because of the huge electric conductivity in NS, , at the temperature . Note that the conductivity depends on the temperature  [20], where . NS cools down as [21],

 TT9=(tt9)−1/6, (5.14)

where is the time after the SN explosion when the NS temperature drops down to . Hence the conductivity changes over time as where and . Using figure 7 in ref. [32] for NS with the mass , one finds that NS cools down to the temperature during .

On this way, multiplying consistently eq. (5.10) by and , substituting given by eq. (5) without diffusion term, we can rewrite eq. (5.10) as

 ∂A∂τ=F−1{−μ2A+∂2A∂θ2−Asin2θ+cotθ∂A∂θ −3×107αemπ[∂A∂θ+cotθA](BMeV2)},

where, in the last term, we substitute the number since in units . Note that this nonlinear term in eq. (5.3) is produced by the mean spin term . It is new compared to that we deal with in the standard MHD.

Substituting the decomposition in eq. (5.11), multiplying by or , and then integrating over the latitude angle , one obtains the system in low mode dynamo approach,

 ˙a1= −1F{(μ2+2)a1+2a2 +2κ[(a1−a2)b1+2a2b2]1+K[(b21+b22)+μ2(a21+a22)+4(a21+5a22+2a1a2)]/F}, ˙a2= −1F{(μ2+12)a2 (5.16)

where , , and the azimuthal field components are normalized on . Note that, in eq. (5.3), we account for the quenching factor in eq. (5.1).

Analogously, multiplying consistently the Faraday eq. (5.2) by and , substituting , given by eq. (5), we rewrite eq. (5.2) as

 (1MeV2)∂Bdτ= F−1{(1MeV2)[−μ2B+∂2B∂θ2+cotθ∂B∂θ−Bsin2θ] +2κ(∂A∂θ+cotθA)[−μ2A+2(∂2A∂θ2+cotθ∂A∂θ−Asin2θ)]}, (5.17)

where the azimuthal potential is normalized on . Thus the components in the decomposition in eq. (5.11) are dimensionless as there. Let us stress again that the second line results from the mean spin contribution .

Multiplying eq.(5.3) by or , integrating over latitude angle , we complete the system in eq. (5.3) by

 ˙b1= −1F{(μ2+6)b1+4b2 +2κ[a21(μ2+4)+24a1a2+a22(μ2+20)]1+K[(b21+b22)+μ2(a21+a22)+4(a21+5a22+2a1a2)]/F}, ˙b2= −1F{(μ2+20)b2 +2κ[a1a2(3μ2+32)+a22(μ2+32)]1+K[(b21+b22)+μ2(a21+a22)+4(a21+5a22+2a1a2)]/F}, (5.18)

where take into account the the quenching factor analogously to eq. (5.3). Thus, we have the four ordinary differential eqs. (5.3) and (5.3) for the four amplitudes and .

#### 5.3.1 Initial conditions

To get the initial conditions for eqs. (5.3) and (5.3), we equate the initial normalized components

 Bθ(RNS,θ,t=0)MeV2= −[a1(0)sinθ+a2(0)sin3θ], Br(RNS,θ,t=0)MeV2= 2cosθ[a1(0)+a2(0)(4cos2θ−1)], (5.19)

that result from eq. (5.9) in the low mode approximation in eq. (5.11), . Then one can find at the same force line of the poloidal field , while at different latitudes when substituting corresponding for where and