UNTWISTING MAGNETOSPHERES OF NEUTRON STARS

The evolution of magnetic field is related to electric field  E according to the induction equation,

$$\begin{equation} \frac{1}{c}\,\frac{\partial {\,{\mathbf B}}}{\partial t}=-\nabla \times {\,{\mathbf E}}. \end{equation} \tag{ 7 }$$

We will express ∇ ×  E in curvilinear coordinates qⁱ defined as follows. Let us label magnetic field lines on a flux surface f by the azimuthal angle ϕ₀ of their northern footpoints at r = r_c (northern footpoints have B_r > 0). Thus, the set of all field lines is parameterized by two coordinates f and ϕ₀. Let s be a parameter running along the field line (increasing in the direction of  B). The parameter s may be chosen, e.g., equal to length l measured along the field line from its northern footpoint (then |e_s| = 1). The coordinate system qⁱ = (s, f, ϕ₀) covers the entire magnetic field that passes through the sphere r = r_c and emerges in the force-free region.

The electric field can be written in components in the new coordinate basis,

$$\begin{equation} {\,\mathbf E}=\sum _{i=1}^3 E^i{\mathbf e}_i, \qquad {\mathbf e}_i=\frac{\partial {\mathbf r}}{\partial q^i}. \end{equation} \tag{ 8 }$$

The length of the basis vectors will be denoted by h_i = |e_i|. Note that

$$\begin{equation} {\mathbf e}_3\equiv \left.\frac{\partial {\mathbf r}}{\partial \phi _0}\right|_{f,s} =\left.\frac{\partial {\mathbf r}}{\partial \phi }\right|_{r,\theta }={\mathbf e}_\phi =h_\phi \hat{{\mathbf e}}_\phi. \end{equation} \tag{ 9 }$$

We will need the determinant g of metric g_ik = e_i · e_k, which is given by $\sqrt{g}={{\mathbf e}}_s\cdot ({\mathbf e}_{f}\times {\mathbf e}_\phi)$. Consider an infinitesimal axisymmetric ring perpendicular to the poloidal component of  B. Its surface element is e_f × e_ϕ df dϕ, and the area of the ring is 2π |e_f × e_ϕ|df. The magnetic flux through the ring, df, is

$$\begin{equation} df=2\pi {\,\mathbf B}\cdot ({\mathbf e}_{f}\times {\mathbf e}_\phi)df, \end{equation} \tag{ 10 }$$

which implies the identity 2π B · (e_f × e_ϕ) = 1. We substitute  B = Be_s/h_s and find

$$\begin{equation} \sqrt{g}=\frac{h_s}{2\pi B}. \end{equation} \tag{ 11 }$$

The general expression for ∇ ×  E in curvilinear coordinates is

$$\begin{equation} \left(\nabla \times {\,\mathbf E}\right)^i=\sum _{j,k=1}^{3} \frac{\epsilon ^{ijk}}{\sqrt{g}}\,\frac{\partial {E_k}}{\partial q^j}, \end{equation} \tag{ 12 }$$

where ^ijk is the Levi-Civita symbol and E_i = e_i ·  E are the covariant components of the electric field in the coordinate system qⁱ. Using ${\mathbf e}_i=h_i\hat{{\mathbf e}}_i$, we obtain the azimuthal component of ∇ ×  E in the normalized basis $\hat{{\mathbf e}}_i$ and find

$$\begin{equation} \frac{1}{c}\,\frac{\partial B_\phi }{\partial t}= h_\phi \,\frac{2\pi B}{h_s}\, \left(\frac{\partial E_s}{\partial f}-\frac{\partial E_f}{\partial s}\right). \end{equation} \tag{ 13 }$$

Let us divide both sides of this equation by 2πBh_ϕ/h_s and integrate it over s along an entire closed field line (including its part at r < r_c). Then the second term on the right-hand side disappears, and we get

$$\begin{equation} \frac{1}{2\pi c}\,\oint \frac{\partial B_\phi }{\partial t}\, \frac{h_s\,ds}{h_\phi B}=\frac{\partial }{\partial f}\oint E_s\,ds. \end{equation} \tag{ 14 }$$

The integral on the right-hand side is the net voltage induced along the magnetic field line,⁸

$$\begin{equation} \Phi _e=\oint E_s\,ds=\oint E_\parallel \,dl, \end{equation} \tag{ 15 }$$

where $E_\parallel ={\,{\mathbf E}}\cdot \hat{{\mathbf e}}_s$ and dl = h_sds is the length element along the field line. Equation (14) shows that the twist evolution is controlled by the longitudinal voltage, as expected. Using B_ϕ = 2I(f, t)/ch_ϕ (Equation (4)), we rewrite Equation (14) as

$$\begin{equation} \frac{1}{2\pi c}\,\oint \frac{1}{I}\,\frac{\partial I}{\partial t}\, \frac{B_\phi \,dl}{B\,h_\phi }=\frac{\partial \Phi _e}{\partial f}. \end{equation} \tag{ 16 }$$

Note that ∂B_ϕ/∂t = 0 and ∂I/∂t = 0 at r < r_c since the magnetic field remains static inside the static ideally conducting inner region. ∂I/∂t jumps from 0 to its value in the force-free region at r ≈ r_c, in a transition layer of thickness Δ ≪ r_c. The contribution from this layer to the integral on the left-hand side of Equation (16) is small and we neglect it. Then, effectively, the integral is taken only along the force-free part of the field line at r > r_c. The current I and its time derivative ∂I/∂t are constant along this part of the field line. Then, using the expression for the twist angle (Equation (3)), we find

$$\begin{equation} \frac{\psi }{2\pi cI}\,\frac{\partial I}{\partial t}\, =\frac{\partial \Phi _e}{\partial f}. \end{equation} \tag{ 17 }$$

Equation (17) describes the evolution of I(f, t) and the corresponding ψ(f, t). The twist is changing with time as the magnetic field lines gradually slip in the resistive magnetosphere and connect new footpoints with different ϕ_s − ϕ_n. Thus, ψ is changing despite the fact that the magnetosphere remains anchored in the static deep crust.

For small twists ψ < 1, the magnetic field can be written as  B ≈  B₀ +  B_ϕ(t) (Section 2), where poloidal field  B_p =  B₀ remains static and B_ϕ < B₀. Then $B=B_0[1+{\cal O}(B_\phi ^2/B^2)]\approx {\rm const}$ and the twist evolution equation (16) simplifies to

$$\begin{equation} \frac{1}{2\pi c}\,\frac{\partial \psi }{\partial t} =\frac{\partial \Phi _e}{\partial f}. \end{equation} \tag{ 18 }$$

This approximate equation quickly becomes accurate for B_ϕ ≪ B: its error is decreasing as (B_ϕ/B)². In the following sections we will use Equation (18) to study the evolution of the twist amplitude ψ. The exact evolution equation (17) would have to be used when the effect of the twist on  B_p is of interest (e.g., for calculations of the spindown rate of the star).

