A Self-regulated Stochastic Acceleration Model of Pulsar Wind Nebulae

Pulsar wind nebulae (PWNe) are clouds of the magnetized relativistic electron/positron plasma supplied from the central pulsar. However, the number of radio-emitting particles inside a PWN is larger than the expectation from the study of pulsar magnetospheres and then their origin is still unclear. A stochastic acceleration of externally injected particles by a turbulence inside the PWN is proposed by our previous studies. In this paper, the previous stochastic acceleration model of the PWN broadband spectra is improved by taking into account the time evolution of the turbulent energy and then the total energy balance inside a PWN is maintained. The turbulent energy supplied from the central pulsar is wasted by the backreaction from the stochastic particle acceleration and the adiabatic cooling according the PWN expansion. The model is applied to the Crab Nebula and reproduce the current broadband emission spectrum, especially the flat radio spectrum although time evolution of the turbulent energy (diffusion coefficient) is a bit complicated compared with our previous studies, where we assumed an exponential behavior of the diffusion coefficient.


Introduction
Pulsars are powerful wind-blowing objects that form relativistic magnetized plasma bubbles called pulsar wind nebulae (PWNe).PWNe are expanding with time and a typical size of a few pc at an age of a few kyr (c.f.Ref. (1)).Broadband emission from radio through PeV gamma-rays is produced by high-energy particles accelerated at the termination shock of the relativistic pulsar wind (2; 3).Dynamics of relativistic magnetized plasmas, particle acceleration mechanism, and the properties of their central pulsars can be studied from the PWNe (4; 5).
One of the well-known problems for the dynamics of the PWN is the σ-problem, where the magnetization of the pulsar wind at the termination shock is required to be much smaller ).The three-dimensional magnetohydrodynamic simulation (8) pointed out that a possible answer to the σ-problem is 'turbulence' inside the PWN (c.f.Refs.(9; 10)).The turbulence would also play a role in particle 'spatial' diffusion (11; 12; 13; 14), where the observed extent of PWNe in X-rays (15) and in γ-rays (16) can be explained by spatial diffusion of the accelerated particles.Polarization observations (17; 18; 19) also require the turbulent component of the magnetic field comparable with the ordered magnetic field (20; 21; 22).
In this paper, we study the stochastic particle acceleration by the turbulence inside the PWN (particle diffusion in momentum space).Tanaka and Asano (2017, hereafter TA17) (23) developed a stochastic acceleration model of PWNe in order to explain the origin of the low-energy (γ 10 6 ) radio-emitting particles which are thought as a distinct population from the high-energy (γ 10 6 ) X-ray-emitting particles injected from the central pulsar (c.f.Refs.(24; 25)).The observed amount of the radio-emitting particles are much more than the central pulsar can produce in its life time (26; 27; 28; 29).The radio-emitting particles would be served externally, such as the surrounding supernova ejecta for example, and they are accelerated stochastically to form the observed hard radio spectra (23).However, in their model, the total energy was not conserved because the origin of the turbulence was not taken into account (23; 30).
Here, we develop a self-regulated turbulent acceleration model of PWNe.We solve the evolution of the turbulent energy E T (t) regarding it as a component of the PWN.The turbulent energy is continuously supplied from the central pulsar and is subsequently converted into the energy of the radio-emitting particles by the stochastic acceleration process (e.g.Ref. (31)).In Section 2, we describe our self-regulated acceleration model, and a relation with the diffusion coefficient in the momentum space and the energy of the turbulence is introduced.In Section 3, the results for the present model applying to the Crab Nebula are shown.We summarize the present paper in Section 4.

