© ESO, 2014

1. Introduction

The pulsar magnetospheres consist mainly of electron-positron plasma. This plasma can affect the radiation produced in the inner region of the magnetosphere or at the stellar surface and, owing to this, the pulsar emission can provide information regarding the physical conditions in the magnetosphere. This might be a powerful method for the diagnostics of the magnetosphere. For instance, fluctuations of the pulsar emission can be caused by non-stationary phenomena in the magnetospheric plasma (such as instabilities, waves, etc.), which are determined by the conditions in the magnetosphere. Therefore, the spectrum and characteristic timescale of the detected fluctuations can provide important information regarding the state of plasma in the magnetosphere. That is why understanding non-stationary properties of the magnetospheric plasma is of crucial importance for the interpretation of observational data.

The physical processes in the pulsar magnetosphere are very particular. The mean free path of particles is typically short compared to the characteristic lengthscale and, hence, the magnetohydrodynamic (MHD) description is justified. For typical values of the magnetic field, the electromagnetic energy density is greater than the kinetic energy density. This suggests that the force-free equation is a good approximation for determining the magnetic field structure over much of the magnetosphere. The growing observational data on spectra and pulse profiles of isolated pulsars prompt continued improvement of theoretical models of such force-free magnetospheres (see, e.g. Goodwin et al. 2004; Contopoulos et al. 1999; Komissarov 2006; McKinney 2006; Petrova 2013). In the axisymmetric case, the equations governing the structure of a pulsar magnetosphere can be reduced to the well-known Grad-Shafranov equation (see, e.g. Michel 1973; Mestel 1973; Mestel & Shibata 1994; see also Beskin 1997, for general overview). Most models based on this equation have a “dead zone” with field lines that are close within the light-cylinder and a “wind zone” with poloidal field lines that cross the light-cylinder. Poloidal currents in the “wind zone” maintain a toroidal field component, whereas currents are vanishing in the “dead zone”.

Many MHD phenomena in the force-free magnetosphere, however, are still poorly understood. Particularly, this concerns the non-stationary processes, such as the various types of instability, that can occur in the magnetosphere. This question is of particular interest because the existence of a stationary force-free configuration even rises doubts (see, e.g. Timokhin 2006). One of the MHD instabilities that can arise in the pulsar magnetosphere is the so-called diocotron instability, which is the non-neutral plasma analog of the Kelvin-Helmholtz instability. This instability has been studied extensively in the context of laboratory plasma devices (see, e.g. Levy 1965; Davidson 1990; Davidson & Felice 1998). The possible existence of pulsars having a differentially rotating equatorial disk with a non-vanishing charge density could trigger instability of diocotron modes (Petri et al. 2002). In the non-linear regime, the diocotron instability might cause diffusion outwards of the charged particles across the magnetic field lines outwards (Petri et al. 2003). The role of a diocotron instability in causing drifting subpulses in radio pulsar emission has been discussed by Fung et al. (2006).

Recently, a new mode of the magnetospheric oscillations has been considered by Urpin (2011). This mode is closely related to the Alfvénic waves from standard magnetohydrodynamics, which have been modified by the force-free condition and non-vanishing electric charge density. This type of magnetospheric waves can be unstable because there is a number of destabilising factors in the magnetosphere (such as differential rotation, electric currents, non-zero charge density, etc.). For example, many models of the magnetosphere predict that rotation should be differential (see, e.g. Mestel & Shibata 1994; Contopoulos et al. 1999) but it is known that differential rotation in plasma with the magnetic field leads to the so-called magnetorotational instability (Velikhov 1959). In the axisymmetric model of a magnetosphere suggested by Countopoulos et al. (1999), the angular velocity decreases inversely proportional to the cylindrical radius beyond the light cylinder and even stronger in front of it. For such rotation, the growth time of unstable magnetospheric waves is of the order of the rotation period (Urpin 2012). Numerical modelling by Komissarov (2006) showed that plasma rotates differentially basically near the equator and poles within the light cylinder. Such strong differential rotation should lead to instability that also arises on a timescale of the order of a rotation period. Note that the magnetorotational instability in the pulsar magnetosphere differs essentially from the standard magnetorotational instability because of a non-vanishing charge density and the force-free condition (Urpin 2012).

The electric currents flowing in plasma also provide a destabilising influence that leads to the so-called Tayler instability (see, e.g. Tayler 1973a,b). This instability is well studied in both laboratory and stellar conditions. It arises basically on the Alfvén time scale and is particularly efficient if the strengths of the toroidal and poloidal field components differ significantly (see, e.g. Bonanno & Urpin 2008a,b). This condition is satisfied in many magnetospheric models (see, e.g. Contopoulos et al. 1999), and these models can be unstable.

However, this instability might have a number of qualitative features in the pulsar magnetosphere because of the force-free condition and non-zero charge density. In the present paper, we consider the instability of the pulsar magnetosphere relevant to magnetospheric waves considered by Urpin (2011). We show that these waves can be unstable in the regions where the magnetic pressure gradient is non-vanishing. The considered instability arises on a short timescale and can be responsible for a short-term variability of the pulsar emission.

2. Basic equations

Despite uncertainties in estimates of many parameters, plasma in the pulsar magnetosphere is likely collisional and the Coulomb mean free path of electrons and positrons is small compared to the characteristic length scale. Therefore, the MHD description can be applied to such highly magnetized plasma (see Urpin 2012, for more details).