The linearized description of twisted configurations is useful for the first calculations of the resistive untwisting. However, the limitations of this approximation should be kept in mind. The linearized description may fail at large distances from the star (Low 1986). Besides, the force-free configuration may be unable to smoothly adjust to the growing twist. Sudden relaxation to a new topological configuration is possible, with partial opening of the field lines and the loss of connectivity between the footpoints. The twist evolution equation derived above (linear or nonlinear) does not describe such transitions.

The energy of a linear twist B_ϕ < B₀ is given by

$$\begin{eqnarray} E_{\rm tw}\equiv \int _{r>r_c} \frac{B_\phi ^2}{8\pi }\, dV &=&\frac{1}{8\pi } \int \int \int \frac{2I}{ch_\phi }\,B_\phi \,\sqrt{g}\,ds\,df d\phi \nonumber\\ &=&\frac{1}{4\pi c}\int _0^{f_c} I(f)\,\psi (f)\,df. \end{eqnarray} \tag{ 19 }$$

Differentiating with respect to time and using Equations (17) and (18), we get

$$\begin{equation} \frac{dE_{\rm tw}}{dt}=\int _0^{f_c} I\,\frac{\partial \Phi _e}{\partial f}\,df. \end{equation} \tag{ 20 }$$

Integrating by parts and taking into account that I(0) = 0 and Φ_e(f_c) = 0 (E_∥ = 0 on flux surfaces confined to the perfect conductor), we find

$$\begin{equation} \frac{dE_{\rm tw}}{dt} =-\int \Phi _e \,dI. \end{equation} \tag{ 21 }$$

The right-hand side represents the net Ohmic losses. E_tw is the free energy of the twist that is gradually dissipated as the magnetic field evolves according to Equation (18).

The evolution equation derived above assumes that the magnetosphere is anchored in the deep static crust following a starquake. It is straightforward to generalize Equation (18) for magnetospheres with moving footpoints,

$$\begin{equation} \frac{\partial \psi }{\partial t} =2\pi c\,\frac{\partial \Phi _e}{\partial f}+\omega (f,t), \end{equation} \tag{ 22 }$$

where ω = dϕ_s/dt − dϕ_n/dt is the differential angular velocity of the northern and southern footpoints of the magnetic field lines. A crust that remains motionless at all times except a sudden starquake at t = 0 is described by

$$\begin{equation} \omega (f,t)=\psi _0(f)\,\delta (t), \end{equation} \tag{ 23 }$$

where δ(t) is the Dirac function and ψ₀ is the amplitude of the twist imparted by the starquake. The impulsive twisting is a good approximation for recurring starquakes if the time between subsequent starquakes is longer than the timescale of untwisting. In the opposite limit, the crust is frequently deformed my mini-starquakes or moves plastically. Then the footpoint motion can be described by a continuous function ω(f, t). In this paper, we focus on the case of impulsive twisting described by Equation (23).

For the study of the untwisting mechanism in the next section, it is useful to consider a concrete simple magnetic configuration. We will consider a dipole field with the symmetry axis passing through the center of the star. Let $\vec{\mu }$ be the dipole moment. The poloidal flux function for the dipole is given by (Appendix A),

$$\begin{equation} f(r,\theta)=2\pi \mu \,\frac{\sin ^2\theta }{r} =\frac{2\pi \mu }{R_{\rm max}}, \end{equation} \tag{ 24 }$$

where R_max = r/sin²θ is the maximum radius reached by the flux surface passing through given r, θ. Hereafter, instead of f or R_max we will label flux surfaces by the dimensionless coordinate

$$\begin{equation} u\equiv \frac{f}{f_R}=\frac{R}{R_{\rm max}}. \end{equation} \tag{ 25 }$$

The last flux surface in the force-free region (marginally emerging from the inner crust: R_max = r_c) has u = u_c = R/r_c ≈ 1.1. The region 0 < u < 1 corresponds to the magnetospheric flux surfaces. In this region, u = sin²θ₁, where θ₁ is the polar angle of the northern footprint of the field line on the star surface r = R.

Suppose now that the field has been twisted by a starquake that was symmetric about the dipole axis and resulted in a differential rotation of the crust through angle ψ(u) (Section 2). The poloidal current I and twist angle ψ will be viewed below as functions of u and time t. The following relation holds between ψ and I,

$$\begin{equation} \psi =\frac{4R^2I}{c\mu u^2}\sqrt{1-\frac{u}{u_c}}, \end{equation} \tag{ 26 }$$

(see Appendix A). The twist evolution equation (18) becomes

$$\begin{equation} \frac{\partial \psi }{\partial t} =\frac{cR}{\mu }\,\frac{\partial {\Phi _e}}{{\partial u}}, \end{equation} \tag{ 27 }$$

or

$$\begin{equation} \frac{\partial I}{\partial t} =\frac{c^2u^2}{4R\sqrt{1-u/u_c}}\,\frac{\partial {\Phi _e}}{{\partial u}}. \end{equation} \tag{ 28 }$$

The current density in the twisted dipole magnetosphere is given by (cf. Equation (5)),

$$\begin{equation} j=\frac{RB}{2\pi \mu }\frac{\partial I}{\partial u}. \end{equation} \tag{ 29 }$$

Consider first what would happen if the magnetosphere was filled by a medium with a fixed conductivity σ; we will assume in this toy model σ(r) = const. Then the current density is related to E_∥ by Ohm's law j = σE_∥. The star will be modeled as a perfect conductor, so E_∥ = 0 is assumed inside the star. Then the voltage Φ_e along a magnetospheric field line is given by

$$\begin{equation} \Phi _e=\int _{r>R} \frac{j}{\sigma }\,dl =\frac{\sqrt{1-u}}{\pi R\sigma }\frac{\partial I}{\partial u}. \end{equation} \tag{ 30 }$$

(The integral has been calculated using Equations (29) and (A10).) Substitution of this result to the twist evolution equation (28) yields the equation for I(u, t),

$$\begin{equation} \frac{\partial I}{\partial t} =\frac{c^2u^2}{4\pi \sigma R^2}\left(1-\frac{u}{u_c}\right)^{-1/2} \frac{\partial }{\partial u} \left(\sqrt{1-u}\,\frac{\partial I}{\partial u}\right). \end{equation} \tag{ 31 }$$

Given an initial twist with current function I(u, 0), one can calculate the evolution of I(u, t) by solving this differential equation with two boundary conditions: I(0) = 0 and I(1) = const. The latter condition is valid as long as the perfect-conductor approximation is used for the star.

