A Shear-density instability in neutron star crusts

Magnetar activity mediated by plastic deformations of neutron star crust


We advance a ”Solar flare” model of magnetar activity, whereas a slow evolution of the magnetic field in the upper crust, driven by electron MHD (EMHD) flows, twists the external magnetic flux tubes, producing persistent emission, bursts and flares. At the same time the neutron star crust plastically relieves the imposed magnetic field stress, limiting the strain to values well below the critical strain of a brittle fracture, .

Magnetar-like behavior, occurring near the magnetic equator, takes place in all neutron stars, but to a different extent. The persistent luminosity is proportional to cubic power of the magnetic field (at a given age), and hence is hardly observable in most rotationally powered neutron stars. Giant flares can occur only if the magnetic field exceeds some threshold value, while smaller bursts and flares may take place in relatively small magnetic fields.

Bursts and flares are magnetospheric reconnection events that launch Alfvén shocks which convert into high frequency whistlers upon hitting the neutron star surface. The resulting whistler pulse induces a strain that increases with depth both due to the increasing electron density (and the resulting slowing of the waves), and due to the increasing coherence of a whistler pulse with depth. The whistler pulse is dissipated on a time scale of approximately a day at shallow depths corresponding to ; this energy is detected as enhanced post-flare surface emission.

1 Two competing models of magnetar bursts and flares

Two closely related classes of young neutron stars – Anomalous X-ray Pulsars (AXPs) and the Soft Gamma-ray Repeaters (SGRs) – both show X-ray flares and quiescent X-ray emission (see Woods & Thompson, 2006, for review). Their high energy emission is powered by dissipation of super-strong magnetic fields, G (Thompson & Duncan, 1995).

Two models of magnetar bursts and flares have been proposed. In the first, a star-quake model, a flare relies on a sudden fracture of the neutron star crust that lead to fast, on time scale of hundred milliseconds, untwisting of the internal magnetic field and a resulting twisting-up of the external magnetic field (Thompson & Duncan, 1995). Alternatively, a Solar paradigm postulates that slow evolution of the internal magnetic field leads to a gradual twisting of magnetospheric field lines on a time scale much longer than that of the giant flare (Lyutikov, 2003, 2006). Eventually, the increasing twist associated with the current-carrying magnetic field in the magnetosphere lead to a sudden relaxation of the external magnetic field, accompanied by associated dissipation and change of magnetic topology.

Perhaps the best observational argument in favor of external magnetic dissipation during the burst is the observed sharp rise of -ray flux during giant flares, on the time scale similar to the Alfvén crossing time of the inner magnetosphere, msec (Palmer et al., 2005). This points to the magnetospheric origin of giant flares (Lyutikov, 2006).

The efficient crack formation, needed in the star-quake model, requires the solid’s ability to form a small void at the crack’s location. Since in the neutron-star crust the pressure is greater than the shear modulus by two orders of magnitude, the conventional crack that relies on a formation of the void cannot occur (Jones, 2003) (though so called deep Earth quakes Frohlich, 2006, apparently violate this condition). Secondly, even if elastic properties of the crust did allow brittle cracking, the energy release would be strongly suppressed by the magnetic tension (even though the magnetic field itself leads to cracking!). Levin & Lyutikov (2012) demonstrated that the energy release from a thin crack would be strongly suppressed, since the magnetic field provides mechanical connection between the two slipping sides of the crack and rapidly suppresses the slippage. The key point is that magnetically induced cracks form orthogonally to the magnetic field; then the tension strongly suppresses the slippage and the energy release. Thin cracks release energy via magnetic-field diffusion too slowly to be able to contribute to the energetics of magnetar flares. In addition, the waves generated within neutron star experience many internal reflections due to impedance mismatch. As a result, the rise time of the emission is minutes to hours, too long account for millisecond rise time of bursts and flares (Link, 2013). Thus, even if mechanical properties of the crust allowed brittle fracture, the resulting energy release on the time scale of flares is very small.

2 Evolution of magnetic field in neutron stars

2.1 Generation and initial relaxation of magnetic field

We envision the following paradigm for evolution of magnetic field in neutron star crusts, and specifically in magnetars. Strong magnetic fields, of the order of G in case of magnetars, may be created by a dynamo mechanism, e.g.,  of the type, operating at birth of neutron stars (Thompson & Duncan, 1993). It is expected that the dynamo operates most efficiently in the outer layers of the proto-neutron star, where neutrino-driven turbulence is most efficient (Colgate & Fryer, 1995; Herant et al., 1994).

Typically, neutrino driven turbulence dies out on time scales of seconds, well before the crust solidifies (though some fluid motion, like a shearing rotation, may persist for longer times). As a result, during the time that the star remains fluid the magnetic field reaches an MHD-type equilibrium. In MHD equilibria Lorentz forces are balanced by gradients of pressure and gravitational forces,


where is the gravitational potential. By dividing Eq. (1) by and taking a curl, we find (e.g.,  Reisenegger et al., 2007)


Stability of fluid stars requires that the configuration involves both toroidal and poloidal magnetic fields (Flowers & Ruderman, 1977; Braithwaite & Spruit, 2004). In addition, it appear that purely barotropic stars are unstable (Lander & Jones, 2012). In stable configurations the energy of the toroidal magnetic field is at least comparable to the energy of the poloidal field, and may greatly exceed it (Akgün et al., 2013). In addition, toroidal magnetic field are confined to relatively small volume near the surface at the magnetic equator. As a result, the value of the toroidal magnetic field may exceed the poloidal by as much as two orders of magnitude, reaching as high as G (Akgün et al., 2013).

2.2 Freezing and reviving the magnetic dynamics: electron MHD

