MAGNETAR FIELD EVOLUTION AND CRUSTAL PLASTICITY

The initial rise timescale for giant flares is 1 ms, suggestive of some explosive reconnection process in the charge-filled corona surrounding the star (Lyutikov 2003; Huang & Yu 2014; Elenbaas et al. 2016), in analogy with solar flares (Masada et al. 2010). For the timescale to be sufficiently short, this reconnection must occur within a few stellar radii (∼10 km) of the surface. Direct release of energy from the crust is unlikely: the shortest characteristic timescale is ∼0.2 s, from shear-wave propagation—and the release of energy may be slower still (Lyutikov 2006; Link 2014).

The coronal-reconnection scenario for giant flares requires a huge amount of magnetic energy to be stored in strongly twisted exterior field lines. The long periods of magnetars, ∼1–10 s, mean that the rotation-powered mechanism invoked for pulsar magnetospheres would only twist a small region of field lines near the polar cap (Glampedakis et al. 2014). Instead, shearing of the crust could move magnetospheric footpoints (Thompson et al. 2002), inducing a twist and forming an equatorial lobe of current, with charges stripped off the surface (Beloborodov 2009).

This corona, in turn, requires a mechanism for build-up of crustal stresses. If the star were born rotating rapidly, its crust would freeze into an oblate spheroid shape, and the star's later spindown would induce stresses as the preferred shape of the star became more spherical (Ruderman 1969)—but this axisymmetric process would displace footpoints vertically, without inducing any coronal twist. The only credible candidate to generate the required azimuthal crustal motion is the star's magnetic field. Magnetar-corona simulations show how crustal motion (albeit added by hand) produces exterior twist (Parfrey et al. 2013), potentially leading to an overtwisting instability (Wolfson 1995) and flare. The build-up of stresses due to a changing global equilibrium can produce the right kind of crustal displacement (Lander et al. 2015), and power giant flares, but this has not been verified with full core-crust field evolutions. On the other hand, magneto-thermal evolutions of the crust show how the crustal field can evolve to become locally intense and fail in small events (Pons & Rea 2012). While this provides a credible explanation of the phenomena of magnetar bursts and their persistent luminosity, it is less clear whether it can be applied to giant flares, with their shorter rise timescales and greater energy release.

The relevant region for energy storage is whatever portion of the interior stellar magnetic field can be tapped for giant flare energy; a deeply buried core field may be irrelevant over magnetar lifetimes. The most energetic giant flare, from SGR 1806–20, was estimated to be $3\times {10}^{46}$ erg (Hurley et al. 2005); the other two giant flares, from SGR 1900+14 (Feroci et al. 2001) and SGR 0526–66 (Mazets et al. 1979), were both $\sim {10}^{45}$ erg. In order to have a mechanism that can credibly provide enough energy to power any giant flare, we will adopt a magnetar model harboring 10 times the energy of the largest one observed: $3\times {10}^{47}$ erg. Note that the more detailed discussion in Thompson & Duncan (2001) suggests a magnetar energy budget $\gtrsim {10}^{47}$ erg; this is probably a conservative value since the paper predated the largest giant flare.

We assume that an NS crust responds in an elastic, reversible manner to stresses below some yield value ${\tau }_{\mathrm{el}}$ and plastically above it. Many terrestrial media can also be approximated as having elastic or plastic responses, depending on the magnitude of the applied stress; they include toothpaste, crude oil, mud, and concrete.

Let us first consider a magnetized crust below ${\tau }_{\mathrm{el}}$. ${\boldsymbol{B}}$ need not be in a global fluid equilibrium with the star since the crust can absorb magnetic stresses:

$$\begin{eqnarray}&&0=-{\rm{\nabla }}P-\rho {\rm{\nabla }}{\rm{\Phi }}+\displaystyle \frac{1}{4\pi }({\rm{\nabla }}\times {\boldsymbol{B}})\times {\boldsymbol{B}}+{\rm{\nabla }}\;\cdot \;{\boldsymbol{\tau }},\end{eqnarray} \tag{ 2 }$$

where P is pressure, ρ mass density, Φ the gravitational potential, and ${\boldsymbol{\tau }}$ the elastic stress tensor. We will think of the crust reaching ${\tau }_{\mathrm{el}}$ for some Lorentz force $({\rm{\nabla }}\times {{\boldsymbol{B}}}_{\mathrm{el}})\times {{\boldsymbol{B}}}_{\mathrm{el}}$, which defines the "yield magnetic field" ${{\boldsymbol{B}}}_{\mathrm{el}}$. For ${\boldsymbol{B}}\gt {{\boldsymbol{B}}}_{\mathrm{el}}$ the crust's response to the force