Equation (31) is of diffusion type. It has a special feature: the effective diffusion coefficient is proportional to $\sqrt{1-u}$ and vanishes at u = 1. This fact allows one to qualitatively understand the evolution of I(u, t) before solving the equation numerically. It is instructive to consider the analogous problem of particle diffusion with a position-dependent diffusion coefficient D(u) that vanishes at the boundary u = 1. The vanishing of D means that particles "stick" to the boundary. With time, more and more particles get stuck, and eventually the particle density vanishes everywhere except at the boundary. The magnetospheric current behaves in a similar way. Far from the boundary, it simply spreads diffusively: ∂I/∂t = const u² ∂²I/∂u² at u ≪ 1. At the same time, near the boundary u = 1, the current is sucked toward u = 1. The current keeps accumulating at u = 1 until I = 0 at all u < 1.

For illustration, we solved numerically Equation (31) for a twist that initially has a uniform amplitude ψ₀(u) = 0.2. The corresponding initial current function I(u, 0) and ∂I/∂u ∝ j/B are shown in Figure 3. The evolution of this twist is shown in Figure 4. The characteristic diffusion timescale is t_σ = R²c²/4πσ, and we express time in units of t_σ. As expected, the current tends to spread to smaller u and, at the same time, it is quickly drawn into the current sheet at u = 1. The sum of the currents flowing in the magnetopshere and in the current sheet, I_tot, remains constant. Our numerical model assumes Φ_e = 0 at u > 1 (inside the star), which allows the current sheet to persist at u = 1. A real star has a finite conductivity and the current sheet will slowly spread into the star (toward larger u > 1).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Note that the twist amplitude ψ grows near u = 0, because ∂Φ_e/∂u > 0 in this region. A strongly twisted bundle develops near the magnetic axis; however, its thickness shrinks with time, so its net current decreases. At late times, the twist growth near the axis enters a self-similar regime. The amplitude peak ψ_peak would grow indefinitely in the limit t → ∞, if the twist remained stable. In fact, the growth must be stopped by the MHD instability expected at ψ ≳ 1. An upper limit on ψ is set by the field-line opening, which may give a partially opened configuration with a lower energy (Wolfson & Low 1992). Overtwisted configurations are prone to kink instability, which is likely to reconnect part of the twisted region away from the star. These dynamic processes are not described by our model. We shall assume that they keep ψ below $\psi _{\rm max}={\cal O}(1)$.

The toy model discussed in Section 5.1 is deficient: the real magnetospheres of neutron stars are not filled by a medium of a fixed conductivity. Instead, the magnetospheric currents are maintained through a discharge with a threshold voltage that is approximately the same for any current density j ≠ 0 (or, more precisely, for any j exceeding a small value j_⋆). The threshold nature of the magnetospheric voltage changes the evolution of the twist. Remarkably, it simplifies the solution: the partial differential equation (27) will be reduced to an ordinary differential equation.

We denote the threshold voltage by ${\cal V}$. Once Φ_e reaches ${\cal V}$, copious particle supply is available to carry any large current. The particle supply and voltage regulation in magnetars was studied in BT07. In principle, there are two sources of particles: (1) the surface of the star and (2) e^± creation in the magnetosphere. The current is carried by charges of both signs, since the net charge density must be nearly zero to avoid huge voltages. If no e^± are created, electrons and ions must be lifted from the star's surface (if ions are available in an atmospheric layer atop the solid crust). Maintaining a flow of ions along a magnetic loop requires a minimum voltage,

$$\begin{equation} e{\cal V}_{ei}=\frac{GMm_i}{R}-\frac{GMm_i}{R_{\rm max}}, \end{equation} \tag{ 32 }$$

where R_max is the maximum radius reached by the loop and m_i is the ion mass. This voltage is ∼0.2m_ic² for loops with R_max ≫ R. The plasma is lifted into the loop by the self-induction electric field (BT07). The exact solution was obtained for this process in the one-dimensional circuit model. It shows that voltage keeps growing even after reaching ${\cal V}_{ei}$, and the circuit evolves toward a global double-layer configuration. This result may not hold in the complete three-dimensional model that includes the excitation of transverse waves in the magnetosphere. Nevertheless, it suggests that lifting plasma from the surface is not the ultimate regulator of the voltage.

A robust mechanism for limiting the voltage is the e^± discharge. If Φ_e exceeds a certain threshold ${\cal V}_{\pm }$, the exponential runaway of pair creation occurs (cf. Figure 5 in BT07), and the e^± pairs screen E_∥. In a magnetar magnetosphere, the e^± avalanche is triggered when the accelerated electrons resonantly scatter X-rays streaming from the star, and the scattered photons convert to e^± off the magnetic field. The threshold for the discharge is given by

$$\begin{equation} e{\cal V}_{\pm }\sim \gamma _{\rm res} m_ec^2\sim \frac{ceB}{\omega _X} \approx 1\left(\frac{B}{10^{14}\,{\rm G}}\right) \left(\frac{\omega _X}{10^{18}\,{\rm Hz}}\right)\;{\rm GeV}, \end{equation} \tag{ 33 }$$

where ω_X is the typical frequency of target X-rays and γ_res ∼ (B/B_Q)(m_ec²/ℏω_X) is the electron Lorentz factor at which the electron begins to scatter >1 X-rays as it travels along the magnetic loop. Here B_Q = m²_ec³/eℏ ≈ 4.4 × 10¹³ G. Numerical experiments in BT07 show that the e^± discharge is intermittent on a timescale ∼r/c and proceeds in the regime of self-organized criticality. The time-averaged current equals the current imposed by the magnetospheric twist, (c/4π)∇ ×  B, and the time-averaged voltage Φ_e is close to ${\cal V}_{\pm }$.

The discharge voltage remains almost independent of the imposed current j unless j is reduced below j_⋆. The value of j_⋆ is unknown but small. It may be comparable to cρ_GJ, where ρ_GJ is the corotation charge density (Goldreich & Julian 1969) that should be maintained in the magnetosphere in the absence of electric currents. For the purposes of the present paper, the following description for Φ_e(j) will be sufficient:

$$\begin{equation} \Phi _e(j)=\left\lbrace \begin{array}{@{}ll@{}} {\cal V}& \mbox{if $j\gg j_\star $} \\ 0 & \mbox{if $j\ll j_\star $} \end{array} \right. \end{equation} \tag{ 34 }$$

where j_⋆ is much smaller than the characteristic currents induced by the starquake. As will be shown below, neither the value of j_⋆ nor the behavior of Φ_e(j) in the transition region j ∼ j_⋆ matters. In essence, $\Phi _e={\cal V}(u)\Theta (j)$, where Θ is the Heaviside step function. In computer simulations, we use a smoothed step function,

$$\begin{eqnarray} \Phi _e &=&{\cal V}(u)\;W\left(\frac{j}{j_\star }\right), \nonumber\\ W(x)&=&\frac{\arctan (1/\Delta)+\arctan \left[(x-1)/\Delta \right]}{\arctan (1/\Delta)+\pi /2}, \end{eqnarray} \tag{ 35 }$$