Let us define the hydrodynamic velocity and electric current of an electron-positron plasma as ${V}\mathrm{=}\frac{\mathrm{1}}{\mathit{n}}\mathrm{(}{\mathit{n}}_{\mathrm{e}}{{V}}_{\mathrm{e}}\mathrm{+}{\mathit{n}}_{\mathrm{p}}{{V}}_{\mathrm{p}}\mathrm{)}\mathit{,}{j}\mathrm{=}\mathit{e}\mathrm{(}{\mathit{n}}_{\mathrm{p}}{{V}}_{\mathrm{p}}\mathrm{-}{\mathit{n}}_{\mathrm{e}}{{V}}_{\mathrm{e}}\mathrm{)}\mathit{,}$(1)where (V_e,n_e) and (V_p.n_p) are the partial velocities and number densities of electrons and positrons, respectively; n = n_e + n_p. Then, partial velocities of the electrons and positrons can be expressed in terms of V and j: ${{V}}_{\mathrm{e}}\mathrm{=}\frac{\mathrm{1}}{\mathrm{2}{\mathit{n}}_{\mathrm{e}}}\left(\mathit{n}{V}\mathrm{-}\frac{{j}}{\mathit{e}}\right)\mathit{,}{{V}}_{\mathrm{p}}\mathrm{=}\frac{\mathrm{1}}{\mathrm{2}{\mathit{n}}_{\mathrm{p}}}\left(\mathit{n}{V}\mathrm{+}\frac{{j}}{\mathit{e}}\right)\mathrm{·}$(2)If the number density of plasma, n, is much greater than the charge number density, | n_p − n_e |, then V ≫ j/en. In the general case, the hydrodynamic and current velocities can be comparable in the electron-positron plasma.

MHD equations governing the electron-positron plasma can be obtained from the partial momentum equations for the electrons and positrons in the standard way (see Urpin 2012). Assuming that plasma is non-relativistic, the momentum equation for particles of the sort α (α = e,p) reads $\begin{array}{ccc}{\mathit{m}}_{\mathit{\alpha }}{\mathit{n}}_{\mathit{\alpha }}& & \left[\dot{{V}}\mathit{\alpha }\mathrm{+}\mathrm{(}{{V}}_{\mathit{\alpha }}\mathrm{·}\mathrm{\nabla }\mathrm{)}{{V}}_{\mathit{\alpha }}\right]\mathrm{=}\mathrm{-}\mathrm{\nabla }{\mathit{p}}_{\mathit{\alpha }}\mathrm{+}{\mathit{n}}_{\mathit{\alpha }}{{F}}_{{\alpha }}\\ & & \mathrm{+}\hspace{0.17em}{\mathit{e}}_{\mathit{\alpha }}{\mathit{n}}_{\mathit{\alpha }}\left({E}\mathrm{+}\frac{{{V}}_{\mathit{\alpha }}}{\mathit{c}}\mathrm{\times }{B}\right)\mathrm{+}{{R}}_{\mathit{\alpha }}\end{array}$(3)(see, e.g. Braginskii 1965 where the general plasma formalism is considered); the dot denotes the partial time derivative. Here, V_α is the mean velocity of particles α; n_α and p_α are their number density and pressure, respectively; F_α is an external force acting on the particles α (in our case F_α is the gravitational force); E is the electric field; and R_α is the internal friction force caused by collisions of the particles α with other sorts of particles. Since R_α is the internal force, the sum of R_α over α is zero in accordance with Newton’s Third Law. Therefore, we have in the electron-positron plasma R_e = −R_p.

The inertial terms on the l.h.s. of Eq. (3) give a small contribution to the force balance because of a small mass of both electrons and positrons. Gravitational force can also be neglected because of the same reason. A gas pressure is much smaller than the magnetic pressure in the force-free pulsar magnetosphere. Therefore, the momentum Eq. (3) reads in the electron-positron plasma ${\mathit{e}}_{\mathit{\alpha }}{\mathit{n}}_{\mathit{\alpha }}\left({E}\mathrm{+}\frac{{{V}}_{\mathit{\alpha }}}{\mathit{c}}\mathrm{\times }{B}\right)\mathrm{+}{{R}}_{\mathit{\alpha }}\mathrm{=}\mathrm{0.}$(4)Generally, the friction force, R_α contains two terms: one proportional to the difference of partial velocities (V_e − V_p) and another proportional to the temperature gradient (see, e.g. Braginskii 1965). We will neglect the thermal contribution to R_α and take into account only friction caused by a difference in the partial velocities. This is equivalent to neglecting the thermal diffusion of particles compared to their hydrodynamic velocities. The friction force is related to the velocity difference, (V_e − V_p), by a tensor that generally has components along and across the magnetic field and the so-called Hall component, which is perpendicular to the both magnetic field and velocity difference. In a strong magnetic field, parallel and perpendicular components are comparable but the Hall component is small (see Braginskii 1965). Therefore, we mimic the friction force between electrons and positrons by ${{R}}_{\mathrm{e}}\mathrm{=}\mathrm{-}\frac{{\mathit{m}}_{\mathrm{e}}{\mathit{n}}_{\mathrm{e}}}{{\mathit{\tau }}_{\mathrm{e}}}\mathrm{(}{{V}}_{\mathrm{e}}\mathrm{-}{{V}}_{\mathrm{p}}\mathrm{)}\mathit{,}$(5)where τ_e is the relaxation time of electrons. Note that this simple model for the friction force is often used even for a magnetized plasma in laboratory conditions (Braginskii 1965) and yields qualitatively correct results. We assume that accuracy of Eq. (5) is sufficient to study the magnetosphere of pulsars.

It is usually more convenient to use linear combinations of Eq. (4) than to solve partial equations. The sum of electron and positron Eq. (4) yields the equation of hydrostatic equilibrium in the magnetosphere $\begin{array}{ccc}{\mathit{\rho }}_{\mathrm{e}}{E}\mathrm{+}\frac{\mathrm{1}}{\mathit{c}}{j}\mathrm{\times }{B}\mathrm{=}\mathrm{0}\mathit{,}& & \end{array}$(6)where ρ_e = e(n_p − n_e) = eδn is the charge density. Taking the difference between electron and positron Eq. (4), we obtain the Ohm’s law in the form ${j}\mathrm{=}{\mathit{\rho }}_{\mathrm{e}}{V}\mathrm{+}\mathit{\sigma }\left({E}\mathrm{+}\frac{{V}}{\mathit{c}}\mathrm{\times }{B}\right)$(7)where σ = e²n_pτ_e/m_e is the conductivity of plasma.