$$\begin{eqnarray}&&({\rm{\nabla }}\times {\rm{\Delta }}{\boldsymbol{B}})\times {\rm{\Delta }}{\boldsymbol{B}}\equiv ({\rm{\nabla }}\times {\boldsymbol{B}})\times {\boldsymbol{B}}-({\rm{\nabla }}\times {{\boldsymbol{B}}}_{\mathrm{el}})\times {{\boldsymbol{B}}}_{\mathrm{el}}\gt 0\end{eqnarray} \tag{ 3 }$$

will be plastic, unless the crust is in hydromagnetic equilibrium (for which ${\boldsymbol{\tau }}=0$). Equation (3) defines the field ΔB that sources the plastic flow; note that in general $\Delta {\boldsymbol{B}}\ne {\boldsymbol{B}}-{{\boldsymbol{B}}}_{{\rm{el}}}$.

Next, we wish to make a depth-dependent estimate of ${B}_{\mathrm{el}}$, rather than using the standard assumption that ${B}_{\mathrm{el}}\sim {10}^{15}$ G. As in Lander et al. (2015), we make fits to NS equation-of-state parameters from Douchin & Haensel (2001) and a magnetar temperature profile from Kaminker et al. (2009). Using these, we calculate ${\tau }_{\mathrm{el}}$ from the formula of Chugunov & Horowitz (2010). Now, Equation (2) evaluated at the yield stress shows us that a natural definition for ${B}_{\mathrm{el}}$ is

$$\begin{eqnarray}&&{B}_{\mathrm{el}}\equiv \sqrt{4\pi {\tau }_{\mathrm{el}}}.\end{eqnarray} \tag{ 4 }$$

Plotting this in Figure 1, we see that ${B}_{\mathrm{el}}$ may be approximated by the linear relation

$$\begin{eqnarray}&&{B}_{\mathrm{el}}\approx 1.8\times {10}^{16}\left(1-\displaystyle \frac{r}{{R}_{\ast }}\right)\;{\rm{G}},\end{eqnarray} \tag{ 5 }$$

where r is the radial coordinate and R_* the surface radius.

Zoom In Zoom Out Reset image size

Beyond the elastic yield stress, the crust must be in motion, with some velocity ${\boldsymbol{v}}$. Since we expect this flow to be approximately incompressible, we take ${\rm{\nabla }}\cdot {\boldsymbol{v}}=0$. General models of viscoplastic flow are based on a relation between two tensorial quantities: the rate of strain $\dot{{\boldsymbol{\varepsilon }}}$ and the stress. Let us however assume that the crust is under simple shear stress, for which each tensor has a single independent component, $\dot{\varepsilon }$ and τ, and the problem reduces to a scalar one. We may now adopt the Bingham model to relate $\dot{\varepsilon }$ and τ when $\tau \gt {\tau }_{\mathrm{el}}$:

$$\begin{eqnarray}&&\dot{\varepsilon }=\displaystyle \frac{1}{2\nu }(\tau -{\tau }_{\mathrm{el}})H(\tau -{\tau }_{\mathrm{el}}),\end{eqnarray} \tag{ 6 }$$

where ν is the dynamical viscosity of the plastic flow and $H(\cdot )$ the Heaviside function. Under Equation (6), one can derive an equation of motion valid above ${\tau }_{\mathrm{el}}$ (Prager 1961); for our magnetar crust model this is

$$\begin{eqnarray}\rho \dot{{\boldsymbol{v}}}+\rho ({\boldsymbol{v}}\cdot {\rm{\nabla }}){\boldsymbol{v}} & = & -{\rm{\nabla }}P-\rho {\rm{\nabla }}{\rm{\Phi }}+\nu {{\rm{\nabla }}}^{2}{\boldsymbol{v}}\\ & & +\displaystyle \frac{1}{4\pi }({\rm{\nabla }}\times {\boldsymbol{B}})\times {\boldsymbol{B}}+{\rm{\nabla }}\;\cdot \;{{\boldsymbol{\tau }}}_{\mathrm{el}}.\end{eqnarray} \tag{ 7 }$$