At approximately 100 seconds the crust cools sufficiently and ions in the neutron star crusts form a fixed lattice, while for slow perturbations the electrons behave as an inertialess fluid. The resulting dynamical system is called electron MHD (EMHD) (Lighthill, 1960; Kingsep et al., 1987; Gordeev et al., 1994; Goldreich & Reisenegger, 1992). Within the framework of EMHD, the magnetic field is frozen into electron fluid. In the limit of infinite conductivity, the electric field then satisfies the condition


and the induction equation becomes


For infinitely rigid lattice the the electric current is produced by the flow of electrons, while ions provide a neutralizing background. In case of plastic deformations the flowing ions also provide current density, so that , where , , and and are the local number density, charge, and average velocity of the electrons. As we demonstrate below, see Eq. (24), for the chosen model of plastic deformations the ion velocity is typically much smaller than the electron velocity. Neglecting the contribution of ions to the current density, , the only dynamical variable in this equation is the magnetic field (Kingsep et al., 1987):


(density is an externally prescribed function). So, given an initial field configuration and appropriate boundary conditions, the induction equation uniquely determines its evolution.

After the crust freezes, it no longer obeys MHD, but EMHD equations. The initial field, satisfying the MHD equilibrium condition, Eq. (2), starts evolving according to Eq. (5), which initially (replacing Eq. [1] in it) takes the form


In the initial MHD equilibrium, due to the Lorentz forces and the stable, compositional stratification of the neutron-star matter (Goldreich & Reisenegger, 1992), the gradients in this equation are not exactly radial and not exactly parallel to each other, although their relative angles are generally small, (Reisenegger et al., 2007). Thus, an MHD-equilibrium magnetic field structure is, generally, not an equilibrium configuration of EMHD. Two terms will drive evolution: the non-barotropic condition and the stratification of the neutron star material, constant.

There are exception to the above (that freezing of MHD equilibrium results in non-equilibrium EMHD state): (i) if MHD configuration is force-free, (ii) if the fluid is barotropic, when and has constant mean molecular weight (or constant proton fraction in case of neutron stars), const, (iii) if the system is effectively one-dimensional, will all the quantities depending only on one particular coordinate. In these cases MHD equilibrium corresponds to EMHD equilibrium as well. We do not expect that any of the three condition listed above are satisfied in proto-neutron stars. Thus, after the crust freezes the resulting magnetic field state is, generally, not in Hall equilibrium, and a system start dynamical evolution. Unfortunately, we cannot give an explicit example of this very fundamental statement, since no analytical examples of at least two-dimensional stable MHD equilibria in stars exist.

The normal modes of EMHD plasma are whistlers. Harmonic whistlers of arbitrary amplitude are exact solution of EMHD (Lyutikov, 2013a),


where and are electron cyclotron and plasma frequencies, and is the eigenvector.

A whistler wave with fluctuating magnetic field  induces stresses in the crust of the order of . This stress is initially zero and builds up to this value on the Hall time scale (Goldreich & Reisenegger, 1992)


where is a typical scale of magnetic field fluctuation (crust thickness) and is typical density at the neutron drip point, . The Hall time scale has a steep dependence on the scale, density and magnetic field, and can vary considerably within the crust.

Thus, freezing of MHD equilibrium results in an non-equilibrium EMHD state, which revives the evolution of the magnetic field and leads to generation of whistler waves. At the moment of crust solidification, there are no shear stresses, but as the magnetic field evolves on the Hall time scale, shear stresses build-up, reaching a maximum value on the time scale of order of the crust Hall time.

The initial MHD-stable state is dominated by low multipoles (e.g.,  axisymmetric dipole-like configuration of Braithwaite & Spruit (2004)). The evolution of the magnetic field in the crusts then proceeds to a number of processes: (i) turbulent EMHD cascade, that creates higher multipole (Vainshtein, 1973; Biskamp et al., 1999; Galtier & Bhattacharjee, 2003; Cho & Lazarian, 2004; Howes et al., 2008; Lyutikov, 2013b); (ii) large scale Hall drifts that may lead to the formation of current sheets approximately in one Hall time (Kingsep et al., 1987; Gordeev et al., 1994); (iii) crustal yields that dissipate the energy of whistler waves. All this effects proceed, roughly on a Hall time scale. It is not clear at the moment which of the above is the dominant mechanism of the magnetic field evolution.

2.3 Crustal plastic properties do not determine magnetar activity

Importantly, the evolution of the magnetic field in the crust proceeds in the electron MHD (EMHD) regime, where the magnetic field is frozen in the electron fluid. This may lead to the evolution of the external fields and externally-generated flares even without any crustal yielding, be it ductile or brittle. The crust may still respond to the imposed stress (e.g.,  relieving the stress plastically, as we discuss in this paper), but it does not need to in order to produce a flare: even if the crust were absolutely rigid the crustal magnetic field would evolve (and twist the external magnetic field) due to the EMHD drift.

On the other hand, the post-burst enhanced surface thermal emission does indicate substantial energy deposition deep within the crust during the bursts and flares (Eichler et al., 2006; Scholz et al., 2012). Scholz et al. (2012) inferred the deposition density of g cm, close to the neutron drip line, the outer core boundary located approximately hundred meters below the surface. One of the points of the paper is to reconcile the external energy deposition and observation of a long-term (weeks to months) post-burst surface cooling.

3 A model problem: plastic deformations of the crust

Evolution of the magnetic field due to Hall drift of the electron fluid induces shear stresses in the ion lattice, which may yield. Thus, the crustal dynamics is described by two separate processes: Hall electron drifts and the response of the crust (which can take various forms: plastic, ductile or brittle).

In the NS crust the shear modulus increases monotonically with depth, and then quickly goes to zero once the crust dissolves. The shear modulus changes according to (Negele & Vautherin, 1973; Thompson & Duncan, 2001)


We expect that when the magnetic field energy density exceeds the shear modulus, for fields