It was shown by Urpin (2012) that Eqs. (6), (7) are equivalent to two equations ${j}\mathrm{=}{\mathit{\rho }}_{\mathrm{e}}{V}\mathit{,}{E}\mathrm{=}\mathrm{-}\frac{{V}}{\mathit{c}}\mathrm{\times }{B}\mathit{.}$(8)These equations imply that the force-free condition and the Ohm’s law (Eqs. (6), (7)) are equivalent to the conditions of a frozen-in magnetic field and the presence of only advective currents in the magnetosphere. Departures from this equivalence can be caused, for instance, by general relativistic corrections (see Palenzuela 2013) but they are very small in the magnetosphere. Note that the cross-production of the frozen-in condition and B yields the well-known expression for the transverse to B component of the velocity: V_⊥ = c(E × B/B²). Certainly, the electric current should be non-vanishing in the magnetosphere, j ≠ 0, because it maintains the magnetic configuration. Hence, the hydrodynamic velocity should be non-zero as well since the current is advective. Therefore, the force-free magnetosphere can only exist if hydrodynamic motions are non-vanishing (Urpin 2012).

3. Equation governing the magnetospheric waves

Equation (8) should be complemented by the Maxwell equations. Then, the set of equations, governing MHD processes in the force-free pulsar magnetosphere reads $\begin{array}{ccc}\mathrm{\nabla }\mathrm{·}{E}& \mathrm{=}& \mathrm{4}\mathit{\pi }{\mathit{\rho }}_{\mathrm{e}}\mathit{,}\mathrm{\nabla }\mathrm{\times }{E}\mathrm{=}\mathrm{-}\frac{\mathrm{1}}{\mathit{c}}\frac{\mathit{\partial }{B}}{\mathit{\partial t}}\mathit{,}\\ \mathrm{\nabla }\mathrm{·}{B}& \mathrm{=}& \mathrm{0}\mathit{,}\mathrm{\nabla }\mathrm{\times }{B}\mathrm{=}\frac{\mathrm{1}}{\mathit{c}}\frac{\mathit{\partial }{E}}{\mathit{\partial t}}\mathrm{+}\frac{\mathrm{4}\mathit{\pi }}{\mathit{c}}{j}\mathit{,}\\ {j}& \mathrm{\approx }& {\mathit{\rho }}_{\mathrm{e}}{V}\mathit{,}{E}\mathrm{\approx }\mathrm{-}\frac{{V}}{\mathit{c}}\mathrm{\times }{B}\mathit{.}\end{array}$(9)Consider the properties of MHD waves with small amplitudes as descibed by these equations. We assune that the electric and magnetic fields are equal to E₀ and B₀ in the unperturbed magnetosphere. The corresponding electric current, charge density, and velocity are j₀, ρ_e0. and V₀, respectively. For the sake of simplicity, we assume that motions in the magnetosphere are non-relativistic (V₀ ≪ c). Linearising Eq. (9), we obtain the set of equations that describes the behaviour of modes with a small amplitude. Small perturbations are indicated by subscript 1. We consider waves with a short wavelength. The space-time dependence of such waves can be taken in the form ∝exp(iωt − ik·r), where ω and k are the frequency and wave vector, respectively. Such waves exist if their wavelength λ = 2π/k is short compared to the characteristic length scale of the magnetosphere, L. Typically, L is greater than the stellar radius. We search in MHD modes with the frequency to satisfy the condition ω < 1/τ_e, since we use the MHD approach.

Substituting the frozen-in condition, E = −V × B/c, into the equation c∇ × E = −∂B/∂t and linearising the obtained induction equation, we have $\mathrm{i}\mathit{\omega }{{B}}_{\mathrm{1}}\mathrm{=}\mathrm{\nabla }\mathrm{\times }\mathrm{(}{{V}}_{\mathrm{1}}\mathrm{\times }{{B}}_{\mathrm{0}}\mathrm{+}{{V}}_{\mathrm{0}}\mathrm{\times }{{B}}_{\mathrm{1}}\mathrm{)}\mathit{.}$(10)A destabilising effect of shear already has been studied by Urpin (2012). In the present paper, we concentrate on the instability caused by the presence of electric currents in the magnetosphere. Therefore, we assume that shear is small and neglect terms proportional to | ∂V_0i/∂x_j |. Then, Eq. (10) reads $\mathrm{i}\begin{array}{c}\mathrm{˜}\\ \mathit{\omega }\end{array}{{B}}_{\mathrm{1}}\mathrm{=}\mathrm{i}{{B}}_{\mathrm{0}}\mathrm{(}{k}\mathrm{·}{{V}}_{\mathrm{1}}\mathrm{)}\mathrm{-}\mathrm{i}{{V}}_{\mathrm{1}}\mathrm{(}{k}\mathrm{·}{{B}}_{\mathrm{0}}\mathrm{)}\mathrm{-}\mathrm{(}{{V}}_{\mathrm{1}}\mathrm{·}\mathrm{\nabla }\mathrm{)}{{B}}_{\mathrm{0}}\mathit{,}$(11)where $\begin{array}{c}{}_{\mathrm{˜}}\\ \mathit{\omega }\end{array}\mathrm{=}\mathit{\omega }\mathrm{-}{k}\mathrm{·}{{V}}_{\mathrm{0}}$. The last term on the r.h.s. is usually small compared to the third term (~λ/L) in a short wavelength approximation. However, it becomes crucially important if the wavevector of perturbations is almost perpendicular to B₀.