The fluid terms in this equation have the same form as in the standard Navier–Stokes equation for a viscous medium. In general, however, viscoplastic dynamics can be richer than those of viscous fluids—if, for example, the medium yields in a more complicated manner than through simple shearing, a tensor generalization of Equation (6) is required. This in turn results in the appearance of a new, uniquely plastic force term in the equation of motion, ${\tau }_{\mathrm{el}}{\rm{\nabla }}\cdot (\dot{{\boldsymbol{\varepsilon }}}/| | \dot{{\boldsymbol{\varepsilon }}}| | )$, where $| | \cdot | |$ is a tensor norm (Prager 1961). This may be important for more sophisticated modeling of NS crustal failure.

Returning to Equation (7), we anticipate the plastic flow to be slow and steady, and hence neglect the advective term $({\boldsymbol{v}}\cdot {\rm{\nabla }}){\boldsymbol{v}}$ for being quadratic in ${\boldsymbol{v}}$ and the acceleration term $\dot{{\boldsymbol{v}}}=0$. The resulting viscous crustal motion is Stokes flow. Comparing this limiting case of Equation (7) with (2), and dropping the tiny differences between the hydrostatic terms ${\rm{\nabla }}P+\rho {\rm{\nabla }}{\rm{\Phi }}$ in the elastic and plastic regimes, we arrive at the intuitive result that the unbalanced piece of the Lorentz force sources a viscous flow:

$$\begin{eqnarray}&&\nu {{\rm{\nabla }}}^{2}{\boldsymbol{v}}=-\displaystyle \frac{1}{4\pi }({\rm{\nabla }}\times {\rm{\Delta }}{\boldsymbol{B}})\times {\rm{\Delta }}{\boldsymbol{B}},{\rm{\nabla }}\cdot {\boldsymbol{v}}=0.\end{eqnarray} \tag{ 8 }$$

If this unbalanced Lorentz force is curl-free, Equation (8) takes an extremely compact form. Since ${\rm{\nabla }}\cdot {\boldsymbol{v}}=0$, we may write ${\boldsymbol{v}}={\rm{\nabla }}\times {\boldsymbol{\psi }}$ for some potential ${\boldsymbol{\psi }}$. This ${\boldsymbol{\psi }}$ is only fixed up to transformations of the form ${\boldsymbol{\psi }}\to {\boldsymbol{\psi }}+{\rm{\nabla }}\phi$, a freedom which allows us to choose ${\rm{\nabla }}\cdot {\boldsymbol{\psi }}=0$ (the Coulomb gauge). Now taking the curl of Equation (8) and using these results, together with various vector identities (including those for a triple curl), we find that the flow is governed by the biharmonic equation:

$$\begin{eqnarray}&&{{\rm{\nabla }}}^{4}{\boldsymbol{\psi }}=0.\end{eqnarray} \tag{ 9 }$$

The overwhelming majority of studies of magnetized NS crusts assume electron magnetohydrodynamics: the crustal ion lattice is strictly static, magnetic forces are balanced by elastic stresses (rendering the equation of motion irrelevant), and only the electrons are mobile (see, e.g., Gourgouliatos & Cumming 2014; Wood & Hollerbach 2015). Field evolution is then governed by the interplay of two secular terms (Cumming et al. 2004): the conservative Hall drift and Ohmic decay (the second and third terms, respectively, in Equation (10)). One situation where this is clearly no longer valid is for $\tau \gt {\tau }_{\mathrm{el}}$, when we need to add a third term involving the plastic-flow velocity ${\boldsymbol{v}}$:

$$\begin{eqnarray}\displaystyle \frac{\partial {\boldsymbol{B}}}{\partial t} & = & {\rm{\nabla }}\times ({\boldsymbol{v}}\times {\boldsymbol{B}})-{\rm{\nabla }}\times \left[\displaystyle \frac{({\rm{\nabla }}\times {\boldsymbol{B}})\times {\boldsymbol{B}}}{4\pi {\rho }_{e}}\right]\\ & & -{\rm{\nabla }}\times \left[\displaystyle \frac{{\rm{\nabla }}\times {\boldsymbol{B}}}{4\pi {\sigma }_{0}}\right],\end{eqnarray} \tag{ 10 }$$

where ${\rho }_{e}$ is charge density and ${\sigma }_{0}$ electrical conductivity. Next, we check when this plastic-flow term may dominate the crust's evolution. Dimensional analysis of Equation (8) shows that

