On the Escape of Low-frequency Waves from Magnetospheres of Neutron Stars

We study the nonlinear decay of the fast magnetosonic (fms) into the Alfvén waves in relativistic force-free magnetohydrodynamics. The work has been motivated by models of pulsar radio emission and fast radio bursts (FRBs), in which the emission is generated in neutron star magnetospheres at conditions when not only the Larmor but also the plasma frequencies significantly exceed the radiation frequency. The decay process places limits on the source luminosity in these models. We estimated the decay rate and showed that the phase volume of Alfvén waves available for the decay of an fms wave is infinite. Therefore, the energy of fms waves could be completely transferred to the small-scale Alfvén waves not via a cascade, as in the Kolmogorov turbulence, but directly. Our results explain the anomalously low radio efficiency of the Crab pulsar and show that FRBs could not be produced well within magnetar magnetospheres.


INTRODUCTION
Nonlinear effects play an important role in powerful compact sources of radio emission, such as pulsars and fast radio bursts (FRBs). For example, in nonmagnetized plasmas, the induced Compton and Raman scattering could even prevent the escape of the waves. In the strongly magnetized magnetosphere of neutron stars, the electromagnetic waves propagate in two orthogonally polarized modes: the so-called O-mode is polarized in the plane set by the background magnetic field and the propagation direction, whereas the X-mode is polarized perpendicularly to the magnetic field and the propagation direction. Only the O-mode is subject to induced scattering because the electric field of this mode has a component along the background magnetic field. In the field of the X-mode, the particles oscillate only due to the weak E × B drift; therefore, the scattering is suppressed.
The O-mode could propagate only if, in the rest frame of the plasma, the wave's frequency exceeds the plasma frequency. If the density of the plasma is high enough so that not only the Larmor but also the plasma frequency is well above the wave frequency, the two magnetohydrodynamic (MHD) waves could propagate: the fast magnetosonic (fms) and the Alfvén waves. The fms wave is polarized perpendicularly to the background magnetic field and the propagation direction. When this wave propagates towards decreasing plasma density, it is smoothly converted into the X-mode and could escape from the system. The Alfvén wave does not escape; it follows the curved magnetic field lines and eventually decays via the Landau damping (Arons & Barnard 1986).
Thus, in dense magnetospheres (such that both the Larmor and the plasma frequencies are above the emission frequency), the radiation propagates in the form of fms waves, independently of the emission mechanism. In this case, the nonlinear decay of fms into Alfvén waves could strongly affect the outgoing radiation. The goal of this paper is to study the decay process. In the magnetospheres of neutron stars, the magnetic energy significantly exceeds the plasma energy; therefore, the wave interaction could be considered in the scope of relativistic force-free MHD. The wave energy is well below the energy of the background field; therefore, we could employ the methods of the weak turbulence theory (see, e.g., Zakharov et al. 1992). Namely, we write down and solve the kinetic equation for the waves. Note that in the non-relativistic case, solutions to the kinetic equations for fms and Alfvén waves were investigated both analytically (Kuznetsov 2001) and numerically (Chandran 2005(Chandran , 2008.
The paper is organized as follows. In sect. 2, we write down the kinetic equations for MHD waves in the relativistic force-free regime. In sect. 3, we analyze the equations, estimate the nonlinear decay rate of fms waves, and qualitatively describe the kinetics of the decay of fms into Alfvén waves. In sect. 4, we solve the kinetic equations numerically, confirming our qualitative analysis. The implications of our findings for FRBs are outlined in sect. 5. Conclusions are presented in sect. 6.

FORCE-FREE REGIME
In this paper, we address weakly nonlinear interactions of MHD waves in the relativistic force-free regime when the plasma energy density, including the rest mass energy, is negligible compared to the magnetic energy density. This implies that the widely used magnetization parameter, σ = B 2 4πρc 2 , is infinite. Here B is the background magnetic field, and ρ is the plasma density.
In the force-free limit, there are two MHD waves: the fms wave, polarized perpendicularly to the background magnetic field and to the propagation direction, and the Alfvén wave, polarized perpendicularly to the background magnetic field in the plane set by the field and the propagation direction. The dispersion equations in the force-free limit are very simple: for the fms waves and ω = ck cos θ for the Alfvén waves. Here ω is the frequency, k the wave vector, and θ the angle between the wave vector and the background magnetic field. The interaction of weakly nonlinear MHD waves in the force-free limit was studied by Thompson & Blaes (1998) and Lyubarsky (2019). The strongest is the interaction of three waves satisfying the resonance conditions which in fact represent energy and momentum conservation. It follows from the dispersion relations that the conservation laws are satisfied for the decay of an fms wave into an fms and an Alfven waves (S ↔ S + A) and into two Alfven waves (S ↔ A + A). Of course, reverse merging processes are also possible. The conservation laws also permit the process S ↔ S + S for the aligned fms waves, but in the force-free limit, the probability of the process is zero (Lyubarsky 2019). At a finite σ, the weakly nonlinear interaction of aligned fms waves becomes possible; it leads to the steepening of the waves and formation of shocks (Levinson & van Putten 1997;Lyubarsky 2003). The three-wave interaction of Alfvén waves, A + A → A, is possible only if one of two waves has zero frequency. Then an anisotropic Alfvén cascade develops, transferring the energy to waves with large components of the wave vector perpendicular to the background field, k ⊥ ≫ k (Montgomery & Matthaeus 1995;Ng & Bhattacharjee 1996;Goldreich & Sridhar 1997). The important point is that the zero-frequency Alfvén waves cannot be treated as linear waves. Such interaction occurs when the field lines of the background magnetic field wander away. The Alfvén waves are stretched when propagating along diverging field lines, so the wave vector component perpendicular to the background field increases. Thus a cascade is formed, redistributing the Alfvén waves towards the high-k ⊥ domain, where they eventually decay. Goldreich & Sridhar (1997) presented a qualitative explanation of how small turbulent fluctuations, δB ≪ B, could lead to divergent field lines. Assume that the mean magnetic field is directed along z-axis and describe the turbulence as an ensemble of localized wave packets with the longitudinal and transverse scales l and l ⊥ , correspondingly, and the amplitude δB. The local magnetic field line turns within a wave packet by the angle θ ∼ δB B, so the field line deviates from the initial position by θl . Adding random deviations, one finds that the average displacement of the field line grows with the distance, s ∼ θ l z. Therefore, the field lines that were initially separated by l ⊥ diverge. This picture implicitly assumes a long wavelength tail in the fluctuation spectrum. Namely, the amplitude of the fluctuations at the scale z is ∆B ∼ Bs z ∼ δB l z. This implies ∆B 2 ∝ k so that the spectral power of the turbulence, ∆B 2 k = d∆B 2 dk , goes to a constant at k → 0. Such a spectrum is obtained if the turbulence is presented as an ensemble of bell-shaped fluctuations. We here deal with Alfvén waves produced by the decay of fms waves. The fms waves with the frequency ω produce Alfvén waves with the wave vector k ∼ ω c. The production rate is proportional to ω and the energy density of fms waves at this frequency, U (see eq. 21). Both these quantities decrease towards smaller frequencies, so the spectrum of the produced waves is cut off at long wavelengths. Therefore, in the case of interest, the three-wave Alfvén cascade is suppressed and could be neglected.
Of course, the Alfvén cascade could develop via fourwave interactions, A + A → A + A (Sridhar & Goldreich 1994;Goldreich & Sridhar 1995). However, we assume that the turbulence is weak, i.e., the system's evolution is governed by the lowest-order processes. At k ⊥ ∼ k , the rate of four-wave processes is lower than that of threewave processes by the ratio of the wave energy to the energy of the background field, 8πU B 2 0 . We will show that the fms waves with the frequency ω decay first into Alfvén waves with k ⊥ ∼ k ∼ ω c, so if the above ratio is small, the fms-to-Alfvén transformation initially occurs in the weak turbulence regime. However, higher k ⊥ Alfvén waves are produced in the course of time. It is well known that the role of nonlinearity grows at higher k ⊥ so that turbulence ceases to be weak at k ⊥ k ∼ c δv, where δv = cδB B 0 is the velocity of turbulent motions, δB the fluctuating magnetic field (Goldreich & Sridhar 1995). However, we will see that the spectrum of Alfvén waves reaches high k ⊥ only when the fms energy decreases ∼ (ck ⊥ ω) 2 times. Therefore, most of the transformation process occurs in the weak turbulence regime.
In astrophysical applications, we typically deal with wide spectra and random phases of waves. Then the wave field is conveniently described in the quantum language via the occupation numbers, n k , which are related to the wave energy density of a mode with the wave vector k as E k = ω k n k . Since we deal with two types of waves, we denote the occupation numbers of Alfvén waves by n k and of fms waves by N k . The evolution of the system is described by the kinetic equations for the waves ( e.g., Zakharov et al. 1992). For fms waves, these equations may be written as where R S↔S+A k,k1,k2 and R S↔A+A k,k1,k2 are the rates of the S ↔ S + A and S ↔ A + A processes for the given set of wave vectors k, k 1 , k 2 . Here, the first term describes the decay of the wave k into an fms and an Alfvén waves and the reverse process; the second is for the production of the fms k-wave via decay of an fms wave with the frequency ω 1 > ω and the reverse process, and the third term is for the decay of the fms k-waves into two Alfvén waves and the reverse process. The factor 1 2 in the third term takes into account double counting in the case of the decay into two waves of the same type. Similarly, the kinetic equation for the Alfvén waves is written as Now the factor 1 2 does not appear in the corresponding term because one of the two Alfvén quanta is fixed. Denoting the probabilities of the spontaneous decay processes as W S→A+S k,k1,k2 and W S→A+A k,k1,k2 , correspondingly, taking into account the induced processes and using the detailed balance principle, one can write In all cases of interest, N k ≫ 1, therefore we neglect the linear in N k terms, which describe spontaneous processes, and get The interaction probability is expressed via the amplitudes of the processes by the golden rule: 2 , and analogously for W S→A+A k,k1,k2 . The corresponding amplitudes were calculated by Lyubarsky (2019): Here the indexes and ⊥ describe the components of the k vector parallel and perpendicular to the background magnetic field. In equation (10), the index 1 is for the Alfvén wave and 2 for the fms wave. In equation (11), the Heaviside step function, H(x), expresses the well-known fact that Alfvén waves do not interact if they propagate in the same direction along the magnetic field.

CONSIDERATIONS
We study the possible decay of fms radiation produced by a strong enough source in a highly magnetized medium. Let us consider for simplicity the time evolution of a spatially homogeneous, isotropic fms radiation with the characteristic frequency of the order of ω 0 = ck 0 . Assume that the spectrum is moderately wide, ∆ω ∼ ω, and the total radiation energy is The fms waves decay into fms and Alfven waves of smaller frequency. The induced decay is possible only into states which are not empty; therefore, we assume that a weak background of fms and Alfven waves is present in the whole phase space. The initial pulse decays into waves satisfying the conservation laws (3); therefore, the occupation numbers in the new states grow until the reverse merging process balances the decay. This happens when the occupation numbers in all three states become comparable. The decay S → S +A occurs into states with k 1 and k 2 comparable with k 0 , and each initial fms quantum could produce only one quantum in the state k 1 and one quantum in the state k 2 . Therefore the equilibrium is achieved, and the decay stops after the energy of the initial peak decreases roughly two times.
The important point is that the phase volume available for decay S → A + A is in fact infinite because the Alfven waves could have the perpendicular component of their wave vector arbitrarily large. To demonstrate this, let us write the conservation laws (3) with the account of the dispersion laws (1) and (2) explicitly: One sees that for a given k, one finds k 1 and k 2 from the first two equations, whereas the third equation is satisfied with an arbitrary k 1⊥ by choosing k 2⊥ = k ⊥ − k 1⊥ . In particular, the fms wave could decay into two Alfven waves with arbitrarily large but nearly oppositely directed perpendicular components of the wave vectors. This implies that the phase volume available for the decay of an fms wave into a pair of Alfvén waves is infinite, so the fms pulse could decay significantly, practically to the background level. Substituting k 2⊥ = k ⊥ − k 1⊥ into equation (11) and expanding in small k k 1⊥ , one finds the interaction amplitude in the limit k 1⊥ ≫ k: (16) One sees that the rate of the fms decay into Alfvén waves with large perpendicular components of the wave vector does not go to zero.
The population of such states is described by equation (5), in which only the second term in the rhs should be retained. Taking into account that initially the occupation numbers for these states are small, one can write (note that in this equation, k 2 is referred to fms waves whereas k to Alfvén waves) Here we could take n k ≈ n k2−k because adding a quantum into one of these states is accompanied by adding a quantum into another state. Then one finally finds that population of the Alfvén waves with k ⊥ ≫ k 0 grows initially (as soon as n k ≪ N k0 ) exponentially, with the rate Here we used equation (16) with an appropriate permutation of indexes k, k 1 and k. Taking into account that k ∼ k 2 ∼ k 2⊥ ∼ k 2 = ω c (according to the conservation laws, the longitudinal components of the wave vectors of the Alfvén and fms waves are comparable, whereas all components of the fms wave vector are typically of the order of ω c), one can find a rough estimate: The coefficient in the definition of q 0 is chosen to show explicitly that the interaction rate is proportional to the small ratio of the wave energy and the energy of the background magnetic field. The total number of the produced Alfvén quanta, N = ∫ n k dk, grows as Here the integration over k is limited, due to the conservation laws, by the region k ∼ k 0 , whereas the integral over k ⊥ is unlimited. Taking into account that q is independent of k ⊥ at large k ⊥ , one sees that the total production rate of the Alfvén quanta diverges unless n k decreases with k ⊥ faster than k −2 ⊥ . Of course, Alfvén waves are erased by dissipation processes at small enough wavelengths, but the corresponding k is typically too high so that the Alfvén waves would be produced at an inappropriately high rate unless the background spectrum is steep enough. Here we adopt a natural assumption that the background waves have the spectrum n k ∝ k −α with α > 2.
In this case, the evolution of the system could be qualitatively described as follows. At the initial stage, the background Alfvén spectrum exponentially grows with the rate of eq. (21) at all k ⊥ . The waves with k ⊥ ∼ k 0 reach the saturation level, n k ∼ N k , first. At this stage, the fms pulse weakens a few times. The fms waves keep decaying into Alfvén waves with larger k ⊥ so that the energy of fms waves decreases monotonically. Then the Alfvén waves with smaller k ⊥ begin to merge into fms wave to maintain the equilibrium n k ∼ N k . Thus the Alfvén population of the states with larger k ⊥ gradually increases until the saturation level, n k ∼ N k . At smaller k ⊥ , a plateau is formed, which gradually decreases with time. In such a way, the fms pulse eventually could decay to the background level. The important point is that the decay rate is determined by the instant energy of the fms pulse. Since the energy of fms waves decreases, the decay process slows down with time. Therefore in real systems, the energy of fms waves decays until the decay rate becomes equal to the wave escape rate. Therefore the process places limits on the luminosity of the source. In the next section, we present numerical simulations that confirm the picture described.

SIMULATIONS
We study the time evolution of a spatially homogeneous, isotropic system of waves described by the kinetic equations (4) and (5). We assume that the initial spectral distribution of fms waves has a bell-shaped form, In addition, we assume that the whole phase space is filled by the background fms and Alvén waves with the spectrum In simulations, we used the power law index κ from 2.5 to 3.5 and found that the results differ insignificantly. Below, we present results for κ = 3. We solve the system of equations (4) and (5) using the Runge-Kutta method of the second order with correction (see Appendix). The system is closed; therefore, the total energy of the system, U 0 = ∫ ω(N k + n k )d 3 k, is conserved. Figure 1 shows the time evolution of the energies of fms and Alfvén waves normalized to the system's total energy, which is chosen to be equal to unity. One sees that the fms waves are efficiently converted to the Alfvén waves with the rate given by equation (21). The evolution of spectra is shown in figs. 2 and 3. We plotted the spectral distributions integrated over the longitudinal wave vectors because the evolution of the perpendicular component of the wave vector is the most interesting. The presented results confirm the qualitative picture described in the previous section. Namely, the population of the fms gradually decreases. The population of the Alfvén waves grows across the whole spectrum until both populations become comparable, n k ∼ N k , at k ⊥ ∼ k 0 . This occurs at q 0 t ∼ 1 − 2. At larger times, only the Alfvén population of the high k ⊥ states grows, whereas at smaller k ⊥ , a plateau is formed, the level of which slowly decreases.

ASTROPHYSICAL APPLICATIONS
The nonlinear interactions of low-frequency waves play an important role in compact sources of powerful radio emission, such as pulsars and FRBs. The results of this paper could be applied if the waves may be described in the MHD limit, i.e., when the emission frequency is well below both the Larmor and the plasma frequencies.

Fast radio bursts
FRBs are radio pulses of millisecond duration coming from cosmological distances and having isotropic luminosities L iso ∼ 10 42 − 10 45 erg⋅s −1 . The origin of these pulses is still not known; however, there is some evidence that they are associated with magnetar flares (see, e.g., review by Zhang 2022 and references therein). In the magnetar magnetospheres, both the Larmor and the plasma frequencies are well above the radio band; therefore, independently of the emission mechanism, radio waves propagate as fms waves. Then, the nonlinear decay of fms into the Alfvén waves places severe limits on radio emission power if FRBs are produced well within the magnetar magnetospheres.
The transformation efficiency is determined by the product of the transformation rate (21) and the propagation time, r c. Assuming the dipole magnetic field in the magnetosphere, B = µ r 3 , and expressing the fms energy density via the isotropic FRB luminosity, U = L FRB 4πcr 2 , one finds Here f = ω 2π is the radiation frequency, and we employ the standard short-hand notation, s = 10 x s x in cgs units. One sees that at any FRB power, the radiation is absorbed at a distance from a few to a few dozen stellar radii. This estimate assumes that the turbulence is weak, i.e., the ratio of the wave energy to the energy of the background field is small. A simple estimate shows that this condition is fulfilled in the region of interest: This justifies the neglect of four-wave processes 1 . Note that the wave amplitude exceeds the background field when this ratio exceeds unity. In this case, the fms radiation is heavily absorbed because the MHD condition, E < B, is violated (Beloborodov 2021(Beloborodov , 2023. If the background field is a dipole, this happens at distances of a few hundred stellar radii. However, the magnetic disturbance from the magnetar flare propagates away as a large-scale electro-magnetic pulse, whose amplitude decreases as r −1 . In this case, the high-frequency waves propagate at the top of the pulse. Then, the wave amplitude remains smaller than the background field, so Beloborodov's mechanism does not work. On the other hand, the fms-to-Alfvén transformation provides an effective absorption in any case.
The important point is that the wave transformation is a stimulated process, i.e., the transformation rate is proportional to the occupation number of the waves in the final state. The high rate is obtained because the initially small density of Alfvén waves grows exponentially and rapidly becomes comparable with the density of fms waves. The transformation could be suppressed if an efficient absorption mechanism does not permit the growth of Alfvén waves.
The newly produced Alfvén waves could be absorbed because of the current starvation, i.e., when the parameter ξ = j enc (27) exceeds unity (Thompson & Gill 2014;Thompson 2023).
Here n is the plasma density, j = ck ⊥ δB 4π the current density in the Alfvén wave with the amplitude δB. Let us consider this process.
In the magnetosphere of an active magnetar, the electron-positron pairs are produced by slow untwisting of magnetospheric magnetic field lines (Beloborodov & Thompson 2007;Beloborodov 2013). The plasma density is estimated as (Beloborodov 2020) where M ∼ 10 3 is the pair multiplicity, r ± = 5⋅10 6 µ 1 3 33 cm the distance from the star where the magnetic field falls to 10 13 G, so the pair production stops. Substituting this estimate into eq. (27) and assuming that the FRB energy is completely transferred to Alfvén waves, L FRB = δB 2 r 2 c, one finds One sees that the current starvation sets in only when a significant fraction of the FRB energy is transferred to Alfvén waves. In this case, the FRB is absorbed via the transformation to Alfvén waves, which decay because of current starvation. However, the important point is that the pairs are heated to high Lorentz factors and produce new pairs. Let us assume that a fraction ζ < 1 of the FRB energy is absorbed. Then the acquired Lorentz factor is estimated as γ = ζL 4πr 2 m e c 3 n = 10 7 r 7 ζL FRB,43 M 3 µ 2 3 33 .
The particles with this Lorentz factor emit curvature photons with the energy MeV.
These photons produce pairs just as in pulsars. The condition for single photon pair production is Here θ is the angle between the photon direction and the magnetic field, B q = m 2 e c 3 ̵ he 2 = 4.4 ⋅ 10 13 G the quantum magnetic field. The photon is emitted along the magnetic field and, after passing the distance x, acquires the angle θ = x r. Now one finds The condition (32) is fulfilled, so all the emitted photons are converted to pairs. One can easily check that the power of the curvature emission is sufficient for the particle to lose the whole energy to radiation. Therefore each particle emits More pairs are produced from the synchrotron photons emitted by the newly produced pairs. In any case, one sees that absorption of a fraction of the FRB energy produces enough pairs to provide conditions for fms-to-Alfvén decay in the MHD regime. The above consideration shows that FRBs could not be generated well within the magnetar magnetosphere. However, one could not directly extrapolate the obtained conclusion to the outer magnetosphere or the magnetar wind. First of all, the plasma density rapidly decreases with the distance, so eventually, waves in the radio band could not be described in the MHD approximation. Moreover, the magnetic perturbation from the magnetar flare propagates away as a large-scale MHD pulse, which amplitude decreases as 1 r, so that in the outer magnetosphere, the magnetic field of the pulse exceeds the dipole field and the plasma is pushed away with relativistic velocities. Therefore the ratio of the radiation energy density to the energy density of the background field stops growing, and moreover, one has to take into account the relativistic slowing down of time. This implies that magnetar flares could produce FRBs only far enough from the magnetar (see, e.g., review by Lyubarsky 2021 and references therein).

Radio emission of the Crab pulsar
A typical pulsar produces a pencil beam of radio emission, presumably generated in the electron-positron plasma flowing along the magnetic axis of the neutron star within a narrow open field line tube. The emission is generally attributed to plasma oscillations in the flow (see the recent review by Philippov & Kramer 2022 and references therein). The frequency of these waves is comparable with the plasma frequency in the plasma rest frame; therefore, they could not be considered MHD waves. The nonlinear process under consideration is irrelevant to this emission.
In pulsars with large magnetic fields at the light cylinder, such as the Crab and millisecond pulsars, there is another emission site, namely, the current sheet separating, beyond the light cylinder, the oppositely directed magnetic fields. The energy release due to the magnetic reconnection in the current sheet feeds the powerful synchrotron emission in the gamma-ray, and sometimes also in the X-ray and optical, band (Lyubarskii 1996;Bai & Spitkovsky 2010;Cerutti et al. 2016). The fan beam thus formed rotates with the neutron star so that the observer typically sees two peaks per pulsar period. Some of these pulsars also exhibit radio pulses in phase with high-energy pulses. This radio emission could be produced because magnetic islands in the reconnecting current sheet continuously merge, giving rise to magnetic perturbations that propagate away in the form of fms waves, which further away are transformed into radio waves (Uzdensky & Spitkovsky 2014;Lyubarsky 2019;Philippov et al. 2019). The nonlinear interaction of fms waves places limits on the radio luminosity of these pulsars.
According to simulations by Philippov et al. (2019), about 0.5% of the total energy release in the current sheet is radiated away in the form of low-frequency waves. The luminosity of the Crab pulsar in the X-and γ-ray bands is roughly L hard = 10 36 erg⋅s −1 . This quantity could be considered a proxy for the energy release rate in the current sheet. Then, one would expect the radio luminosity of the Crab to be of the order of L ′ = 5⋅10 33 erg⋅s −1 . However, the observed radio luminosity is two orders of magnitude smaller, L radio = 7 ⋅ 10 31 erg⋅s −1 (e.g., Malov et al. 1994). This may be attributed to the decay of the fms into the Alfvén waves on the way out of the magnetosphere.
The decay rate (21) is calculated in the zero electric frame of the plasma because the non-linear interactions in the force-free regime are not affected by plasma moving along the magnetic field lines. Just beyond the light cylinder, the magnetospheric electric and magnetic fields are of the same order but not too close to each other so that the velocity of the zero electric field frame, v = cE × B B 2 , is only mildly relativistic. Therefore we use the parameters in the lab frame. We find the radiation energy density, U , from the condition that the decay rate (21) is comparable with the wave escape rate, q 0 ∼ c r. The radiation is produced near the light cylinder; therefore, r ∼ c Ω, where Ω = 2π P is the angular velocity of the neutron star, and P is the pulsar period. Now this condition yields The magnetic field in the equatorial zone at the distance of the light cylinder is B = 2µ(Ω c) 3 , where µ is the magnetic moment of the neutron star. The last is related to the pulsar spin-down power: where ψ is the angle between the magnetic and rotational axes (Spitkovsky 2006). The radio emission forms a fan beam with the opening angle α ∼ 0.1, so that the radio luminosity may be presented as Now one finds where f = ω 2π is the radiation frequency. The period of the Crab pulsar is P = 0.033 s, and the spectrum is very steep without the low-frequency cutoff down to the decameter band. Therefore, we take a low frequency f = 30 MHz (Malov et al. 1994). Then one gets L radio L sd ∼ 10 −7 . Taking into account that the observed slowing down rate of the Crab corresponds to the spindown power L sd = 5 ⋅ 10 38 erg⋅s −1 , one sees that the obtained estimate is compatible with the observed radio luminosity. Now let us check that the MHD conditions are fulfilled, i.e., the plasma density is sufficient to maintain Alfvén waves. The plasma density in pulsars may be presented as where κ is the multiplicity. In young pulsars, κ > 10 5 (Timokhin & Harding 2015), which is compatible with the observations of the Crab Nebula (de Jager et al. 1996). Assume that the emitted power, L ′ , (which is larger than the observed luminosity, L radio , see above) is transformed into Alfvén waves so that the relation between L ′ and the Alfvén energy density is the same as the relation (37) between the observed luminosity and the radiation energy density. Now we estimate the parameter ξ (see eq. 27) as = 0.2 f 7.5 P −1.5 L ′1 2 33.7 α −1 κ 5 (1 + sin 2 ψ)L 1 2 sd,38.7 k ⊥ c ω .
One sees that the density is marginally sufficient to maintain fms-to-Alfvén transformation at the expected luminosity of the current sheet. The produced Alfvén waves decay because of current starvation when the energy is transferred to higher k ⊥ , as is described in previous sections. This justifies our explanation of the low luminosity of the Crab pulsar.

CONCLUSIONS
In this paper, we addressed the nonlinear decay of fms waves in relativistic force-free MHD. There are models of pulsars and FRBs, in which the radio emission is generated in dense magnetospheres such that not only the Larmor but also the plasma frequency is well below the radiation frequencies (see, e.g., reviews by Zhang 2022; Lyubarsky 2021; Philippov & Kramer 2022). Within these sources, the radiation propagates in the form of fms waves, and the nonlinear decay of fms into Alfvén waves could strongly affect the properties of the outgoing radiation.
Using the kinetic equations for the waves, we estimated the decay rate and studied the kinetics of the decay process. We have shown that an fms wave could decay into two Alfvén waves with arbitrary large wave vectors if these wave vectors are nearly perpendicular to the background magnetic field and nearly oppositely directed. Therefore the phase volume available for the decay of an fms wave is in fact infinite. In this case, the energy of fms waves could be completely transferred to the small-scale Alfvén waves not via a cascade, as in the Kolmogorov turbulence, but directly. Numerical solutions of the kinetic equations confirmed these conclusions. Our results explain the anomalously low radio efficiency (the ratio of the radio to the spin-down power) of the Crab pulsar and demonstrate that FRBs could not be produced well within magnetar magnetospheres.
The set of equations (4) and (5) was solved by the Runge-Kutta method. It is known that the determination of the error and stability of the Runge-Kutta method is quite difficult Butcher (2015). As a rule, a very small time step is necessary to minimize the error and ensure the stability of the method. However, if there are conserved parameters in the system (e.g., total energy, total mass, helicity, etc.), one can achieve better convergency and stability by making use of a correction procedure that explicitly exploits the conservation laws (see, e.g., Christlieb et al. 2011;Palha & Gerritsma 2017;Coppola et al. 2019 and references there).
In our case, the conserved parameter is the total energy of the system, U = ω(N k + n k )d 3 k.
This means that summing up the rhs of equations (4) and (5), multiplying the obtained expression by ω and integrating over the whole phase space yields zero. However, the conservation law is violated when the integrals in the rhs of these equations are evaluated numerically.
To keep the total energy exactly conserved, we introduce a correction to the obtained values of N (k) and n(k) at each time step. Let U 1 be the energy of the system obtained after the i-th time step, and U 0 be the initial energy of the system. We multiply N i+1 (k) and n i+1 (k) by ξ = 2 − U1 U0 . Taking into account the identity one sees that after such a correction, the total energy of the system remains equal to U 0 to within (1 − U1 U0 ) 2 . In our calculations the value of 1 − U1 U0 does not exceed 10 −3 . Consequently, the energy is conserved to within ∼ 10 −6 . In figure 4, we show how the above correction procedure affects the stabilization of the calculations. The calculations were carried out with the time step q 0 ∆t = 10 −4 and values ∆k ⊥ = ∆k ∥ = 0.1k 0 . One sees that without correction (Fig.4a), the system quickly becomes unstable. When using the correction procedure (Fig.4b), the system remains stable, and the total energy of the system is conserved. Figure 5 shows the spectra of Alfvén waves calculated with different input parameters. The spectra are presented at q 0 t = 1, when the energies of the fms and Alfvén waves are comparable. Figure 5a shows the spectra calculated with different time steps, q 0 ∆t, at value ∆k ⊥ k 0 = 0.1. One sees that when the time resolution is improved, the "waves" in the spectra disappear, whereas the overall spectral shape remains intact. Figure 5b shows the spectra obtained for different values ∆k ⊥ k 0 at q 0 ∆t = 10 −4 . One sees that improving the resolution in the wave vector space does not significantly affect the shape of the spectrum.