Theoretical and simulation studies of relativistic ion holes in astrophysical plasmas

Theoretical and numerical studies of relativistic ion holes in a relativistically hot electron–ion plasma are presented. Previous particle-in-cell (PIC) simulations have shown that the ion holes are formed as a result of relativistic beam-plasma instabilities in the foreshock region of internal shocks of gamma-ray bursts and the relativistic jets of active galactic nuclei. In this process, the electrons are heated to ultra-relativistic temperatures so that their relativistic mass becomes comparable to the proton mass, and relativistic ion holes are formed by a secondary ion beam instability. The electrostatic potentials associated with the ion holes are large enough to accelerate particles to GeV energies. We use a semi-analytical model to construct relativistic ion holes and investigate their stability by means of fully relativistic Vlasov simulations. This investigation is relevant for astrophysical settings where the ion holes may work as efficient particle accelerators.

Theoretical and numerical studies of relativistic ion holes in a relativistically hot electron-ion plasma are presented. Previous particle-in-cell (PIC) simulations have shown that the ion holes are formed as a result of relativistic beam-plasma instabilities in the foreshock region of internal shocks of gammaray bursts and the relativistic jets of active galactic nuclei. In this process, the electrons are heated to ultra-relativistic temperatures so that their relativistic mass becomes comparable to the proton mass, and relativistic ion holes are formed by a secondary ion beam instability. The electrostatic potentials associated with the ion holes are large enough to accelerate particles to GeV energies. We use a semi-analytical model to construct relativistic ion holes and investigate their stability by means of fully relativistic Vlasov simulations. This investigation is relevant for astrophysical settings where the ion holes may work as efficient particle accelerators. The space x is divided by X u = 2πv b /ω pe , where v b ≈ c is the initial beam speed. (After [25].) In this paper, we present an analytical description of relativistic ion holes that are relevant for extreme environments of astrophysical objects. The investigation is motivated by the fact that the present relativistic ion holes are associated with electrostatic potentials of several GeV over short length-scales, and therefore they can be regarded as one of the strongest particle accelerators in the Universe. The dynamics and stability of the relativistic ion holes and their interaction with the electrons are investigated by means of a fully relativistic Vlasov simulation.

The excitation of relativistic ion holes
Recent PIC simulations [24,25] reveal the formation of relativistic ion (proton) holes by a beam instability relevant for the foreshock region of internal shocks of GRBs and the relativistic jets of active galactic nuclei. In this process, large-amplitude electron waves are excited through a relativistic Buneman instability which is saturated by electron trapping in the electrostatic potential of the wave. The trapped electrons are accelerated to ultra-relativistic (GeV) energies by a process similar to surfing acceleration [29]. As the relativistic mass of the electrons become comparable to the proton mass, a secondary instability is triggered in which the ion beam becomes unstable. Via this instability, the ions are partially thermalized to relativistic temperatures, and large-scale electrostatic structures are formed on the proton beam [24,25], identified as relativistic ion phase-space vortices/ion holes in which ions are trapped [25]. One of these ion holes, taken from the simulation of [25], is shown in figure 1. In this simulation, two counter-propagating proton beams are launched a speed of 99% of the speed of light, yielding a gamma factor of 7.1 for the protons relative to the bulk plasma. We see in figure 1 that the ions have been accelerated to gamma factors of ∼10 by the relativistically large potential of the ion hole. In [25], where the simulation box length was much larger than the ion hole, the large-scale electrostatic structures were relatively short-lived, while in [25] where the box length was the same as the width of the ion hole, the ion hole was stabilized and long-lived. As a result of the surfing acceleration and the following collapse of the electrostatic structures in [24], the electron distribution function became ultra-relativistic and flat-topped in momentum space up to a gamma factor of ∼10 3 -10 4 [24].