$$\begin{eqnarray}&&v\sim \displaystyle \frac{{L}_{\mathrm{pl}}^{2}{({\rm{\Delta }}B)}^{2}}{{L}_{{\rm{B}}}\nu }\end{eqnarray} \tag{ 11 }$$

where ${L}_{\mathrm{pl}},{L}_{{\rm{B}}}$ are the lengthscales of the plastic flow and the field, respectively. Assuming the plastic flow to be more localized than the magnetic field, i.e., ${L}_{\mathrm{pl}}\leqslant {L}_{{\rm{B}}}$, and also that ${\rm{\Delta }}B\sim B$ (reasonable unless the field is just above ${B}_{\mathrm{el}}$), Equations (11) and (10) combine to give the following plastic-flow evolution timescale:

$$\begin{eqnarray}&&{t}_{\mathrm{pl}}\sim \displaystyle \frac{{L}_{{\rm{B}}}}{{L}_{\mathrm{pl}}}\displaystyle \frac{\nu }{{B}^{2}},\quad \mathrm{when}\ B\gt {B}_{\mathrm{el}}.\end{eqnarray} \tag{ 12 }$$

The plastic viscosity ν is an unknown crustal parameter, whose calculation might require molecular-dynamics simulations, but we may make an estimate by arguing from dimensional analysis that $\nu \sim {t}_{\mathrm{char}}{\tau }_{\mathrm{el}}$ for some characteristic plastic timescale ${t}_{\mathrm{char}}$. We already know ${\tau }_{\mathrm{el}}$, from the calculation described before (4). For ${t}_{\mathrm{char}}$, let us demand it be short enough to allow for a persistent corona. Since coronal dissipation occurs over ∼1–10 years (Beloborodov & Thompson 2007), we need to take ${t}_{\mathrm{char}}=10\;{\rm{years}}$ to allow sufficiently fast regeneration of the exterior current. ν will increase with crustal depth and decrease with temperature, but rather than exploring these dependences here we will simply evaluate ${t}_{\mathrm{char}}{\tau }_{\mathrm{el}}$ at the crust–core boundary for our (hot) magnetar model and use this as a reference value: $\nu \sim {10}^{38}$ poise. If this estimate is reliable, then the $\nu \lesssim {10}^{35}$ poise condition for thermo-plastic instability (Beloborodov & Levin 2014) would only be attained in the outer crust. For comparison, a famously viscous terrestrial material is pitch; in an experiment begun in 1930, nine drops of pitch have fallen to date, yielding a viscosity estimate of $2\times {10}^{9}$ poise (Edgeworth et al. 1984).

Using our viscosity estimate, we may compare ${t}_{\mathrm{pl}}$ with the timescale ${t}_{\mathrm{Hall}}$ for Hall drift (which is faster than Ohmic decay for high B), using the results of Cumming et al. (2004). We forsee three regimes:

• 1.  
  If $B\gtrsim {10}^{15}$ G, plastic flow always dominates, with a characteristic timescale $\lesssim 3{L}_{{\rm{B}}}/{L}_{\mathrm{pl}}$ years.
• 2.  
  In the range $B\sim (2\mbox{--}7)\times {10}^{14}$ G, ${B}_{\mathrm{el}}$ is exceeded for crustal depths ≲100–400 m, and ${t}_{\mathrm{pl}}\sim {t}_{\mathrm{Hall}}\;\sim 10\mbox{--}100\;{\rm{years}}$ there, so the two field-evolution mechanisms compete. Deeper into the crust there is no plastic flow.
• 3.  
  If $B\lesssim {10}^{13}$ G, plastic flow shuts off essentially everywhere.

We now consider a simple illustrative example of magnetically induced plastic flow that can be solved analytically. Taking a segment of the crust and neglecting its curvature, we are left with a Cartesian slab of material, as shown in Figure 2: the x and z coordinates represent azimuthal and radial directions in the crust, respectively.

Zoom In Zoom Out Reset image size

Assume a Lorentz force that stresses the crust beyond ${\tau }_{\mathrm{el}}$ in some band $-{y}_{0}\lt y\lt {y}_{0}$. To produce a constant, non-radial force, we choose

$$\begin{eqnarray}&&{\rm{\Delta }}{\boldsymbol{B}}=b\sqrt{x}{{\boldsymbol{e}}}_{z},\end{eqnarray} \tag{ 13 }$$