where Δ ≪ 1. Note that the discharge voltage ${\cal V}$ can be different for different flux surfaces, i.e., ${\cal V}$ in general depends on u. There is a sharp drop in ${\cal V}(u)$ at u = 1 (voltage is small inside the highly conducting star). In computer simulations, we model this drop by introducing the factor exp{[/(1 − u)]¹⁰} with ≪ 1. The exact form of this factor plays no role for the twist dynamics.

Consider again the twist shown in Figure 3, and let us calculate its evolution with the new threshold relation between j and Φ_e (Equation (34)). The numerical solution of Equations (28), (29), and (35) is shown in Figure 5 for the simplest model that assumes ${\cal V}(u)={\rm const}$ in the magnetosphere and ${\cal V}(u)=0$ inside the star.

Zoom In Zoom Out Reset image size

Initially, the drop in Φ_e near u = 1 is very steep, i.e., ∂Φ_e/∂u is large and negative, and hence the current density here is quickly reduced (cf. Equation (27)). The reduction continues until j ∼ j_⋆, which permits $\Phi _e<{\cal V}$. Then a smooth profile of $0<\Phi _e(u)<{\cal V}$ is established in a region u_⋆(t) < u < 1. The threshold nature of the discharge leads to the formation of two distinct regions in the magnetosphere:

• 1.  
  "Cavity" u_⋆ < u < 1 where $0<\Phi _e<{\cal V}$ and j ∼ j_⋆ (essentially zero current). The current originally injected in this region is sucked into the star and flows in the current sheet at u = 1.
• 2.  
  Region u < u_⋆ where $\Phi _e={\cal V}$. Here the current remains equal to its initial value at t = 0, i.e., the twist remains static.⁹

The voltage profile in the cavity has ∂Φ_e/∂u < 0 which implies ∂I/∂t < 0 (Equation (28)). Hence the cavity must grow, as indeed seen in the simulation. The structure of the untwisting magnetosphere is shown in Figure 6.

Zoom In Zoom Out Reset image size

The boundary of the cavity forms a sharp front, which resembles a shock wave. The front starts as a tiny arc emerging from the star at the magnetic equator, and propagates toward smaller u (i.e., outward and poleward, to flux surfaces with larger R_max). We denote its instantaneous position by u_⋆(t). The profile of the front—the shape of j(u) near u_⋆—is controlled by the behavior of Φ_e at j ∼ j_⋆.¹⁰ However, we are not interested in the exact profile of the front. It is sufficient to know that it is steep and can be treated as a step function. The quantity of interest is the speed of the front propagation du_⋆/dt. In the limit of small j_⋆ (steep front), one can use

$$\begin{equation} \Phi _e={\cal V}(u)\Theta (j), \end{equation} \tag{ 36 }$$

and derive an explicit expression for du_⋆/dt (see Appendix B),

$$\begin{equation} \displaystyle \frac{du_{\star }}{dt}=-\frac{\displaystyle {\cal V}(u_{\star })\,\frac{u_{\star }(1+1.4u_{\star })}{1-u_{\star }} +\frac{u_{\star }^2{\cal V}^\prime (u_{\star })}{\sqrt{1-u_{\star }/u_c}}}{\displaystyle \frac{4R}{c^2}\left.\frac{dI_0}{du}\right|_{u_{\star }} +\left.\frac{d}{du}\frac{u^2{\cal V}^\prime (u)}{\sqrt{1-u/u_c}}\right|_{u_{\star }} t}, \end{equation} \tag{ 37 }$$

where ${\cal V}^\prime =d{\cal V}/du$ and I₀(u) ≡ I(u, 0) is the initial current function. This ordinary differential equation can be solved for u_⋆(t). If ${\cal V}^\prime =0$ (as in the model in Figure 5) the front equation simplifies to

$$\begin{equation} \frac{du_{\star }}{dt}=-\frac{c^2{\cal V}\,u_{\star }(1+1.4u_{\star })}{4R\,(1-u_{\star })(dI_0/du)_{u_{\star }}}. \end{equation} \tag{ 38 }$$

The history of the front propagation is shown in Figure 7. The front starts with extremely high speed near u = 1, then decelerates and approaches the axis u = 0 with $du_{\star }/dt=-cR{\cal V}/2\mu \psi _0$. It reaches u = 0 (and erases all of the twist) at $t_{\rm end}=\mu \psi _0/cR{\cal V}$.

Zoom In Zoom Out Reset image size

The current function I of the evolving twist is given by Equation (B8) in Appendix B. The corresponding twist amplitude is given by

$$\begin{equation} \psi (u,t)=\left\lbrace \begin{array}{@{}ll@{}} \displaystyle \psi _0(u)+\frac{cR}{\mu }\,{\cal V}^\prime (u)\,t & 0<u<u_{\star }\\ \displaystyle \frac{4R^2}{c\mu u^2}\sqrt{1-\frac{u}{u_c}}\,I_{\star }(t) & u_{\star }<u<1 \\ \end{array} \right. \end{equation} \tag{ 39 }$$

where I_⋆(t) ≡ I[u_⋆(t), t] is the net current that flows through the magnetosphere at time t. Equation (39) together with u_⋆(t) gives a complete analytical description for untwisting magnetospheres (in the linear-twist approximation; see Section 4.2).

Note that the numerical simulation shown in Figure 5 gives $\Phi _e/{\cal V}<1$ for u ≪ 1, i.e., for flux surfaces with footpoints near the magnetic axis. This is not predicted by the analytical model, which assumes j_⋆ = 0. The drop appears in the numerical simulation at a finite u ≪ 1 because the numerical model assumes a finite j_⋆. The twist with ψ₀(u) = const has the current density near the axis j(u) ∝ u² and hence, in a small polar region where j ≲ j_⋆, Φ_e must drop. In the limit j_⋆ → 0, this effect disappears in the sense that the polar region with j ∼ j_⋆ shrinks to one point u = 0.

The twist behavior on the axis becomes important when the spindown of the star is of interest. The model in Figure 5 assumes a finite j_⋆; however, it does not take into account rotation of the star; therefore j(0) = 0 and Φ_e(0) = 0. Rotation with angular velocity Ω implies the additional twisting and opening of the field lines that extend to the light cylinder R_lc = c/Ω. This persistent "external" twisting induces small but finite currents on the magnetic dipole axis, along the bundle of open field lines. The open bundle has the parameter u_lc = R/R_lc ≪ 1. In particular, observed magnetars have u_lc ≲ 10⁻⁴. As long as the behavior of the closed magnetosphere u > u_lc is concerned, the rotation can be neglected, and j_⋆ → 0 is a good approximation.

As the cavity expands, the current-carrying region u < u_⋆ shrinks. We will call the current-carrying bundle of magnetic field lines "j-bundle," for brevity. The numerical model of Section 5.3 (Figures 5 and 7) assumed that the discharge voltage is the same for all magnetospheric flux surfaces, ${\cal V}(u)={\rm const}$. Then the twist remains static inside the j-bundle. It freezes and waits while it is eaten by the expanding front u_⋆(t).