Substituting the expression j = ρ_eV into Ampere’s law (second line of Eq. (9)) and linearising the obtained equation, we have ${{V}}_{\mathrm{1}}\mathrm{=}\mathrm{-}\frac{\mathrm{i}}{\mathrm{4}\mathit{\pi }{\mathit{\rho }}_{\mathrm{e}\mathrm{0}}}\mathrm{(}\mathit{c}{k}\mathrm{\times }{{B}}_{\mathrm{1}}\mathrm{+}\mathit{\omega }{{E}}_{\mathrm{1}}\mathrm{)}\mathrm{-}\frac{{\mathit{\rho }}_{\mathrm{e}\mathrm{1}}}{{\mathit{\rho }}_{\mathrm{e}\mathrm{0}}}{{V}}_{\mathrm{0}}\mathit{.}$(12)We search in relatively low-frequency MHD modes with the frequency ω < ck. Note that the frequency of MHD modes must satisfy the condition ω < 1/τ_e because of the MHD approach used. The relaxation time can be estimated as τ_e ~ ℓ_e/c_e, where c_e and ℓ_e are the thermal velocity and mean free path of particles, respectively. The frequency ck can be greater or smaller than 1/τ_e depending on a wavelength λ. If λ > 2πℓ_e(c/c_e), then we have ck < 1/τ_e. If λ < 2πℓ_e(c/c_e), then we have ck > 1/τ_e.

By eliminating E₁ from Eq. (12) by making use of the linearised frozen-in condition and neglecting terms of the order of (ω/ck)(V₀/c), we obtain the following for such modes: ${{V}}_{\mathrm{1}}\mathrm{+}\frac{\mathrm{i}\mathit{\omega }}{\mathrm{4}\mathit{\pi c}{\mathit{\rho }}_{\mathrm{e}\mathrm{0}}}{{B}}_{\mathrm{0}}\mathrm{\times }{{V}}_{\mathrm{1}}\mathrm{=}\mathrm{-}\frac{\mathrm{i}\mathit{c}}{\mathrm{4}\mathit{\pi }{\mathit{\rho }}_{\mathrm{e}\mathrm{0}}}{k}\mathrm{\times }{{B}}_{\mathrm{1}}\mathrm{-}\frac{{\mathit{\rho }}_{\mathrm{e}\mathrm{1}}}{{\mathit{\rho }}_{\mathrm{e}\mathrm{0}}}{{V}}_{\mathrm{0}}\mathit{.}$(13)The perturbation of the charge density can be calculated from the equation ρ_e1 = ∇·E₁/4π. We have then with accuracy in terms of the lowest order in (λ/L)${\mathit{\rho }}_{\mathrm{e}\mathrm{1}}\mathrm{=}\frac{\mathrm{1}}{\mathrm{4}\mathit{\pi c}}\left[\mathrm{i}{{B}}_{\mathrm{0}}\mathrm{·}\mathrm{(}{k}\mathrm{\times }{{V}}_{\mathrm{1}}\mathrm{)}\mathrm{-}\mathrm{i}{{V}}_{\mathrm{0}}\mathrm{·}\mathrm{(}{k}\mathrm{\times }{{B}}_{\mathrm{1}}\mathrm{)}\right]\mathit{.}$(14)Substituting Eq. (14) into Eq. (13) and neglecting terms of the order of ${\mathit{V}}_{\mathrm{0}}^{\mathrm{2}}\mathit{/}{\mathit{c}}^{\mathrm{2}}$, we obtain the second equation, which couples B₁ and V₁, $\begin{array}{ccc}\mathrm{4}\mathit{\pi c}{\mathit{\rho }}_{\mathrm{e}\mathrm{0}}{{V}}_{\mathrm{1}}\mathrm{+}\mathrm{i}\mathit{\omega }{{B}}_{\mathrm{0}}\mathrm{\times }{{V}}_{\mathrm{1}}\mathrm{=}\mathrm{-}\hspace{0.17em}\mathrm{i}{\mathit{c}}^{\mathrm{2}}{k}\mathrm{\times }{{B}}_{\mathrm{1}}\mathrm{-}\mathrm{i}{{V}}_{\mathrm{0}}\left[{{B}}_{\mathrm{0}}\mathrm{·}\mathrm{(}{k}\mathrm{\times }{{V}}_{\mathrm{1}}\mathrm{)}\right]\mathit{.}& & \end{array}$(15)Eliminating B₁ from Eqs. (11) and (15) in favor of V₁ and again neglecting terms of the order of (ω/ck)(V₀/c) and (ω/ck)², we obtain the equation for V₁ in the form $\mathrm{4}\mathit{\pi c}{\mathit{\rho }}_{\mathrm{e}\mathrm{0}}{{V}}_{\mathrm{1}}\mathrm{-}\mathrm{i}\frac{{\mathit{c}}^{\mathrm{2}}}{\begin{array}{c}\mathrm{˜}\\ \mathit{\omega }\end{array}}\mathrm{(}{k}\mathrm{·}{{B}}_{\mathrm{0}}\mathrm{)}{k}\mathrm{\times }{{V}}_{\mathrm{1}}\mathrm{=}\frac{{\mathit{c}}^{\mathrm{2}}}{\begin{array}{c}\mathrm{˜}\\ \mathit{\omega }\end{array}}{k}\mathrm{\times }\mathrm{\left[}\mathrm{\right(}{{V}}_{\mathrm{1}}\mathrm{·}\mathrm{\nabla }\mathrm{)}{{B}}_{\mathrm{0}}\mathrm{-}\mathrm{i}{{B}}_{\mathrm{0}}\mathrm{(}{k}\mathrm{·}{{V}}_{\mathrm{1}}\mathrm{\left)}\mathrm{\right]}\mathit{.}$(16)It immediately follows from this equation that (k·V₁) = 0 and, hence, the magnetospheric waves are transverse. Equation (16) is simplified to $\mathit{\alpha }{{V}}_{\mathrm{1}}\mathrm{-}\mathrm{i}\mathrm{(}{k}\mathrm{·}{{B}}_{\mathrm{0}}\mathrm{)}{k}\mathrm{\times }{{V}}_{\mathrm{1}}\mathrm{=}{k}\mathrm{\times }\mathrm{(}{{V}}_{\mathrm{1}}\mathrm{·}\mathrm{\nabla }\mathrm{)}{{B}}_{\mathrm{0}}\mathit{,}$(17)where $\mathit{\alpha }\mathrm{=}\mathrm{4}\mathit{\pi }{\mathit{\rho }}_{\mathrm{e}\mathrm{0}}\frac{\begin{array}{c}\mathrm{˜}\\ \mathit{\omega }\end{array}}{\mathit{c}}$(18)In the case of a uniform magnetic field, Eq. (17) reduces to the equation considered by Urpin (2011).