where b is a constant. Note that ${\rm{\Delta }}{\boldsymbol{B}}$ is divergence-free, and its corresponding force curl-free, as required. Our assumptions, however, force us to have a depth-independent ν; the intention is to relax this restriction in future work. Now, returning to Equation (8),

$$\begin{eqnarray}&&{{\rm{\nabla }}}^{2}{\boldsymbol{v}}=({{\rm{\nabla }}}^{2}{v}_{x}){{\boldsymbol{e}}}_{x}+({{\rm{\nabla }}}^{2}{v}_{y}){{\boldsymbol{e}}}_{y}=\displaystyle \frac{{b}^{2}}{8\pi \nu }{{\boldsymbol{e}}}_{x}.\end{eqnarray} \tag{ 14 }$$

Here, v_z = 0 to ensure incompressibility of the flow. Equation (14) implies ${{\rm{\nabla }}}^{2}{v}_{y}=0;$ we will satisfy this by taking v_y = 0 (this is not actually a simplification; one can show that keeping a ${v}_{y}\ne 0$ term is incompatible with confining the flow into a channel $-{y}_{0}\lt y\lt {y}_{0}$).

In terms of the velocity potential, ${\boldsymbol{v}}={{\boldsymbol{e}}}_{x}\partial \psi /\partial y-{{\boldsymbol{e}}}_{y}\partial \psi /\partial x$, so v_y = 0 implies $\psi =\psi (y)$, and Equation (9) becomes simply ${d}^{4}\psi /{{dy}}^{4}=0$. Integrating once and comparing with Equation (14) gives

$$\begin{eqnarray}&&\displaystyle \frac{{d}^{3}\psi }{{{dy}}^{3}}=\displaystyle \frac{{d}^{2}{v}_{x}}{{{dy}}^{2}}=\displaystyle \frac{{b}^{2}}{8\pi \nu }.\end{eqnarray} \tag{ 15 }$$

Integrating twice more and imposing the boundary conditions that the flow must stop when there is no unbalanced Lorentz force, i.e., ${v}_{x}(-{y}_{0})={v}_{x}({y}_{0})=0$, we find that the plastic-flow velocity is

$$\begin{eqnarray}&&{\boldsymbol{v}}=\displaystyle \frac{{b}^{2}}{16\pi \nu }({y}^{2}-{y}_{0}^{2}){{\boldsymbol{e}}}_{x}.\end{eqnarray} \tag{ 16 }$$

The evolution of a crustal field ${\boldsymbol{B}}={B}_{x}(x,y,t){{\boldsymbol{e}}}_{x}\;+{B}_{y}(x,y,t){{\boldsymbol{e}}}_{y}+{B}_{z}(x,y,t){{\boldsymbol{e}}}_{z}$ under this flow is given by ${\rm{\nabla }}\times ({\boldsymbol{v}}\times {\boldsymbol{B}})$, i.e.,

$$\begin{eqnarray}\displaystyle \frac{\partial {\boldsymbol{B}}}{\partial t} & = & \displaystyle \frac{{b}^{2}}{16\pi \nu }\left\{\displaystyle \frac{\partial }{\partial y}[({y}^{2}-{y}_{0}^{2}){B}_{y}]{{\boldsymbol{e}}}_{x}\right.\\ & & +\left.({y}_{0}^{2}-{y}^{2})\left(\displaystyle \frac{\partial {B}_{y}}{\partial x}{{\boldsymbol{e}}}_{y}+\displaystyle \frac{\partial {B}_{z}}{\partial x}{{\boldsymbol{e}}}_{z}\right)\right\}.\end{eqnarray} \tag{ 17 }$$

By separation of variables, writing each field component in the form $X(x)Y(y)T(t)$, we find that we require ${B}_{y}={B}_{y}(y)$ and ${B}_{z}={B}_{z}(y)$ for consistency, so that only B_x is time-varying, with

$$\begin{eqnarray}&&{B}_{x}=\displaystyle \frac{{b}^{2}}{16\pi \nu }\displaystyle \frac{d}{{dy}}[({y}^{2}-{y}_{0}^{2}){B}_{y}]t.\end{eqnarray} \tag{ 18 }$$

This toy model suggests that a magnetic field with an unbalanced radial component induces an azimuthal plastic flow, whose velocity depends (through b) on how much one exceeds ${\tau }_{\mathrm{el}};$ field evolution along this flow provides an ever-increasing twist to the corona. More realistically, the twist would be relieved gradually through viscoplastic or Ohmic dissipation, or suddenly in coronal flare events.