Table Of ContentThree-dimensional simulations of the magnetic stress in a neutron star crust
∗
T. S. Wood
School of Mathematics and Statistics, Newcastle University,
Newcastle upon Tyne, NE1 7RU, United Kingdom
R. Hollerbach
Department of Applied Mathematics, University of Leeds, Leeds, LS2 9JT, United Kingdom
(Dated: June15, 2015)
Wepresentthefirstfullyself-consistentthree-dimensionalmodelofaneutronstar’smagneticfield,
5 generated by electric currents in the star’s crust via the Hall effect. We find that the global-scale
1 field converges to a dipolar Hall-attractor state, as seen in recent axisymmetric models, but that
0 small-scalefeaturesinthemagneticfieldsurviveevenonmuchlongertimescales. Thesesmall-scale
2 featurespropagatetowardthedipoleequator,wherethecrustalelectriccurrentsorganizethemselves
into a strong equatorial jet. By calculating the distribution of magnetic stresses in the crust, we
n predict that neutron stars with fields stronger than 1014G can still be subject to starquakes more
u than 105yr after their formation.
J
2
1 I. INTRODUCTION However,because ofdifficulties solvingthe EMHDequa-
tions computationally in 3D, in all of these studies the
]
R Neutron stars are of interest not only for the exotic magnetic field was assumed to be axisymmetric. Al-
S statesofmattertheycontain,butalsofortheirmagnetic though there have been 3D studies of EMHD in simpli-
. fields, which are the strongest in the universe. In the fiedperiodic-boxgeometry[20–22],theresultscannotbe
h
case of pulsars, the rotation of the magnetic field pro- usedtopredicttheglobal-scalefieldmorphologyinareal
p
- duces beams of non-thermal radiation that can be de- neutron star. Currently, the only predictions regarding
o tected thousands of light-years away. The field can also the global-scale field come from axisymmetric models.
tr affect the dynamics of the star itself, through the mag- Using one such axisymmetric model, Gourgouliatos
s netic stress it exerts on the star’s solid outer layer, or and Cumming [18] found that an initially dipolar mag-
a
[ “crust”. Inneutronstarswithespeciallystrongmagnetic netic field evolves towards a quasi-steady configuration
fields (known as magnetars) the magnetic stress can be that they called a “Hall attractor”. After ∼ 105yr the
2 strong enough to fracture the crust, producing a star- electriccurrentsinthecrustareconcentratedinanarrow
v
quake [1]. This is believed to be the mechanism behind jetaroundthe dipole equator,producingastrongbeltof
9
4 the gamma-ray flares and X-ray outbursts detected in poloidal field. However, it is unknown whether this ax-
1 these objects [2]. isymmetricattractorwouldbestableto3Dperturbations
5 Themagneticfields ofpulsars,inferredfromspindown in a real neutron star.
0 measurements,canbeanywhereintherange108–1015G. Here,wepresentforthefirsttimeglobal3Dnumerical
1. However, the spindown rate depends only on the large- simulations of the magnetic field in a neutron star crust.
0 scale component of the field at the magnetic poles, so We use a pseudo-spectralcode that allowsus to perform
5 these measurements may underestimate the actual field simulations that are not only fully 3D, but also higher
1 strength in the star. This could explain why some neu- resolutionthananyofthosepreviouslypresentedevenin
v: tronstarsproducemagnetar-likeemissionsyethaveslow 2D.WedescribehowtheHallattractorismodifiedin3D,
i spin-downrates[3–5]. Thereisalsoevidencethatatleast and discuss the magnitude and distribution of magnetic
X
some neutron stars have stronger magnetic features on stresses within the crust.
r
smaller scales [6–9]. To interpret these observations it is
a
necessary to develop a self-consistent, three-dimensional
II. THE MODEL
(3D) model of neutron star magnetic fields.
The external field of a neutron star is generated by
electric currents flowing within its crust and core. The We work in a reference frame that corotates with the
ions in the crust form a rigid lattice, and the currents neutron star’s crust. Because the ions are fixed within
there arise purely through the flow of electrons, whose the crust, in this frame the electric current, J, depends
dynamics depend primarily on the Hall effect [10]; this only onthe electronfluid velocity, v. The currentis also
situation is commonly referred to as electron magneto- directly related to the magnetic field, B, via Amp`ere’s
hydrodynamics (EMHD). Recently there have been nu- Law, and we have (in Gaussian cgs units)
merous studies of EMHD in neutron star crusts [11–19].
c
J=−env= ∇×B, (1)
4π
where n is the electron number density, e is the elemen-
∗ [email protected] tary charge, and c is the speed of light. If the crust is
2
either sufficiently cool, or close to isothermal, then the where we have assumed that the characteristic length-
magnetic field is frozen to the electron fluid [10], and scale for the magnetic field is R=10km. If the field has
evolves according to the equation strong, small-scale features, then these will evolve on a
shorter timescale.
∂B cB
=∇× ×(∇×B) −∇×(η∇×B), (2) To solve the EMHD equation(2) we haveadapted the
∂t (cid:18)4πen (cid:19) 3DMHDcodePARODY,developedbyDormyet al.[24]
and Aubert et al. [25]. The code is pseudo-spectral, and
whereηisthemagneticdiffusivity. Thetwotermsonthe
uses spherical harmonic expansions in latitude and lon-
right-hand side of Equation (2) represent the Hall effect
gitude, and a discrete grid in radius, making it perfectly
and Ohmic dissipation, respectively.
suited to solving problems in spherical-shell geometry.
Inthisstudyweapproximatethestructureofthecrust
The resultsthatwe presenthere havea resolutionof128
as fixed and spherically symmetric, and so the electron
grid-points in radius, and spherical harmonics up to de-
density n and magnetic diffusivity η are fixed functions
gree l = 100. We have benchmarked the code against
of radius r. We do not therefore take account of any
previouslypublishedaxisymmetricresults[12,23,26]and
deformationsin the crustarisingfrom magnetic stresses.
find excellent agreement in all cases.
However,wecanuseourresultstocalculatethemagnetic
stress and determine whether it would be large enough
to induce fractures in the crust. We take the top and
bottom of the crust to be the spherical surfaces r = R III. RESULTS
andr=0.9Rrespectively,whereR=10kmistheradius
of a typical neutron star. For n(r) and η(r) we adopt A. Robustness of the Hall attractor
the same analytical profiles used by Gourgouliatos and
Cumming [23]. Specifically, Todeterminewhetherthe“Hallattractor”seeninear-
4 −8/3 lier axisymmetric simulations is robust against 3D per-
1−r/R 1−r/R
n=n0 1+ andη =η0 1+ , turbations,wehaverepeatedonesimulationofGourgou-
(cid:18) 0.0463 (cid:19) (cid:18) 0.0463 (cid:19) liatos and Cumming [23], using the same boundary con-
where n0 =2.5×1034cm−3 and η0 =4.0×10−4cm2s−1 ditions anddipolar “Ohmic eigenmode” initial condition
fortheaxisymmetriccomponentofthemagneticfield. To
are the values at r = R. These profiles are only rough
this initial axisymmetric poloidal field we add a low am-
approximations to the (highly uncertain) profiles in real
plitude, small-scale 3D perturbation with both poloidal
neutron stars, but fortunately our results are not sensi-
and toroidal components, comprising spherical harmon-
tive to the specific profiles used.
icsofdegree206l 640. AsshowninFigure1,themag-
Finally, we must impose boundary conditions on the
netic field evolves towards a state broadly similar to the
magnetic field at the top and bottom boundaries of the
axisymmetric Hall attractor. However, significant non-
crust. At the top we impose vacuum boundary condi-
axisymmetric features persist in the simulation even on
tions, i.e., we match to a current-free field outside the
long timescales (comparable to the global Hall timescale
star. At the bottom we impose either the idealized
tHall). These features would rapidly decay by the action
boundary conditions used by Gourgouliatos and Cum-
of Ohmic dissipation alone, so their longevity can be at-
ming [23], which are B = 0 and J = 0, or the more
r r tributed to the Hall effect. In factit has previously been
realistic boundary conditions used by Hollerbach and
shownthat, in the presence of a strong density gradient,
Ru¨diger [12], which model the star’s core as a type-I su-
theHalleffectcansustainorevenamplifysmall-scalefea-
perconductor. From here on we refer to these two sets
tures in the magnetic field [27]. The length-scale of the
of conditions as GC and HR respectively. In order to
longest-livednon-axisymmetricfeaturesiscomparableto
implement the HR boundary conditions, it is convenient
the largest scales present in the initial 3D perturbations
tomakeasmallmodificationtotheelectrondensitypro-
(degreel=20),suggestingthattheyformbyanup-scale
file n(r), asdescribedby HollerbachandRu¨diger[12], to
transfer of magnetic energy.
make 1/n vanish at the bottom of the crust. We there-
The convergence of the global-scale magnetic field to
fore use a slightly modified density profile n˜(r), defined
theattractorstateoccursontheHalltimescale(3),after
as 1/n˜(r)=1/n(r)−1/n(0.9R).
which the evolution proceeds on the slower timescale of
The relative importance of the Hall effect and Ohmic
Ohmic diffusion. This is in agreement with results from
dissipation terms in Equation (2) depends on the
axisymmetricsimulations [12, 18]. (In the simulations of
strength of the magnetic field, B0 say, and is quantified Marchantet al. [19], by contrast,the fieldwasinitialized
by the Hallparameter,H ≡cB0/(4πen0η0). In a typical
magnetar, with B0 = 1014G, we have H ≃ 50, implying alreadyinaHallequilibriumstate,andthereforeevolved
on the Ohmic timescale throughout.) The poloidalmag-
that the Hall effect dominates the dynamics of the mag-
neticfluxbecomesconcentratedintoa“belt”aroundthe
netic field. In that case we expect the field to evolve on
equator of the magnetic dipole, implying a strong equa-
the Hall timescale,
torial jet of electrons within the crust. The angular ve-
4πen0R2 locity in this jet is approximately constant along each
tHall ≡ cB0 ≃1.6Myr, (3) poloidalfieldline(seeFigure1). Wedonotfindthatthe
3
(a) t=0.01Myr t=0.54Myr t=1.77Myr (b) t=0.11Myr Br/1014G
1.2 1.2 1.2 0.6
0.4
1.0 1.0 1.0 0.2
-0.0
-0.2
0.8 0.8 0.8 -0.4
-0.6
0.6 0.6 0.6 t=1.10Myr Br/1 0 14G
0.15
0.4 0.4 0.4 0.10
0.05
0.2 0.2 0.2 -0.00
-0.05
0.0 0.0 0.0 -0.10
-0.15
FIG. 1. (a) Lines of the azimuthally averaged poloidal field and electron angular velocity at successive times, illustrating
convergence to the Hall attractor. The dashed lines indicate the boundaries of the crust. (b) The radial component of the
magnetic field at thesurface at early and late times in the same simulation.
Halleffectenhancesthediffusionofthemagneticfield,as 1.2 Bφ/1014G1.2
0.0
originally suggested by Goldreich and Reisenegger [10]. 1.0 1.0
However, by concentrating the poloidal flux around the -0.2
0.8 0.8
equator,the Halleffect doesleadtoa morerapiddecline
-0.4
in the strength of the radial field at the poles. 0.6 0.6
Themagneticfieldisprimarilypoloidal,withaweaker -0.6
0.4 0.4
toroidalcomponent in the inner crust. However,the GC
-0.8
boundary conditions used in this simulation impose zero 0.2 0.2
toroidal field both at the top and bottom of the crust. 0.0 -1.0 0.0
In a real neutron star, currents can flow between the 0.0 0.2 0.4 0.6 0.8 1.0 1.2 0.0 0.2 0.4 0.6 0.8 1.0 1.2
crust and core, generating significant toroidal fields. We FIG. 2. Poloidal fieldlines and contours of toroidal field at
t = 1.2Myr. The left and right panels show the cases with
have therefore repeated this simulation using the more
GC and HR boundary conditions, respectively.
realistic HR boundary conditions, which allow for finite
toroidal field at the lower boundary. Figure 2 compares
theazimuthallyaveragedmagneticfieldatthesametime
surface. The Hall effect may therefore be even more sig-
in both simulations, revealing that although the mor-
nificantin neutronstarsthanpreviouslythought,andso
phology of the poloidal field is modified, a similar Hall
we haverepeated the same simulationshownin Figure 3
attractor state still exists. The toroidal field in the sec-
with a stronger initial magnetic field (i.e., a larger Hall
ond simulation is much stronger, as expected, but re-
parameter H). Because the evolution of the global-scale
mains weaker than the poloidal field. This toroidal field
magneticfieldoccursontheHalltimescale(3),weexpect
is associated with a meridional flow of electrons that is
that a stronger field will converge more rapidly to the
equatorwardatthe star’ssurface,anddragsthe poloidal
Hall attractor, before the small-scale magnetic features
magnetic field lines. As a result, the equatorial belt of
have been dissipated by resistivity. This is confirmed in
poloidal flux, and the corresponding electron jet, is even
Figure 4, which shows the surface radial field in three
strongerinthesimulationwithHRboundaryconditions,
simulations with increasing magnetic field strengths. In
and is pushed slightly deeper into the crust. Figure 3
each case the field is plotted at time t = 0.6tHall, which
illustrates the development of the equatorial belt, which
correspondsto t≃1, 0.5, and 0.25Myrrespectively, and
manifests as bands of positive and negative B on either
r an equatorial jet of electrons has already formed. In the
side of the equator at the star’s surface.
case with the strongestmagnetic field, the jet is broader
andmorespatiallydisordered;atlater times inthe same
simulation,the jet becomes increasinglylaminar, but re-
B. The effect of field strength mains broader than in the other simulations.
Spin-down measurements of neutron stars tell us only
the strength of the large-scale magnetic field near the C. The crustal magnetic stress
poles. However, it is clear from Figure 3 that the field
strength elsewhere in the crust can greatly exceed this If the magnetic shear stress within the crust exceeds
observedvalue. At the time plotted in Figure 2b, for ex- the breakingstressof the ionic lattice then it will induce
ample, the strength of the surface magnetic field varies a crustal fracture. Molecular dynamics models indicate
between≃5×1012Gatthedipole axisand≃3×1013G that the breaking stress is approximately 5 × 10−3 of
nearthe equator. Near the bottom ofthe crust,the field the electron pressure, Pe [28]. At each point within the
is stronger by a further order of magnitude than at the crust, the strongest magnetic shear stress is exerted on
4
t=0.21Myr t=0.67Myr t=1.40Myr
FIG. 3. Evolution of the magnetic field in the simulation with HR boundary conditions. Top row: B at r=R. Bottom row:
r
B at r=0.9R. Colorbars are adjusted to themaximum values in each plot, and use a linear scale.
φ
Br/B0
0.3
0.2
0.1
0.0
-0.1
-0.2
-0.3
FIG. 4. The surface radial field in threesimulations with B0 =1014G, 2×1014G, 4×1014G (i.e., H =50,100,200).
a surface that makes an angle of 45◦ to the local mag- 1.2
netic field direction, and this stress is exactly equal to
the magnetic pressure, Pm. The condition for fracturing 1.0
is therefore Pm & 5 × 10−3Pe, where Pm = |B|2/(8π)
and Pe ≃ (3π2n)4/3ℏc/(12π2). Fractures are most likely
nearthesurfaceofthecrust,wherethedensityandpres- 0.8
surearelowest,andtheenergyreleasedinanear-surface
fracture can directly power flares and outbursts. In Fig-
0.6
ure 5 we plot the ratio of magnetic pressure to breaking
stress at the surface of the crust in the same simulation
showninthelastpanelofFigure4,revealingthatpatches 0.4
aroundthedipole equatoraresusceptibletocrustalfrac-
turing. These patches are localized in both latitude and
0.2
longitude, as a consequence of the three-dimensionality
ofthe magnetic field. They persistonlong timescales,of
ordertHall,andpropagateazimuthallyinthedirectionof 0.0
the electron jet.
FIG. 5. Lines of the magnetic field generated by currents
in the crust. The coloring indicates the ratio of magnetic
pressure to breaking stress at the surface of the crust. A
IV. DISCUSSION
portion of the surface has been cut away to show fieldlines
inside thecrust, whose lower boundary is the gray sphere.
Our results indicate that neutron starquakes are most
common in the vicinity of the electron jet within the
crust,whichtypicallyformsaringaroundthedipoleaxis. equatorward as a result of the toroidal magnetic field in
This jet forms on the Hall timescale, which is dependent thecrustandtheassociatedmeridionalflowofelectrons.
on the overall strength of the magnetic field. At earlier This migration could potentially be reversed if a strong
times the structureofthe magnetic fieldis dependent on toroidal field of the opposite sign were present, corre-
the initial conditions for the proto-neutron star, which sponding to a poleward flow of electrons at the star’s
unfortunately are not well known. The surface field is surface [29]. There is no obvious process that can gen-
strongestinlocalizedpatchesaroundthe magneticequa- erate such a strong toroidal field within the crust, but
tor, where the local field strength typically exceeds that it could be a remnant from the star’s formation [30] or
at the poles by an order of magnitude. the result of toroidal flux expulsion from the core [31].
In our simulations, surface magnetic features migrate To fully describe the processes that can impart a strong
5
toroidal field to the crust, it will be necessary to couple ACKNOWLEDGMENTS
our crust model to a model of the superconducting core.
The flow of heat from a cooling neutron star is inhib-
ited across magnetic field lines, so the strong concentra-
tionofpoloidalmagneticfluxintheequatorialbeltcould T.S.W. and R.H. were supported by STFC grant
trapheatwithintheelectronjet[32]. Thebreakingstress ST/K000853/1. ComputationswereperformedonARC1
is very sensitive to temperature [28], and so the equato- andARC2,partofthe HighPerformanceComputing fa-
rial region of the crust may be even more susceptible to cilities at the University of Leeds. We thank Konstanti-
fracturing than we have found here. We are currently nos Gourgouliatos, C´eline Guervilly, David Ian Jones,
extending our model to describe the thermal–magnetic MaximLyutikov,andthe anonymousrefereesforhelpful
evolution of the crust, and its coupling to the core. comments and suggestions.
[1] M. Ruderman,Astrophys.J. 382, 587 (1991). [17] J. A. Pons, D. Vigan`o, and N. Rea,
[2] C. Thompson and R. C. Duncan, Nature Physics 9, 431 (2013).
Astrophys.J. 561, 980 (2001). [18] K. N. Gourgouliatos and A. Cumming,
[3] F. P. Gavriil, M. E. Gonzalez, E. V. Gotthelf, Phys. Rev.Lett. 112, 171101 (2014).
V. M. Kaspi, M. A. Livingstone, and P. M. Woods, [19] P. Marchant,A.Reisenegger, J.Alejandro Valdivia,and
Science 319, 1802 (2008). J. H.Hoyos, Astrophys.J. 796, 94 (2014).
[4] N. Rea, P. Esposito, R. Turolla, G. L. Israel, S. Zane, [20] D. Biskamp, E. Schwarz, A. Zeiler, A. Celani, and J. F.
L. Stella, S. Mereghetti, A. Tiengo, D. G¨otz, E. G¨o˘gu¨¸s, Drake, Phys. Plasmas 6, 751 (1999).
and C. Kouveliotou, Science330, 944 (2010). [21] J. Cho and A. Lazarian, Astrophys.J. 615, L41 (2004).
[5] P. Scholz, V. M. Kaspi, and A. Cumming, [22] C. J. Wareing and R. Hollerbach,
Astrophys.J. 786, 62 (2014). J. Plasma Phys. 76, 117 (2010).
[6] T. Gu¨ver, E. G¨oˇgu¨¸s, and F. O¨zel, [23] K. N. Gourgouliatos and A. Cumming,
Mon. Not. R.Astron. Soc. 418, 2773 (2011). Mon. Not. R.Astron. Soc. 438, 1618 (2014).
[7] E. V. Gotthelf, J. P. Halpern, and J. Alford, [24] E. Dormy, P. Cardin, and D. Jault,
Astrophys.J. 765, 58 (2013). Earth Planet. Sci. Lett. 160, 15 (1998).
[8] A. Tiengo, P. Esposito, S. Mereghetti, R. Turolla, [25] J. Aubert, J. Aurnou, and J. Wicht,
L. Nobili, F. Gastaldello, D. G¨otz, G. L. Israel, Geophys. J. Int.172, 945 (2008).
N. Rea, L. Stella, S. Zane, and G. F. Bignami, [26] R. Hollerbach and G. Ru¨diger,
Nature500, 312 (2013). Mon. Not. R.Astron. Soc. 337, 216 (2002).
[9] U. Geppert, J. Gil, and G. Melikidze, [27] T. S. Wood, R. Hollerbach, and M. Lyutikov,
Mon. Not. R.Astron. Soc. 435, 3262 (2013). Phys. Plasmas 21, 052110 (2014).
[10] P. Goldreich and A. Reisenegger, [28] A. I. Chugunov and C. J. Horowitz,
Astrophys.J. 395, 250 (1992). Mon. Not. R.Astron. Soc. 407, L54 (2010).
[11] D.A.ShalybkovandV.A.Urpin,A&A321,685(1997). [29] U. Geppert and D. Vigan`o,
[12] R. Hollerbach and G. Ru¨diger, Mon. Not. R.Astron. Soc. 444, 3198 (2014).
Mon. Not. R.Astron. Soc. 347, 1273 (2004). [30] C. Thompson and R. C. Duncan,
[13] J. A. Pons and U.Geppert, A&A 470, 303 (2007). Astrophys. J. 408, 194 (1993).
[14] R.PernaandJ.A.Pons,Astrophys.J. 727, L51 (2011). [31] S.K.Lander,Mon. Not. R. Astron.Soc. 437, 424 (2014).
[15] Y. Kojima and S. Kisaka, [32] D. Vigan`o, N. Rea, J. A. Pons, R. Perna,
Mon. Not. R.Astron. Soc. 421, 2722 (2012). D. N. Aguilera, and J. A. Miralles,
[16] D. Vigan`o, J. A. Pons, and J. A. Miralles, Mon. Not. R.Astron. Soc. 434, 123 (2013).
Computer PhysicsCommunications 183, 2042 (2012).