Theory of solitary ion holes in a relativistic plasma
In this section, we extend the non-relativistic ion hole theory of [1]- [3] to include relativistic mass variations of the plasma particles. The dynamics is governed by the relativistic Vlasov-Poisson system where the ion and electron number densities (normalized by n 0 ) are obtained as N i = ∞ −∞ f i dp and N e = ∞ −∞ f e dp, respectively. Here, M = m i /m e is the ratio between the ion and electron rest masses, γ i = 1 + α 2 e p 2 /M 2 and γ e = 1 + α 2 e p 2 are the relativistic ion and electron gamma factors, respectively, and α e = V Te /c is the ratio between the electron thermal speed V Te = (T e /m e ) 1/2 and the speed of light c. We will take M = 1836 in our numerical work. The distribution functions f i and f e have been normalized by n 0 /m e V Te , the time t by ω −1 pe , the space x by V Te /ω pe , the momentum p by m e V Te and the electrostatic potential φ by T e /e. Here, T e is the electron temperature, e is the magnitude of the electron charge, and ω pe = (4πn 0 e 2 /m e ) 1/2 is the non-relativistic electron plasma frequency, where n 0 is the unperturbed electron number density.
We now seek steady-state structures moving with velocity v 0 relative to the observer. With the ansatz f(p, ξ) and φ(ξ), where ξ = x − v 0 t and v 0 is the constant propagation speed of the ion hole, the system of equations (1)-(3) takes the form In the ion hole, a population of ions are trapped in a localized negative potential φ which goes to zero far away from the ion hole, while a population of the electrons are reflected on each side Institute of Physics ⌽ DEUTSCHE PHYSIKALISCHE GESELLSCHAFT of the potential minimum. The equilibrium of the ion hole can be constructed from the general solutions of equations (2) and (5), which are of the form f i = f i0 (E i ) and f e = f e0 (E e ), where and are the conserved energy integrals along the particle trajectories, and where trapped ions and reflected electrons have a negative energy, while the untrapped and unreflected (free) particles have a positive energy. Here, we have denoted γ 0 = 1/ 1 − α 2 e v 2 0 . Far away from the ion hole, where the potential is zero, we assume a relativistic Jüttner-Synge distribution function [30] for the ions, while the trapped ions inside the ion hole are taken to have a Jüttner-Synge distribution with negative 'temperature', leading to an excavated vortex distribution for the trapped ions. From these assumptions, the ion distribution function for the free ions can be written as 2 − M 2 are the separatrices between the free and trapped ions in momentum space, which are found by setting E i = 0 and solving for p. The modified momenta 2 have the property that they equal the momentum p when φ = 0, and thus the ion distribution function will converge to a relativistic Jüttner-Synge distribution function far away from the hole where the potential goes to zero. The trapped ion population, for where β i is a trapping parameter and τ = T i /T e . Negative values of β i , which we will consider here, leads to a negative 'temperature' of the trapped ions and a depletion in the ion distribution function. This choice of the distribution function for the trapped ion leads to an overall continuous distribution function for the ions. The normalization constant A i is chosen so that the ion number density is normalized to unity, For the electrons, we choose a relativistic Maxwellian function where the potential takes its minimum φ min , while the distribution function for the reflected electrons on both sides of the ion hole are taken to be flat-topped. Flat-topped distributions for the electrons have been observed to arise self-consistently in simulation studies of nonrelativistic colliding ion holes [31], as well as in the saturated stage of both weakly relativistic Buneman instabilities [32], relativistic Raman scattering instabilities [33] and in strongly relativistic ion beam instabilities [24]. The electron distribution functions for the free electrons take the form 1 are the separatrices between the free and reflected electrons, which are found by setting E e = 0 and solving for p. The modified momenta p e± (E e where φ = φ min , and there the electron distribution function will be set to a relativistic Jüttner-Synge function. The distribution function for the reflected electrons, in the momentum interval p e− (φ) p p e+ (φ), is simply taken to be the constant f e = A e exp[α −2 e (1 − γ 0 )] which equals the electron distribution for the free electrons for p = p e± (φ) and thus the electron distribution function is also continuous. The value of the normalization constant A e is determined by requiring that N e = f e dp = 1 at |ξ| = ∞ where φ = 0 far away from the ion hole. There we have the integration limits p e± = γ 0 v 0 (1 − γ 0 α 2 e φ min ) ± (γ 0 /α e ) (1 − γ 0 α 2 e φ min ) 2 − 1 between the free and reflected electrons and the electron energy integral is E e (0) = −v 0 p + α −2 e (γ e − γ −1 0 ) + φ min . By integrating the electron distribution function over momentum space, we obtain the normalization constant where we have used (in the second term inside curly brackets) that f e is independent of p in the interval p e− < p < p e+ . Integrating the distribution functions over momentum space we obtain the total ion and electron number densities as a function of the potential, which is calculated self-consistently by means of the Poisson equation (6). The integration of equation (6) can be done analytically with the potential method [2] to obtain the equation where the classical potential is For solitary ion holes it is required that ∂V/∂φ equals zero at φ = 0, and that V is negative between φ = 0 and φ = −ψ, where ψ = |φ| max denotes the amplitude of the ion hole potential. Equation (6) is useful for finding existence criteria [2] for solitary relativistic ion holes. However, in our numerical treatment below we solve equation (6) directly as a nonlinear boundary value problem, where the potential is set to zero far away from the ion hole, by an iterative numerical method (similar to Newton's method) to obtain profiles of the ion hole.

Numerical results
In order to study the influence of relativistic effects on the ion holes, we have solved the Poisson equation (6) numerically for the potential. In figure 2, we have plotted the amplitude of the potential as a function of the speed (normalized by the non-relativistic ion thermal speed) of the ion hole. We see that the amplitude of the ion hole decreases with increasing hole speed, which is similar to the non-relativistic case [2] for solitary ion holes. The solitary ion hole speed has an upper limit of ≈1.4 V Ti relative to the bulk plasma beyond which the ion hole ceases to exist. The potential increases strongly with increasing values on the relativistic factor α e . Here the potential is normalized by the electron thermal speed, so that in order to obtain the potential in volts, the numerical value should be multiplied by the factor α 2 e × 0.5 × 10 6 . For α e = 40, the amplitudes of the potential is of the order 0.5 GV, which can accelerate protons to relativistic gamma factors of a few and the electrons to gamma factors of ∼ 10 3 . As a comparison with [2], the dotted line in figure 2 shows the case with a non-relativistic Boltzmann distribution, N e = exp(φ), instead of the flat-topped distribution for the weakly relativistic electrons (the solid line). In figure 3, we have presented the existence diagrams for the ion hole, for one strongly relativistic case with α e = 40 and for one weakly relativistic case with α e = 0.1. We see that the amplitudes increase for smaller negative values of β i , and that there is a smallest value of β i below which we could not find ion hole solutions. For the non-relativistic case with Boltzmann-distributed electrons [2], the minimum value of β i (the symbol α used in [2]) was given as 0.71 for small values of ψ. This agrees quite well with figure 3, where the smallest β i is allowed to be ≈ −0.5 to −0.7 for small values of ψ. We also found that the ion hole has a maximum speed of v i0 ≈ 1.38 above which there are no ion hole solutions. In the non-relativistic case with Boltzmann-distributed electrons presented in [2], the upper limit on the ion hole speed was v i0 = 1.307(1 + τ) ≈ 1.44 for τ = 0.1. This result is thus very similar to the relativistic case presented here. This upper limit on the phase speed of the ion hole indicates that it rests on the slow ion acoustic mode [1,2,26], which in the non-relativistic case is a nonlinear and marginally stable mode that is not subject to linear Landau damping [9]. On the other hand, there are some relativistic modifications of the existence diagram for large values of ψ (and small values of β i ) which is seen by comparing the left-and right-hand panels of figure 3, in that the scaled amplitude of the potential ψ becomes larger in the non-relativistic case. In comparing the weakly and strongly relativistic cases, one has to keep in mind that the scaling encompasses the absolute amplitude of the potential (in volts), These diagrams have a similar appearance as for the non-relativistic case with Boltzmann-distributed electrons, see figure 7 in [2]. To obtain the potential in volts, ψ should be multiplied by the factor α 2 e × 0.5 × 10 6 .
which equals 0.5 × 10 6 α 2 e ψ. In the weakly relativistic case (α e = 0.1), we thus have amplitudes of the order of 10 kV, while in the relativistic case (α e = 40) we have amplitudes of the order of 0.5 GV. There is also a restriction T e /T i > 3.5, given for the existence of the non-relativistic ion hole in [2], and there should be some relativistic corrections to this existence criterion. A more thorough study of this matter is postponed to a future careful parametric study of relativistic ion holes. In figure 4, we display the profiles of the potential and the electron and ion densities for a few values of the relativistic factor α e . The amplitude of the localized potential increases for larger values of α e , and the ion hole width increases: to have the space ξ in physical units (cm), its numerical values should be multiplied by the factor α e c/ω pe . The electron and ion distribution functions are visualized in figure 5, for the same values of α e as in figure 4. We see that the ions are accelerated to relativistic momenta of p ≈ 1.2 m i c in the centre of the ion hole with α e = 40, while the electrons on both sides of the ion hole have a flat-topped distribution up to relativistic momenta of p ∼ 10 3 , which makes the relativistic mass of the electrons comparable to that of the ions. The dynamics of two interacting ion holes is studied by means of direct simulation of the coupled equations (1)-(3), and the results are displayed in figure 6. The temperature ratio τ = 0.1 and trapping parameter β i = −1 are chosen to show some typical behaviour of large-amplitude ion holes. Usually in simulations (e.g. [24,25]) a Buneman instability leads to a strong heating of the electrons, while the ions are heated by wave turbulence to a lesser extent. Hence, we have taken the electron temperature one order of magnitude larger than the ion temperature. The choice of β i = −1 was done to demonstrate the dynamics of large-amplitude ion holes, while  larger negative values of β i gives smaller amplitude ion holes (cf figure 3). We see in figure 6 that the ion holes, which are stable before the interaction, merge and form a new and slightly larger ion hole. We continued the simulation up to the time t = 15 000 and the single ion hole remained stable and essentially unchanged up to this point. This is in line with the stability analysis of [34]- [36], which shows that non-relativistic electron and ion holes are stable in one spatial dimension, but unstable in two spatial dimensions. In the small-amplitude limit, the nonrelativistic ion hole is governed by a modified Korteweg-de Vries equation, being referred to as the Schamel equation [1,3,26], [37]- [39]. The latter is not integrable, and hence it may turn out that small-amplitude ion holes collide inelastically. The ions that are trapped by the ion holes are accelerated to relativistic momenta so that their relativistic gamma factor is of the order ∼ 2. The ultra-relativistic electrons form a neutralizing background for the ion holes and are strongly depleted at the ion holes. Note that p should be multiplied by α e in order to have the momentum in units of m e c and by α e /M to have it in units of m i c, as has been done in figures 5 and 6.
In the simulation, we used a pseudo-spectral (Fourier transform) method to approximate the x derivatives in equations (1) and (2), and for integrating the Poisson equation (3). The p derivatives were approximated with a fourth order compact 'classical' Padé scheme [40], where f i and f e were set to zero at the boundary at the largest p. The time-stepping was performed with a fourth-order Runge-Kutta scheme. The domain size was 0 < x < 80 with periodic boundary conditions, while the momentum space was −300 p 300 for the electrons and −200 p 200 for the ions. We used 500 grid intervals in x space and 400 intervals in p space and a time step of t = 1.

Summary
We have presented a theoretical model for relativistic ion holes in hot plasmas that are relevant for astrophysical environments in the vicinity of GRBs and the relativistic jets of active galactic nuclei. Previous PIC-simulations [24,25] have demonstrated that relativistic ion beams can accelerate electrons to ultra-relativistic energies via an initial Buneman instability that excites electrostatic plasma waves that trap the electrons and accelerate them to ultra-relativistic energies by surfing acceleration across a magnetic field. The ion beam becomes unstable and is partially thermalized, and relativistic ion holes are formed in the saturation of the ion beam instability. The electrostatic potential of these ion holes can be in the GeV range, and trapped ions are accelerated to relativistic energies. By a relativistic Vlasov simulation, we have shown that solitary ion holes are robust structures that can accelerate ions to gamma factors of a few, while electrons are reflected and accelerated up to gamma factors of ∼10 3 by the ion hole potential. Future studies of the relativistic ion holes should include a more thorough parametric study to explore the existence and amplitude of the ion holes depending on the ion hole speed, the ion-toelectron temperature ratio and mass ratio etc, as well as theories for periodic ion holes that arise naturally by relativistic ion beam instabilities. Although our theoretical work excludes relative drift motions between the electrons and ions, they can easily be incorporated following [12]. In such circumstances, there would appear a new parameter of the differential equilibrium speed which can control the formation of relativistic ion holes.