4. Dispersion equation and instability of magnetospheric modes

Generally, the stability properties of perturbations are complicated even in the local approximation. Equation (17) can be transformed to a more convenient form that does not contain a cross production of k and V₁. Calculating the cross production of k and Eq. (17) and taking into account k·V₁ = 0, we have ${k}\mathrm{\times }{{V}}_{\mathrm{1}}\mathrm{=}\mathrm{-}\frac{\mathrm{1}}{\mathit{\alpha }}\mathrm{\left\{}\mathrm{i}{\mathit{k}}^{\mathrm{2}}\mathrm{\right(}{k}\mathrm{·}{{B}}_{\mathrm{0}}\mathrm{)}{{V}}_{\mathrm{1}}\mathrm{-}\mathrm{(}{{V}}_{\mathrm{1}}\mathrm{·}\mathrm{\nabla }\mathrm{\left)}\mathrm{\right[}{k}\mathrm{\times }\mathrm{(}{k}\mathrm{\times }{{B}}_{\mathrm{0}}\mathrm{)}{\mathrm{\left]}}^{\mathrm{\right\}}}\mathit{.}$(19)Substituting this expression into Eq. (17), we obtain $\begin{array}{ccc}\mathrm{[}{\mathit{\alpha }}^{\mathrm{2}}\mathrm{-}{\mathit{k}}^{\mathrm{2}}\mathrm{(}{k}\mathrm{·}{{B}}_{\mathrm{0}}{\mathrm{)}}^{\mathrm{2}}\mathrm{]}{{V}}_{\mathrm{1}}\mathrm{=}\mathit{\alpha }\mathrm{(}{{V}}_{\mathrm{1}}\mathrm{·}\mathrm{\nabla }\mathrm{)}{k}\mathrm{\times }{{B}}_{\mathrm{0}}& & \\ \mathrm{+}\hspace{0.17em}\mathrm{i}\mathrm{(}{k}\mathrm{·}{{B}}_{\mathrm{0}}\mathrm{)}\mathrm{(}{{V}}_{\mathrm{1}}\mathrm{·}\mathrm{\nabla }\mathrm{)}\mathrm{\left[}{k}\mathrm{\right(}{k}\mathrm{·}{{B}}_{\mathrm{0}}\mathrm{)}\mathrm{-}{\mathit{k}}^{\mathrm{2}}{{B}}_{\mathrm{0}}\mathrm{]}\mathit{.}& & \end{array}$(20)The magnetospheric waves can exist in the force-free pulsar magnetosphere only if the wavevector k and the unperturbed magnetic field B₀ are almost (but not exactly) perpendicular and the scalar production (k·B₀) is small but non-vanishing (see Urpin 2011, 2012). The reason for this is clear from simple qualitative arguments. The magnetospheric waves are transverse (k·V₁ = 0), and the velocity of plasma is perpendicular to the wave vector. However, wave motions across the magnetic field in a strong field are suppressed and the velocity component along the magnetic field should be much greater than in the transverse (see, e.g. Mestel & Shibata 1994). Therefore, the direction of a wavevector k should be close to the plane perpendicular to B₀. That is why we treat Eq. (20) only in the case of small (k·B₀).

Consider Eq. (20) in the neighbourhood of a point, r₀, using local Cartesian coordinates. We assume that the z-axis is parallel to the local direction of the unperturbed magnetic field and the corresponding unit vector is b = B₀(r₀)/B₀(r₀). The wavevector can be represented as k = k_∥b + k_⊥, where k_∥ and k_⊥ are parallel and perpendicular to the magnetic field components of k, respectively. Then, we have from the continuity equation ${\mathit{V}}_{\mathrm{1}\mathit{z}}\mathrm{=}\mathrm{-}\frac{\mathrm{1}}{{\mathit{k}}_{\mathrm{\parallel }}}\mathrm{(}{{k}}_{\mathrm{\perp }}\mathrm{·}{{V}}_{\mathrm{1}\mathrm{\perp }}\mathrm{)}\mathit{.}$(21)Since k_⊥ ≫ k_∥, we have V_1z ≫ V_1 ⊥ and, hence, (V₁·∇) ≈ V_1z∂/∂z. Therefore, the z-component of Eq. (20) yields the following dispersion relation ${\mathit{\alpha }}^{\mathrm{2}}\mathrm{+}\mathit{A\alpha }\mathrm{+}\mathrm{i}\mathit{D}\mathrm{=}\mathrm{0}\mathit{,}$(22)where $\begin{array}{ccc}\mathit{A}& \mathrm{=}& \mathrm{(}{k}\mathrm{\times }{b}\mathrm{)}\mathrm{·}\frac{\mathit{\partial }{{B}}_{\mathrm{0}}}{\mathit{\partial z}}\mathit{,}\\ \mathit{D}& \mathrm{=}& {\mathit{k}}^{\mathrm{2}}\mathrm{(}{k}\mathrm{·}{{B}}_{\mathrm{0}}\mathrm{)}\left[{b}\mathrm{·}\frac{\mathit{\partial }{{B}}_{\mathrm{0}}}{\mathit{\partial z}}\mathrm{+}\mathrm{i}\mathrm{(}{k}\mathrm{·}{{B}}_{\mathrm{0}}\mathrm{)}\right]\mathit{.}\end{array}$(23)We neglect in D corrections of the order ~λ/L to k²(k·B₀)². The roots of Eq. (22) correspond to two modes with the frequencies given by ${\mathit{\alpha }}_{\mathrm{1}\mathit{,}\mathrm{2}}\mathrm{=}\mathrm{-}\frac{\mathit{A}}{\mathrm{2}}\mathrm{\pm }{\left(\frac{{\mathit{A}}^{\mathrm{2}}}{\mathrm{4}}\mathrm{-}\mathrm{i}\mathit{D}\right)}^{\mathrm{1}\mathit{/}\mathrm{2}}\mathit{.}$(24)If the magnetic field is approximtely uniform along the field lines, then ∂B₀/∂z ≈ 0, and, hence, A ≈ 0 and D ≈ ik²(k·B₀)². In this case, the magnetospheric modes are stable and ${\mathit{\alpha }}_{\mathrm{1}\mathit{,}\mathrm{2}}\mathrm{\approx }\mathrm{\pm }\sqrt{\mathrm{-}\mathrm{i}\mathit{D}}$. The corresponding frequency is $\begin{array}{c}\mathrm{˜}\\ \mathit{\omega }\end{array}\mathrm{\approx }\mathrm{\pm }\frac{\mathit{ck}}{\mathrm{4}\mathit{\pi }{\mathit{\rho }}_{\mathrm{e}\mathrm{0}}}\mathrm{(}{k}\mathrm{·}{{B}}_{\mathrm{0}}\mathrm{)}\mathit{.}$(25)These waves have been first studied by Urpin (2011). Deriving the dispersion Eq. (25), it was assumed that ω ≪ ck. Therefore, the magnetospheric waves can exist only if (k·B₀) is small, as was discussed above: (k·B₀) ≪ 4πρ_e0. Sometimes, it is convenient to measure the true charge density, ρ_e0, in units of the Goldreich-Julian charge density, ρ_GJ = ΩB₀/2πc, where Ω is the angular velocity of a pulsar. Then, we can suppose ρ_e0 = ξρ_GJ, where ξ is a dimensionless parameter and, hence, the condition ω ≪ ck transforms into $\mathrm{2}\mathit{\xi }\mathrm{\Omega }\mathrm{\gg }\mathit{c}\mathrm{|}{k}\mathrm{·}{b}\mathrm{|}\mathit{.}$(26)Obviously, this condition can be satisfied only for waves with the wavevector almost perpendicular to B₀. Note, however, that if k is exactly perpendicular to B₀ the magnetospheric waves do not exist.