the crust must respond plastically (but, it can also respond plastically at much smaller magnetic fields). In fact, plastic response may dominate the lattice response even for stresses much smaller than the shear modus, especially for slow strain rates.

In this paper we consider plastic response of the crust given an externally imposed Hall magnetic stresses, from the magnetic field frozen in the evolving electron fluid. The typical time scale of the electron drift, the Hall time (8), is much longer compared to the observed duration of the magnetars’ bursts and flares, so we can assume that the stress builds slowly, quasi-statically. The normal modes of the Hall plasma - whistlers - carry a shear stress, assumed to be balanced by the lattice stress. As a simplified model of this whistler-carried stress on the ion lattice, we assumed that whistlers create a simple two-dimensional current layer where the -displacement of the electron fluid evolves according to


where is the typical amplitude of the lateral displacement of the electron fluid, is a typical thickness of the current layer (realistically, ) and is a typical time of the growth of the electron stress on the ion lattice. The strain in the electron fluid then is growing linearly in time according to the law


where dot denotes time derivative. The parallel magnetic field, frozen in the electron fluid, then evolves according to


The magnetic field exerts a stress on the ion lattice, that results in an elastic strain rate


where is Alfvén velocity and is shear speed. We assume that the ion lattice responds in a viscoelastic way to the stress induced by the Hall motions of electrons (12). Thus, in addition to the elastic strain rate there is also a plastic relaxation rate.

Plastic deformations in crystals are typically divided into two categories. First, there is plasticity that occurs beyond some critical strain; below this critical stain there are no irreversible deformations. Secondly, plastic deformations can occur at small strains. This type of plastic deformations is called creep (or, related, a plastic relaxation) - in this case even small stress leads to irreversible deformations (like in a pure unseeded aluminum). Generally, all materials exhibit some viscoelastic response, while the existence of a well-defined yield stress is the exception rather than the rule (Lubliner, 2008).

Plastic creep/relaxation in crystals is induced by the motions of crystal deformations (Orowan, 1934; Gilman, 1969). Shear deformations are due to the motions of the line defects. The plastic rate of strain depends on the surface density of defects (dimensions cm), the rate of the defects’s motion and a typical distance a defect travels in one jump (the Burgers vector) via the Orowan equation


(the minus sign reflects the fact that the rate of the plastic stress, like viscosity, works as a damper.) The length of the Burgers vector is of the order of the inter-particle distance, , , where is the mass atomic number of the ion lattice.

The rates of strain of the ion lattice due to the EMHD dynamics and plastic relaxation add:


Eq. (16) describes viscoelastic deformation due to externally imposed, time-dependent strain and the plastic relaxation .

The density of the dislocations and the defect drift velocity are, generally, functions of the stress, activation energy of the defect’s movement, temperature and history of stress (e.g.,  Gilman, 1969) – determining the density of dislocation and their drift velocity is a complicated problem. Typically, at higher temperatures the density and especially the mobility of dislocations obey an Arrhenius-type scaling, , where is an activation energy. On the other hand, at small temperatures/high stresses the mobility and the density of dislocations are determined by the applied stress. In this case the density of dislocations may be estimated as (Gilman, 1969)


(this scaling uses the fact that a stress from a dislocation line varies with distance according to , so a mean stress is ).

We parametrize the dislocations’ velocity as (Gilman, 1960)


where is the critical fracture strain. (So that at the critical strain the dislocation velocity is a fraction of the shear speed).

The total rate of strain of the lattice is given by the EMHD strain (12) and the plastic strain (15)


Equation (19) represents a balance of the EMHD stress and the plastic stress on the lattice, it determines the evolution of the strain as a function of time and coordinate (as a parameter) and properties of the medium. The lattice strain induced by the motion of the electron fluid is relieved via creep.

Beyond the Maxwell time, Eq. (22), the flow reaches a steady state, 1 with the terminal strain (up to logarithmic accuracy),


The argument of the logarithm in (20) is very large, the shear velocity times the Hall time divided by the inter-particle distance . There is also a weak, logarithmic dependence on the parameter , the ratio of the lateral deformation of the electron fluid to the thickness of the current layer and on the ratio of shear speed to Alfvén speed . Typically, for parameters in the NS crust,


where we adopted for the critical strain (Hoffman & Heyl, 2012), (see also Horowitz & Kadau, 2009, who argue for larger ).

A time to reach the steady state (Maxwell time) is


(we assumed that the lateral displacement of the plastically flowing layer is of the order of its thickness .)

Thus, after a fairly short time a plastic creep relieves the crustal strain at a level of a few percent of the critical strain. At the same time the assumed deformations of the electron fluid proceed linearly in time, increasing on the Hall time scale. Under the EMHD approximation the magnetic field is frozen in the electron fluid, so that the deformations of the magnetic field increase linearly with time (recall, we assume that the plasma is ideal, so there is no slipping of the magnetic field through the electron fluid). The deformations of the magnetic field induced by the motion of the electron fluid within the crust will lead to the similar deformations of the magnetic field outside the neutron star – twisting of the external magnetic field.

Figure 1: Evolution of the strain in the plastically relaxing medium subject to the externally-imposed stress (linearly increasing with time). Left Panel: strain as function of at different times; Right Panel: strain as function of time at different (top to bottom curves, increases by ). Initial strain is zero; at early times the response of the medium is mostly elastic, balancing the externally-imposed stress. As time progresses, the system reaches a steady state, balancing the rates of the elastic strain and plastic relaxation. At smaller , where the strain grows quicker, the steady state is reached earlier in time. For the chosen driving, the terminal shape is, approximately, , Eq. (20).

Lorentz stresses lead to the distortion of the crustal lattice and, thus, do work on the crystal (this is a deviation from the EMHD, where the lattice is assumed to be infinitely ridged). At the steady state all the work done by the magnetic field is dissipated plastically, so that the plastically dissipated energy is