In contrast, if ${\cal V}(u)\ne {\rm const}$, the twist in the j-bundle will change linearly with time as it waits for the front to come. This change is described by Equation (27) (note that $\Phi _e={\cal V}(u)$ inside the bundle u < u_⋆) or Equation (39). The twist amplitude ψ decreases if ${\cal V}^\prime <0$ and grows if ${\cal V}^\prime >0$. Observational data (discussed below) suggest ${\cal V}^\prime >0$ and the growth of ψ near the axis u = 0. Despite the twist growth at small u, its total energy E_tw is decreasing with time as the j-bundle shrinks. This evolution is consistent with the energy conservation law (Equation (21)).

For illustration, we calculated the same model as in Figure 5 but with new threshold voltage ${\cal V}(u)=(0.04+2u)^{1/2}\bar{{\cal V}}$, where $\bar{{\cal V}}$ approximately equals the average of ${\cal V}(u)$ in the magnetosphere 0 < u < 1. Figure 8 shows the evolution of the twist amplitude ψ in this case and compares it with the case of ${\cal V}(u)={\rm const}$.

Zoom In Zoom Out Reset image size

The twist evolution is described by simple analytical formulae when the j-bundle is narrow, u_⋆ ≪ 1. In the leading order of u_⋆ ≪ 1 the front equation becomes

$$\begin{equation} \frac{du_{\star }}{dt}=-\frac{u_{\star }{\cal V}(u_{\star })}{(4R/c^2)\,I_0^\prime (u_{\star }) +[u^2{\cal V}^\prime ]^\prime |_{u_{\star }} t} [1+{\cal O}(u_{\star })], \quad u_{\star }\ll 1, \end{equation} \tag{ 40 }$$

where prime denotes the differentiation with respect to u. The relation between I and ψ (Equation (26)) gives $I_0(u)=(c\mu \psi _0/4R^2)u^2+{\cal O}(u^3)$ where ψ₀(u) = ψ(u, 0) is the initial amplitude of the twist. Then we get,

$$\begin{equation} \frac{du_{\star }}{dt}=-\frac{{\cal V}(0)}{[2\mu \psi _0(0)/cR]+2{\cal V}^\prime (0)\,t}\,\left[1+{\cal O}(u_{\star })\right]. \end{equation} \tag{ 41 }$$

If ${\cal V}^\prime (0)\ne 0$, integration of this equation yields

$$\begin{equation} u_{\star }(t)=\frac{{\cal V}(0)}{2{\cal V}^\prime (0)}\,\ln \frac{t_0+t_{\rm end}}{t_0+t}, \qquad t_0\equiv \frac{\mu \psi _0(0)}{cR{\cal V}^\prime (0)}, \end{equation} \tag{ 42 }$$

where t_end is the time when the front reaches u = 0; it is finite in all cases.

The twist amplitude ψ(u, t) inside the j-bundle grows according to Equation (39) (unless ψ reaches ψ_max ≳ 1),

$$\begin{equation} \psi (0,t)=\psi (0,0)\left(1+\frac{t}{t_0}\right), \qquad 0<t<t_{\rm end}. \end{equation} \tag{ 43 }$$

It grows by a large factor if t_end ≫ t₀.

For example, consider a model with linear ${\cal V}(u)={\cal V}(0)+{\cal V}^\prime (0)u$ at $u<\hat{u}$ and constant ${\cal V}(u)={\cal V}(\hat{u})$ at $u>\hat{u}$. Suppose ${\cal V}(\hat{u})\gg {\cal V}(0)$. Then we find,

$$\begin{equation} \frac{\psi (0,t_{\rm end})}{\psi (0,0)} \sim \exp \left[2\left(\frac{{\cal V}(\hat{u})}{{\cal V}(0)}-1\right)\right]\gg 1. \end{equation} \tag{ 44 }$$

The drop in ${\cal V}(u_{\star })$ at small u_⋆ delays the arrival of the front to u = 0 and gives an exponentially longer time for the twist growth at u = 0. Thus, even a small twist with ψ₀ ≪ 1 can grow to ψ_max inside the j-bundle. Further growth is impeded by the MHD instability.

Starquakes may rotate part of the crust and leave the rest of it untouched. Suppose that a ring u₂ < u < u₁ has been rotated. Then the twist ψ ≠ 0 is created only in the region u₂ < u < u₁. This implies I(u > u₁) = 0, i.e. the net ejected current is zero, and hence the current density must change sign at some u_m in the region u₂ < u < u₁. The current function I₀(u) reaches a maximum at u_m. An example of such a localized twist and its evolution are shown in Figure 9. Voltage ${\cal V}(u)=(\frac{1}{2}+u)\bar{{\cal V}}\ne {\rm const}$ is assumed in the model.

Zoom In Zoom Out Reset image size

Three stages may be noted in the evolution of the ring twist.

• 1.  
  One current front is launched from u = 1. It immediately jumps to the twist boundary u₁, and continues to propagate toward smaller u. In addition, two divergent current fronts are immediately launched from u = u_m. Thus, two cavities form in the magnetosphere. Both cavities grow until they merge: the two fronts erasing the spike of negative current (Figure 9) eventually meet and "annihilate." Stage 1 ends at this point (at time $t\approx 0.005t_{\cal V}$ where $t_{\cal V}=\mu /cR\bar{{\cal V}}$).
• 2.  
  The merged cavity continues to expand poleward and erase the remaining spike of positive current at smaller u. All of the initially injected current is erased (and stage 2 ends) at $t\approx 0.05t_{\cal V}$.
• 3.  
  The current front proceeds toward the axis, erasing the currents that have grown there (from zero) since the beginning of the twist evolution.

Voltage ${\cal V}(u)\ne {\rm const}$ was chosen in the model to allow the twist growth near the axis. The simulation shows, however, that the growth is slow compared to the expansion of the cavity. Only at the last stage (stage 3) does the grown current create a significant twist ψ near the axis and somewhat decelerate the expansion of the cavity.

Figure 10 shows the evolution of the twist luminosity L(t). Its decrease is quickest during stage 1, as the spike of negative current is quickly erased. Figure 10 also shows the propagation of the current front u_⋆(t) that starts at u_m and moves to smaller u, erasing the ejected positive current.

Zoom In Zoom Out Reset image size

Instead of a ring u₂ < u < u₁, an axisymmetric localized starquake may rotate a cap u < u₁. The cap-twist and ring-twist evolve in a similar way. In particular, stages I and II are similar. However, stage 3 is absent for the cap-twist. It has u₂ = 0, i.e., the initially ejected positive currents occupy the entire region around the polar axis. Erasing these currents takes a longer time. As a result, the twist grows to a large amplitude at u ≪ 1 before the cavity reaches the axis (Section 5.4). Then a maximally twisted narrow bundle with ψ = ψ_max forms.