Self-regulated Stochastic Acceleration Model
Here, the self-regulated stochastic acceleration model is described focusing on improvements from the previous studies (23; 30).The turbulence inside the PWN (E T (t)) is explicitly set to an element in addition to the magnetic field (E B (t)) and the relativistic electrons/positrons (E e (t)).We introduce the fraction parameters, η T , η B and η e , where they control the injection rate of the turbulent ĖT , magnetic ĖB and particle (electron/positron) Ėe energy from the 2/16 central pulsar.The spin-down power of the central pulsar (L spin (t)) is divided into the three components, such as Ėi = η i L spin (i = T, B, e), and they satisfy the condition, where we expect η T is the same order as η e , while η B ≪ 1.In the present model, we do not require the sum of a fraction parameters (Eq.( 1)) to be unity, because η T is attributed to acceleration only for electrons while the ions injected externally (see Eq. ( 7)) would also be accelerated (see discussion in Section 4).The turbulent energy interacts with ions are ignored because the emission from the accelerated ions are inefficient.First, one of the roles of the turbulence within the one-zone spectral evolution model is a particle diffusion in the momentum space.Introducing the time-scale of the stochastic acceleration t acc , the diffusion coefficient is written as where γ is the Lorentz factor of electrons/positrons and t acc depends both on γ and t in general.Based on the quasilinear approximation of the gyroresonant scattering off the particles by turbulence, this acceleration time is inversely proportional to the spectral energy of the turbulence whose wavenumber corresponds the gyro-radius of an electron of γ (e.g.Refs.(31; 32)).Nevertheless, here, we use the hardsphere formula assuming the non-resonant scattering between the turbulence, i.e., t acc is independent from γ according to Tanaka and Kashiyama (2023, hereafter TK23 (30)).The acceleration time is still expected to be inversely proportional to the energy of the spectrally integrated turbulent energy E T (t), and then we introduce where E rot (t) = L spin dt is the total injected rotational energy from the central pulsar.In the present model, the acceleration time-scale t acc (t) changes with time because of the decay mechanisms of the turbulence discussed below, i.e., we consider the case E T (t) = η T E rot (t).τ acc is the acceleration time-scale without any decay processes of the turbulence and is an only parameter that controls the stochastic acceleration.Eq. ( 3) is a major difference from the previous studies, where an artificial parameter of the turbulence decay time-scale was implemented (see Eq. ( 7) of TK23).We wiil discuss validity of Eq. ( 3) in section 4.
The evolution of the turbulent energy E T (t) would be expressed as where t adi (t) ≡ R PWN (t)/ ṘPWN (t) is the adiabatic cooling time-scale and R PWN is the radius of the PWN.The first term in the right-hand side of Eq. ( 4) is injection from the pulsar as a random bulk motion of electron/positron plasma as a magnetohydodynamic fluid.The second term represents that the turbulence is attenuated by the backreaction of the stochastic acceleration of electrons/positrons (31) and will be defined in Eq. ( 6).Finally, the third term is the adiabatic cooling by the expansion of the PWN (33).Note that the turbulence should have the pressure, which we assume one-third of the turbulent energy density and the 3/16 turbulent pressure is taken into account for the expansion of the PWN R PWN (t) in addition to the particle and magnetic pressure (c.f.Refs.(30; 34; 35)).Eq. ( 4) is solved with the energy distribution of the accelerated particles N (γ, t), where γcool includes the adiabatic, synchrotron and inverse Compton coolings, Q ext represents an external source of the radio-emitting particles and Q PSR is the (single power-law) particle injection from the pulsar.The total particle energy inside the PWN is given by E e (t) = dγγm e c 2 N (γ, t), where m e and c are the electron mass and the speed of light, respectively.The third term on the left-hand side of Eq. ( 5) is the interaction between the turbulence and particles, and relates with the damping term in Eq. ( 4) as The backreaction of the stochastic particle acceleration to the turbulence is formulated self-consistently and the acceleration is suppressed by the decay of turbulence.
The external particle injection is described as where γ inj ∼ 1, m p is the proton mass and Ṁsh (t) is the rate of the mass of the supernova (SN) ejecta swept-up by the PWN expansion (see section 2.2 of TK23).f imp and f cont are the parameters which gives the particle injection rate to the stochastic acceleration process and should be smaller than unity.TA17 found that the externally injected low energy particles can be accelerated stochastically to be the radio-emitting particles.For the impulsive injection case Q ext (γ, t) ∝ f imp , we consider that a fraction f imp of the rotational energy at an early phase t init (≪ 10 yr) is almost instantaneously converted into the low energy e ± of γ = γ inj .For the continuous injection case Q ext (γ, t) ∝ f cont , we consider that a fraction f cont of the swept-up ejecta neutrals is continuously injected to the stochastic acceleration process by photoionization inside the PWN and then a similar number of hadrons can also be injected and accelerated inside the nebula.We ignore the acceleration and emission processes of the hadrons here, while we allow a fraction of the spin-down energy converted to ion acceleration, i.e., we set η T + η B + η e < 1 for the continuous injection case.Note that the upper limit of the pressure (energy) of the accelerated hadrons is at most the sum of the other components inside the PWN, because they can affect the PWN expansion if the hadrons become dominant in pressure.This upper limit is a requirement from the past studies of the dynamical and spectral evolution of PWNe (34; 36).The other part of the model is the same as our previous studies, where synchrotron radiation and inverse Compton scattering off synchrotron radiation (SSC) and off the cosmic microwave background radiation (IC/CMB) from the expanding PWN inside the parent supernova remnant are calculated (23; 26; 27; 28; 30; 37).The parameters of the model are those for the evolution of the supernova remnant (E SN , M SN , n ISM ), those for the evolution of the central pulsar (P 0 , B 0 , M NS , R NS , n), those for the injection of particles, magnetic fields, and turbulence from the pulsar (p, γ min , γ max , η e , η B , η T ), and those for external particle injection (f imp or f cont ).