where is the magnetic energy density.

Next, let us verify that the velocity of ions due to plastic deformations is much smaller than the Hall drift velocity of electrons (which is ). The typical velocity of the ion component can be estimated as, see Eq. (19)


where we used Eq. (20) to estimate the typical strain. Since the term is logarithmically large, and is raised to large negative power, the ion velocity can be neglected.

Finally, let us comment on the possibility of crust melting during propagation of crustal faults (Beloborodov & Levin, 2014). The plastic flowing layer layer is heated viscously, while cools by conduction. Can the resulting dissipation lead to local melting of the layer? The answer critically depends on whether an initially distributed stress localizes into the narrow bands. At higher temperatures the effective viscosity is smaller; this could affect the overall evolution of the crust via so called shear thinning. If there is no localization of stress, the answer to the above question is clearly not. Balancing energy generation rate, , and the heat flow per unit area of the crack, of the order of , where is thermal conductivity and is melting temperature (Flowers & Itoh, 1976, 1979, 1981; Potekhin et al., 1999; Page & Reddy, 2012), the required size of the plastically flowing layer,


is much larger than the size of the neutron star.

On the other hand, some material do show stress localization - formation of so-called shear bands (Wright, 2002). This dynamic localization of stress (as opposed to stress localization due to preexisting defects) depends critically on the details of the shape of the stress-strain curve (Rudnicki & Rice, 1975; Montési & Zuber, 2002), which are not known for the neutron star crusts. Typically, stress localization occurs for compressive deformations (as opposed to shear deformations Rice, 1976); compressive deformations are not likely to be important in neutron star crusts.

4 Evolution of the external twist

4.1 Persistent emission

Motion of the electron fluid inside the crust twists the outside magnetic field - the electric current is pushed outside of the neutron star. As the current-carrying particles hit the surface, they are stopped by the interaction with the crust, heating the surface - this plays a role of resistivity, leading to continuos untwisting of the magnetic flux. In order to produce flares through instabilities of the twisted flux tube the rate of twist, of the order of the Hall time scale, should exceed the rate of this resistive relaxation.

Consider a magnetic flux tube that at the stellar surface has a radius and carries total current (assuming that the charge carries move relativistically). By flux conservation, the radius of the flux tube changes with height according to . Typical toroidal magnetic field within the flux tube then changes according to


The field line twist changes with height according to


The maximal twist of the magnetic field line is reached at the highest point :


For a stable flux tube it is required that .

We can use Eq. (28) to parametrize the current in term of the maximal height reached by the flux tube and the maximal twist




Note, that the energy stored in the toroidal magnetic field per unit is independent of radius


(For a given maximal twist and given radius at the neutron star surface the flux tubes that extend to larger heights carry less toroidal energy.) The total energy in the toroidal magnetic field can be estimated as


If current-carrying particles are accelerated though a typical potential (Beloborodov, 2009), the dissipated power is


where we assumed that . This value is close to the observed persistent luminosity of magnetars.

For those flux tubes that satisfy , Eq. (37), the twisting due to Hall drift is relieved gradually via resistive untwisting. Estimating , the dissipation rate is


Note the strong dependence of the dissipated power on the magnetic field, (available energy , dissipation time scale ). This explains why magnetar activity, especially giant flares, is mostly seen in the high field neutron stars: in our model all neutron stars show magnetar-like activity, but the continuously dissipated power has a cubic dependence on the magnetic field. The scaling is consistent with data (Fig. 13 of Olausen & Kaspi, 2014). Also major flares have a threshold magnetic field, Eq. (38), while smaller flares can occur in lower field/older neutron stars.

4.2 Flares: Hall twisting versus resistive untwisting of flux tubes