Energy dissipation in the untwisting magnetosphere is confined to the bundle of current-carrying field lines u < u_⋆ (j-bundle). It can be a bright source of radiation, with luminosity equal to the rate of Ohmic dissipation L. The luminosity generally decreases as the magnetosphere untwists. Examples of the evolution of L(t) are shown in Figures 7 and 10.

A large fraction of the dissipated power may be radiated quasi-thermally at the footprints of the j-bundle as the accelerated magnetospheric particles run into the star (BT07), creating a hot spot on the surface θ < θ_⋆. The area of this spot is given by

$$\begin{equation} A\approx \pi (R\sin \theta _\star)^2=\pi R^2u_{\star }. \end{equation} \tag{ 45 }$$

As the cavity expands in an untwisting magnetosphere (Figure 6), the spot shrinks.

The evolution of a narrow j-bundle (u_⋆ ≪ 1) with a uniform twist ψ (e.g., ψ ≈ ψ_max ∼ 1) is described by simple formulae. The free energy of the twist (Equation (19)) is then given by¹¹

$$\begin{eqnarray} E_{\rm tw}&=&\frac{\mu }{2cR}\int _0^1 I(u)\,\psi (u)\,du \approx \frac{\mu ^2\psi ^2u_{\star }^3}{6R^3} \nonumber\\ &\approx & 4\times 10^{44}\,B_{14}^2\,R_6^3\,\psi ^2\,u_{\star }^3\;{\rm erg}, \end{eqnarray} \tag{ 46 }$$

where R₆ ≡ R/10⁶ cm, B₁₄ ≡ B_pole/10¹⁴ G, and B_pole ≡ 2μ/R³. The luminosity of the j-bundle with ψ ≈ const can be immediately calculated for a given voltage ${\cal V}(u)$,

$$\begin{equation} L=\int _0^{I_{\star }}{\cal V}\,dI =\frac{c\mu \psi }{4R^2}\int _0^{u_{\star }} {\cal V}(u)\, d\left(\frac{u^2}{\sqrt{1-u/u_c}}\right). \end{equation} \tag{ 47 }$$

The simplest model with ${\cal V}(u)={\rm const}$ gives

$$\begin{equation} L={\cal V}I_{\star }\approx \frac{c\mu }{4R^2}\,\psi \,{\cal V}\,u_{\star }^2 \approx 1.3\times 10^{36}\,B_{14}\,R_6\,\psi \,{\cal V}_9\,u_{\star }^2\;{\rm erg\; s}^{-1}. \end{equation} \tag{ 48 }$$

where ${\cal V}_9\equiv {\cal V}/10^9$ V. The j-bundle with ${\cal V}(u)={\rm const}$ shrinks according to Equation (38), which yields at u_⋆ ≪ 1

$$\begin{equation} \frac{du_{\star }}{dt}\approx -\frac{cR{\cal V}}{2\mu \psi } \quad \Rightarrow \quad u_{\star }(t)\approx \frac{cR{\cal V}}{2\mu \psi }\left(t_{\rm end}-t\right). \end{equation} \tag{ 49 }$$

This equation is also easy to derive from dE_tw/dt = −L, using Equations (46) and (48). The evolution timescale of the luminosity is given by

$$\begin{equation} t_{\rm ev}=-\frac{L}{dL/dt}\approx \frac{\mu u_{\star }}{cR{\cal V}} \approx 15\,{\cal V}_9^{-1}\,B_{14}\,R_6^2\,\psi \,u_{\star }\;{\rm yr}. \end{equation} \tag{ 50 }$$

Approximate formulae can also be derived for ${\cal V}(u)\ne {\rm const}$. For example, for ${\cal V}(u)={\cal V}(0)+{\cal V}^\prime (0)u$ we find

$$\begin{equation} L=\frac{c\mu {\cal V}(0)}{4R^2}\,\psi \,u_{\star }^2\left[1+ \left(\frac{2{\cal V}^\prime (0)}{3{\cal V}(0)}+\frac{1}{2u_c}\right)u_{\star }+{\cal O}\left(u_{\star }^2\right)\right]. \end{equation} \tag{ 51 }$$

This equation again assumes ψ(u < u_⋆) ≈ const. It may approximately describe e.g., the maximally twisted j-bundle, ψ ≈ ψ_max.

The decay of the twist luminosity L(t) was previously estimated assuming that the current is decaying uniformly in the twisted region, which gave a linear L(t) ∝ t − t_end (BT07). The electrodynamic theory developed in this paper shows that the untwisting is strongly non-uniform: the twist is erased by the propagating front that resembles a shock wave. The speed of this front depends on the initial twist configuration ψ₀(u). Similar to the simple estimate L(t) ∝ t − t_end, we find that the twist is erased in a finite time and L(t) vanishes at t_end (unless new starquakes occur). However, no universal linear shape of L(t) is predicted. In some cases L(t) may resemble a linear decay (e.g., segment II in Figure 10 is almost linear in a linear plot). In observed sources, L(t) was close to linear in AXP 1E 1048.1−5937 (Dib et al. 2009) and nonlinear in XTE J1810−197 (Gotthelf & Halpern 2007).

The energy released in the j-bundle can power nonthermal magnetospheric emission. Such emission is observed in most magnetars. Two distinct nonthermal components are detected in their spectra: (1) soft X-ray tail that extends from 1 keV to ∼10–20 keV with a photon index Γ ∼ 2–4 (e.g., Woods & Thompson 2006 and references therein) and (2) hard X-ray component that extends to ∼300 keV with Γ ∼ 0.8–1.5 (e.g., Kuiper et al. 2008 and references therein).

The 1–20 keV tail is usually explained by resonant scattering of thermal radiation by the magnetospheric plasma (Thompson et al. 2002; Lyutikov & Gavriil 2006; Fernández & Thompson 2007; Nobili et al. 2008; Rea et al. 2008). Ions resonantly scatter thermal photons near the star where ℏeB/m_ic∼ keV, and e^± scatter at radii r ∼ 10R where ℏeB/m_ec∼ keV. The growth of cavity in the untwisting magnetosphere implies that the scattering in the inner magnetosphere is suppressed, because the dense plasma is confined to the narrow j-bundle. At larger radii r ∼ R/u_⋆, the j-bundle broadens and forms an outer corona that subtends a large solid angle as viewed from the star. The cyclotron energy in this region is

$$\begin{equation} \hbar \,\frac{eB}{m_ec}\sim 1 \left(\frac{B_{\rm pole}}{10^{14}\,{\rm G}}\right) \left(\frac{u_{\star }}{0.1}\right)^3 {\rm keV}, \end{equation} \tag{ 52 }$$