If ∂B₀/∂z ≠ 0, the magnetospheric waves turn out to be unstable. The instability is especially pronounced if | k·B₀ |  < B₀/L. In this case, the second term in the brackets of Eq. (24) is smaller than the first one and, therefore, the roots are ${\mathit{\alpha }}_{\mathrm{1}}\mathrm{\approx }\mathrm{-}\mathit{A}\mathrm{+}\mathrm{i}\frac{\mathit{D}}{\mathit{A}}\mathit{,}{\mathit{\alpha }}_{\mathrm{2}}\mathrm{=}\mathrm{-}\mathrm{i}\frac{\mathit{D}}{\mathit{A}}\mathrm{·}$(27)The coefficient D is approximately equal to $\mathit{D}\mathrm{\approx }{\mathit{k}}^{\mathrm{2}}\mathrm{(}{k}\mathrm{·}{{B}}_{\mathrm{0}}\mathrm{)}\left({b}\mathrm{·}\frac{\mathit{\partial }{{B}}_{\mathrm{0}}}{\mathit{\partial z}}\right)\mathrm{·}$(28)The first and second roots of Eq. (27) correspond to oscillatory and non-oscillatory modes, respectively. The occurence of instability is determined by the sign of the ratio D/A. If this ratio is positive for some wavevector k, then the non-oscillatory mode is unstable but the oscillatory one is stable for such k. If D/A < 0, then the oscillatory mode is unstable but the non-oscillatory one is stable for corresponding k. Note, however, that the frequency of oscillatory modes often can be very high and ω ≫ ck. Our consideration does not apply in this case. Indeed, we have α₁ ~ A and, hence, $\begin{array}{c}{}_{\mathrm{˜}}\\ {\mathit{\omega }}_{\mathrm{1}}\end{array}\mathrm{~}\mathit{ck}\mathrm{(}{\mathit{B}}_{\mathrm{0}}\mathit{/}\mathrm{4}\mathit{\pi }{\mathit{\rho }}_{\mathrm{e}\mathrm{0}}\mathit{L}\mathrm{)}$. The condition ω ≪ ck implies that B₀/4πρ_e0L < 1. Expressing the charge density in units of the Goldreich-Julian density, ρ_e0 = ξρ_GJ, we transform this inequality into $\frac{\mathrm{1}}{\mathrm{2}\mathit{\xi }}\frac{\mathit{c}}{\mathrm{\Omega }\mathit{L}}\mathrm{\ll }\mathrm{1.}$(29)This condition can be satisfied only in regions where ξ ≫ 1 and the charge density is much greater than the Goldreich-Julian density. If inequality (29) is not fulfilled, then Eq. (27) for the oscillatory mode α₁ does not apply, and only the non-oscillatory modes exist. For example, the charge density is large in the region where the electron-positron plasma is created. Therefore, condition (29) can be satisfied there, and, hence, the oscillatory instability can occur in this region.

The non-oscillatory modes have a lower growth rate and usually can occur in the pulsar magnetosphere. For any magnetic configuration, it is easy to verify that one can choose the wavevector of perturbations, k, in such a way that the ratio D/A becomes positive, and, hence, the non-oscillatory mode is unstable. Indeed, we can represent k as the sum of components parallel and perpendicular to the magnetic field, k = k_∥ + k_⊥. Obviously, A ∝ k_⊥ and D ∝ k_∥ and, hence, A/D ∝ k_∥/k_⊥. Therefore, if A/D < 0 for a value of k = (k_∥,k_⊥), this ratio changes the sign for k = (−k_∥,k_⊥) and k = (k_∥, −k_⊥), and waves with such the wavevectors are unstable. It turns out that there always exists the range of wavevectors for which the non-oscillatory modes are unstable and, hence, the force-free magnetosphere is always the subject of instability.