4/16
We apply the present model to the Crab Nebula in the next section 3. The parameters of the parent supernova are set to E SN = 10 51 erg, M ej = 9.5 M ⊙ , and n ISM = 0.1 cm −3 (36; 38) and those of the Crab pulsar are also set to P 0 = 18.7 msec, B 0 = 3.35 × 10 12 G, n = 2.51, M NS = 1.4 M ⊙ and R NS = 12 km (L 0 = 2.6 × 10 39 erg s −1 , τ 0 = 1.1 kyr).The values of γ min , γ max and p are almost unique from previous studies (c.f.Refs.(24; 25)).The main parameter of the present study is only τ acc and we study the cases with τ acc = 1, 3, 10, 30, 100 yr.The other parameters, (f imp , f cont ), η B , η T , η e , are allowed to change while they have already been constrained by our previous studies.Both impulsive and continuous injection cases are studied.

Results
Figs. 1 and 2 show evolution of spectra of (a) photons and (b) electrons for the impulsive and continuous injection cases, respectively.We found that, for the both cases, the observed broadband spectrum of the Crab Nebula at an age of 1 kyr is reproduced with τ acc = 10 yr and then it seems difficult to distinguish the impulsive injection from the continuous one with the one-zone spectral evolution model as the previous studies (23 ; 30).Below, we just discuss the impulsive injection case in detail.The used parameters are summarized in Table 1.Note that we obtained the parameters fitting by eye because the calculation is too heavy to make a statistical fit by doing a lot of runs.
The top panel of Fig. 1(c) shows time evolution of the typical time-scales: the acceleration time t acc (t) (red solid), the synchrotron cooling time t syn (γ, t) (yellow dashed) and the adiabatic cooling time t adi (t) (blue dotted).The acceleration time-scale adopted by TK23 t TK23 (gray solid) is also plotted for reference, where t TK23 (t) = τacc exp(t/τ turb ) with (τ acc , τturb ) = (10, 70) yr.Only t syn (γ, t) ∝ γ −1 B −2 PWN (t) depends on γ and rapidly increases with time because the magnetic field decreases with the PWN expansion.The particles injected from the pulsar Q PSR whose Lorentz factor is γ ≥ γ min = 10 6 are not stochastically accelerated because t acc (t) is always longer than the other time-scales for γ γ min .The particles of γ 10 6 are stochastically accelerated at a time when t acc (t) < t adi (t) is satisfied, i.e., the stochastic acceleration starts at t age ≈ 25 yr and ends at t age ≈ 0.2 kyr.The complicated behaviors of t acc (t), such as an initial increase of a factor two around t age ∼ yr, an exponential increase phase in 50 t age 130 yr and a power-law increase phase after that, are not essential to reproduce the current observation, because the relatively simple t TK23 (t) still reproduce the observed flat radio spectrum (30).Evolution of the turbulent energy E T (t) reflects the complicated behaviors of t acc (t) (Eq.( 4)).
The middle and bottom panels of Fig. 1(c) are time evolution of each energy content per E rot (t) and its time derivative (luminosity) per L spin (t), respectively.The adiabatic loss of the turbulent energy plays a role for the initial (∼ yr) factor two increase of t acc (t) (factor two decrease of E T /E rot , see dot-dashed green line in the bottom panel of Fig. 1(c)).t acc (t) ∝ E rot (t)/E T (t) is almost constant, i.e., dE T /dt ∼ η T L spin , until the backreaction of the stochastic acceleration to the turbulence becomes significant ≤ 0.1 kyr and after that dE T /dt ∼ (δE T /δt) damp (see Eq. ( 6) and red solid line in the bottom panel of Fig. 1(c)).Finally, t acc (t) has a power-law dependence on time because η T L spin ∼ (δE T /δt) damp .
We find that the current broadband observations can also be reproduced by taking τ acc = 1 and 3 yr with very similar sets of the parameters, while we do not find the solution with τ acc = 30 and 100 yr.The results for τ acc = 3 and 30 yr are briefly summarized in Appendix 5/16  (c) Evolution of the typical time-scales (top), the fractional energies (middle) and luminosities (bottom) Fig. 1: The results for the impulsive injection case.(a) The spectral energy distribution and (b) the particle energy spectra (bottom) of the Crab Nebula are plotted for t age = 30 yr (yellow dashed), 300 yr (blue dotted), 1 kyr (red solid) and 3 kyr (green dot-dashed), respectively.The black points on (a) are the observational data in radio (39; 40; 41; 42; 43; 44), in IR and optical (45; 46; 47; 48; 49; 50), in X-rays (51), and in γ-rays (52; 53; 54; 55; 56; 57; 58).(c) Top panel: t acc (t) (Eq.( 3), red solid) and t adi (t) (Eq.( 4), blue dotted) evolve with time t, while the synchrotron cooling time-scales t syn (γ, t) (yellow dashed) depends on both t and γ, where γ = 10 4 , 10 6 , 10 8 are plotted.t TK23 (t) (definition in the text) is the acceleration time-scale adopted by TK23 (gray solid).(c) Middle panel: the fractional energy (per E rot (t)) of the particles E e (t) (black solid), of the turbulence E T (t) (red solid), of the magnetic field E B (t) = η B E rot (t) (purple dot-dot dashed), of the particle coolings by synchrotron (yellow dashed) and adiabatic loss (blue dotted) t 0 dt ′ dγ| γcool (γ, t ′ )|m e c 2 N (γ, t ′ ), of the adiabatic loss of the turbulent energy (green dot-dashed) t 0 E e,T (t ′ )/t adi (t ′ )dt ′ and of the sum of all of them (gray solid).(c) Bottom panel: the fractional luminosity (per L spin (t)) of stochastic acceleration backreaction (Eq.( 6), red solid), synchrotron (yellow dashed) and adiabatic cooling rates for the particle (blue dotted) and the turbulence (green dot-dashed) energy.The energy injection rates η e,T,B (black solid, gray solid, and purple dot-dot dashed) are also shown.The adopted parameters are summarized in Table 1.
6/16  (c) Evolution of the typical time-scales (top), the fractional energies (middle) and luminosities (bottom) Fig. 2: The results for the continuous injection case.The adopted parameters are summarized in Table 1.
A. We also studied the dependence on t init by taking an order of magnitude small value of t init = 0.1 yr and then the behavior of t age a few years are almost the same as the results with t init = 1 yr shown in this section.