and resonant scattering by e^± can give 1–20 keV photons. The luminosity expected from a strongly twisted j-bundle with u_⋆ ∼ 0.1 is consistent with the typical nonthermal luminosity of magnetars, L ∼ 10³⁵ erg s⁻¹ (see Equation (48) and substitute the typical μ ∼ 3 × 10³² G cm³ and ψ = 1–2). If no new starquakes occur, the j-bundle must shrink toward the magnetic axis with time, and u_⋆ will be reduced below 0.1. Then the resonant scattering must be suppressed. This suppression is caused by two reasons: ℏeB/m_ec in the outer corona r ∼ R/u_⋆ decreases below keV (it is proportional to u³_⋆), and the power dissipated in the j-bundle becomes small as L ∝ u²_⋆.

The j-bundle with u_⋆ ∼ 0.1–0.2 dissipates sufficient energy to explain also the hard X-ray component. This component was detected in three AXPs and two soft gamma-ray repeaters (SGRs; see Kuiper et al. 2008 for a recent review). In AXPs, the hard X-ray emission has a huge pulsed fraction, approaching 100% at high energies (Kuiper et al. 2006; den Hartog et al. 2008a, 2008b). This may be explained if the emission is produced in the narrow j-bundle near the star.

Observations indicate that the nonthermal emission can be stable on timescales as long as a decade (den Hartog et al. 2008a). There may be two reasons for this stability. (1) The discharge voltage ${\cal V}$ is relatively low (below 1 GeV) and the j-bundle is relatively thick, u_⋆ ∼ 0.2. Then the untwisting timescale becomes long (see Equation (50)). (2) The j-bundle is kept in a quasi-steady, maximally twisted state ψ ∼ ψ_max by frequently repeating (possibly continual) shearing motion of the crust in a fixed region u < u_⋆.

The plasma filling the j-bundle may also produce optical and infrared radiation by mechanisms discussed in BT07. Besides, it can be a bright source of radio waves, as suggested in Section 7.2.

The magnetospheric twist is expected to impact the spindown rate when the twist amplitude is large, ψ ≳ 1. This impact may occur in two ways.

• 1.  
  The strong twist inflates the poloidal field lines and increases the magnetic field at the light cylinder, which leads to stronger spindown torque acting on the star (Thompson et al. 2002). The poloidal inflation is common for twisted configurations. It is seen, e.g., in the self-similarly twisted dipole (Wolfson 1995). A similar inflation must occur when the currents are confined to the j-bundle. Its calculation will require a full nonlinear model. Here we limit our consideration to simple estimates.One can think of poloidal inflation as an increase of magnetic dipole moment with radius. This effect is easiest to evaluate for moderate twists ψ < 1. The current dI flowing along a bundle of twisted field lines du creates a toroidal current dI_ϕ ≈ dI ψ(u)/2π. The field lines carrying this current extend to R_max = R/u, and the dipole moment created by dI_ϕ is dμ ∼ dI_ϕR²_max/c. Integrating over u, we obtain the net change of dipole moment of the star due to the twist

  $$\begin{equation} \Delta \mu \sim \frac{R^2}{2\pi c}\int \frac{\psi (u)}{u^2}\,dI. \end{equation} \tag{ 53 }$$

  The j-bundle u_lc < u < u_⋆ with a uniform twist ψ increases the dipole moment of the star by

  $$\begin{equation} \frac{\Delta \mu }{\mu }\sim \frac{\psi ^2}{4\pi }\,\ln \frac{u_{\star }}{u_{\rm lc}}. \end{equation} \tag{ 54 }$$

  This effect is quadratic in ψ and quickly becomes small for ψ < 1. For strong twists ψ ∼ ψ_max, the estimate (54) must be replaced by a full nonlinear calculation. Qualitatively, it suggests that Δμ/μ is reduced as the j-bundle shrinks (u_⋆ decreases). Therefore, the spindown torque is expected to be reduced with time and gradually come back to the standard dipole torque as u_⋆ → u_lc.
• 2.  
  When the twist has grown to ψ_max, it will begin to "boil over" through a repeated instability. The energy that would be stored in the magnetosphere if ψ kept growing above ψ_max is then carried away by an intermittent magnetic outflow. The outflow may also carry away a significant angular momentum.The outflow can be generated where the overtwisted field lines open up. This opening occurs because the twist is constantly pumped near the axis by ${\cal V}^\prime >0$. A possible structure of untwisting magnetospheres with ${\cal V}^\prime >0$ is schematically shown in Figure 11. It resembles the picture of a rotationally powered pulsar (see, e.g., Arons 2008 for a review); however, the opening is caused by the internal resistive dynamics of the magnetosphere and depends on the profile of ${\cal V}(u)$. By contrast, in ordinary pulsars the twist is pumped by the star rotation.

Zoom In Zoom Out Reset image size

Mechanisms (1) and (2) start immediately after the starquake if it implants a strong initial twist ψ₀ > 1 near the magnetic axis. If ψ₀ < 1, the spindown torque may not be affected until ψ grows to ∼1, which takes time (cf. Equation (27) or Equation (39))

$$\begin{equation} t_{\rm delay}\approx \frac{\mu }{cR{\cal V}^\prime }, \end{equation} \tag{ 55 }$$

where ${\cal V}^\prime =d{\cal V}/du$ is evaluated near the axis u ≈ 0; ${\cal V}^\prime$ may be large, because u is small. For example, a change in ${\cal V}$ from ${\cal V}=10^9$ V at u = 0 to 2 × 10⁹ V at u = 10⁻² corresponds to ${\cal V}^\prime =10^{11}$ V. The delay is observed in some objects, with a characteristic t_delay ∼ 10⁷ s (e.g., Gavriil & Kaspi 2004). This requires ${\cal V}^\prime \sim 10^{11}$ V.

In contrast to luminosity L(t) (which always tends to decrease after the starquake), the torque behavior is generally non-monotonic. If ψ₀ < 1, the torque is expected to grow as ψ grows in the shrinking j-bundle. (Such an anti-correlation between the torque and the X-ray luminosity was observed in 1E 1048.1−5937; see Gavriil & Kaspi 2004). Once ψ has reached ψ_max ∼ 1 the torque should start to decrease, because the shrinking of the j-bundle at constant ψ ≈ ψ_max leads to the reduction of Δμ. The power of the outflow from the magnetosphere is also expected to decrease. The simultaneous decrease in torque and luminosity was observed, e.g., in XTE J1810−197 (Camilo et al. 2007).

Useful preliminary estimates can be made if we assume a uniform voltage ${\cal V}(u)\approx {\rm const}$ and uniform twist ψ(u) ≈ const across the j-bundle u < u_⋆ (Section 6.1). Then the produced luminosity L(t) is given by Equation (48). Suppose a large fraction of L is emitted thermally at the footprint of the j-bundle. The estimates for L (Equation (48)) and its evolution timescale t_ev (Equation (50)) for a given spot area A (Equation (45)) may be compared with observations. For example, in the fall of 2004, the spot area was A ≈ 10¹¹ cm (and hence u_⋆ ≈ 0.03). One then finds that a strongly twisted j-bundle (ψ ∼ 1–1.5) with voltage ${\cal V}_9\sim 3$ explains both observed L ≈ 2 × 10³⁴ erg s⁻¹ and t_ev ≈ 0.6 yr.¹²