The necessary condition of instability is D ≠ 0. As it was mentioned, the magnetospheric waves exist only if the wavevector k is close to the plane perpendicular to the unperturbed magnetic field, B₀, and the scalar production (k·B₀) is small (but non-vanishing). Therefore, the necessary condition D ≠ 0 is equivalent to b·(∂B₀/∂z) ≠ 0. Since b = B₀/B₀, we can rewrite this condition as ${{B}}_{\mathrm{0}}\mathrm{·}\frac{\mathit{\partial }{{B}}_{\mathrm{0}}}{\mathit{\partial z}}\ne \mathrm{0.}$(30)This condition is satisfied if the magnetic pressure gradient along the magnetic field is non-zero. The topology of the magnetic field can be fairly complicated in the magnetosphere, particularly in a region close to the neutron star. This may happen because the field geometry at the neutron star surface should be very complex (see, e.g. Bonanno et al. 2005, 2006). Therefore, Eq. (30) can be satisfied in different regions of the magnetosphere. However, this condition can be fulfilled even if the magnetic configuration is relatively simple. As a possible example, we consider a region near the magnetic pole of a neutron star. It was shown by Urpin (2012) that a special type of cylindrical waves can exist there with m ≠ 0, where m is the azimuthal wavenumber. The criterion of instability (30) is certainly satisfied in this region and, hence, the instability can occur. Indeed, one can mimic the magnetic field by a vacuum dipole near the axis. The radial and polar components of the dipole field in the spherical coordinates (r,θ,ϕ) are ${\mathit{B}}_{\mathit{r}}\mathrm{=}{\mathit{B}}_{\mathrm{p}}{\left(\frac{\mathit{a}}{\mathit{r}}\right)}^{\mathrm{3}}\cos \mathit{\theta ,}{\mathit{B}}_{\mathit{\theta }}\mathrm{=}\frac{\mathrm{1}}{\mathrm{2}}{\mathit{B}}_{\mathrm{p}}{\left(\frac{\mathit{a}}{\mathit{r}}\right)}^{\mathrm{3}}\sin \mathit{\theta ,}$(31)where B_p is the polar strength of the magnetic field at the neutron star surface and a is the stellar radius (see, e.g. Urpin et al. 1994). The radial field is much greater than the polar one near the axis and, therefore, it is easy to check that the criterion of instability (30) is fulfilled in the polar gap. Hence, filament-like structures can be formed there. In some models, note that the force-free field at the top of the polar gap can differ from that of a vacuum dipole (see, e.g. Petrova 2012) but this cannot change our conclusion qualitatively. We will consider the instability in the polar gap in more detail elsewhere.

From Eq. (30), it follows that the instability in pulsar magnetospheres is driven by a non-uniformity of the magnetic pressure and, hence, it can be called “the magnetic pressure-driven instability”. Note that this instability can occur only in plasma with a non-zero charge density, ρ_e0 ≠ 0, and does not arise in a neutral plasma with ρ_e0 = 0.

5. Discussion

We have considered stability of the electron-positron plasma in the magnetosphere of pulsars. The pair plasma in the magnetosphere is likely created in a two-stage process: primary particles are accelerated by an electric field parallel to the magnetic field near the poles up to extremely high energy, and these produce a secondary, denser pair plasma via a cascade process (see, e.g. Michel 1982). The number density of this secondary plasma greatly exceeds the Goldreich-Julian number density, n_GJ =  | ρ_GJ | /e, required to maintain a corotation electric field and, hence, the multiplicity factor ξ can be very large. Unfortunately, this factor is model dependent and rather uncertain with estimates in a wide range from 10² to 10⁶ (see, e.g. Gedalin et al. 1998). For example, the model of bound-pair creation above polar caps (Usov & Melrose 1996) results in a relatively low value of ξ. This model postulates that the photons emitted by primary particles as a result of curvature emission create bound pairs (positronium) rather than free pairs. While the pairs remain bound, the screening of the component E parallel to the magnetic field is inefficient because screening is attributed to free pairs that can be charge-separated as a result of acceleration by E_∥. In the absence of screening, the height of the polar gap and the maximum energy of primary particles increase. The main part of the energy that primary particles gain during their motion through the polar gap is transformed into the energy of curvature photons and then into the energy of secondary pairs. Assuming that the electric field in the gap has more or less standard value (~10¹⁰ V cm^-1), Usov & Melrose (1996) obtain the following estimate for the multiplicity factor above polar caps: $\mathit{\xi }\mathrm{~}\mathrm{4}\mathrm{\times }{\mathrm{10}}^{\mathrm{2}}{\left(\frac{\mathit{P}}{\mathrm{0.1}\mathrm{s}}\right)}^{\mathrm{-}\mathrm{3}\mathit{/}\mathrm{4}}\mathit{,}$(32)where P is the pulsar period. Note that this value can be higher if dissociation of bound pairs is taken into account. However, the instability can operate in the region of pair creation even if ξ is such small. It can be efficient in regions with ξ ~ 1 as well.

The geometry of motions in the unstable magnetospheric waves is simple. Since these waves are transverse (k·V₁ = 0) and the wavevector of such waves should be close to the plane perpendicular to B₀, plasma motions are almost parallel or anti-parallel to the magnetic field. The velocity across B₀ is small. In our model, we have only considered the stability of plane waves using a local approximation. In this model, the instability should lead to formation of filament-like structures with filaments alongside the magnetic field lines. Note that plasma can move in the opposite directions in different filaments. The characteristic timescale of formation of such structures is ~1/Imω. Since the necessary condition (30) is likely satisfied in a major fraction of a magnetosphere, one can expect that filament-like structures can appear (and disappear) in different magnetospheric regions. We used the hydrodynamic approach in our consideration, which certainly does not apply to a large distance from the pulsar where the number density of plasma is small. Therefore, the considered instability is most likely efficient in the inner part of a magnetosphere where filament-like structures can be especially pronounced. The example of a region where the instability can occur is the so-called dead zone. Most likely, the hydrodynamic approximation is valid in this region and hydrodynamic motions are non-relativistic, as it was assumed in our analysis. Note that a particular geometry of motions in the basic (unperturbed) state is not crucial for the instability and cannot suppress the formation of filament-like structures. These structures can be responsible for fluctuations of plasma and, hence, the magnetospheric emission can fluctuate with the same characteristic timescale.