External twist is driven by the Hall evolution of the electron fluid within the crust and is relieved via resistive dissipation of the external currents (which is linear in the twist angle, Eq. (33). Thus, the twist angle evolves according to


where the resistive time scale is


In order to reach instability, , it is required that the resistive time be longer than the Hall time scale,


where we assumed .

Note strong dependence of this ratio on the magnetic field, , on the size of evolving region and on the electron density (crustal depth where the driving occurs). Since the size of the region is related to the energy dissipated in the flare, we expect that smaller magnetic fields can induce smaller size flares, at shallower depths. This explains why magnetar activity is observed in low field pulsars, but only producing relatively small burst and flares. Major flares (like giant flares) that involve large fraction of the crust do require high magnetic fields.

The condition for a magnetar to produce flares, , can be expressed as


The last terms being of the order of unity, the most important terms are , which evaluates to an order of unity. For G magnetic field, the corresponding density is


Summarizing, smaller scale bursts and flares can be produced in smaller fields neutron stars and/or at smaller crustal depth. (So that the twisting of external field by Hall drift occurs faster than the resistive untwisting of the external currents). Giant flares, which require large scale re-configuration of the magnetosphere require large magnetic fields and shallower depths, where crustal density is .

5 Internal dissipation of magnetospheric Alfvén pulse during magnetar flares

Observations of the post-flare evolution of the surface thermal emission indicate that in addition to external dissipation of energy, a considerable amount is also dissipated inside the neutron star (Eichler et al., 2006; Scholz et al., 2012). In the Section we discuss how external trigger can in addition lead to internal dissipation. Note, that this is opposite to what is assumed in the star-crack model, where internal trigger leads to external dissipation.

In the “Solar flare” model of magnetospheric energy release the flares are associated with instabilities developing in the twisted current-carrying magnetosphere, e.g.,  when a magnetic flux tube becomes unstable to kinks and reconnection. As a result Alfvén waves are launched from the reconnection site down to the stellar surface. As the Alfvén shock hits the surface it launched a whistler waves in the crust. On the Hall time scale of the crust such an impulse is nearly instantaneous. Below we demonstrate that the plastic dissipation of the elastic strain produced by such a pulse within the crust is broadly consistent with observations of the post-burst magnetar activity.

5.1 Interaction of Alfvén wave with the EMHD surface

Consider an Alfvén wave propagating along magnetic field with amplitude ( is a wave vector, is a displacement) and hitting a neutron star surface with plasma frequency . First, for typical plasma parameters most of the Alfvén pulse is reflected. This can be seen from the Fresnel’s reflection coefficients for the reflection of normal electromagnetic at refractive index jump,


For highly magnetized magnetospheric plasma ( refers to the Alfvén velocity in the magnetosphere), while inside the neutron star the whistler modes have . In this limit the reflection coefficient is .

Reflection of Alfvén pulse by the surface of a neutron star launches a whistlers pulse in the crust with the same frequency as the incoming Alfvén wave, so that


where is a typical frequency of the Alfvén waves in the magnetosphere. Since , the whistler waves launched in the curst have wavelength much larger that the skin depth, , consistent with the EMHD assumption. Note also, that the transmission coefficient


Thus, as the Alfvén wave launched at the reconnection region propagates toward the star, it is mostly reflected at the surface, launching a whistler pulse in the crust. In the next section we consider the propagation of a whistler pulse in the crust.

5.2 Green’s function for whistlers

Consider a half-space of plasma described by constant density electron MHD, while the region obeys MHD equations. Let the unperturbed magnetic field  be in the direction. Introducing the electron displacement we can write the magnetic field perturbations in the EMHD regime as


where the prime signifies the derivative with respect to and is the strain in the electron fluid. The EMHD Eq. (5) then gives the non-relativistic Schrodinger-type equation for the electron displacement


The solution to Eq. (44) is given by the Green’s function of a free quantum-mechanical particle, or, equivalently, by the Green’s function of the diffusion equation with a complex diffusion coefficient:


(We assumed that the displacement is non-vanishing on the boundary; otherwise there is also Green’ function .) This is the Green’s function for whistler modes; it is normalized to ; the minus sign chosen so that the strain, Eq. (47), is positive. The key feature of Eq. (45) is the wave vector dependence of whistler modes: initial -function has a broad spectral content, that gets spatially separated at later times/larger distances. Thus the whistler pulse becomes more coherent with distance.

For a given initial displacement the general solution is


The corresponding Green’s function for the strain is


Note, that for a fixed time the amplitude of the strain fluctuations increases with depth.

5.3 Rotational discontinuity in EMHD

A magnetar flare lasts typically msec; this time is much shorter that the period; the resulting field disturbance can be treated as Alfvén rotational discontinuity. Next we consider how such a structure interacts with the neutron star surface. In an Alfvén shock magnetic field in the plane of the discontinuity experiences a rotation by some angle. Thus, there is a surface current density associated with it:


where is the jump in the transverse magnetic field. The current pulse can be characterized by the transverse magnetic field  and the rotation angle ,


In our representation of the magnetic field in terms of the displacement vector, such initial condition corresponds to the discontinuous second derivative of , or, equivalently, the discontinuous first derivative of . The corresponding Green’s function for the strain is


where is the Fresnel integral of the second kind. (The superscript indicates that this is the Green’s function for the Alfven shock). In this case the strain in the electron fluid is given by


see Fig. 2.

Figure 2: Strain induced by Alfven shock in the crust. The plot shows a strain at three different times, illustrating evolution of the point of maximal strain to larger depths.

Note, that on the surface . The maximum value of the strain is reached at and equals . Note that this value is independent of time. Thus a whistler pulse propagates diffusively into the crust, conserving its maximal strain (in constant density and assuming no dissipation - see discussion below where these constraints are relaxed).

The amplitude of the electron strain depends only on the properties of the initial Alfvén pulse, while the propagation depends on the properties on the medium. Thus, in a varying density, in the WKB approximation, the strain is given by (64) with and amplitude depending on the depth, see §5.4.

In time the whistler pulse reaches typical depth


This is approximately, the depth and the density at the outer-inner crust boundary, close to the neutron drip point. Importantly, this is the depth where Scholz et al. (2012) infer deposition of the energy during the outburst in Swift J1822.3-1606.

5.4 Propagation of whistlers along density gradient

In §5.3 we considered propagation of whistlers in constant density, let us next take variations of density into account. Consider magnetic field along in plasma with density changing along the field, . Let a whistler with frequency propagate along direction. Then


where is defined in terms of and prime denotes derivative with respect to coordinate. The equation for becomes


which can be rewritten for ,


We can solve this equation in a WKB-like approximation. Note, that the classical WKB is developed for the second order differential equations, while Eq. (55) is fourth order. Rewriting (55) as


we seek solution in the form . Expanding for , the zeroth order gives






like in a classical WKB approach. Note that there are no resonances/reflection points. Qualitatively,


Solution (59) implies that . Thus, the amplitude of the whistlers increases with depth. This relation can be understood from the conservation of the Poynting flux carried by the wave:


from which it follows that .

5.5 Propagation of whistlers across density gradient

Next consider magnetic field along in a plasma with density changing in a transverse direction, . Let a whistler with frequency propagate obliquely to the magnetic field, so that . We find


Equation for has a form of the Schrodinger’s equation with effective potential . Solutions of (62) are well known: they include propagating and evanescent waves. For given and the reflection occurs at given by


Thus, oblique whistlers are reflected from the low density regions in a direction of the component of the density gradient perpendicular to the magnetic field, Fig. 3.

Figure 3: Propagation of oblique whistlers in density gradient. If magnetic field is not aligned with the density gradient whistlers are reflected from the high density regions in a direction of the component of the density gradient perpendicular to the magnetic field. The shade indicates increasing electron density.

5.6 Whistler pulse with plastic response

In the previous sections we considered propagation of the whistler pulse within the crust neglecting possible plastic response of the crust. If we take plasticity into account, the evolution will be approximately described by (20) with the driving term


This expression for incorporates both the Green’s fucntion in the constant density, Eq. (51), as well as WKB amplitude changes with depth, Eq. (61).

Similarly, the strain rate,


sharply increases with depth due to: (i) increasing density and the correspondingly increasing amplitude of the magnetic field fluctuations (the factor in Eq. (65); (ii) dispersive effects, so that from an initial broad pulse the higher-frequency/higher wave number components spatially separate (the factor in Eq. (65)).

Overall, at time equal to the local Hall time the flow will reach a plasticity-mediated terminal value given approximately by Eq. (20) with ,


The increasing density, the term under logarithm, makes the terminal larger, contributing to higher rate of dissipation at larger depths, Eq. (23).

We conclude that an Alfvén pulse generated at the magnetospheric reconnection event launches a whistler pulse in the crust, that is effectively dissipated around densities g cm.

6 Shear-density EMHD instability in neutron star crusts

In addition to linear evolution on Hall time scale (that builds crustal stress and leads to plastic deformation of the crust) and to the non-linear interaction of whistlers (that produces turbulent cascade to smaller scales, e.g.,  Lyutikov, 2013a), the EMHD plasma is also subject to a shear-density EMHD instability (Wood et al., 2014) that requires (i) magnetic field gradients on scales ; (ii) electron density gradient across the field ; (iii) density gradient should be on scales smaller than magnetic field gradients, . All these conditions are satisfied in the equatorial regions of neutron star’s crusts, where the magnetic field is mostly in the direction, varying in the radial direction with electron density varying in radial direction as well. Scales of density variations are typically smaller than those of the magnetic field: below a total mass density of , the density drops to zero over about 100 meters.

The unstable modes have a length scale longer than the transverse density scale, and a growth-rate of the order of the inverse Hall timescale. Qualitatively, the instability is driven by the current (velocity) shear that stretchers the perturbations and the rotation back by the whistler mode; the instability also requires density gradient, making it somewhat analogous to the magneto-buoyancy instability in regular MHD.

Next we consider shear-density EMHD instability in cylindrical geometry. This approximates the equatorial regions of the NS crust. Consider cylindrical EMHD configuration with radially-dependent axial magnetic field , . For axially-symmetric perturbations the EMHD equation (5) gives


Multiplying by and partially integrating this gives


The instability can be driven by the last term in the numerator. The instability requires that - both magnetic field and density decreasing with radius. From Eq. (68) it follows that the constant density configurations are stable, since the last term is zero, while the first two are positively defined. Also, comparing the density and magnetic field gradients in the second and third term in Eq. (68), we conclude that for the third term to dominate it is required that density changes on scales smaller than magnetic field, (the exact conditions on that leads to instability are mathematically challenging to derive.)

The typical growth of the instability can be estimated from the last term in Eq. (68) as the whistler travel time over the scale . For further discussion of the shear-density instability in neutron star crusts, see Appendix A.

7 Discussion and Summary

In this paper we advance a model of magnetar activity whereby twisting of the external magnetic fields by the electron Hall drift in the crust leads to both persistent emission as well as bursts and flares. In our model the elastic/plastic properties of the crust are not important for the production of flares: since the magnetic field in the crust is frozen in the electron fluid, and not the ion lattice, magnetic field changes can proceed without crustal faults, plastic or brittle. Evolution of the internal magnetic fields drives the electric currents into the magnetosphere, which becomes susceptible to plasma instabilities.

We suggest that at low strain rates the response of the crust is mostly plastic and present a model of slow plastic deformations. Employing a simple model of plastic deformations in metals, due to Gilman (1960), we found that the strain saturates at values well below the critical strain, Eq. (20). Importantly, the saturation level depends on the parameters of the problem only logarithmically, and thus remains insensitive to the variations of these badly constrained parameters.

The critical strain of metals is still an open equation. Ideal metals are expected to have a critical stress Gilman (1969). Taking accounts of defects reduces the critical strain by approximately two orders of magnitudes, e.g., Peierls-Nabarro critical stress Peierls (1940); Nabarro (1947). In practice, critical stresses are another two orders of magnitude lower, . Our result, that the plastically relaxing strain saturates at values is, qualitatively, consistent with this general scaling.

The scaling (20) also incorporates the fact that the behavior of a material depends on the rate of strain. For higher rates of strain (smaller ), the saturation strain is higher, reaching, in the limit of very fast driving , the critical strain. Modern molecular dynamic numerical experiments (e.g.,  Horowitz et al., 2011), that show a very high strain before development of a crustal faults, occur in this vary fast regime due to computational limitations. As a results, we suggest, they cannot capture slow plastic deformation regime.

While plastic recovery releases the crustal stress, the magnetic field, frozen in the electron fluid, undergoes twisting deformations, which are pushed outside the neutron star and are dissipated both continuously and in a Solar flare-like reconnection events. These magnetospheric reconnection events launch Alfvén shocks that convert into fast high frequency whistlers upon hitting the neutron star surface. The resulting whistler pulse can propagate to the outer core boundary and dissipate the energy of the initial Alfvén pulse on a time scale of approximately a day. This dissipated energy is detected as an enhanced post-flare surface emission.

We expect that all neutron stars show some kind of magnetar activity. Plastic deformations dissipate magnetic energy, the rate of dissipation, Eq. (23), is a cubic function of the magnetic field (two powers come from the energy density of the magnetic field and one from the inverse scaling of the Hall time). For crustal magnetic fields of the order of G the total dissipated power may reach erg s.

The persistent magnetar-like emission is unobservable in most rotationally powered neutron stars. Giant flares, which involve rearrangement of the magnetic field on the scale of a whole neutron star magnetosphere, can occur only if the magnetic field exceeds some threshold value, while smaller flares and persistent emission occurs in smaller fields. Two competing effects are at play: larger scales EMHD deformation lead to increase twist of the external magnetic field. On the other hand, non-linear interaction of whistlers lead to generation of smaller scales. Both processes occur, approximately, on Hall time scales .

In addition to large scale whistler motions, the non-linear interaction of whistlers will lead to the development of a turbulent cascade, generating smaller scales ion older pulsars. Adopting quasi-isotropic turbulent cascade, advocated by Lyutikov (2013a), the magnetic field power at scales , so that ; where is the outer scale of the turbulence; this gives - in a turbulent cascade more energy is released on smaller scales. Also, the non-linear energy transfer accelerates toward small scales, . This explains weak magnetar-like activity in an older, low magnetic field pulsar (Rea et al., 2012).

We argued that magnetic fields in the equatorial regions of neutron star crusts are unstable to shear-density instability, which develops on whistler time scale of the upper crust. In a stark contrast to MHD, this is an axially-symmetric instability of purely poloidal magnetic field – there is no analogue in the MHD case. The regions near the magnetic equator are preferred both due to the development of the above mentioned shear-density instability, as well as higher magnetic fields (and thus shorter time scales) in the equatorial regions (Lyutikov, 2010; Reisenegger, 2013). The fact that the magnetar activity occurs mostly near the magnetic equator, while the spin down rate measured the current flowing near the magnetic poles is related to the highly complicated timing behavior during flares (Dib & Kaspi, 2014).

Magnetars share many properties (strength of magnetic field first of all) with normal radio pulsars. Detection of radio emission from magnetars (Halpern et al., 2005) and magnetar-like emission from the rotation-powered pulsar (Gavriil et al., 2008) further blends this distinction. What makes magnetars different? We see two possibilities: first, the initial configuration of the magnetic field at the moment of crust freezing may be an important factor affecting magnetar activity. Numerically, stable magnetic configuration of fluid stars were fond to vary from nearly dipolar to complicated non-axisymmtrical shapes Braithwaite (2008). During magnetar phase, the internal structure of magnetic fields evolve toward the dipolar equilibrium; no toroidal field is required to remain inside a star, though we expect that some toroidal flux is trapped inside flux surfaces fully enclosed within a star. Since no analytical approximation to the structure of magnetic field in fluid stars exists, we cannot estimate, e.g., the difference in magnetic energies between the initial and final states.

Secondly, for a given poloidal fields the crustal magnetic fields can vary in a wide range of values (Akgün et al., 2013). Since magnetar activity is proportional to the thrid power of magnetic field, mild variations of the magnetic field within the crust can lead to large variations in the dissipated power.

The paradigm of magnetar evolution described above assumes that magnetic field is confined mostly to the crust of a neutron star. This is consistent with the possibility that strong magnetic fields are created by neutrino-driven turbulence in proto-neutron stars, which is mostly effective in the outer layers Colgate & Fryer (1995); Herant et al. (1994). An alternative possibility is that magnetar activity is driven by the core of the neutron star, which expels magnetic field after becoming a type-I superconductor.

I would like to thank Andrei Beloborodov, Andrew Cumming, Jeremy Heyl, Konstantinos Gourgouliatos, Victoria Kaspi, Yuri Levin, Mikhail Medvedev and Jay Melosh for discussions.

Appendix A Shear-density instability in neutron star crusts

To illustrate more clearly the shear-density instability in neutron star crusts here we first discuss a special case of the density and magnetic field that analytically shows the instability and then comment on the special case when the poloidal magnetic field becomes zero at some points.

Eq. (67) has a particularly simple solution, that demonstrates the instability, for the case of similar scaling of the magnetic field and density


In this case


showing instability for .

Two other special cases include , ,


and neutrally stable case , . These solution give analytical examples of the shearing density instability in the equatorial regions of the neutron star crusts.

One mathematical subtlety concerns the convergence of integrals in (68): it requires that the initial magnetic field does not pass through zero. If at some radius, , then expanding near the point gives


The requirement that perturbation vanishes at then leads to - neutral stability. We expect that surfaces where will be unstable to resistive tearing instability (Wood et al., 2014).


  1. The steady state is reached only for sufficiently slow driving since the plastic strain rate, , has a maximum value of reached at . So, for plastic relaxation cannot compensate for the externally imposed shear.


  1. Akgün, T., Reisenegger, A., Mastrano, A., & Marchant, P. 2013, MNRAS, 433, 2445
  2. Beloborodov, A. M. 2009, ApJ, 703, 1044
  3. Beloborodov, A. M., & Levin, Y. 2014, ArXiv e-prints
  4. Biskamp, D., Schwarz, E., Zeiler, A., Celani, A., & Drake, J. F. 1999, Physics of Plasmas, 6, 751
  5. Braithwaite, J. 2008, MNRAS, 386, 1947
  6. Braithwaite, J., & Spruit, H. C. 2004, Nature, 431, 819
  7. Cho, J., & Lazarian, A. 2004, ApJ, 615, L41
  8. Colgate, S. A., & Fryer, C. 1995, Phys. Rep., 256, 5
  9. Dib, R., & Kaspi, V. M. 2014, ApJ, 784, 37
  10. Eichler, D., Lyubarsky, Y., Kouveliotou, C., & Wilson, C. A. 2006, ArXiv Astrophysics e-prints
  11. Flowers, E., & Itoh, N. 1976, ApJ, 206, 218
  12. —. 1979, ApJ, 230, 847
  13. —. 1981, ApJ, 250, 750
  14. Flowers, E., & Ruderman, M. A. 1977, ApJ, 215, 302
  15. Frohlich, C. 2006, Deep Earthquakes (Cambridge U. Press, New York)
  16. Galtier, S., & Bhattacharjee, A. 2003, Physics of Plasmas, 10, 3065
  17. Gavriil, F. P., Gonzalez, M. E., Gotthelf, E. V., Kaspi, V. M., Livingstone, M. A., & Woods, P. M. 2008, Science, 319, 1802
  18. Gilman, J. J. 1960, Australian Journal of Physics, 13, 327
  19. —. 1969, Micromechanics of Flow in Solids (McGraw-Hill Book Company; First Edition edition (1969))
  20. Goldreich, P., & Reisenegger, A. 1992, ApJ, 395, 250
  21. Gordeev, A. V., Kingsep, A. S., & Rudakov, L. I. 1994, Phys. Rep., 243, 215
  22. Halpern, J. P., Gotthelf, E. V., Becker, R. H., Helfand, D. J., & White, R. L. 2005, ApJ, 632, L29
  23. Herant, M., Benz, W., Hix, W. R., Fryer, C. L., & Colgate, S. A. 1994, ApJ, 435, 339
  24. Hoffman, K., & Heyl, J. 2012, MNRAS, 426, 2404
  25. Horowitz, C. J., Hughto, J., Schneider, A., & Berry, D. K. 2011, ArXiv e-prints
  26. Horowitz, C. J., & Kadau, K. 2009, Physical Review Letters, 102, 191102
  27. Howes, G. G., Dorland, W., Cowley, S. C., Hammett, G. W., Quataert, E., Schekochihin, A. A., & Tatsuno, T. 2008, Physical Review Letters, 100, 065004
  28. Jones, P. B. 2003, ApJ, 595, 342
  29. Kingsep, A. S., Chukbar, K. V., & Ian’kov, V. V. 1987, Voprosy Teorii Plazmy, 16, 209
  30. Lander, S. K., & Jones, D. I. 2012, MNRAS, 424, 482
  31. Levin, Y., & Lyutikov, M. 2012, MNRAS, 427, 1574
  32. Lighthill, M. J. 1960, Royal Society of London Philosophical Transactions Series A, 252, 397
  33. Link, B. 2013, ArXiv e-prints
  34. Lubliner, J. 2008, Plasticity Theory (Dover Publications)
  35. Lyutikov, M. 2003, MNRAS, 346, 540
  36. —. 2006, MNRAS, 367, 1594
  37. —. 2010, MNRAS, 402, 345
  38. —. 2013a, Phys. Rev. E, 88, 053103
  39. —. 2013b, ArXiv e-prints
  40. Montési, L. G. J., & Zuber, M. T. 2002, Journal of Geophysical Research (Solid Earth), 107, 2045
  41. Nabarro, F. R. N. 1947, Proceedings of the Physical Society, 59, 256
  42. Negele, J. W., & Vautherin, D. 1973, Nuclear Physics A, 207, 298
  43. Olausen, S. A., & Kaspi, V. M. 2014, ApJS, 212, 6
  44. Orowan, E. 1934, Zeitschrift fur Physik, 89, 605
  45. Page, D., & Reddy, S. 2012, ArXiv e-prints
  46. Palmer, D. M., Barthelmy, S., Gehrels, N., Kippen, R. M., Cayton, T., Kouveliotou, C., Eichler, D., Wijers, R. A. M. J., Woods, P. M., Granot, J., Lyubarsky, Y. E., Ramirez-Ruiz, E., Barbier, L., Chester, M., Cummings, J., Fenimore, E. E., Finger, M. H., Gaensler, B. M., Hullinger, D., Krimm, H., Markwardt, C. B., Nousek, J. A., Parsons, A., Patel, S., Sakamoto, T., Sato, G., Suzuki, M., & Tueller, J. 2005, Nature, 434, 1107
  47. Peierls, R. 1940, Proceedings of the Physical Society, 52, 34
  48. Potekhin, A. Y., Baiko, D. A., Haensel, P., & Yakovlev, D. G. 1999, A&A, 346, 345
  49. Rea, N., Israel, G. L., Esposito, P., Pons, J. A., Camero-Arranz, A., Mignani, R. P., Turolla, R., Zane, S., Burgay, M., Possenti, A., Campana, S., Enoto, T., Gehrels, N., Gövgüş, E., Götz, D., Kouveliotou, C., Makishima, K., Mereghetti, S., Oates, S. R., Palmer, D. M., Perna, R., Stella, L., & Tiengo, A. 2012, ApJ, 754, 27
  50. Reisenegger, A. 2013, ArXiv e-prints
  51. Reisenegger, A., Benguria, R., Prieto, J. P., Araya, P. A., & Lai, D. 2007, ArXiv Astrophysics e-prints
  52. Rice, J. R. 1976, in in: W.T. Koiter (Ed.), Theoretical and Applied Mechanics (North-Holland Publishing Company), 207–220
  53. Rudnicki, J. W., & Rice, J. R. 1975, Journal of Mechanics Physics of Solids, 23, 371
  54. Scholz, P., Ng, C.-Y., Livingstone, M. A., Kaspi, V. M., Cumming, A., & Archibald, R. F. 2012, ApJ, 761, 66
  55. Thompson, C., & Duncan, R. C. 1993, ApJ, 408, 194
  56. —. 1995, MNRAS, 275, 255
  57. —. 2001, ApJ, 561, 980
  58. Vainshtein, S. I. 1973, Soviet Journal of Experimental and Theoretical Physics, 37, 73
  59. Wood, T. S., Hollerbach, R., & Lyutikov, M. 2014, ArXiv e-prints
  60. Woods, P. M., & Thompson, C. 2006, Soft gamma repeaters and anomalous X-ray pulsars: magnetar candidates, ed. W. H. G. Lewin & M. van der Klis, 547–586
  61. Wright, T. W. 2002, The Physics and Mathematics of Adiabatic Shear Bands