The model ${\cal V}(u)\approx {\rm const}$ gives good estimates for L, A, and t_ev; however, it has a drawback: it is unable to describe the possible growth of ψ from a smaller ψ₀ before the j-bundle became maximally twisted. The growth occurs if ${\cal V}^\prime >0$ (Section 5.4). Therefore, we adopt a slightly more general model that includes the next (linear) term in the expansion of ${\cal V}(u)$ near u = 0,

$$\begin{equation} {\cal V}(u)\approx {\cal V}_0+{\cal V}^\prime \,u, \qquad u\ll 1. \end{equation} \tag{ 56 }$$

${\cal V}^\prime$ is unknown and probably large near the axis (see the text after Equation (55)).

A simplest model of a twisted magnetosphere has four parameters: magnetic dipole moment of the star μ, radius of the star R, the initial size of the j-bundle created by the starquake u₀ = sin²θ₀, and the initial amplitude of the twist ψ₀ (it equals the angular displacement of the crustal cap rotated by the starquake; in our fiducial model, the starquake imparts a uniform twist ψ₀(u < u₀) = const). The twist evolution after the starquake is controlled by two parameters ${\cal V}_0$ and ${\cal V}^\prime$ that specify voltage ${\cal V}(u)$ in Equation (56). The evolution is described by Equations (39) and (37) that we solve numerically. If ψ reaches $\psi _{\rm max}={\cal O}(1)$, it is assumed to stall at ψ_max.¹³ The luminosity of a maximally twisted j-bundle, L(u_⋆), is approximately given by Equation (51).

Figure 12 compares the model with observations of XTE J1810−197. The star is assumed to have μ = 1.5 × 10³² G cm³ and R = 9 km. The parameters of the starquake are u₀ = 0.15 and ψ₀ = 0.5. The discharge voltage is assumed to drop linearly from 5.5 GeV at u = 0.15 to 1 GeV at u = 0 (which corresponds to ${\cal V}_0=10^9$ V and ${\cal V}^\prime =3\times 10^{10}$ V). Although we did not attempt a formal fitting of the data, the figure suggests that the model is successful in explaining the evolution of L(t) and A(t). The data firmly constrain the voltage to be in 1–6 GeV range, which is close to the theoretical estimate (Equation (33)).

The accuracy of the model is limited by several idealizing assumptions: uniform initial twist in the starquake region, the linear form of ${\cal V}(u)$, and axial symmetry. The actual value of μ is probably smaller than inferred from spindown measurements (Section 6.3). The linear twist-evolution equation becomes approximate when ψ ∼ 1. Note also that our simplest fiducial model assumes that all energy dissipated in the j-bundle is emitted thermally at one (e.g., anode) footprint. A more realistic model can predict a more complicated spectrum with both footprints emitting. Finally, the hot-spot emission may be significantly anisotropic (Perna & Gotthelf 2008), which may increase or reduce its apparent luminosity depending of the average inclination of the spot to the line of sight.

Note that the observed spectrum was well fitted by a two-temperature blackbody. The second blackbody component was cooler and had emission area comparable to that of the star. No nonthermal component from magnetospheric scattering was required by the data. The weakness of the scattered component may be explained by the small size u_⋆ of the j-bundle: the star's radiation does not scatter in the cavity around the j-bundle, and the probability of scattering in the narrow bundle is low (Section 6.2).

In ordinary radio pulsars, radio emission is believed to be produced by the open field lines passing through the light cylinder, because they are the only part of the magnetosphere that carries electric currents. In magnetars, these field lines have a tiny u_lc = R/R_lc ≲ 10⁻⁴ where R_lc = c/Ω is the light-cylinder radius for a star rotating with angular velocity Ω. Radio emission from the open bundles of magnetars may be undetectable for two reasons. (1) e^± discharge has a low threshold in magnetars, ${\cal V}\sim 10^9$ V. When this voltage is multiplied with the current in the bundle, $I_{\rm lc}\approx \hat{I}u_{\rm lc}^2=c\mu /4R_{\rm lc}^2$, one obtains the dissipated power $L_{\rm lc}={\cal V}I_{\rm lc}\sim 10^{28}{\cal V}_9$ erg s⁻¹. This power appears to be too small to feed the observed radio luminosity of XTE J1810−197, L_radio ∼ 10³⁰ erg s⁻¹. (2) The radio beam may be narrow (because of small u_lc). Then the probability of its passing through our line of sight is small. This suggests that radio pulsations are hardly detectable when the magnetar is in quiescence, i.e., when its magnetosphere is untwisted.

In contrast, after the starquake, the j-bundle forms. It is much thicker and more energetic than the bundle passing through the light cylinder and can produce much brighter radio emission with a much broader pulse. The untwisting magnetosphere in XTE J1810−197 had u_⋆/u_lc ≳ 3 × 10². The net current flowing in the j-bundle is ∼10⁵ times larger than I_lc,

$$\begin{equation} \frac{I_{\star }}{I_{\rm lc}}\approx \left(\frac{u_{\star }}{u_{\rm lc}}\right)^2\sim 10^5. \end{equation} \tag{ 57 }$$

This gives L/L_lc ∼ 10⁵ (assuming a comparable discharge voltage at u_⋆ and u_lc). A small fraction _radio of this luminosity may escape as radio waves; _radio ∼ 10⁻³ is consistent with observations.

Radio waves are efficiently absorbed by the magnetospheric plasma. Only waves produced close to the magnetic axis, u < u_esc, are likely to escape, while waves emitted at u > u_esc are trapped in the magnetosphere. The value of u_esc depends on the frequency of the wave and the state of the plasma (particle density and velocity distribution); u_esc might be inferred from the observed opening angle of the radio beam. The radio luminosity should quickly decrease when the j-bundle shrinks to u_⋆ < u_esc.

Radio pulsar XTE J1810−197 is distinguished from ordinary radio pulsars by its hard spectrum, very strong linear polarization, and variable pulse profile. If its emission is produced on the bundle of closed field lines with u ≫ u_lc, it may be expected to be different from ordinary radio pulsars. The spectrum of radio waves may form as the sum of emissions from different (frequency-dependent) u_esc inside the j-bundle. The sporadic changes in the pulse profile may be caused by the instabilities in the maximally twisted outer magnetosphere (Section 6.3).

Radio observations provided accurate measurements of the spindown torque acting on the star (Camilo et al. 2007). The torque was decreasing together with the X-ray luminosity 2–3 years after the starquake. Its history at earlier times was not observed; the torque is believed to have increased after the starquake (Camilo et al. 2007). Such a non-monotonic evolution of the torque would be consistent with the theoretical expectations (Section 6.3).