It should be also noted that the considered instability is basically electromagnetic in origin as followed in our treatment. Hydrodynamic motions in the basic state play no important role in the instability. For instance, the unperturbed velocity does even not enter the expression for the growth rate. Therefore, one can expect that the same type of instability arises in the regions where velocities are relativistic. This case will be considered in detail elsewhere.

Hydrodynamic motions accompanying the instability can be the reason of turbulent diffusion in the magnetosphere. This diffusion should be highly anisotropic because both the criteria of instability and its growth rate are sensitive to the direction of the wave vector. However, the turbulent diffusion caused by motions may be efficient in the transport of angular momentum and mixing plasma with a much higher enhancement of the diffusion coefficient in the direction of the magnetic field since the velocity of motions across the field is much slower than along it.

The characteristic growth rate of unstable waves, Im ω, can be estimated from Eq. (27) as Im ω ~ (c/4πρ_e0)(D/A). Since k and B₀ should be close to orthogonality in magnetospheric waves, we have A ~ kB₀/L and $\mathit{D}\mathrm{~}{\mathit{k}}^{\mathrm{2}}\mathrm{(}{k}\mathrm{·}{b}\mathrm{)}{\mathit{B}}_{\mathrm{0}}^{\mathrm{2}}\mathit{/}\mathit{L}$, where we estimate b·(∂B₀/∂z) as B₀/L. Then, $\mathrm{Im}\mathit{\omega }\mathrm{~}\mathit{ck}\frac{\mathrm{(}{k}\mathrm{·}{{B}}_{\mathrm{0}}\mathrm{)}}{\mathrm{4}\mathit{\pi }{\mathit{\rho }}_{\mathrm{e}\mathrm{0}}}\mathrm{~}\frac{\mathrm{1}}{\mathit{\xi }}\mathit{ck}\frac{\mathit{c}\mathrm{(}{k}\mathrm{·}{b}\mathrm{)}}{\mathrm{\Omega }}\mathrm{·}$(33)Like stable magnetospheric modes, the unstable ones can occur in the magnetosphere if Eq. (26) is satisfied. Generally, this condition requires a position close to orthogonality of (but not orthogonal) k and B₀. Under certain conditions, however, the instability can arise even if departures from orthogonality are not very small but ξ ≫ 1. As it was mentioned, this can happen in regions where the electron-positron plasma is created. The parameter ξ can also be greater than 1 in those regions where plasma moves with the velocity greater ΩL. Indeed, we have ρ_e0 = (1/4π)∇·E₀ for the unperturbed charge density. Since E₀ is determined by the frozen-in condition (8), we obtain ρ_e0 ~ (1/4πcL)V₀B₀. If the velocity of plasma in a magnetosphere is greater than the rotation velocity, then ξ ~ V₀/ΩL. Some models predict that the velocity in the magnetisphere can reach a fraction of c. Obviously, in such regions, condition (26) is satisfied even if departures from orthogonality of k and B₀ are relatively large.

The growth rate of instability (31) is sufficiently high and can reach a fraction of ck. For example, if a pulsar rotates with the period 0.01 s and ξ ~ 1, magnetospheric waves with the wavelength ~10⁵−10⁶ cm grow on a timescale ~10^-4−10^-5 s if a departure from orthogonality between k and B₀ is of the order of 10^-4. The considered instability can occur everywhere in the magnetosphere except regions close to the surfaces where B₀·(∂B₀/∂z) = 0 and instability criterion (30) is not satisfied.

The instability considered is caused by a combined action of non-uniform magnetic field and non-zero charge density. Certainly, this is not the only instability that can occur in the pulsar magnetosphere. There are many factors that can destabilize a highly magnetized magnetospheric plasma and lead to various instabilities with substantially different properties. For instance, differential rotation predicted by many models of the magnetosphere can be the reason of instability as was shown by Urpin (2012). This instability is closely related to the magnetorotational instability (Velikhov 1959) which is well-studied in the context of accretion disks (see, e.g. Balbus & Hawley 1991). Generally, the regions, where rotation is differential and the magnetic field is non-unifom, can overlap. Thus, the criteria of both instabilities can be fulfilled in the same region. However, these instabilities usually have substantially different growth rates. The instability caused by differential rotation arises usually on a time-scale comparable to the rotation period of a pulsar. The growth rate of the magnetic pressure-driven instability is given by Eq. (30) and can even reach a fraction of ck in accordance with our results. Therefore, this instability occurs

typically on a shorter time-scale than the instability caused by differential rotation. If two different instabilities can occur in the same region, then, the instability with a shorter growth time usually turns out to be more efficient and determines plasma fluctuations. It is likely, therefore, that the magnetic pressure-driven instability is more efficient everywhere in the magnetosphere except surfaces where criterion (30) is not satisfied. In the neighbourhood of these surfaces, howevere, the instability associated with differential rotation can occur despite it arises on a longer time-scale. Therefore, it appears that the whole pulsar magnetosphere should be unstable.

Instabilities can lead to a short-term variability of plasma and, hence, to modulate the magnetospheric emission of pulsars. The unstable plasma can also modulate the radiation produced at the stellar surface and propagating through the magnetosphere. Since the growth time of magnetospheric waves can be substantially different in different regions, the instability can lead to a generation of fluctuations over a wide range of timescales, including those yet to be detected in the present and future pulsar searches (Liu et al. 2011; Stappers et al. 2011). Detection of such fluctuations would uncover the physical conditions in the magnetosphere and enable one to construct a relevant model of the pulsar magnetosphere and its observational manifestations beyond the framework of the classical concept (see, e.g. Kaspi 2010).

Acknowledgments

The author thanks the Russian Academy of Sciencs for financial support under the programme OFN-15.

References