Discussion & Conclusions
In this paper, the previous stochastic acceleration model is improved by considering the time evolution of the turbulent energy and then the total energy balance inside a PWN is maintained.The turbulent energy supplied from the central pulsar is wasted by the backreaction from the stochastic particle acceleration and the adiabatic cooling according to the PWN expansion.The model is applied to the Crab Nebula and reproduces the current broadband emission spectrum, especially the flat radio spectrum, although time evolution of the turbulent energy (diffusion coefficient) is different from our previous studies, where we assumed an exponential behavior of the diffusion coefficient (see the top panels of Figs.1(c) and 2(c)).
The obtained acceleration time-scale of τ acc = 10 yr or shorter (see Appendix A) is consistent with TK23.According to the discussion by TK23, the acceleration time-scale can be constrained by observing an early ( a few decades) evolution of PWNe in radio band.It is also difficult to distinguish the impulsive and continuous injection cases only from the present 7/16 observational data (23; 30).However, the continuous injection model would also accelerate ions which are simultaneously injected as electrons and then a characteristic neutrino signal would be expected.We leave the study of the neutrino emission as a future study.
The impulsive and continuous injection models would have different spatial distributions of the radio-emitting particles.For the impulsive injection case, the generated low energy e ± are distributed throughout the dense PWN, and then the radio-emitting particles are expected to be almost uniformly distributed according to the PWN expansion.For the continuous injection case, the particles are injected from the outer edge of the nebula, or more precisely, they are also injected from the supernova ejecta filaments, which penetrate the PWN due to the Rayleigh-Taylor instability (59; 60), so that we expect the radio-emitting particles to be distributed more in the outer part than near the center of the nebula.In both cases, the spatial distribution of the radio-emitting particles should also be different from that of the X-ray-emitting particles injected from the central pulsar.The spatial brightness profile also depends on the magnetic field structure and the spatial diffusion of the particles would also be important (11; 12; 13; 14).Future observations and modeling of the spatial structure of PWNe will distinguish the origin of the radio-emitting particles.
The essence of the present study is the modeling of the acceleration time, i.e., Eq. 3. We consider the non-resonant acceleration with large-scale eddies of a compressible turbulence, which is characterized by the fact that the acceleration time 'does not depend on the energy of the particles' (D γγ ∝ γ 2 ) and 'depends on the (spectrally integrated) energy density of the turbulence' (61; 62).Since the present model is a one-zone model, we use the spatially integrated energy rather than the (local) energy density of the turbulence (D γγ ∝ E T (t)).The energy of the turbulence injected so far by the pulsar (η T E rot (t)) is introduced only to normalize the present turbulent energy inside the PWN.In practice, the local physical quantities of the PWN plasma should be considered to determine the acceleration time, but it is not easy to define 'waves' parameterized by mass density, temperature, etc. in the PWN.The problem arises from the fact that the PWN consists of a relativistic electron-positron plasma with a 'pure' non-thermal distribution, and then the accelerated particles scattered by the turbulence and the background plasma responsible for the turbulence are not easily distinguishable as usual.The important conclusion obtained from the present study is that the requirement for t acc (γ, t) in order to satisfy the radio observations is quite simple, i.e., the stochastic acceleration stops as the turbulent energy decreases (see the top panels of Figs.1(c) and 2(c)).
The turbulence introduced to reproduce the broadband spectrum in this paper is not exactly the same as the turbulent magnetic field introduced to resolve the dynamical σproblem (9; 10).The energy of the turbulence E T (t) is the kinetic energy and is not the magnetic one because we do not include the synchrotron cooling by the turbulent magnetic field.In addition, E T (t) is directly supplied from the central pulsar but Tanaka et al. (2018) (10) considered that the turbulent magnetic field is converted from the ordered magnetic field although the kinetic turbulence would convert the ordered magnetic field into the tangled one.It is interesting that the turbulent energy E T (t) is close to E B (t) at an age of ∼ 1 kyr (see the middle panels of Figs.1(c) and 2(c)), where the ratio of the turbulent to ordered magnetic field is expected to be order unity from the polarization observation (63).The 8/16 turbulence inside PWNe would be a plausible solution to the long-standing problems of the PWN physics and then we should relate and combine them in the future studies.Acknowledgment S. J. T. would like to thank S. Kisaka, K. Murase, Y. Ohira, and K. Nishiwaki for useful discussion.The authors would like to thank the anonymous referee for helpful comments.This work is supported by JPJS Bilateral Program, Grant No. JPJSBP120229940 (SJT) and by JSPS KAKENHI 24H01816 (SJT), 21J01450 (WI).S. J. T. would like to thank the Sumitomo Foundation, and the Research Foundation For Opto-Science and Technology for support.

A. Results for different acceleration time-scales
Figs. A1 and A2 show evolution of the spectral energy distribution and the typical timescales, the fractional energies (middle) and luminosities (bottom) for τ acc = 3 yr (left) and 30 yr (right), respectively.The results for the continuous injection case are shown.The parameters are almost the same as the case of τ acc = 10 yr (Fig. 2).We do not obtain the flat radio spectrum for τ acc = 30 yr while the radio observations can be reproduced for τ acc = 3 yr.The turbulent energy is not fully converted into the energy of the radio-emitting particles for τ acc = 30 yr.We conclude that τ acc ≤ 10 yr is required to reproduce the current radio observations of the Crab Nebula.

Table 1 :
Summary of the model parameters for Figs.1 and 2.