© The Authors 2024

1. Introduction

Fast radio bursts (FRBs) are millisecond-duration and energetic radio pulses (see Cordes & Chatterjee 2019; Petroff et al. 2019, for reviews). They are similar in some respects to single pulses from radio pulsars, however, the bright temperature of FRBs is many orders of magnitude higher than that of normal pulses from Galactic pulsars. The physical origin(s) of FRB are still unknown and the only certainty is that their radiation mechanisms must be coherent.

Magnetars have been invoked as the most likely candidate to interpret the emission of FRBs, especially repeating ones. The models can be generally divided into two categories by considering the distance of the emission region from the magnetar (Zhang 2020): emission within the magnetosphere (pulsar-like model, e.g., Katz 2014; Kumar et al. 2017; Yang & Zhang 2018; Wadiasingh & Timokhin 2019; Lu et al. 2020; Cooper & Wijers 2021; Wang et al. 2022a; Zhang 2022; Liu et al. 2023a; Qu & Zhang 2023) and emission from relativistic shocks region far outside the magnetosphere (GRB-like model, e.g., Metzger et al. 2019; Beloborodov 2020; Margalit et al. 2020; Chen et al. 2022; Khangulyan et al. 2022). A magnetar origin for at least some FRBs was established after the detection of FRB 20200428D¹, an FRB-like burst from SGR J1935+2154, which is a Galactic magnetar (Bochenek et al. 2020; CHIME/FRB Collaboration 2020a). Nevertheless, the localization of FRB 20200120E to an old globular cluster in M 81 (Bhardwaj et al. 2021; Majid et al. 2021; Kirsten et al. 2022) poses a challenge to the assumption of a young magnetar born via a massive stellar collapse. This would mean that the magnetar formed in accretion-induced collapse or a system associated with a compact binary merger (Kremer et al. 2021; Lu et al. 2022).

The spectral-time characteristics of FRBs are significant tools that are useful for gaining a deeper understanding of emission mechanisms. The burst morphologies of FRBs are varied. Some FRBs have a narrow bandwidth, for instance, with the sub-bursts of FRB 20220912A showing that the relative spectral bandwidth of the radio bursts was distributed near at Δν/ν ≃ 0.2 (Zhang et al. 2023). Similar phenomena of narrowband can appear at FRB 20121102A, in which spectra are modeled by simple power laws with indices ranging from −10 to +14 (Spitler et al. 2016). Simple narrowband and other three features (simple broadband, temporally complex, and downward drifting) are summarized as four observed archetypes of burst morphology among the first CHIME/FRB catalog (Pleunis et al. 2021). The narrow spectra may result from small deflection angle radiation or coherent processes. Alternatively, they could be modulated by scintillation and plasma lensing in the FRB source environment (Yang 2023).

Downward drift in the central frequency of consecutive sub-bursts with later-arriving time, has been discovered in at least some FRBs (CHIME/FRB Collaboration 2019a; Hessels et al. 2019; Fonseca et al. 2020). In addition, a variety of drifting patterns including upward drifting, complex structure, and no drifting subpulse, was found in a sample of more than 600 bursts detected from the repeater source FRB 20201124A (Zhou et al. 2022). Such a downward drifting structure was unprecedentedly seen in pulsar PSR B0950+08 (Bilous et al. 2022), suggesting that FRBs originate from the magnetosphere of a pulsar. Downward drifting patterns can be aptly understood by invoking magnetospheric curvature radiation (Wang et al. 2019), however, it is unclear why the spectral bandwidth can be extremely narrow, namely, Δν/ν ≃ 0.1.

Polarization properties carry significant information about radiation mechanisms to shed light on the possible origin of FRBs. High levels of linear polarization are dominant for most sources and most bursts for individual repeating sources (e.g., Michilli et al. 2018; Day et al. 2020; Luo et al. 2020; Sherman et al. 2023), while some bursts have significant circular polarization fractions, which can be even up to 90% (e.g., Masui et al. 2015; Xu et al. 2022; Jiang et al. 2022; Kumar et al. 2022). The properties are reminiscent of pulsars, which show a wide variety of polarization fractions between sources. The emission is highly linearly polarized when the line of sight (LOS) is confined into a beam angle, whereas highly circular polarization is present when the LOS is outside, when considering bunched charges moving at ultra-relativistic speed in the magnetosphere (Wang et al. 2022b; Liu et al. 2023b).

In this work, we attempt to understand the origin(s) of the observed narrow spectra of FRBs from the perspective of intrinsic radiation mechanisms and absorption effects in the magnetosphere. The paper is organized as follows. The intrinsic radiation mechanisms by charged bunches in the magnetosphere are discussed in Sect. 2. We investigate two scenarios of the spectrum, including the bunch mechanism (Sect. 2.1) and perturbations on the bunch moving at curved trajectories (Sect. 2.2). The polarization properties of the intrinsic radiation mechanisms are summarized in Sect. 2.3. In Sect. 3, we discuss several possible mechanisms as triggers that may generate emitting charged bunches (Sect. 3.1), as well as oscillations and Alfvén waves. Some absorption effects during the wave propagating in the magnetosphere are investigated in Sect. 4. We discuss the bunch’s evolution in the magnetosphere in Sect. 5. The results are discussed and summarized in Sects. 6 and 7. The convention Q_x = Q/10^x in cgs units and spherical coordinates (r,  θ,  φ) relating to the magnetic axis are used throughout the paper.

2. Intrinsic radiation mechanism

A sudden trigger may occur on the stellar surface and last at least a few milliseconds to create free charges. Such a trigger event is a sudden and violent process in contrast to a consecutive “sparking” process in the polar cap region of a pulsar (Ruderman & Sutherland 1975). The charges can form bunches and the total emission of charges is coherently enhanced significantly if the bunch size is smaller than the half wavelength. We note that electrons as the emitting charges are the point of interest in the following discussion, leading to negatively net-charged bunches. We attempt to study how to trigger the emitting charges and how to form bunches in Sect. 3.

Within the magnetosphere, the charged particles are suddenly accelerated to ultra-relativistic velocities and stream outward along the magnetic field lines. The charges hardly move across the magnetic field lines because the cyclotron cooling is very fast in a strong magnetic field and charges stay in the lowest Landau level so that charges’ trajectories essentially track with the field lines.

Curvature radiation can be produced via the perpendicular acceleration for a charge moving at a curved trajectory. It is a natural consequence that charges moving at curved trajectories, which has been widely discussed in pulsars (Ruderman & Sutherland 1975; Cheng & Ruderman 1977; Melikidze et al. 2000; Gil et al. 2004; Gangadhara et al. 2021) and FRBs (Kumar et al. 2017; Katz 2018; Lu & Kumar 2018; Ghisellini & Locatelli 2018; Yang & Zhang 2018; Wang et al. 2022a; Cui et al. 2023). However, the noticeable difference is that there should be an electric field parallel (E_∥) to the local magnetic field B in the FRB emission region because the enormous emission power of FRB lets the bunches cooling extremely fast. The E_∥ sustains the energy of bunches so that they can radiate for a long enough time to power FRBs. The equation for motion in the FRB emission region in the lab frame can be written as:

$$\begin{array}{c}\hfill N_{\mathrm{e}}{E}_{\Vert }e\mathrm{d}s-{\mathcal{L}}_{\mathrm{b}}\mathrm{d}t=N_{\mathrm{e}}{m}_{\mathrm{e}}{c}^2\mathrm{d}\gamma ,\end{array}$$(1)

where N_e is the number of net charges in one bunch that are radiating coherently, c is the speed of light, e is the elementary charge, m_e is the positron mass, γ is the Lorentz factor of the bunch, and ℒ_b is the luminosity of a bunch. The timescale for the E_∥ that can be balanced by radiation damping is essentially the cooling timescale of the bunches, which is much shorter than the FRB duration. As a result, bunches stay balanced throughout the FRB emission process.

Bunches in the magnetosphere are accelerated to ultra-relativistic velocity. A stable E_∥ (∂E/∂t = 0) may sustain the Lorentz factor of bunches to be a constant. The emission of bunches is confined mainly in a narrow conal region due to the relativistic beaming effect, so the line of sight (LOS) sweeps across the cone within a very short time. According to the properties of the Fourier transform, the corresponding spectrum of the beaming emission within a short time is supposed to be wide-band. This is a common problem among the emission from relativistic charges. In the following discussion, we strive to understand the narrowband FRB emission by considering the bunching mechanism (Sect. 2.1) and perturbation on the curved trajectories (Sect. 2.2).

2.1. Bunching mechanism

The waves are coherently enhanced significantly when the emitting electrons form charged bunches, which is the fundamental unit of coherent emission. Different photon arrival delays are attributed to charges moving in different trajectories, so that charges in the horizontal plane (a light solid blue line denotes a bunch shown in Fig. 1) have roughly the same phases. The bunch requires a longitude length that is smaller than the halfwavelength, which is ∼10 cm for the GHz wave.

Inside a bunch, the total emission intensity is the summation of single charges that are coherently added. We considered a single charge that moves along a trajectory r_s(t). The energy radiated per unit solid angle per unit frequency interval is interested in considering spectrum, which is given by (Rybicki & Lightman 1979; Jackson 1998):

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill \frac{{\mathrm{d}}^{2}W}{\mathrm{d}\omega \mathrm{d}\mathrm{\Omega }}=& \frac{e^2{\omega }^2}{4{\pi }^{2}c}\hfill \\ \hfill & \times {|{\displaystyle {\int }_{-\infty }^{+\infty }\sum _{\mathrm{s}}^{N_{\mathrm{s}}}}-{\mathit{\beta }}_{\mathrm{s}\perp }exp[i\omega (t-\mathit{n}·{\mathit{r}}_{\mathrm{s}}/c)]\mathrm{d}t|}^2,\hfill \end{array}\end{array}$$(2)

where β_s⊥ is the component of β_s in the plane that is perpendicular to the LOS:

$$\begin{array}{c}\hfill {\mathit{\beta }}_{\mathrm{s}\perp }=-\mathit{n}\times (\mathit{n}\times {\mathit{\beta }}_{\mathrm{s}}).\end{array}$$(3)

ϵ₁ is the unit vector pointing to the center of the instantaneous circle, ϵ₂ = n × ϵ₁ is defined. If there is more than one charged particle, one can use subscripts (i,  j,  k) to describe any charge in the three-dimensional (3D) bunch. The subscripts contain all the information about the location. We assume that bunches are uniformly distributed in all directions of ($\widehat{\mathit{r}},\widehat{\mathit{\theta }},\widehat{\mathit{\phi }}$). Equation (3) can be written as:

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill \frac{{\mathrm{d}}^{2}W}{\mathrm{d}\omega \mathrm{d}\mathrm{\Omega }}=& \frac{e^2{\omega }^2}{4{\pi }^{2}c}\hfill \\ \hfill & \times {|{\displaystyle {\int }_{-\infty }^{+\infty }\sum _i^{N_i}\sum _j^{N_j}\sum _k^{N_k}}-{\mathit{\beta }}_{ijk,\perp }exp\left[i\omega (t-\mathit{n}·{\mathit{r}}_{\mathit{ijk}}(t)/c)\right]\mathrm{d}t|}^2,\hfill \end{array}\end{array}$$(4)

where N_e = N_iN_jN_k. The emission power from charges with the same phase is $N_{\text{e}}^2$ enhanced. Detailed calculations are referred to Appendix A.

The bunches thus act like single macro charges in some respects. We define a critical frequency of curvature radiation, ω_c = 3cγ³/(2ρ), where ρ is the curvature radius. Spectra of curvature radiation from a single charge can be described as a power law for ω ≪ ω_c and exponentially drop when ω ≫ ω_c. The spectrum and polarization properties of a bunch with a longitude size smaller than the half-wavelength and half-opening angle φ_t ≲ 1/γ, are similar to the case of a single charge. In general, the spectra can be characterized as multi-segmented broken power laws by considering a variety of geometric conditions of bunch (Yang & Zhang 2018; Wang et al. 2022a).

2.1.1. Bunch structure

There could be a lot of bunches contributing to instantaneous radiation in a moving bulk and their emissions are added incoherently. The spread angle of relativistic charge is:

$$\begin{array}{c}\hfill {\theta }_{\mathrm{c}}(\omega )\simeq \{\begin{array}{ccc}& \frac{1}{\gamma }{\left(\frac{2{\omega }_{\mathrm{c}}}{\omega }\right)}^{1/3}={\left(\frac{3c}{\omega \rho }\right)}^{1/3},\hfill & \hfill \omega \ll {\omega }_{\mathrm{c}}\\ \hfill & \frac{1}{\gamma }{\left(\frac{2{\omega }_{\mathrm{c}}}{3\omega }\right)}^{1/2},\hfill & \hfill \omega \gg {\omega }_{\mathrm{c}}.\end{array}\end{array}$$(5)

For ω ∼ ω_c, bunches within a layer of ρ/γ can contribute to the observed luminosity at an epoch. The number of that contributing charges is estimated as $N_{\text{b}}~2\rho \nu /(c\gamma )\simeq {10}^4{\rho }_7{\nu }_9{\gamma }_2^{-1}$.

We assume that the number density in a bunch has ∂n_e/∂t = 0. The case of fluctuating bunch density can give rise to a suppressed spectrum, and there is a quasi-white noise in a wider band in the frequency domain (Yang & Zhang 2023). The fluctuating density may cause the integral of Eq. (A.5) to not converge.

Different bunches may have different bulk Lorentz factors. We note that emissions from relativistic charges may have a characteristic frequency determined by the Lorentz factor. The scattering of the Lorentz factor from different emitting units likely enables the characteristic frequency to be polychromatic. Therefore, non-monoenergetic distributed charges give rise to a wide spectrum. However, for the bunch case discussed in Sect. 2, the layer contributing to the observed luminosity simultaneously is ∼ρ/γ much smaller than the curvature radius. Emission power in the layers is thus roughly the same. The Lorentz factors between different bunch layers can be regarded as a constant.

The frequency structure of a burst is essentially the Fourier transform of the spatial structure of the radiating charge density (Katz 2018). Therefore, a relatively narrow spectrum may be produced via the emission from a quasi-periodically structured bulk of bunches. The possible mechanisms that generate such a quasi-periodically structure are discussed in Sect. 3.

2.1.2. Spectrum of the bunch structure

We focus on the radiation properties due to the quasi-periodical structure. The charged bunch distribution is thought to be quasi-period in the time domain, which may be attributed to some quasi-periodic distribution in space between bunches. Each emitting bunch has the same E(t), but with different arrival times. The total received electric field from the multiple bunches simultaneously is given by:

$$\begin{array}{c}\hfill E(t)={\displaystyle \sum _n^{N_{\mathrm{b}}}}{E}_n(t-t_n),\end{array}$$(6)

and its time-shifting property of the Fourier transform is:

$$\begin{array}{c}\hfill E(\omega )={\displaystyle \sum _n^{N_{\mathrm{b}}}}{E}_n(\omega )exp(i\omega t_n).\end{array}$$(7)

A quasi-periodically distributed bulk of bunches means that ωt_n = nω/ω_M + δϕ_n, where 1/ω_M denotes the period of the bunch distribution and ϕ_n is the relative random phase.

A strictly periodically distributed bulk means that the relative random phase is zero strictly. The energy per unit time per unit solid angle is (Yang 2023):

$$\begin{array}{c}\hfill \frac{{\mathrm{d}}^{2}W}{\mathrm{d}\omega \mathrm{d}\mathrm{\Omega }}=c^2{\mathcal{R}}^2{\left|E_n(\omega )\right|}^2\frac{{sin}^2[N_{\mathrm{b}}\omega /(2{\omega }_{\mathrm{M}})]}{{sin}^2[\omega /(2{\omega }_{\mathrm{M}})]},\end{array}$$(8)

where ℛ is the distance from the emitting source to the observer. The coherence properties of the radiation by the multiple bulks of bunches are determined by ω_M. The radiation is coherently enhanced when ω is 2πω_M of integer multiples so that the energy is radiated into multiples of 2πω_M with narrow bandwidths. Even when ϕ_n is not zero, the Monte Carlo simulation shows that the coherent radiation energy is still radiated into multiples of 2πω_M with narrow bandwidths (Yang 2023).

Figure 2 shows the simulation of radiation spectra from quasi-periodically distributed bulk of bunches with ω_M = 10⁸ Hz and ω_M = 5 × 10⁸ Hz. We assume that N_b = 10⁴, γ = 10², δϕ_n = 0 and ρ = 10⁷ cm. Since the flux drops to a small number rapidly when both 1/γ < φ′ and 1/γ < χ′, we fix χ′ = 10⁻³ in the simulations for simplicity. Spectra from uniformly distributed bulk of bunches are plotted as comparisons. By assuming there is an observation threshold exceeding the peak flux of the pure curvature radiation, the spectrum is characterized by multiple narrowband emissions with a bandwidth smaller than 2πω_M. As ω_M becomes larger, the number of spikes is smaller. Even though there are multi-spikes in the spectra, there is a main maximum at near ν = ν_c ordered larger than other spikes, so that multiple narrowband emissions are difficult to detect with an ultra-wideband receiver.

2.2. Perturbation on curved trajectories

We considered two general cases that the perturbations exist both parallel (E_1, ∥) and perpendicular (E_1, ⊥) to the local B-field independently. The perturbations can produce extra parallel $c{\dot{\beta }}_{\Vert }$ and perpendicular accelerations $c{\dot{\beta }}_{\perp }$ to the local B-field, as shown in Fig. 1. The perturbation electric fields are much smaller than the E_∥, which is ∼10⁷ esu, sustaining a constant Lorentz factor (Wang et al. 2019). When they are affected by these accelerations, electrons not only emit curvature radiation, but also emissions caused by the perturbations on the curved trajectory. As we have limited this issue to general discussion, we do not consider the perturbation in this section; instead, we simply require for the velocity to change slightly, namely, δβ ≪ β.

2.2.1. Parallel acceleration

The emission power by the acceleration parallel to the velocity, namely, longitudinal acceleration, is given by (Rybicki & Lightman 1979):

$$\begin{array}{c}\hfill P=\frac{2e^2}{3m_{\mathrm{e}}^2{c}^3}{\left(\frac{\mathrm{d}p}{\mathrm{d}t}\right)}^2.\end{array}$$(9)

According to Eq. (9), if the ratio of two longitudinal forces is N, the ratio of their emission power to N², so that we can ignore the radiation damping caused by E_1, ∥. The energy loss by radiation is mainly determined by the coherent curvature radiation. The balance between the E_∥ and radiation damping is still satisfied approximately, leading to a constant Lorentz factor. Thus, the velocity change caused by E_1, ∥ of a single charge can be described as:

$$\begin{array}{c}\hfill E_{1,\Vert }e\mathrm{d}t\simeq m_{\mathrm{e}}{\gamma }_{\mathrm{s}}^3\mathrm{d}{v}_{\mathrm{s}}.\end{array}$$(10)

The velocity change is so slight that the partial integration of Eq. (A.5) converges, thus Eq. (A.6) can be obtained from Eq. (A.5).

2.2.2. Spectrum of parallel acceleration

The dimensionless velocity of a single charge reads β_s = β_s, 0 + β_1, ∥, where β_s, 0 is the mean dimensionless velocity and β_1, ∥ is its perturbation parallel to the field lines, so that the emitting electric field can be calculated by the curvature radiation and perturbation-induced emission independently. Let us consider a charge with a dentifier of i, j, k in a bunch. The polarization components of the parallel perturbation-induced wave are calculated as:

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill {\stackrel{\sim }{A}}_{ijk,\Vert }\simeq & {\displaystyle {\int }_{-\infty }^{+\infty }}{\beta }_{1,\Vert }(t)sin(\frac{\mathit{ct}}{\rho }+{\chi }_{\mathit{ij}})\hfill \\ \hfill & \times exp\left(i\frac{\omega }{2}[(\frac{1}{{\gamma }^2}+{\phi }_k^2+{\chi }_{\mathit{ij}}^2)t+\frac{c^2{t}^3}{3{\rho }^2}+\frac{c{t}^2}{\rho }{\chi }_{\mathit{ij}}]\right)\mathrm{d}t,\hfill \\ \hfill {\stackrel{\sim }{A}}_{ijk,\perp }\simeq & {\displaystyle {\int }_{-\infty }^{+\infty }}{\beta }_{1,\Vert }(t)sin{\phi }_{k}cos(\frac{\mathit{ct}}{\rho }+{\chi }_{\mathit{ij}})\hfill \\ \hfill & \times exp\left(i\frac{\omega }{2}[(\frac{1}{{\gamma }^2}+{\phi }_k^2+{\chi }_{\mathit{ij}}^2)t+\frac{c^2{t}^3}{3{\rho }^2}+\frac{c{t}^2}{\rho }{\chi }_{\mathit{ij}}]\right)\mathrm{d}t.\hfill \\ \hfill \end{array}\end{array}$$(11)

A complex spectrum can be caused by a complexly varying β_1, ∥(t). A simple case is that the integrals of Eq. (11) can be described by the modified Bessel function K_ν(ξ) if β₁ is a constant. The function is approximated to be power-law-like for ξ ≪ 1 and exponential-like for ξ ≫ 1, when ν ≠ 0. Here, it is worth noting that K_ν(ξ = 0)→+∞ when ξ = 0, exhibiting a very narrow spectrum.

In particular, if a waveform E(t) is quasi-sinusoid periodic, the spectrum would be quasi-monochromatic, equivalently to the narrowband. The quasi-sinusoid periodic waveform can be produced by the intrinsic charge radiation if the charge is under a periodic acceleration during its radiation beam pointing to the observer. For a perturbation E_1.∥, which is monochromatic oscillating with time of interest, according to Eq. (A.1), the velocity is also monochromatic oscillating. Therefore, we introduce a velocity form with monochromatic oscillating in the lab frame, that is, β_1, ∥(t) = β₁exp(−iω_at). The two parameters introduced in the curvature radiation in Eq. (A.11) are replaced by:

$$\begin{array}{c}\hfill \begin{array}{cc}& u_{\mathrm{a}}=\frac{\mathit{ct}}{\rho }{(\frac{1}{{\gamma }^2}+{\phi }_k^2+{\chi }_{\mathit{ij}}^2-\frac{2{\omega }_{\mathrm{a}}}{\omega })}^{-1/2},\hfill \\ \hfill & {\xi }_{\mathrm{a}}=\frac{\omega \rho }{3c}{(\frac{1}{{\gamma }^2}+{\phi }_k^2+{\chi }_{\mathit{ij}}^2-\frac{2{\omega }_{\mathrm{a}}}{\omega })}^{3/2}.\hfill \end{array}\end{array}$$(12)

The radiation can be narrowband when γ is a constant and the center frequency is:

$$\begin{array}{c}\hfill \omega =\frac{2{\omega }_{\mathrm{a}}{\gamma }^2}{1+{\gamma }^2({\chi }_{\mathit{ij}}^2+{\phi }_k^2)}·\end{array}$$(13)

For the bunching emission, charges are located with different χ and φ, which results in the bandwidth getting wider. Now we consider the change in the value of E_∥ is much smaller than l∂E_∥/∂r ≪ E_∥, where l is the longitude bunch size. The charges in the bunch are then monoenergic. However, for an ultra-relativistic particle, a small change in velocity can create a violent change in the Lorentz factor, namely, δγ ≃ γ³βδβ_∥. To keep the electrons monoenergic, it is required for the amplitude to be β₁ ≪ 1/γ².

Following the Appendix A, one can obtain the spectrum of bunching charges. Adopting γ = 100, ρ = 10⁷ cm, φ = 0, and ω_a = 2π × 10⁵ Hz, we simulate the spectra in four cases: φ_t = 0.1/γ, β₁ = 10⁻⁶; φ_t = 0.1/γ, β₁ = 10⁻⁵; φ_t = 1/γ, β₁ = 10⁻⁶; and φ_t = 1/γ, β₁ = 10⁻⁵, shown as Fig. 3. The spectra are normalized to the F_ν of curvature radiation with φ_t = 10⁻³ at ω = ω_c. For φ_t = 0.1/γ, the spectra show spikes at 2ν_aγ². The amplitude of spikes decreases as β₁ decreases. For φ_t = 1/γ, the scattering of φ′ leads the spikes to become wider and shorter. The small amplitude β₁ allows for the spike disappear and the spectrum is then determined by the curvature radiation. As φ becomes larger, more different trajectories may have more different azimuths so that the coherence of emission decreases. The spectrum has a flat component at ω < ω_c when 3c/ρ < ω_c(χ^′2 + φ^′2)^3/2 (Yang & Zhang 2018).

2.2.3. Perpendicular acceleration

For a comparable parallel and perpendicular forces, dp/dt, the radiation power from the parallel component is on the order of 1/γ² compared to that from the perpendicular component. The energy loss by radiation is also mainly determined by the coherent curvature radiation. Similar to the parallel scenario, the velocity change caused by E_1, ⊥ of a single charge can be described as:

$$\begin{array}{c}\hfill E_{1,\perp }e\mathrm{d}t\simeq m_{\mathrm{e}}{\gamma }_{\mathrm{s}}c\mathrm{d}{\beta }_{1,\perp },\end{array}$$(14)

when β_1, ⊥ ≪ β_∥. The velocity change is also slight, therefore, we can use Eq. (A.6) to calculate the spectrum (see Appendix A).

2.2.4. Spectrum of perpendicular acceleration

The polarization components of the perpendicular perturbation-induced wave of a charge with an identifier of ijk are calculated as:

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill {\stackrel{\sim }{A}}_{ijk,\Vert }\simeq & {\displaystyle {\int }_{-\infty }^{+\infty }}\frac{1}{2}{\beta }_{1,\perp }(t)[{\phi }_k^2-{({\chi }_{\mathit{ij}}+ct/\rho )}^2]\hfill \\ \hfill & \times exp\left(i\frac{\omega }{2}[(\frac{1}{{\gamma }^2}+{\phi }_k^2+{\chi }_{\mathit{ij}}^2)t+\frac{c^2{t}^3}{3{\rho }^2}+\frac{c{t}^2}{\rho }{\chi }_{\mathit{ij}}]\right)\mathrm{d}t,\hfill \\ \hfill {\stackrel{\sim }{A}}_{ijk,\perp }\simeq & {\displaystyle {\int }_{-\infty }^{+\infty }}{\beta }_{1,\perp }(t){\phi }_k({\chi }_{\mathit{ij}}+ct/\rho )\hfill \\ \hfill & \times exp\left(i\frac{\omega }{2}[(\frac{1}{{\gamma }^2}+{\phi }_k^2+{\chi }_{\mathit{ij}}^2)t+\frac{c^2{t}^3}{3{\rho }^2}+\frac{c{t}^2}{\rho }{\chi }_{\mathit{ij}}]\right)\mathrm{d}t,\hfill \\ \hfill \end{array}\end{array}$$(15)

where χ_ij and φ_k are small angles. The amplitudes are much smaller than those of parallel perturbations with comparable β₁. Similar to the parallel scenario, the total amplitudes of a bunch for the perpendicular perturbation case can be calculated by summation of Eq. (A.17).

A narrow spectrum is expected when γ is a constant with the scattering of χ and φ is small. The spike is at ω = 2γ²ω_a resembling that of the parallel perturbation case. A small change in velocity can also make a violent change in the Lorentz factor but with a noticeable difference from that of longitude motion. Under the energy supply from the E_∥, charges have β_∥ ≈ 1. The change in the Lorentz factor is δγ ≃ γ³β_⊥δβ_⊥ when β_1, ⊥ ≪ 1. To keep electrons being monoenergic, the condition of β₁ ≪ 1/γ should be satisfied at least.

Following the calculation steps in the Appendix A, we simulate the spectra in four cases: φ_t = 0.1/γ, β₁ = 10⁻⁴; φ_t = 0.1/γ, β₁ = 10⁻³; φ_t = 1/γ, β₁ = 10⁻⁴; and φ_t = 1/γ, β₁ = 10⁻³, shown as Fig. 4. Here, we take γ = 100, ρ = 10⁷ cm, φ = 0 and ω_a = 2π × 10⁵ Hz. The spectra are normalized to the F_ν of curvature radiation with φ_t = 10⁻³ at ω = ω_c. The spikes for the four cases are very small so the spectra are dominated by curvature radiations.

2.3. Polarization of the intrinsic radiation mechanisms

Four Stokes parameters can be used to describe the polarization properties of a quasi-monochromatic electromagnetic wave (see Appendix B). Based on Eq. (A.6), the electric vectors of the emitting wave reflect the motion of charge in a plane normal to the LOS. For instance, if the LOS is located within the trajectory plane of the single charge, the observer will see electrons moving in a straight line so that the emission is 100% linearly polarized. Elliptical polarization would be seen if the LOS is not confined to the trajectory plane. We define that φ < θ_c as on-beam and φ > θ_c as off-beam. For the bunch case, similar on-beam and off-beam cases can also be defined in which θ_c is replaced by φ_t + θ_c (Wang et al. 2022b). The general conclusion of the bunching curvature radiation is that high linear polarization appears for the on-beam geometry, whereas high circular polarization is present for the off-beam geometry.

We consider the polarization properties of the bunch structure case (Sect. 2.1.2), and the perturbation scenarios including the parallel and the perpendicular cases (Sect. 2.2). For the quasi-periodically distributed bulk of bunches scenario, both A_∥ and A_⊥ are multiplied by the same factor, so that the polarization properties are similar with a bulk of bunches, which bunches are uniformly distributed (e.g., Wang et al. 2022a). For the perturbation scenario, the perturbation on the curved trajectory can bring new polarization components, which generally obey the on(off)-beam geometry. However, the intensity of the perturbation-induced components is too small then the polarization profiles are mainly determined by the bunching curvature radiation. Consequently, the polarization properties for both the bunching mechanism and the perturbation scenarios can be described by the uniformly distributed bulk of bunches.

We took γ = 10², φ_t = 10⁻² and ρ = 10⁷ cm. The beam angle at ω = ω_c is 0.02. According to the calculations in Appendix B, we simulate the I, L, V, and PA profiles, shown in Fig. 5. These four parameters can completely describe all the Stokes parameters. The total polarization fraction is 100% because the electric vectors are coherently added. The polarization position angle (PA) for the linear polarization exhibits some variations but within 30°. The differential of V is the largest at φ = 0 and becomes smaller when |φ| gets larger. A high circular polarization with flat fractions envelope would appear at off-beam geometry. In this case, the flux is smaller but may still be up to the same magnitude as the peak of the profile (e.g., Liu et al. 2023b). The spectrum at the side of the beam becomes flatter due to the Doppler effect (Zhang 2021). The sign change in the V would be seen when φ = 0 is confined into the observation window (Wang et al. 2022b).

3. Trigger of bunches

The trigger mechanism plays an important role in producing bunching charges and creating the E_∥. Intriguingly, a giant glitch was measured 3.1 ± 2.5 days before FRB 20200428A associated with a quasi-periodic oscillations (QPOs) of ∼40 Hz in the X-ray burst from the magnetar SGR J1935+2154 (Mereghetti et al. 2020; Li et al. 2021, 2022; Ridnaia et al. 2021; Tavani et al. 2021; Ge et al. 2022). Another glitch event was found on October 5, 2020 and the magnetar emitted three FRB-like radio bursts in subsequent days (Younes et al. 2023; Zhu et al. 2023). Both the glitches and QPOs are most likely the smoking guns of starquakes. Therefore, magnetar quakes are considered a promising trigger model of FRB in the following.

As discussed in Sect. 2.1.2, we are particularly interested in the mechanisms that can generate quasi-periodically distributed bunches. Two possible scenarios are discussed independently in Sects. 3.1.1 and 3.1.2. Motivated by the fact that Alfvén waves associated with the quake process have been proposed by lots of models (e.g., Wang et al. 2018; Kumar & Bošnjak 2020; Yang & Zhang 2021; Yuan et al. 2022), we discuss the Alfvén waves as the perturbations in Sect. 2.2.

3.1. Structured bunch generation

3.1.1. Pair cascades from the “gap”

A magnetized neutron star’s rotation can create a strong electric field that extracts electrons from the stellar surface and accelerates them to high speeds (Ruderman & Sutherland 1975). The induced electric field region lacks adequate plasma, namely, charge-starvation, referred to as a “gap”. The extracted electrons move along the curved magnetic field lines near the surface, emitting gamma rays absorbed by the pulsar’s magnetic field, resulting in electron-positron pairs. Such charged particles then emit photons that exceed the threshold for effective pair creation, resulting in pair cascades that populate the magnetosphere with secondary pairs. The electric field in the charge-starvation region is screened once the pair plasma is much enough and the plasma generation stops (Levinson et al. 2005). When the plasma leaves the region, a gap with almost no particles inside is formed forms again and the vacuum electric field is no longer screened.

The numerical simulations show that such pair creation is quasi-periodic (Timokhin 2010; Tolman et al. 2022), leading to quasi-periodically distributed bunches. The periodicity of cascades is ∼3 h_gap/c, where h_gap is the height of the gap. The gap’s height is expressed as (Timokhin & Harding 2015):

$$\begin{array}{c}\hfill h_{\mathrm{gap}}=1.1\times {10}^4{\rho }_7^{2/7}{P}_0^{3/7}{B}_{12}^{-4/7}{|cos\alpha |}^{-3/7}\mathrm{cm},\end{array}$$(16)

where P is the spin period and α is the inclination angle of magnetic axis. For a typical radio pulsar, we take B = 10¹² G and P = 1 s. The corresponding ω_M is 10⁷ Hz. According to Eq. (8), radiation energy is emitted into multiples of ∼10 MHz narrow bandwidths, which is consistent with the narrowband drifting structure of PSR B0950+08 (Bilous et al. 2022).

In contrast to the continuous sparking in the polar cap region of normal pulsars, FRB models prefer a sudden, violent sparking process. Starquake as a promising trigger mechanism can occur when the pressure induced by the internal magnetic field exceeds a threshold stress. The elastic and magnetic energy released from the crust into the magnetosphere are much higher than the spin-down energy for magnetars. Crust shear and motion may create hills on the surface. The height of the gap becomes smaller and the electric field is enhanced due to the so-called point-discharging. Quasi-periodically distributed bunches with ω_M ∼ 10⁸ may be formed from such shorter gaps. Observational consequences of small hills on stellar surfaces could also be known by studying regular pulsars (e.g., Wang et al. 2023).

We note that the sparking process to generate FRB-emitting bunches is powered by magnetic energy, rather than spin-down energy. The stellar oscillations during starquakes can provide additional voltage, making the star active via sparking (Lin et al. 2015). The oscillation-driven sparking process may generate much more charges than the spin-down powered sparking of a normal rotation-powered pulsar.

3.1.2. Two-stream Instability

Alternatively, another possibility to form charged bunches is the two-stream instability between electron-positron pairs. The growth rate of two-stream instability is discussed by invoking linear waves. Following some previous works (e.g., Usov 1987; Gedalin et al. 2002; Yang & Zhang 2023), we consider that there are two plasma components denoted as “1” and “2” with a relative motion along the magnetic field line. The Lorentz factor and the particle number density are γ_j and n_j with j = 1, 2. Each plasma component is assumed to be cold in the rest frame of each component. Then the dispersion relation of the plasma could be written as (Asseo & Melikidze 1998):

$$\begin{array}{c}\hfill 1-\frac{{\omega }_{\mathrm{p},j=1}^2}{{\gamma }_{j=1}^3{(\omega -{\beta }_{j=1}kc)}^2}-\frac{{\omega }_{\mathrm{p},j=2}^2}{{\gamma }_{j=2}^3{(\omega -{\beta }_{j=2}kc)}^2}=0,\end{array}$$(17)

where ω_p, j is the plasma frequency of the component j. For the resonant reactive instability as the bunching mechanism, the growth rate and the characteristic frequency of the Langmuir wave are given by (e.g., Usov 1987; Gedalin et al. 2002; Yang & Zhang 2023):

$$\begin{array}{cc}\hfill \mathrm{\Gamma }& \sim {\gamma }_{j=1}^{-1}{\gamma }_{j=2}^{-1}{\left(\frac{n_{j=1}}{n_{j=2}}\right)}^{1/3}{\omega }_{\mathrm{L}},\hfill \end{array}$$(18)

$$\begin{array}{cc}\hfill {\omega }_{\mathrm{L}}& \sim 2{\gamma }_{j=2}^{1/2}{\omega }_{\mathrm{p},j=2}.\hfill \end{array}$$(19)

The linear amplitude of the Langmuir wave depends on the gain, G ∼ Γr/c. According to Rahaman et al. (2020), a threshold gain, G_th, indicating the breakdown of the linear regime is about G_th ∼ a few. The sufficient amplification would drive the plasma beyond the linear regime once G ∼ G_th, and the bunch formation rate is Γ(G ∼ G_th). The bunch separation corresponds to the wavelength of the Langmuir wave for G ∼ G_th, namely,

$$\begin{array}{c}\hfill l_{\mathrm{sep}}\sim \frac{\pi c}{{\gamma }_{j=2}^{1/2}{\omega }_{\mathrm{p},j=2}}·\end{array}$$(20)

A narrowband spectrum with Δν ∼ 100 MHz would be created if ${\gamma }_{j=2}^{1/2}{\omega }_{\mathrm{p},j=2}\sim {10}^8$ Hz.

Generally, the charged bunches would form and disperse dynamically due to the instability, charge repulsion, velocity dispersion, and so on (e.g., Lyutikov 2021; Yang & Zhang 2023). The spectrum of the fluctuating bunches depends on the bunch formation rate, λ_B, lifetime, τ_B, and bunch Lorentz factor, γ, as proposed by Yang & Zhang (2023). The coherent spectrum by a fluctuating bunch is suppressed by a factor of (λ_Bτ_B)² compared with that of a persistent bunch and there is a quasi-white noise in a wider band in the frequency domain. The observed narrow spectrum implies that, at the minimum, the quasi-white noise should not be dominated in the spectrum, in which case, the condition of 2γ²λ_B ≳ min(ω_peak, 2γ²/τ_B) would be required, where ω_peak is the peak frequency of curvature radiation. We refer to Yang & Zhang (2023) for a detailed discussion.

3.2. Quake-induced oscillation and Alfvén waves

An Alfvén wave packet can be launched from the stellar surface during a sudden crustal motion and starquake. The wave vector is not exactly parallel to the field line so the Alfvén waves have a non-zero electric current along the magnetic field lines. Starquakes as triggers of FRBs are always accompanied by Alfvén wave launching. We investigate whether Alfvén waves as the perturbations discussed in Sect. 2.2 can affect the spectrum.

The angular frequencies of the torsional modes Ω caused by a quake are calculated as follows (also see Appendix C). We chose the equation of state (EOS A; Pandharipande 1971) and EOS WFF3 (Wiringa et al. 1988), which describe the core of neutron stars, as well as to consider two different proposed EOSs for the crusts, including EOSs for NV (Negele & Vautherin 1973) and DH (Douchin & Haensel 2001) models. We match the various core EOS to two different EOSs for the crust. The crust-core boundary for the NV and DH EOSs is defined at ρ_m ≈ 2.4 × 10¹⁴ g cm⁻³, and ρ_m ≈ 1.28 × 10¹⁴ g cm⁻³, respectively.

The frequencies of torsional modes are various with a strong magnetic field (Messios et al. 2001). As already emphasized by Sotani et al. (2007), the shift in the frequencies would be significant when the magnetic field exceeds ∼10¹⁵ G. In the presence of magnetic fields, frequencies are shifted as:

$$\begin{array}{c}\hfill {}_{\ell }{f}_n{=}_{\ell }{f}_n^{(0)}{\left[1{+}_{\ell }{\alpha }_n{\left(\frac{B}{B_{\mu }}\right)}^2\right]}^{1/2},\end{array}$$(21)

where _ℓα_n is a coefficient depending on the structure of the star, ${}_{\mathcal{l}}{f}_n^{(0)}$ is the frequency of the non-magnetized star, and B_μ ≡ (4πμ)^1/2 is the typical magnetic field strength.

In Fig. 6, we show the effects of the magnetic field on the frequencies of the torsional modes. The magnetic field strength is normalized by B_μ = 4 × 10¹⁵ G. Different dashed lines in Fig. 6 are our fits to the calculated numerical data with a high accuracy. For B > B_μ, we find that the frequencies follow a quadratic increase against the magnetic field, and tend to become less sensitive to the NS parameters. The oscillation frequencies exceed 10³ Hz for n = 2 and Ω ∼ 10⁴ Hz for typical magnetars.

Alfvén waves are launched from the surface at a frequency of Ω, which is much faster than the angular frequency of the spin, and propagate along the magnetic field lines. Modifying the Goldreich-Julian density can produce electric fields to accelerate magnetospheric charged particles. The charged particles oscillate at the same frequency as the Alfvén waves in the comoving frame. In the lab frame, the perturbation has ${\omega }_{\text{a}}=\mathrm{\Omega }{\gamma }^2/{\gamma }_{\text{A}}^2$ due to the Doppler effect, where ${\gamma }_{\mathrm{A}}=\sqrt{1+B^2/(4\pi {\rho }_{\mathrm{m}}{c}^2)}$. ω_a is much smaller than ${10}^8{\mathrm{\Omega }}_4{\gamma }_2^2$ Hz because of the strong magnetic field and low mass density in the magnetosphere.

Therefore, by considering that Alfvén waves play a role in the electric perturbations discussed in Sect. 2.2, the spectra have no spike from 100 MHz to 10 GHz, which is the observation bands for FRBs so far. The spectral spike would appear at the observation bands when the perturbations are global oscillations, meaning that all positions in the magnetosphere oscillate simultaneously with Ω rather than a local Alfvén wave packet.

4. Absorption effects in the magnetosphere

In this section, we investigate absorption effects during the wave propagating the magnetosphere. The absorption effects that have been discussed as a possible origin for the narrowband emission are independent of the intrinsic mechanisms in Sect. 2. To create narrowband emission, there must be a steep frequency cut-off on the spectrum that is due to the absorption conditions. Absorption conditions should also change cooperatively with the emission region, that is, the cut-off would appear at a lower frequency as the flux peak occurs at a lower frequency. Curvature radiation spectra can be characterized as a power law for ω < ω_c and exponential function at ω > ω_c (Jackson 1998). In the following discussion, we mainly consider two absorption mechanisms including Landau damping and self-absorption of curvature radiation.

4.1. Wave modes

Here, we consider the propagation of electromagnetic waves in magnetospheric relativistic plasma. The magnetosphere mainly consists of relativistic electron-positron pair plasma, which can affect the radiation created in the inner region. We assume that the wavelength is much smaller than the scale lengths of the plasma and of the magnetic field in which the plasma is immersed. The scale lengths of the magnetic field can be estimated as B/(∂B/∂r)∼r. The plasma is consecutive and its number density is thought to be multiples of the Goldreich-Julian density (Goldreich & Julian 1969), leading to the same scale lengths as that of the magnetic field. Both the scale lengths of magnetic field and plasma are much larger than the wavelength $\lambda =30{\nu }_9^{-1}$ cm, even for the outer magnetosphere at r ∼ 4.8 × 10⁹P₀ cm.

No matter whether it is a neutron star or a magnetar, the magnetospheric plasma is highly magnetized. It is useful to define the (non-relativistic) cyclotron frequency and plasma frequency:

$$\begin{array}{cc}& {\omega }_B=\frac{\mathit{eB}}{m_{\mathrm{e}}c}=1.78\times {10}^{19}{B}_{12}\mathrm{Hz},\hfill \end{array}$$(22)

$$\begin{array}{cc}& {\omega }_{\mathrm{p}}={\left(\frac{4\pi n_{\mathrm{p}}{e}^2}{m_{\mathrm{e}}}\right)}^{1/2}=4.73\times {10}^{11}{\mathcal{M}}_3^{1/2}{B}_{12}^{1/2}{P}_0^{1/2}\mathrm{Hz},\hfill \end{array}$$(23)

where n_p = ℳn_GJ = ℳB/(ceP), in which n_GJ is the Goldreich-Julian density and ℳ is a multiplicity factor. In the polar-cap region of a pulsar, a field-aligned electric field may accelerate primary particles to ultra-relativistic speed. The pair cascading may generate a lot of plasma with lower energy and a multiplicity factor of ℳ ∼ 10² − 10⁵ (Daugherty & Harding 1982). The relationships ω ≪ ω_B and ω_p ≪ ω_B are satisfied at most regions inside the magnetosphere.

The vacuum polarization may occur at the magnetosphere for high energy photons and radio waves by the plasma effect in a strong magnetic field (Adler 1971). In the region of r ≳ 10 R, where FRBs can be generated, in which R is the stellar radius, the magnetic field strength $B<B_Q=m_{\text{e}}^2{c}^3/(e\hslash )$ with assumption of B_s = 10¹⁵ at stellar surface. In the low-frequency limit ℏω ≪ m_ec², the vacuum polarization coefficient can be calculated as ${\delta }_V=1.6\times {10}^{-4}{(B/B_Q)}^2=2.65\times {10}^{-8}{B}_{12}^2$ (Heyl & Hernquist 1997; Wang & Lai 2007). Therefore, the vacuum polarization can be neglected in the following discussion about the dispersion relationship in the magnetosphere.

Waves have two propagation modes: O-mode and X-mode. In the highly magnetized relativistic case, the dispersion of X-mode photon is n ≈ 1. O-mode photons have two branches in the dispersion relationship, which are called subluminous O-mode (n > 1) and superluminous O-mode (n < 1), shown as Fig. D.1. The subluminous O-mode cuts off at 1/cosθ_B, where θ_B is the angle between angle between $\widehat{\mathit{k}}$ and $\widehat{\mathit{B}}$. Following the calculation shown in Appendix D, the dispersion relation of magnetospheric plasma for subluminous O-mode photons can be found following (Arons & Barnard 1986; Beskin et al. 1988):

$$\begin{array}{c}\hfill ({\omega }^2-c^2{k}_{\Vert }^2)(1-{\displaystyle \sum _{\mathrm{s}}}\frac{{\omega }_{\mathrm{p},\mathrm{s}}^2}{{\omega }^2}{g}_{\mathrm{s}})-c^2{k}_{\perp }^2=0,\end{array}$$(24)

with

$$\begin{array}{c}\hfill g_{\mathrm{s}}=\mathcal{P}{\displaystyle {\int }_{-\infty }^{\infty }}\frac{f(u)\mathrm{d}u}{{\gamma }^3{(1-n_{\Vert }\beta )}^2}+i\pi {\beta }_n^2{\gamma }_n^3{\left[\frac{\mathrm{d}f}{\mathrm{d}u}\right]}_{\beta ={\beta }_n},\end{array}$$(25)

where ${\beta }_n=n_{\Vert }^{-1}=\omega /(k_{\Vert }c)$, ${\gamma }_n={(1-{\beta }_n^2)}^{-1/2}$ and f(u) is the distribution function of plasma.

4.2. Landau damping

Wave forces can exist many collective waves in a beam and exchange energy between waves and particles. If the distribution function is df/dv < 0 close to the phase speed of the wave, the number of particles with speeds v < ω/k will be larger than that of the particles having v > ω/k. Then, the transmission of energy to the wave is then smaller than that absorbed from the wave, leading to a Landau damping of the wave. According to the dispersion relationship given in Sect. 4.1, there is no linear damping of the X-mode since wave-particle resonance is absent. The superluminous branch of the O-mode also has no Landau damping because the phase speed exceeds the speed of light. Landau damping only occurs at the subluminous O-mode photons.

The Landau damping decrement can be calculated by the imaginary part of the dispersion relationship (Boyd & Sanderson 2003). We let k = k_r + ik_i, where k_i ≪ k_r ≈ k. If a wave has a form of E ∝ exp(−iωt + ikr), the imaginary k_i can make the amplitude drop exponentially. The imaginary part of Eq. (24) by considering a complex waveform of k is given by:

$$\begin{array}{c}\hfill k_{\mathrm{i}}=\frac{n^2{\omega }_{\mathrm{p}}^2{cos}^2{\theta }_B-{\omega }_{\mathrm{p}}^2}{2k{c}^2(1-g_{\mathrm{r}}{\omega }_{\mathrm{p}}^2{cos}^2{\theta }_B/{\omega }^2)}{g}_{\mathrm{i}},\end{array}$$(26)

where g_r and g_i denote the real and imaginary parts of Eq. (25). We consider the emission at a radius, r_e, along a field line whose footprint intersects the crust at magnetic colatitude Θ_e(R/r_e)^1/2 for Θ_e ≪ 1, one can obtain the angle between k and B at each point:

$$\begin{array}{c}\hfill {\theta }_B=\frac{3}{8}{\mathrm{\Theta }}_{\mathrm{e}}{\left(\frac{r}{r_{\mathrm{e}}}\right)}^{1/2}{(1-\frac{r_{\mathrm{e}}}{r})}^2.\end{array}$$(27)

Consequently, the optical depth for waves emitted at r_e to reach r is:

$$\begin{array}{c}\hfill {\tau }_{\mathrm{LD}}={\displaystyle {\int }_{r_{\mathrm{e}}}^r}2k_i\mathrm{d}r.\end{array}$$(28)

Even the distribution function for the magnetospheric electrons and positrons is not known, a relativistic thermal distribution is one choice. Since the vertical momentum drops to zero rapidly, for instance, a 1D Jüttner distribution is proposed here (e.g., Melrose et al. 1999; Rafat et al. 2019):

$$\begin{array}{c}\hfill f(u)=\frac{1}{2K_1({\rho }_{\mathrm{T}})}exp(-{\rho }_{\mathrm{T}}{\gamma }_{\mathrm{p}}),\end{array}$$(29)

where ρ_T = m_ec²/(k_BT), in which T is the temperature of plasma and k_B is the Boltzmann constant. Alternatively, another widely used distribution is a Gaussian distribution (e.g., Egorenkov et al. 1983; Asseo & Melikidze 1998):

$$\begin{array}{c}\hfill f(u)=\frac{1}{\sqrt{2\pi }{u}_{\mathrm{m}}}exp(-\frac{u^2}{2u_{\mathrm{m}}^2}),\end{array}$$(30)

where u_m can be interpreted as the average u over the distribution function. Besides a thermal distribution, power law distributions $f({\gamma }_{\mathrm{p}})=(p-1){\gamma }_{\mathrm{p}}^{-p}$, where p is the power law index (p > 2), are representative choice to describe the plasma in the magnetosphere (e.g., Kaplan & Tsytovich 1973).

The optical depths as a function of wave frequency with four different distributions are shown in Fig. 7. Here, we assume that r_e = 10⁷ cm, ℳ = 10³, γ_p = 10 and Θ_e = 0.05. For all considered distributions, optical depths decrease with frequencies which means that the Landau damping is more significant for lower-frequency subluminous O-mode photons. The optical depths are much larger than 1 for the FRB emission band so that the subluminous O-mode photon cannot escape from the magnetosphere. Superluminous O-mode photons which frequency higher than ${\omega }_{\mathrm{cut}}={\omega }_{\mathrm{p}}{\gamma }_{\mathrm{p}}^{-3/2}$ can propagate. The bunching curvature radiation is exponentially dropped at ω > ω_c. If the superluminous O-mode cut-off frequency is close to ω_c, the spectrum of the emission would be narrow. A simple condition could be given when a dipole field is assumed:

$$\begin{array}{c}\hfill r\simeq 13.8R{\mathcal{M}}_3{B}_{\mathrm{s},15}{P}_0^{-1}{\gamma }_2^{-6}{\gamma }_{\mathrm{p},1.3}^{-3},\end{array}$$(31)

where the emission region is considered at closed field lines in which ρ ∼ r. The cut-off frequency ω_cut ∝ r^−3/2 which evolves steeper than ω_c ∝ r⁻¹. The bandwidth is suggested to become wider as the burst frequency gets lower.

Coherent curvature radiation can propagate out of the magnetosphere via superluminous O-mode with a frequency cut-off, which likely leads to narrow spectra. However, the condition that the frequency cut-off is close to the characteristic frequency of curvature radiation could not be satisfied at all heights, leading to narrow spectra sometimes. The waves can also propagate via X-mode photon which has no frequency cut-off. The spectrum of X-mode photons ought to be broad. The waves propagate via whether X-mode or O-mode is hard to know because the real magnetic configuration is complex. Narrow spectra have been found in some FRBs, while some FRBs seem to have broad spectra and there are also some not sure due to the limitation of bandwidth of the receiver (Zhang et al. 2023).

In principle, the spectrum of a radio magnetar would be narrow, if the emission is generated from r > 10 R propagating via O-mode. However, radio pulses are rarely detected among magnetars. For normal pulsars, the lower magnetic fields may induce lower number density so that the Landau damping is not as strong as in magnetars.

4.3. Self-absorption by curvature radiation

The self-absorption is also associated with the emission. We assume that the momentum of the particle is still along the magnetic field line even after having absorbed the photon. The cross-section can be obtained through the Einstein coefficients when hν ≪ γ_sm_ec² (e.g., Locatelli & Ghisellini 2018):

$$\begin{array}{c}\hfill \frac{\mathrm{d}{\sigma }_{\mathrm{c}}}{\mathrm{d}\mathrm{\Omega }}=\frac{1}{2m_{\mathrm{e}}{\nu }^2}\frac{\partial }{\partial {\gamma }_{\mathrm{s}}}\left[j_{\mathrm{s}}(\nu ,{\gamma }_{\mathrm{s}},{\phi }_{\mathrm{s}})\right],\end{array}$$(32)

where σ_c is the curvature cross section and j_s is the emissivity.

There are special cases which the stimulated emission becomes larger than the true absorption. As a result, the total cross section becomes negative, and there is a possible maser mechanism. However, such maser action is ineffective when the magnetic energy density is much higher than that for the kinetic energy of charges (Zhelezniakov & Shaposhnikov 1979). The isotropic equivalent luminosity for the curvature radiation is ∼10⁴² erg s⁻¹ if N_e = 10²² at r = 10⁷ cm for γ = 10² (Kumar et al. 2017). The number density in a bunch is estimated as n_e ∼ 10¹¹ cm⁻³ by Wang et al. (2022a), so that n_eγm_ec² ≪ B²/(8π). Thus, the maser process is ignored in the following discussion.

By inserting the curvature power per unit frequency (see e.g., Jackson 1998), the total cross section integrated over angles is given by Ghisellini & Locatelli (2018):

$$\begin{array}{c}\hfill {\sigma }_{\mathrm{c}}\approx 2^{5/3}\sqrt{3}\mathrm{\Gamma }(5/3)\frac{\pi e^2}{{\gamma }_{\mathrm{s}}^3{m}_{\mathrm{e}}c\nu },\nu \ll {\nu }_{\mathrm{c}},\end{array}$$(33)

where Γ(ν) is the gamma function. As explained in the discussion in Sect. 2, the charges are approximately monoenergetic due to the rapid balance between E_∥ and the curvature damping. The absorption optical depth is:

$$\begin{array}{c}\hfill {\tau }_{\mathrm{CR}}\simeq {\sigma }_{\mathrm{c}}{n}_{\mathrm{e}}\frac{\rho }{{\gamma }_{\mathrm{s}}}·\end{array}$$(34)

The absorption is more significant for lower frequencies.

Figure 8 exhibits the comparison between purely bunching curvature radiation and it after the absorption. We adopt that γ_s = 10², n_e = 10¹¹ cm⁻³ and ρ = 10⁷ cm. The plotted spectra are normalized to the F_ν of curvature radiation with φ_t = 10⁻³ at ω = ω_c. The spectral index is 5/3 at lower frequencies due to the self-absorption, which is steeper than pure curvature radiation, but it is still wide-band emission. Alternatively, the E_∥ can let electrons and positrons decouple, which leads to a spectral index of 8/3 at lower frequencies (Yang et al. 2020). The total results of these independent processes can make a steeper spectrum, with an index of 13/3.

5. Bunch expansion

As the bunch moves to a higher altitude, the thickness of the bunch (the size parallel to the field lines) becomes larger due to the repulsive force of the same sign charge. If the bunch travels through the emission region faster than the LOS sweeps across the whole emission beam, the thickness can be estimated as cw, where w is the burst width (Wang et al. 2022b). For simplicity, we consider the motion of electrons at the boundary due to the Coulomb force from bunches in the co-moving rest frame. The acceleration of an electron at the boundary is estimated as:

$$\begin{array}{c}\hfill a_{\mathrm{b}}^{\prime }={\displaystyle \int \frac{e^2{n}_{\mathrm{e}}{l}^{\prime }}{{(R_{\mathrm{b}}^{\prime 2}+l^{\prime 2})}^{3/2}\gamma m_{\mathrm{e}}}}\mathrm{d}{V}^{\prime },\end{array}$$(35)

where $R_{\text{b}}^{\prime }=R_{\text{b}}$ is the transverse size of bunch. We note that the transverse size of a bunch is much larger than the longitude size even in the co-moving rest frame. Equation (35) is then calculated as $a_{\text{b}}^{\prime }\simeq \pi e^2{n}_{\text{e}}{s}^{\prime }/(\gamma m_{\text{e}})$. Charged particles in the bunch cannot move across the B-field due to the fast cyclotron cooling. The transverse size is essentially R_b ∝ r due to the magnetic freezing, thus, we have $a_{\text{b}}^{\prime }\propto r^{-2}$. The distance that the electron travels in the co-moving frame is s′∝r by assuming r ≫ r₀ (r₀ is an initial distance). As a result, the thickness of the bunch is cw ∝ r in the lab frame.

A burst width, w, would be larger if the emission is seen from a higher altitude region. If there is more than one bulk of bunches being observed coincidentally but with different trajectory planes, we would see the drifting pattern (Wang et al. 2019). The geometric time delay always lets a burst from a lower region be seen earlier, that is, a downward drifting pattern is observed. A subpulse structure with such a drifting pattern is a natural consequence caused by the geometric effect when the characteristic frequency of the emission is height-dependent.

Following the previous work (Wang et al. 2022a), we continuously investigate the drifting rate as a function of the burst central frequency and width. The drift pattern is a consequence of narrowband emission. Compared with our previous work, we add the fitting results of FRB 20180916B. Figures 9 and 10 show the best fitting results for FRB 20121102A (blue dashed-dotted lines), FRB 20180916B (red dashed lines) and all FRBs (black solid lines) and their corresponding 1σ regions. The best fitting of $\dot{\nu }-\nu$ relationship for FRB 20121102A, FRB 20180916B, and all FRBs are $(\dot{\nu }/\mathrm{MHz}{\mathrm{ms}}^{-1})=-{10}^{-5.30}\times {(\nu /\mathrm{MHz})}^{2.29}$, $(\dot{\nu }/\mathrm{MHz}{\mathrm{ms}}^{-1})=-{10}^{-4.52}\times {(\nu /\mathrm{MHz})}^{1.93}$, and $(\dot{\nu }/\mathrm{MHz}{\mathrm{ms}}^{-1})=-{10}^{-6.25}\times {(\nu /\mathrm{MHz})}^{2.52}$. The best fitting of $\dot{\nu }$ as a function of ν/width for FRB 20121102A, FRB 20180916B and all FRBs are $(\dot{\nu }/\mathrm{MHz}{\mathrm{ms}}^{-1})=-{10}^{-1.75}\times {(\nu /\mathrm{MHz})}^{1.10}{(w/\mathrm{ms})}^{-1.10}$, $(\dot{\nu }/\mathrm{MHz}{\mathrm{ms}}^{-1})=-{10}^{-1.11}\times {(\nu /\mathrm{MHz})}^{0.88}{(w/\mathrm{ms})}^{-0.88}$, and $(\dot{\nu }/\mathrm{MHz}{\mathrm{ms}}^{-1})=-{10}^{-0.99}\times {(\nu /\mathrm{MHz})}^{0.94}{(w/\mathrm{ms})}^{-0.94}$. The fitting results generally match the case of curvature radiation, which has $\dot{\nu }\propto {\nu }^2$. The relationship $\dot{\nu }\propto \nu w^{-1}$ is attributed to w ∝ r. Data are quoted from references as follows: Michilli et al. (2018), Gajjar et al. (2018), Hessels et al. (2019), Josephy et al. (2019), Caleb et al. (2020), Platts et al. (2021), Kumar et al. (2023), CHIME/FRB Collaboration (2019b, 2020b), Chamma et al. (2021), Chawla et al. (2020), Pastor-Marazuela et al. (2021), Sand et al. (2022), Day et al. (2020), Hilmarsson et al. (2021), Zhou et al. (2022), Gopinath et al. (2024).

6. Discussion

6.1. Curvature radiation

The curvature radiation is formed since the perpendicular acceleration on the curved trajectory and its spectrum is naturally wide-band. Spectra of bunching curvature radiation are more complex due to the geometric conditions and energy distribution (e.g., Yang & Zhang 2018; Wang et al. 2022a). Despite considering the complex electron energy distributions, it remains challenging for the spectrum to appear in the narrowband.

However, quasi-periodically distributed bunches are able to create such narrowband emission. The radiation from the bulk of bunches is quasi-periodically enhanced and the bandwidth depends on ω_M and the threshold of the telescope. As shown in Fig. 2, narrowband emissions with roughly the same bandwidth are expected to be seen from multi-band observations simultaneously. The energy is radiated into several multiples of 2πω_M with narrow bandwidths and the required number of particles in the bunch is smaller than that of pure curvature radiation.

6.2. Inverse Compton scattering

Inverse Compton scattering (ICS) occurs when low-frequency electromagnetic waves enter bunches of relativistic particles (Zhang 2022). The bunched particles are induced to oscillate at the same frequency in the comoving frame. The outgoing ICS frequency in the lab frame is calculated as:

$$\begin{array}{c}\hfill \omega \simeq {\omega }_{\mathrm{in}}{\gamma }^2(1-cos{\theta }_{\mathrm{in}}),\end{array}$$(36)

where ω_in is the angular frequency of the incoming wave and θ_in is the angle between the incident photon momentum and the electron momentum. This ICS process could produce a narrow spectrum for a single electron by giving a low-frequency wave of narrow spectrum (Zhang 2023). In order to generate GHz waves, dozens of kHz incoming wave is required for bunches with γ = 10² and θ_in should be larger than 0.5. This requires the low-frequency waves to be triggered from a location very far from the emission region. The X-mode of these waves could propagate freely in the magnetosphere, but O-mode would not be able to unless the propagation path of the waves approaches a vacuum-like environment.

As discussed in Sect. 2, the Lorentz factor is constant due to the balance between E_∥ and radiation damping. θ_in varies slightly as the LOS sweeps across different field lines. If the incoming wave is monochromatic, the outgoing ICS wave would have a narrow value range matching the observed narrow-band FRB emissions. The polarization of ICS from a single charge is strongly linearly polarized, however, the electric complex vectors of the scattered waves in the same direction are added for all the particles in a bunch. Therefore, circular polarization would be seen for off-beam geometry (Qu & Zhang 2023), similar to the case of bunching curvature radiation but they are essentially different mechanisms (Wang et al. 2022c). Doppler effects may be significant at off-beam geometry so that the burst duration and spectral bandwidth would get larger. The model predicts that a burst with higher circular polarization components tends to have wider spectral bandwidth. This may differ from the quasi-periodically distributed bunch case, in which the mechanism to create narrow spectra is independent of polarization.

The perturbation discussed in Sect. 2.2 can let charges oscillate at the same frequency in the comoving frame. This process is essentially the same as the ICS emission. The acceleration caused by the perturbation is proposed to be proportional to exp(−iω_at), which leads to the wave frequency of ω_a(1 − βcosθ_in) due to the Doppler effect in the lab frame. Taking the wave frequency into Eq. (13), the result is the same as Eq. (36). We let $F(t)={\displaystyle \int E}(t)e{m}_{\mathrm{e}}^{-1}\mathrm{d}t$ and ignore the radiation damping as a simple case. The dimensionless velocity is then given by $\beta =F(t)/\sqrt{1+F{(t)}^2}$. By assuming F(t) has a monochrome oscillating form, according to Eq. (A.6) if F(t)≪1, the result is reminiscent of the perturbation cases and ICS. If F(t)≫1, we have β ∼ 1 which implies a wide-band emission.

6.3. Interference processes

Some large-scale interference processes, for instance, scintillation, gravitational lensing, and plasma lensing, could change the radiation spectra via wave interference (Yang 2023). The observed spectra could be coherently enhanced at some frequencies, leading to the modulated narrow spectra. We consider that the multipath propagation effect of a given interference process (scintillation, gravitational lensing, and plasma lensing) has a delay time of δt, then the phase difference between the rays from the multiple images is δϕ ∼ 2πνδt, where ν is the wave frequency. Thus, for GHz waves, a significant spectral modulation requires a delay time of $\delta t~{10}^{-10}{\nu }_9^{-1}$ s, which could give some constraints on scintillation, gravitational lensing, or plasma lensing – if the observed narrowband is due to these processes.

For scintillation, according to Yang (2023), the plasma screen should be at a distance of ∼10¹⁵ cm from the FRB source, and the medium in the screen should be intermediately dense and turbulent. For plasma lensing, the lensing object needs to have an average electron number density of ∼10 cm⁻³ and a typical scale of ∼10⁻³ AU at a distance of ∼1 pc from the FRB source. The possibility of gravitational lensing could be ruled out because the gravitational lensing event is most likely a one-time occurrence; however, it cannot explain the remarkable burst-to-burst variation of the FRB spectra.

7. Summary

In this work, we discuss the spectral properties from the perspectives of radiation mechanisms and absorption effects in the magnetosphere. We have drawn the following conclusions:

1. A narrow spectrum may reflect the periodicity of the distribution of bunches, which have roughly the same Lorentz factor. The energy is radiated into multiple narrowbands of 2πω_M when the bulk of bunches is quasi-periodic distributed (Yang 2023). Such quasi-periodically distributed bunches may be produced due to quasi-monochromatic Langmuir waves or from a “gap”, which experiences quasi-periodically sparking. The model predicts that narrowband emissions may be seen in other higher or lower frequency bands.

2. The quasi-periodically distributed bunches and the perturbation scenarios share the same polarization properties with the uniformly distributed bunches’ curvature radiation (Wang et al. 2022a). The emission is generally dominated by linear polarization for the on-beam case, whereas it shows highly circular polarization for the off-beam case. The differential of V is the largest and sign change in the circular polarization appears at φ = 0. The PA across the burst envelope may show slight variations due to the existence of circular polarization.

3. Spectra of bunching curvature radiation with normal distribution are wide-band. We investigate the spectra from charged bunches with force perturbations on curved trajectories. If the perturbations are sinusoid with a frequency of ω_a, there would be a spike at 2γ²ω_a on the spectra. However, the perturbations are constrained by the monoenergy of bunches unless the scattering of the Lorentz factor can make the spike wider leading to a wide-band spectrum. Monochromatic Alfvén waves as the perturbations would not give rise to a narrowband spectrum from 100 MHz to 10 GHz.

4. Landau damping as an absorption mechanism may generate a narrowband spectrum. We consider three different forms of charge distribution in the magnetosphere. Subluminous O-mode photons are optically thick, but there is no Landau damping for X-mode and superluminous O-mode photons. The superluminous O-mode photons have a low-frequency cut-off of ${\omega }_{\mathrm{p}}{\gamma }_{\mathrm{p}}^{-3/2}$ which depends on the height in the magnetosphere. If so, the spectrum would be narrow when ${\omega }_{\mathrm{p}}{\gamma }_{\mathrm{p}}^{-3/2}$ is close to ω_c and the bandwidth is getting wider as the frequency becomes lower. The condition of ${\omega }_{\mathrm{p}}{\gamma }_{\mathrm{p}}^{-3/2}\approx {\omega }_{\mathrm{c}}$ can only be satisfied sometimes. For some FRBs, the spectrum of some FRBs can be broader than several hundred MHz. This may be due to the frequency cut-off of the superluminous O-mode being much lower than the characteristic frequency of curvature radiation or the waves propagating via X-mode.

5. Self-absorption by curvature radiation is significant at lower frequencies. By considering the Einstein coefficients, the cross-section is larger for lower frequencies. The spectral index is 5/3 steeper than that of pure curvature radiation at ω ≪ ω_c, however, the spectrum is not as narrow as Δν/ν = 0.1.

6. The two intrinsic mechanisms by invoking quasi-periodic bunch distribution and the perturbations are more likely to generate the observed narrow spectra of FRBs, rather than the absorption effects.

7. As a bulk of bunches moves to a higher region, the longitude size of the bulk becomes larger, suggesting that burst duration is larger for bulks emitting at higher altitudes. By investigating drifting rates as functions of the central frequency and ν/w, we find that drifting rates can be characterized by power law forms in terms of the central frequency and ν/w. The fitting results are generally consistent with curvature radiation.

---

Acknowledgments

We are grateful to the referee for constructive comment and V.S. Beskin, He Gao, Biping Gong, Yongfeng Huang, Kejia Lee, Ze-Nan Liu, Jiguang Lu, Rui Luo, Chen-Hui Niu, Jiarui Niu, Yuanhong Qu, Hao Tong, Shuangqiang Wang, CHengjun Xia, Zhenhui Zhang, Xiaoping Zheng, Enping Zhou, Dejiang Zhou, and Yuanchuan Zou for helpful discussion. W.Y.W. is grateful to Bing Zhang for the discussion of the ICS mechanism and his encouragement when I feel confused about the future. This work is supported by the National SKA Program of China (No. 2020SKA0120100) and the Strategic Priority Research Program of the CAS (No. XDB0550300). Y.P.Y. is supported by the National Natural Science Foundation of China (No.12003028) and the National SKA Program of China (2022SKA0130100). J.F.L. acknowledges support from the NSFC (Nos. 11988101 and 11933004) and from the New Cornerstone Science Foundation through the New Cornerstone Investigator Program and the XPLORER PRIZE.

References

Appendix A: Bunch radiation

In order to obtain the spectrum from bunches, we briefly summarize the physics of radiation from a single charge. As shown in Figure 1, the angle between the electron velocity direction and x-axis at t = 0 is defined as χ_ij; the one between the LOS and the trajectory plane is φ_k. The observation point is far enough away from the region of space where the acceleration occurs. Based on the classical electrodynamics, the radiation electric vector is given by (Jackson 1998):

$$\begin{array}{c}\hfill \mathit{E}(t)=\frac{e}{c}\left\{\frac{\mathit{n}\times [(\mathit{n}-{\mathit{\beta }}_s)\times {\dot{\mathit{\beta }}}_s]}{{(1-{\mathit{\beta }}_s·\mathit{n})}^3\mathcal{R}}\right\},\end{array}$$(A.1)

where the expression in the brackets is evaluated at the retarded time t_ret = t − ℛ/c. By inserting the two accelerations, ${\dot{\beta }}_{\Vert }$ and ${\dot{\beta }}_{\perp }$, the vector component of Equation (A.1) has (the subscript ‘i, j, k’ has been replaced by ‘s’):

$$\begin{array}{c}\hfill \begin{array}{c}\mathit{n}\times [(\mathit{n}-{\mathit{\beta }}_s)\times {\dot{\mathit{\beta }}}_s]\\ =-{\dot{\beta }}_{s\Vert }sin({\chi }_s+v_{s}t/\rho ){\mathit{ϵ}}_1+{\dot{\beta }}_{s\perp }[cos{\phi }_s-cos({\chi }_s+v_{s}t/\rho )]{\mathit{ϵ}}_1\\ +{\dot{\beta }}_{s\Vert }sin{\phi }_{s}cos({\chi }_s+v_{s}t/\rho ){ϵ}_2-{\dot{\beta }}_{s\perp }sin{\phi }_{s}sin({\chi }_s+v_{s}t/\rho ){\mathit{ϵ}}_2.\end{array}\end{array}$$(A.2)

Consider the small value of φ_s, χ_s and v_st/ρ, if ${\dot{\beta }}_{s\Vert }\simeq {\dot{\beta }}_{s\perp }$, the radiation power caused by ${\dot{\beta }}_{s\perp }$ is much smaller than that for ${\dot{\beta }}_{s\Vert }$.

The energy per unit frequency per unit solid angle of a single particle is given by:

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill \frac{{\mathrm{d}}^{2}W}{\mathrm{d}\omega \mathrm{d}\mathrm{\Omega }}& =\frac{c}{4{\pi }^2}{|\mathcal{R}{\mathit{E}}_s(t)exp(i\omega t)\mathrm{d}t|}^2\hfill \\ \hfill & =\frac{e^2}{4{\pi }^{2}c}|{\displaystyle {\int }_{-\infty }^{\infty }}\frac{\mathit{n}\times [(\mathit{n}-{\mathit{\beta }}_s)\times {\dot{\beta }}_s]}{{(1-{\mathit{\beta }}_s·\mathit{n})}^2}\hfill \\ \hfill & {\times exp[i\omega (t^{\prime }-\mathit{n}·{\mathit{r}}_s(t^{\prime })/c)]\mathrm{d}{t}^{\prime }|}^2.\hfill \end{array}\end{array}$$(A.3)

Using the identity:

$$\begin{array}{c}\hfill \frac{\mathit{n}\times [(\mathit{n}-{\mathit{\beta }}_s)\times {\dot{\mathit{\beta }}}_s]}{{(1-{\mathit{\beta }}_s·\mathit{n})}^2}=\frac{\mathrm{d}}{\mathrm{d}{t}^{\prime }}\left[\frac{\mathit{n}\times (\mathit{n}\times {\mathit{\beta }}_s)}{1-{\mathit{\beta }}_s·\mathit{n}}\right],\end{array}$$(A.4)

Equation (A.3) becomes:

$$\begin{array}{c}\hfill \frac{{\mathrm{d}}^{2}W}{\mathrm{d}\omega \mathrm{d}\mathrm{\Omega }}=\frac{e^2}{4{\pi }^{2}c}{|{\displaystyle {\int }_{-\infty }^{\infty }}exp[i\omega (t^{\prime }-\mathit{n}·{\mathit{r}}_s(t^{\prime })/c)]\mathrm{d}\left[\frac{\mathit{n}\times (\mathit{n}\times {\mathit{\beta }}_s)}{1-{\mathit{\beta }}_s·\mathit{n}}\right]|}^2.\end{array}$$(A.5)

It can be proven that by adding or subtracting the integrals over the times, while maintaining a constant velocity from Equation (A.5), one can use partial integration to obtain:

$$\begin{array}{c}\hfill \frac{{\mathrm{d}}^{2}W}{\mathrm{d}\omega \mathrm{d}\mathrm{\Omega }}=\frac{e^2{\omega }^2}{4{\pi }^{2}c}{|{\displaystyle {\int }_{-\infty }^{\infty }}\mathit{n}\times (\mathit{n}\times {\mathit{\beta }}_s)exp[i\omega (t^{\prime }-\mathit{n}·{\mathit{r}}_s(t^{\prime })/c)]\mathrm{d}{t}^{\prime }|}^2,\end{array}$$(A.6)

as long as a convergence factor (e.g., exp(−ϵ|t|), where taking the limit ϵ → 0) is inserted to the integrand function of Equation (A.5) to eliminate instability at t = ±∞.

The energy radiated per unit solid angle per unit frequency interval in terms of the two polarization components can be written as:

$$\begin{array}{c}\hfill \frac{{\mathrm{d}}^{2}W}{\mathrm{d}\omega \mathrm{d}\mathrm{\Omega }}=\frac{e^2{\omega }^2}{4{\pi }^{2}c}{|-{\mathit{ϵ}}_{\Vert }{A}_{\Vert }+{\mathit{ϵ}}_{\perp }{A}_{\perp }|}^2.\end{array}$$(A.7)

As shown in Figure 1, we consider a charge of identifier of ijk. The vector part of Equation (A.6) can be written as:

$$\begin{array}{c}\hfill \begin{array}{cc}& \mathit{n}\times (\mathit{n}\times {\mathit{\beta }}_{ijk,\Vert })\hfill \\ \hfill & ={\beta }_{\Vert }[-{\mathit{ϵ}}_{\Vert }sin(\frac{\mathit{vt}}{\rho }+{\chi }_{\mathit{ij}})+{\mathit{ϵ}}_{\perp }cos(\frac{\mathit{vt}}{\rho }+{\chi }_{\mathit{ij}})sin{\phi }_k],\hfill \end{array}\end{array}$$(A.8)

and the argument of the exponential is:

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill \omega (t-\frac{\mathit{n}·{\mathit{r}}_{\mathit{ijk}}(t)}{c})& =\omega [t-\frac{2\rho }{c}sin\left(\frac{\mathit{vt}}{2\rho }\right)cos(\frac{\mathit{vt}}{2\rho }+{\chi }_{\mathit{ij}})cos{\phi }_k]\hfill \\ \hfill & \simeq \frac{\omega }{2}[(\frac{1}{{\gamma }^2}+{\phi }_k^2+{\chi }_{\mathit{ij}}^2)t+\frac{c^2{t}^3}{3{\rho }^2}+\frac{c{t}^2}{\rho }{\chi }_{\mathit{ij}}],\hfill \end{array}\end{array}$$(A.9)

where t is adopted to replace t′ in the following calculations. Therefore, the two amplitudes for Equation (A.7) are given by:

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill A_{ijk,\Vert }& \simeq {\displaystyle {\int }_{-\infty }^{\infty }}(\frac{\mathit{ct}}{\rho }+{\chi }_{\mathit{ij}})\hfill \\ \hfill & \times exp\left(i\frac{\omega }{2}[(\frac{1}{{\gamma }^2}+{\phi }_k^2+{\chi }_{\mathit{ij}}^2)t+\frac{c^2{t}^3}{3{\rho }^2}+\frac{c{t}^2}{\rho }{\chi }_{\mathit{ij}}]\right)\mathrm{d}t,\hfill \\ \hfill A_{ijk,\perp }& \simeq {\phi }_k{\displaystyle {\int }_{-\infty }^{\infty }}exp\left(i\frac{\omega }{2}[(\frac{1}{{\gamma }^2}+{\phi }_k^2+{\chi }_{\mathit{ij}}^2)t+\frac{c^2{t}^3}{3{\rho }^2}+\frac{c{t}^2}{\rho }{\chi }_{\mathit{ij}}]\right)\mathrm{d}t.\hfill \end{array}\end{array}$$(A.10)

Making the changes of variables:

$$\begin{array}{c}\hfill \begin{array}{cc}& u=\frac{\mathit{ct}}{\rho }{(\frac{1}{{\gamma }^2}+{\phi }_k^2+{\chi }_{\mathit{ij}}^2)}^{-1/2},\hfill \\ \hfill & \xi =\frac{\omega \rho }{3c}{(\frac{1}{{\gamma }^2}+{\phi }_k^2+{\chi }_{\mathit{ij}}^2)}^{3/2},\hfill \end{array}\end{array}$$(A.11)

Equation (A.10) becomes:

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill A_{ijk,\Vert }& \simeq \frac{\rho }{c}(\frac{1}{{\gamma }^2}+{\phi }_k^2+{\chi }_{\mathit{ij}}^2)\hfill \\ \hfill & \times {\displaystyle {\int }_{-\infty }^{\infty }}(u+\frac{{\chi }_{\mathit{ij}}}{\sqrt{1/{\gamma }^2+{\chi }_{\mathit{ij}}^2+{\phi }_k^2}})\hfill \\ \hfill & \times exp\left[i\frac{3}{2}\xi (u+\frac{1}{3}{u}^3+\frac{{\chi }_{\mathit{ij}}}{\sqrt{1/{\gamma }^2+{\phi }_k^2+{\chi }_{\mathit{ij}}^2}}{u}^2)\right]\mathrm{d}u,\hfill \\ \hfill A_{ijk,\perp }& \simeq \frac{\rho }{c}{\phi }_k{(\frac{1}{{\gamma }^2}+{\phi }_k^2+{\chi }_{\mathit{ij}}^2)}^{1/2}\hfill \\ \hfill & \times {\displaystyle {\int }_{-\infty }^{\infty }}exp\left[i\frac{3}{2}\xi (u+\frac{1}{3}{u}^3+\frac{{\chi }_{\mathit{ij}}}{\sqrt{1/{\gamma }^2+{\phi }_k^2+{\chi }_{\mathit{ij}}^2}}{u}^2)\right]\mathrm{d}u.\hfill \\ \hfill \end{array}\end{array}$$(A.12)

We note that:

$$\begin{array}{c}\hfill \begin{array}{cc}& \underset{u\to +0}{lim}(u+\frac{1}{3}{u}^3+\frac{{\chi }_{\mathit{ij}}}{\sqrt{1/{\gamma }^2+{\phi }_k^2+{\chi }_{\mathit{ij}}^2}}{u}^2)=u,\hfill \\ \hfill & \underset{u\to +\infty }{lim}(u+\frac{1}{3}{u}^3+\frac{{\chi }_{\mathit{ij}}}{\sqrt{1/{\gamma }^2+{\phi }_k^2+{\chi }_{\mathit{ij}}^2}}{u}^2)=\frac{1}{3}{u}^3.\hfill \end{array}\end{array}$$(A.13)

The integrals in Equation (A.12) are identifiable as Airy integrals. Consequently, the two amplitudes are:

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill A_{ijk,\Vert }& \simeq \frac{i2}{\sqrt{3}}\frac{\rho }{c}(\frac{1}{{\gamma }^2}+{\phi }_k^2+{\chi }_{\mathit{ij}}^2)K_{\frac{2}{3}}(\xi )\hfill \\ \hfill & +\frac{2}{\sqrt{3}}\frac{\rho }{c}{\chi }_{\mathit{ij}}{(\frac{1}{{\gamma }^2}+{\phi }_k^2+{\chi }_{\mathit{ij}}^2)}^{1/2}{K}_{\frac{1}{3}}(\xi ),\hfill \\ \hfill A_{ijk,\perp }& \simeq \frac{2}{\sqrt{3}}\frac{\rho }{c}{\phi }_k{(\frac{1}{{\gamma }^2}+{\phi }_k^2+{\chi }_{\mathit{ij}}^2)}^{1/2}{K}_{\frac{1}{3}}(\xi ),\hfill \end{array}\end{array}$$(A.14)

where K_ν(ξ) is the modified Bessel function, in which K_ν(ξ)→(Γ(ν)/2)(ξ/2)^−ν for ξ ≪ 1 and ν ≠ 0, and $K_{\nu }(\xi )\to \sqrt{\pi /(2\xi )}exp(-\xi )$ for ξ ≫ 1, leading to a wide-band spectrum.

The curvature radius can be regarded as a constant since the bunch’s size is much smaller than its altitude. We assume that each bunch has roughly the same χ_u and χ_d even with different k so that an average height can be derived to replace the complex boundary condition of the bulk (Wang et al. 2022a). Consider that the bunch longitude size is l < 4ρ/(3γ³). The emission with ω < ω_c is added coherently. For the emission band with ω > ω_l, the emission is added incoherently and the fluxes drop with ω very quickly. The spectrum becomes much deeper than an exponential function. The total amplitudes of a bunch are then given by:

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill A_{\Vert }& \simeq \frac{2}{\sqrt{3}}\frac{\rho }{c}\frac{N_j}{\mathrm{\Delta }\chi }\frac{N_k}{2{\phi }_t}{N}_i\hfill \\ \hfill & \times {\displaystyle {\int }_{{\chi }_{d,i}}^{{\chi }_{u,i}}}d\chi \prime {\displaystyle {\int }_{{\phi }_d}^{{\phi }_u}}[i(\frac{1}{{\gamma }^2}+{\phi }^{\prime 2}+{\chi }^{\prime 2})K_{\frac{2}{3}}(\xi )\hfill \\ \hfill & +\chi \prime {(\frac{1}{{\gamma }^2}+{\phi }^{\prime 2}+{\chi }^{\prime 2})}^{1/2}{K}_{\frac{1}{3}}(\xi )]cos{\phi }^{\prime }\mathrm{d}{\phi }^{\prime },\hfill \\ \hfill A_{\perp }& \simeq \frac{2}{\sqrt{3}}\frac{\rho }{c}\frac{N_j}{\mathrm{\Delta }\chi }\frac{N_k}{2{\phi }_t}{N}_i\hfill \\ \hfill & \times {\displaystyle {\int }_{{\chi }_{d,i}}^{{\chi }_{u,i}}}d\chi \prime {\displaystyle {\int }_{{\phi }_d}^{{\phi }_u}}{(\frac{1}{{\gamma }^2}+{\phi }^{\prime 2}+{\chi }^{\prime 2})}^{1/2}{K}_{\frac{1}{3}}(\xi ){\phi }^{\prime }cos{\phi }^{\prime }\mathrm{d}{\phi }^{\prime },\hfill \end{array}\end{array}$$(A.15)

where χ_u is the upper boundary χ′, χ_d is the lower boundary χ′, and the boundaries of φ′ read φ_u = φ_t + φ and φ_d = −φ_t + φ.

With the particular interest of β_1, ∥(t) = β₁exp(−iω_at), the total amplitudes of a bunch for the parallel perturbation are given by:

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill {\stackrel{\sim }{A}}_{\Vert }& \simeq \frac{2}{\sqrt{3}}\frac{\rho }{c}\frac{N_j}{\mathrm{\Delta }\chi }\frac{N_k}{2{\phi }_t}{N}_i{\beta }_1\hfill \\ \hfill & \times {\displaystyle {\int }_{{\chi }_{d,i}}^{{\chi }_{u,i}}}d\chi \prime {\displaystyle {\int }_{{\phi }_d}^{{\phi }_u}}[i(\frac{1}{{\gamma }^2}+{\phi }^{\prime 2}+{\chi }^{\prime 2}-\frac{2{\omega }_a}{\omega })K_{\frac{2}{3}}({\xi }_a)\hfill \\ \hfill & +{\chi }^{\prime }{(\frac{1}{{\gamma }^2}+{\phi }^{\prime 2}+{\chi }^{\prime 2}-\frac{2{\omega }_a}{\omega })}^{1/2}{K}_{\frac{1}{3}}({\xi }_a)]cos{\phi }^{\prime }\mathrm{d}{\phi }^{\prime },\hfill \\ \hfill {\stackrel{\sim }{A}}_{\perp }& \simeq \frac{2}{\sqrt{3}}\frac{\rho }{c}\frac{N_j}{\mathrm{\Delta }\chi }\frac{N_k}{2{\phi }_t}{N}_i{\beta }_1\hfill \\ \hfill & \times {\displaystyle {\int }_{{\chi }_{d,i}}^{{\chi }_{u,i}}}d\chi \prime {\displaystyle {\int }_{{\phi }_d}^{{\phi }_u}}{(\frac{1}{{\gamma }^2}+{\phi }^{\prime 2}+{\chi }^{\prime 2}-\frac{2{\omega }_a}{\omega })}^{1/2}\hfill \\ \hfill & \times K_{\frac{1}{3}}({\xi }_a){\phi }^{\prime }cos{\phi }^{\prime }\mathrm{d}{\phi }^{\prime }.\hfill \end{array}\end{array}$$(A.16)

The monochromatic perturbation frequency is Doppler boosted to 2γ²ω_a, so that there is a spike in the spectrum at ω = 2γ²ω_a, and the width of the spike is determined by the boundaries of χ′ and φ′.

We consider the case with a similar formula but for the perpendicular perturbation, namely, β_1, ⊥(t) = β₁exp(−iω_at). Following the discussion in Section 2.2.4, we insert Equation (12) into Equation (15) and the amplitudes are calculated as:

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill {\stackrel{\sim }{A}}_{ijk,\Vert }& \simeq \frac{1}{\sqrt{3}}\frac{\rho }{c}{\beta }_1({\phi }_k^2-{\chi }_{\mathit{ij}}^2){(\frac{1}{{\gamma }^2}+{\phi }^{\prime 2}+{\chi }^{\prime 2}-\frac{2{\omega }_a}{\omega })}^{1/2}{K}_{\frac{1}{3}}({\xi }_a)\hfill \\ \hfill & -\frac{i2}{\sqrt{3}}\frac{\rho }{c}{\beta }_1{\chi }_{\mathit{ij}}(\frac{1}{{\gamma }^2}+{\phi }^{\prime 2}+{\chi }^{\prime 2}-\frac{2{\omega }_a}{\omega })K_{\frac{2}{3}}({\xi }_a)\hfill \\ \hfill & +\frac{\rho }{c}{\beta }_1{(\frac{1}{{\gamma }^2}+{\phi }^{\prime 2}+{\chi }^{\prime 2}-\frac{2{\omega }_a}{\omega })}^{3/2}\hfill \\ \hfill & \times [\frac{1}{\sqrt{3}}{K}_{\frac{1}{3}}({\xi }_a)-\pi \delta (3{\xi }_a/2)],\hfill \\ \hfill {\stackrel{\sim }{A}}_{ijk,\perp }& \simeq -\frac{i2}{\sqrt{3}}\frac{\rho }{c}{\phi }_k{\beta }_1(\frac{1}{{\gamma }^2}+{\phi }_k^2+{\chi }_{\mathit{ij}}^2)K_{\frac{2}{3}}({\xi }_a)\hfill \\ \hfill & -\frac{2}{\sqrt{3}}\frac{\rho }{c}{\phi }_k{\beta }_1{\chi }_{\mathit{ij}}{(\frac{1}{{\gamma }^2}+{\phi }_k^2+{\chi }_{\mathit{ij}}^2)}^{1/2}{K}_{\frac{1}{3}}({\xi }_a),\hfill \\ \hfill \end{array}\end{array}$$(A.17)

where δ(x) is the Dirac function. The total amplitudes of a bunch for the perpendicular perturbation case can be given by the integrals of Equation (A.17).

Appendix B: Stokes parameters

Four Stokes parameters are significant tools to describe the polarization properties of a quasi-monochromatic electromagnetic wave. In this section, we consider the monochrome Stokes parameters which can be derived by the Fourier-transformed electric vector. The Fourier-transformed polarized components of electric waves can be written as:

$$\begin{array}{c}\hfill E_{\Vert }(\omega )=\frac{e\omega A_{\Vert }}{2\pi c\mathcal{R}},E_{\perp }(\omega )=\frac{e\omega A_{\perp }}{2\pi c\mathcal{R}}.\end{array}$$(B.1)

By inserting Equation (B.1), four Stokes parameters of observed waves can be calculated as:

$$\begin{array}{c}\hfill \begin{array}{cc}& I=\mu (A_{\Vert }{A}_{\Vert }^{\ast }+A_{\perp }{A}_{\perp }^{\ast }),\hfill \\ \hfill & Q=\mu (A_{\Vert }{A}_{\Vert }^{\ast }-A_{\perp }{A}_{\perp }^{\ast }),\hfill \\ \hfill & U=\mu (A_{\Vert }{A}_{\perp }^{\ast }+A_{\perp }{A}_{\Vert }^{\ast }),\hfill \\ \hfill & V=-i\mu (A_{\Vert }{A}_{\perp }^{\ast }-A_{\perp }{A}_{\Vert }^{\ast }),\hfill \end{array}\end{array}$$(B.2)

where μ = ω²e²/(4π²ℛ²c𝒯) is the proportionality factor, I defines the total intensity, Q and U define the linear polarization and its position angle, and V describes circular polarization. The timescale 𝒯 is chosen as the average timescale the flux density I repeated. The corresponding linear polarization and PA are then given by:

$$\begin{array}{c}\hfill \begin{array}{cc}& L=\sqrt{U^2+Q^2},\hfill \\ \hfill & \mathrm{\Psi }=\frac{1}{2}{tan}^{-1}\left(\frac{U_s}{Q_s}\right),\hfill \end{array}\end{array}$$(B.3)

where

$$\begin{array}{c}\hfill \left(\begin{array}{c}{U}_s\hfill \\ Q_s\hfill \end{array}\right)=\left(\begin{array}{cc}cos2\psi & sin2\psi \\ -sin2\psi & cos2\psi \end{array}\right)\left(\begin{array}{c}U\hfill \\ Q\hfill \end{array}\right),\end{array}$$(B.4)

where ψ is given by the rotation vector model (RVM; Radhakrishnan & Cooke 1969).

Appendix C: Quake-induced oscillation

We calculated the frequencies of the torsional modes for the relativistic star caused by a quake and discuss the effect of the magnetic field on the frequencies of the various torsional modes. Oscillations after the quake have two types: spherical and toroidal, in which only the toroidal oscillations can effectively modify the spin of the star and the Goldreich-Julian density. The axial perturbation equations for the elastic solid star in the Cowling approximation is written as (Samuelsson & Andersson 2007; Sotani et al. 2012):

$$\begin{array}{c}\hfill \begin{array}{cc}& Y[\frac{{\rho }_m+P}{\mu }{\mathrm{\Omega }}^{2}exp(-2\mathrm{\Phi })-\frac{(\ell +2)(\ell -1)}{r^2}]exp(2\mathrm{\Lambda })+Y^{″}\hfill \\ \hfill & +(\frac{4}{r}+{\mathrm{\Phi }}^{\prime }-{\mathrm{\Lambda }}^{\prime }+\frac{{\mu }^{\prime }}{\mu })Y^{\prime }=0,\hfill \end{array}\end{array}$$(C.1)

where μ is the shear modulus, Φ and Λ are functions of r, ρ_m is the mass-energy density, Y(r) describes the radial part of the angular oscillation amplitude, and the integer ℓ is the angular separation constant which enters when Y(r) is expanded in spherical harmonics Y_ℓm(θ, ϕ).

To discuss the possible effects of the magnetic field, we extended our study to the torsional oscillation of a magnetized relativistic star (e.g., Li et al. 2024). Sotani et al. (2007) derived the perturbation equations of the magnetized relativistic star using the relativistic Cowling approximation. The final perturbation equation is:

$$\begin{array}{c}\hfill A_{\ell }(r)Y^{″}+B_{\ell }(r)Y^{\prime }+C_{\ell }(r)Y=0,\end{array}$$(C.2)

where the coefficients are given in terms of the functions describing the equilibrium metric, fluid, and the magnetic field of the star:

$$\begin{array}{cc}\hfill A_{\ell }(r)& =\mu +(1+2{\lambda }_1)\frac{{a_1}^2}{\pi r^4},\hfill \end{array}$$(C.3)

$$\begin{array}{cc}\hfill B_{\ell }(r)& =(\frac{4}{r}+{\mathrm{\Phi }}^{\prime }-{\mathrm{\Lambda }}^{\prime })\mu +{\mu }^{\prime }\hfill \\ \hfill & +(1+2{\lambda }_1)\frac{a_1}{\pi r^4}[({\mathrm{\Phi }}^{\prime }-{\mathrm{\Lambda }}^{\prime })a_1+2a_1^{\prime }],\hfill \end{array}$$(C.4)

$$\begin{array}{cc}\hfill C_{\ell }(r)& =(2+5{\lambda }_1)\frac{a_1}{2\pi r^4}\{({\mathrm{\Phi }}^{\prime }-{\mathrm{\Lambda }}^{\prime })a_1^{\prime }+a_1^{″}\}\hfill \\ \hfill & -(\lambda -2)(\frac{\mu exp(2\mathrm{\Lambda })}{r^2}-\frac{{\lambda }_1{a_1^{\prime }}^2}{2\pi r^4})\hfill \\ \hfill & +[({\rho }_m+P+(1+2{\lambda }_1)\frac{{a_1}^2}{\pi r^4})exp(2\mathrm{\Lambda })-\frac{{\lambda }_1{a_1^{\prime }}^2}{2\pi r^2}]\hfill \\ \hfill & \times {\mathrm{\Omega }}^{2}exp(-2\mathrm{\Phi }),\hfill \end{array}$$(C.5)

where λ = ℓ(ℓ + 1), and λ₁ = −ℓ(ℓ+1)/[(2ℓ−1)(2ℓ+3)]. To solve Eqs. (C.1) and (C.2) and determine the oscillation frequencies, the boundary conditions require that the traction vanishes at the top and the bottom of the crust.

The toroidal oscillations can effectively modify the spin and the electric field in the gap. The velocity components of toroidal oscillations can be written as:

$$\begin{array}{c}\hfill \delta v^{\widehat{i}}=[0,\frac{1}{sin\theta }{\partial }_{\varphi }{Y}_{\ell m}(\theta ,\varphi ),-{\partial }_{\theta }{Y}_{\ell m}(\theta ,\varphi )]\stackrel{\sim }{\eta }(r)e^{-i\mathrm{\Omega }t},\end{array}$$(C.6)

where $\stackrel{\sim }{\eta }(r)$ is a parameter denoting the amplitude of the oscillations. Owing to unipolar induction, we have (Lin et al. 2015):

$$\begin{array}{c}\hfill E_{\theta }=-\frac{1}{\mathit{Nc}}[\frac{2\pi r}{P}sin\theta (1-\frac{R^3\kappa }{r^3})-{\partial }_{\theta }{Y}_{\ell m}\stackrel{\sim }{\eta }(r)]B_{\mathrm{s}}\frac{f(r)}{f(R)}\frac{R^{3}cos\theta }{r^3},\end{array}$$(C.7)

where $N=\sqrt{1-2GM/(r{c}^2)}$, κ = 2GI_m/(R³c²), I_m is the moment of inertia, and f(r) is:

$$\begin{array}{c}\hfill f(r)=-3{\left(\frac{r{c}^2}{2GM}\right)}^3[ln(1-\frac{2GM}{r{c}^2})+\frac{2GM}{r{c}^2}(1+\frac{\mathit{GM}}{r{c}^2})].\end{array}$$(C.8)

Appendix D: Dispersion relationship of ultra-relativistic plasma

In this section, we briefly summarize the physics of waves in ultra-relativistic plasma. Basically, the equation of motion of a given charge species is expressed as:

$$\begin{array}{c}\hfill \frac{\mathrm{d}}{\mathrm{d}t}\left({\gamma }_s{m}_{\mathrm{e}}{\mathit{v}}_s\right)=q_s\mathit{E}+\frac{q_s}{c}{\mathit{v}}_s\times \mathit{B},\end{array}$$(D.1)

where v_s = v_s, 0 + δv_s, E = δE, B = B₀ + δB, in which δv_s, δE, δB are associated with the disturbance in the plasma and have a form of exp[ − i(ωt − k ⋅ r)]. Combined with

$$\begin{array}{c}\hfill \mathrm{d}\delta {\mathit{v}}_s/\mathrm{d}t=-i\omega \delta {\mathit{v}}_s+ic(\mathit{k}·{\mathit{v}}_{s,0})\delta {\mathit{v}}_s,\end{array}$$(D.2)

and

$$\begin{array}{c}\hfill \delta \mathit{B}=(c\mathit{k}/\omega )\times \delta \mathit{E},\end{array}$$(D.3)

δv_s can be then described in terms of δE. The continuity equation is expressed as:

$$\begin{array}{c}\hfill \frac{\partial n_s}{\partial t}+\mathrm{\nabla }·(n_s{\mathit{v}}_{s,0})=0,\end{array}$$(D.4)

so that the density perturbation is given by:

$$\begin{array}{c}\hfill \delta n_s=\frac{n_s(\mathit{k}·\delta {\mathit{v}}_s)}{\omega -\mathit{k}·{\mathit{v}}_{s,0}}.\end{array}$$(D.5)

The conductivity tensor σ is defined by:

$$\begin{array}{c}\hfill \delta \mathit{J}=\mathit{\sigma }·\delta \mathit{E}={\mathrm{\Sigma }}_s(n_s{q}_s\delta {\mathit{v}}_s+\delta n_s{q}_s{\mathit{v}}_{s,0}),\end{array}$$(D.6)

and the dielectric tensor is given by:

$$\begin{array}{c}\hfill \mathit{ϵ}=\mathit{I}+i\frac{4\pi }{\omega }\mathit{\sigma },\end{array}$$(D.7)

where I is the unit tensor. Since the electric displacement vector is defined as D = ϵ ⋅ E, we can obtain the equation for plane waves in a coordinate system (x′, y′, z′) where B₀ goes along z′ as follows:

$$\begin{array}{c}\hfill ϵ=\left(\begin{array}{ccc}S& \mathit{iD}& A\\ -iD& S& -iC\\ A& \mathit{iC}& P\end{array}\right),\end{array}$$(D.8)

with

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill S& =1+{\displaystyle \sum _s}{f}_{s,11},\hfill \\ \hfill D& ={\displaystyle \sum _s}{f}_{s,12},\hfill \\ \hfill A& ={\displaystyle \sum _s}{\zeta }_s{f}_{s,11},\hfill \\ \hfill C& ={\displaystyle \sum _s}{\zeta }_s{f}_{s,12},\hfill \\ \hfill P& =1+{\displaystyle \sum _s}(f_{s,\eta }+{\zeta }_s^2{f}_{s,11}),\hfill \end{array}\end{array}$$(D.9)

where

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill f_{s,11}& =-\frac{{\omega }_{\mathrm{p}}^2{\gamma }_s^{-1}}{{\omega }^2-{\omega }_B^2{\gamma }_s^{-2}{(1-n{\beta }_{s}cos{\theta }_B)}^{-2}},\hfill \\ \hfill f_{s,12}& =-\frac{\text{sign}\left(q_s\right){\omega }_B{\omega }_{\mathrm{p}}^2{\omega }^{-3}{\gamma }_s^{-2}{(1-n{\beta }_{s}cos{\theta }_B)}^{-1}}{1-{\omega }_B^2{\omega }^{-2}{\gamma }_s^{-2}{(1-n{\beta }_{s}cos{\theta }_B)}^{-2}},\hfill \\ \hfill f_{s,\eta }& =-\frac{{\omega }_{\mathrm{p}}^2}{{\omega }^2{\gamma }_s^3{(1-n{\beta }_{s}cos{\theta }_B)}^2},\hfill \\ \hfill {\zeta }_s& =\frac{n{\beta }_{s}sin{\theta }_B}{1-n{\beta }_{s}cos{\theta }_B}.\hfill \end{array}\end{array}$$(D.10)

Here, the subscript ‘0’ has been suppressed. In the coordinate system x₀y₀z₀ with k along the z₀-axis and B in the x₀ − z₀ plane, the components of dielectric tensor are given by:

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill {ϵ}_{\mathit{xx}}& =S{cos}^2{\theta }_B-2Asin{\theta }_{B}cos{\theta }_B+P{sin}^2{\theta }_B,\hfill \\ \hfill {ϵ}_{\mathit{yy}}& =S,\hfill \\ \hfill {ϵ}_{\mathit{zz}}& =S{sin}^2{\theta }_B+2Asin{\theta }_{B}cos{\theta }_B+P{cos}^2{\theta }_B,\hfill \\ \hfill {ϵ}_{\mathit{xy}}& =-{ϵ}_{\mathit{yx}}=i(Dcos{\theta }_B-Csin{\theta }_B),\hfill \\ \hfill {ϵ}_{\mathit{xz}}& ={ϵ}_{\mathit{zx}}=Acos2{\theta }_B+(S-P)sin{\theta }_{B}cos{\theta }_B,\hfill \\ \hfill {ϵ}_{\mathit{yz}}& =-{ϵ}_{\mathit{zy}}=-i(Dsin{\theta }_B+Ccos{\theta }_B).\hfill \end{array}\end{array}$$(D.11)

The three electric components have:

$$\begin{array}{c}\hfill E_z=-{ϵ}_{\mathit{zz}}^{-1}({ϵ}_{\mathit{zx}}{E}_x+{ϵ}_{\mathit{zy}}{E}_y).\end{array}$$(D.12)

Reinserting this back into the definition of electric displacement vector yields:

$$\begin{array}{c}\hfill \left(\begin{array}{cc}{\eta }_{\mathit{xx}}-n^2& {\eta }_{\mathit{xy}}\\ {\eta }_{\mathit{yx}}& {\eta }_{\mathit{yy}}-n^2\end{array}\right)\left(\begin{array}{c}{E}_x\\ E_y\end{array}\right)=0,\end{array}$$(D.13)

where

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill {\eta }_{\mathit{xx}}& =\frac{1}{{ϵ}_{\mathit{zz}}}({ϵ}_{\mathit{zz}}{ϵ}_{\mathit{xx}}-{ϵ}_{\mathit{xz}}{ϵ}_{\mathit{zx}})=\frac{1}{{ϵ}_{\mathit{zz}}}(SP-A^2),\hfill \\ \hfill {\eta }_{\mathit{yy}}& =\frac{1}{{ϵ}_{\mathit{zz}}}({ϵ}_{\mathit{zz}}{ϵ}_{\mathit{yy}}-{ϵ}_{\mathit{yz}}{ϵ}_{\mathit{zy}}),\hfill \\ \hfill & =\frac{1}{{ϵ}_{\mathit{zz}}}\left[(S^2-D^2-SP+C^2){sin}^2{\theta }_B\right]\hfill \\ \hfill & +\frac{1}{{ϵ}_{\mathit{zz}}}[2(AS-CD)sin{\theta }_{B}cos{\theta }_B+SP-C^2]\hfill \\ \hfill {\eta }_{\mathit{yx}}& =-{\eta }_{\mathit{xy}}=\frac{1}{{ϵ}_{\mathit{zz}}}({ϵ}_{\mathit{zz}}{ϵ}_{\mathit{yx}}-{ϵ}_{\mathit{yz}}{ϵ}_{\mathit{zx}}),\hfill \\ \hfill & =\frac{-i}{{ϵ}_{\mathit{zz}}}[PDcos{\theta }_B-SCsin{\theta }_B+A(Dsin{\theta }_B-Ccos{\theta }_B)].\hfill \end{array}\end{array}$$(D.14)

For a normal pulsar or magnetar, one can obtain ω_B ≫ ω and ω_B ≫ ω_p inside the magnetosphere. Under these conditions, some parameters can be reduced to:

$$\begin{array}{c}\hfill \begin{array}{cc}& f_{11}\sim {\omega }_{\mathrm{p},s}^2{\gamma }_s/{\omega }_B^2\sim 0,f_{12}\sim 0,\hfill \\ \hfill & f_{\eta }\simeq -{\omega }_{\mathrm{p},s}^2{\omega }^{-2}{\gamma }_s^{-3}{(1-n{\beta }_{s}cos{\theta }_B)}^{-2},\hfill \\ \hfill & S\simeq 1,P\simeq 1+f_{\eta },D\simeq 0,A\simeq 0,\hfill \\ \hfill & C\simeq 0,{ϵ}_{\mathit{zz}}\simeq 1+f_{\eta }{cos}^2{\theta }_B,\hfill \\ \hfill & {\eta }_{\mathit{xx}}\simeq \frac{1+f_{\eta }}{1+f_{\eta }{cos}^2{\theta }_B},{\eta }_{\mathit{yy}}\simeq 1,{\eta }_{\mathit{yx}}\simeq -{\eta }_{\mathit{xy}}\simeq 0.\hfill \end{array}\end{array}$$(D.15)

Reinserting this back into Equation (D.13) and the solution of refractive index is given by (Arons & Barnard 1986):

$$\begin{array}{c}\hfill \begin{array}{cc}& n^2=1,\hfill \\ \hfill & n^2={\eta }_{\mathit{xx}},\hfill \\ \hfill & ({\omega }^2-c^2{k}_{\Vert }^2)[1-\frac{{\omega }_{\mathrm{p},s}^2}{{\gamma }_s^3{\omega }^2(1-{\beta }_{s}c{k}_{\Vert }/\omega )}]-c^2{k}_{\perp }^2=0,\hfill \end{array}\end{array}$$(D.16)

where k_∥ = kcosθ_B and k_⊥ = ksinθ_B. For the mode with n² = 1, the polarization is given by:

$$\begin{array}{c}\hfill \left|\frac{E_x}{E_y}\right|=0,\left|\frac{E_z}{E_y}\right|=0.\end{array}$$(D.17)

As a result, the mode is a transverse wave with the electric field vector in the k × B direction, that is, the X-mode. For the other two refractive index solutions, we can obtain:

$$\begin{array}{c}\hfill \left|\frac{E_y}{E_x}\right|=0,\left|\frac{E_z}{E_x}\right|=\left|\frac{n^2-1}{tan{\theta }_B}\right|.\end{array}$$(D.18)

The mode is called O-mode which the electric field oscillates at k − B plane.

We plot the dispersion relationship of the X-mode and the O-mode for two cases (γ = 2 and γ = 10) in Figure D.1. For both cases, θ_B = 0.2 is adopted. The O-modes have two branches: the superluminous branch (n < 1) and the subluminous branch (n > 1). The superluminous O-mode has a frequency cut-off at ωγ³ = ω_p. The latter one has a refractive index cut-off at n = 1/cosθ_B. This mode looks to have two branches as shown in Figure D.1, but the branches can get crossed at high frequencies.

Appendix E: Notation List

Subscript i, j, k: The identifier of each charged particle

Subscript s: Species of a single particle c: Speed of light

d²W/dωdΩ: Energy radiated per unit solid angle per unit frequency interval

e: Charge (absolute value) of electron

f(u): Distribution function of plasma

h_gap: Height of the gap

k: Wave vector

k_B: Boltzmann constant

l: Longitude bunch size

ℓ: Angular separation constant

m_e: Mass of the electron

n: Unit vector of the line of sight

n: Refractive index

n_e: Number density of charges

n_GJ: Goldreich-Julian density

n_p: Number density of magnetospheric plasma

p: Power law index

r: Distance to the neutron star center

s: Distance that charges traveled along the field lines

u: Space-like part of the four-velocity

u_m: Average u over the Gaussian distribution

v: Velocity of electron

w: Burst width

A_∥: Parallel component of amplitude

A_⊥: Perpendicular component of amplitude

B: Magnetic field strength

E: Wave electric field

E_1, ∥: Perturbation electric field parallel to the B-field

E_1, ⊥: Perturbation electric field perpendicular to the B-field

E_∥: Electric field parallel to the B-field

I,  Q,  U,  V: Stokes parameters

J: Current density

j: Emissivity

K_ν: Modified Bessel function

ℒ_b: Luminosity of a bunch

M: Stellar mass

ℳ: Multiplicity factor

N_b: Number of bunches that contribute to the observed power at an epoch

N_e: Number of electrons in one bunch

N_i: Maximum number of the subscript of i

N_j: Maximum number of the subscript of j

N_k: Maximum number of the subscript of k

P: Spin period

R: Stellar radius

R_b: Transverse size of bunch

ℛ: Distance from the emitting source to the observer

T: Temperature

𝒯: The average timescale of the flux density I repeated

Y_ℓm(θ, ϕ): Spherical harmonics

α: Angle between spin axis and magnetic axis

β: Dimensionless velocity of one charge

${\dot{\beta }}_{\Vert }$: Dimensionless perturbation acceleration parallel to the B-field

${\dot{\beta }}_{\perp }$: Dimensionless perturbation acceleration perpendicular to the B-field

δ_V: Vacuum polarization coefficient

γ: Lorentz factor

μ: Shear modulus ν: Frequency, index of the modified Bessel function

$\dot{\nu }$: Drift rate

χ: Angle between $\widehat{\mathit{\beta }}$ for different trajectories

χ_d: Lower boundary of angle between $\widehat{\mathit{\beta }}$ for different trajectories

χ_u: Upper boundary of angle between $\widehat{\mathit{\beta }}$ for different trajectories

ϵ: Dielectric tensor

λ: Wavelength

ψ: Polarization angle given by the RVM model

ρ: Curvature radius

σ: Conductivity tensor

σ_c: Curvature cross section

τ: Optical depth

θ: Poloidal angle with respect to the magnetic axis

θ_c: Spread angle of the curvature radiation

θ_in: Angle between the incident photon momentum and the electron momentum

θ_s: Angle of the footpoint for field line at the surface

θ_B: Angle between angle between $\widehat{\mathit{k}}$ and $\widehat{\mathit{B}}$

ω: Angular frequency

ω_a: Frequency of the perturbation-induced acceleration

ω_B: Cyclotron frequency

ω_c: Critical frequency of curvature radiation

ω_cut: Angular frequency cut-off of the superluminous O-mode

ω_in: Angular frequency of incoming wave of ICS

ω_M: Frequency of bulk of bunches distribution

ω_p: Plasma frequency

φ: Azimuth angle to the magnetic axis, the angle between LOS and the trajectory plane

φ_d: Lower boundary of angle between LOS and the trajectory plane

φ_t: Half opening angle of a bulk of bunches

φ_u: Upper boundary of angle between LOS and the trajectory plane

Γ(ν): Gamma function

Ψ: Polarization angle

Θ_e: Magnetic colatitude of the emission point

Ω: Solid angle of radiation, angular frequency of the oscillations.

All Figures