Numerical Simulations of Particle Acceleration at Interplanetary Quasi-perpendicular Shocks

Many energetic particle events in space are associated with the acceleration at collisionless shocks (e.g., Heras et al. 1992; Kallenrode & Wibberenz 1997; Zank et al. 2000; Li et al. 2003; Rice et al. 2003; Wang et al. 2012; Qin et al. 2013; Qi et al. 2017). The most popular theory for charged-particle acceleration is diffusive shock acceleration (DSA; Axford et al. 1977; Krymsky 1977; Bell 1978; Blandford & Ostriker 1978), which is capable of predicting a power-law distribution downstream of the shock with a spectral exponent that depends only on the density compression ratio. It is widely accepted that this theory explains many important energetic particle events related to collisionless shocks. However, the energy spectra observed in space are not always in agreement with that of the prediction of DSA theory. In some cases, the observed energy spectra have exponential-like rollovers at higher energies (Ellison & Ramaty 1985). The reason for the spectra rollover is related to the loss of particles, limited acceleration time, or complicated shock geometries, all of which are not considered in the basic diffusive theory.

The angle between the shock normal and the upstream magnetic field, ${\theta }_{\mathrm{Bn}}$, is one of the most important effects that controls the injection of particles at the shocks. For quasi-parallel shocks, thermal particles can easily move along the magnetic field lines to cross the shock front repeatedly, so that they can be accelerated to higher energies (Ellison 1981; Quest 1988; Scholer 1990; Giacalone et al. 1992). In contrast, the acceleration of thermal particles at quasi-perpendicular shocks is more difficult to understand because particles tied to the magnetic field lines nearly parallel to the shocks can be convected across the shock downstream. Many theoretical works have been devoted to the well-known injection problem at quasi-perpendicular shocks. In earlier studies, Ellison et al. (1995) showed that thermal particles could be accelerated by quasi-perpendicular shocks, but that the acceleration efficiency decreases with the increase of the shock-normal angle. Using hybrid simulations, Giacalone & Ellison (2000) and Giacalone (2003) also found that quasi-perpendicular shocks did not accelerate particles efficiently. Subsequently, Giacalone (2005) performed large-scale self-consistent plasma simulations and demonstrated that thermal particles could be accelerated by perpendicular shocks. It is assumed that with large-scale magnetic fluctuations, the particles tied to meandering magnetic field lines can experience multiple crossings of the shock front and thus, be diffusively accelerated.

Recently, injection mechanisms of shock acceleration have been investigated. The injection energy can be derived from diffusion models assuming small anisotropy in the particle distribution (Giacalone & Jokipii 1999; Giacalone & Kóta 2006; Zank et al. 2006). Neergaard Parker & Zank (2012) and Neergaard Parker et al. (2014) explored the injection energy at quasi-parallel and quasi-perpendicular shocks, respectively, through solving the steady-state transport equation. For quasi-parallel shocks, it is suggested that the injection energy is small, so an upstream Maxwellian distribution is sufficient to provide the particle injection for DSA (Neergaard Parker & Zank 2012). The Maxwellian distribution has the form of

$$\begin{eqnarray}&&{f}_{0}({\boldsymbol{p}})={n}_{1}{\left(\displaystyle \frac{1}{2\pi {mkT}}\right)}^{3/2}\exp \left(-\displaystyle \frac{{({\boldsymbol{p}}-{{\boldsymbol{p}}}_{0})}^{2}}{2{mkT}}\right),\end{eqnarray} \tag{ 1 }$$

where n₁ is the upstream particle (solar wind proton) number density, m is the mass of a proton, k is Boltzmann's constant, T is the temperature of the distribution, p is the proton momentum, and p₀ corresponds to the solar wind momentum in the spacecraft frame. But for quasi-perpendicular shocks, the injection energy is assumed to be large. In addition, it is known that for quasi-perpendicular shocks, all reflected ions can return to the shock front and form a suprathermal population in the upstream region, and that a kappa distribution with an enhanced tail of energetic particles is demonstrated to provide an appropriate particle injection for DSA (Neergaard Parker et al. 2014). The kappa distribution can be written as

$$\begin{eqnarray}{f}_{0}({\boldsymbol{p}}) & = & \displaystyle \frac{{n}_{1}}{{[\pi (2\kappa -3){mkT}]}^{3/2}}\displaystyle \frac{{\rm{\Gamma }}(\kappa +1)}{{\rm{\Gamma }}(\kappa -1/2)}\\ & & \times {\left[1+\displaystyle \frac{{({\boldsymbol{p}}-{{\boldsymbol{p}}}_{0})}^{2}}{(2\kappa -3){mkT}}\right]}^{-\kappa -1}.\end{eqnarray} \tag{ 2 }$$

In this paper, we perform test particle simulations of particle acceleration associated with quasi-perpendicular interplanetary shocks, with modeling magnetic turbulence, by numerically solving the Newton–Lorentz equation with a backward-in-time method to clarify the crucial shock features that are responsible for efficient particle acceleration. In Section 2, we describe the shock model and numerical methods. The numerical results are shown in Section 3. In Section 4, we present the conclusions and discussion.

Our numerical simulations, performed by using the code from Zhang et al. (2017) to solve the Newton–Lorentz equation to obtain the trajectories of particles in the shock frame of reference, are similar to those of Decker & Vlahos (1986a, 1986b) except that here we consider a two-component magnetic turbulence model (described below). The Newton–Lorentz equation for a test particle in the electric field ${\boldsymbol{E}}$ and magnetic field ${\boldsymbol{B}}$ is given by

$$\begin{eqnarray}&&\displaystyle \frac{d{\boldsymbol{p}}}{{dt}}=q[{\boldsymbol{E}}({\boldsymbol{r}},t)+{\boldsymbol{v}}\times {\boldsymbol{B}}({\boldsymbol{r}},t)],\end{eqnarray} \tag{ 3 }$$

where ${\boldsymbol{p}}$ is the particle momentum, ${\boldsymbol{v}}$ is the particle velocity, q is the particle charge, t is time, and the frame of reference is moving with the shock front. The convective electric field is ${\boldsymbol{E}}=-{\boldsymbol{U}}\times {\boldsymbol{B}}$, with ${\boldsymbol{U}}$ as the plasma flow speed. We consider a planar, fast-mode, collisionless shock with the geometry illustrated in Figure 1. The shock is located at z = 0, and the plasma flows in the positive z-direction with the speeds of U₁ and U₂ in the upstream and downstream regions as measured in the shock frame of reference, respectively. In the shock transition the plasma speed is given by

$$\begin{eqnarray}&&U(z)=\displaystyle \frac{{U}_{1}}{2s}\{(s+1)+(s-1)\tanh [\tan (-\pi z/{L}_{\mathrm{diff}})]\}.\end{eqnarray} \tag{ 4 }$$

Here, L_diff is the thickness of the shock transition, which is assumed to be small enough compared to the gyro-radii of particles in the magnetic field, and $s={U}_{1}/{U}_{2}$ is the compression ratio. The magnetic field is taken to be time independent as

$$\begin{eqnarray}&&{\boldsymbol{B}}({x}^{{\prime} },{y}^{{\prime} },{z}^{{\prime} })={{\boldsymbol{B}}}_{0}+{\boldsymbol{b}}({x}^{{\prime} },{y}^{{\prime} },{z}^{{\prime} }),\end{eqnarray} \tag{ 5 }$$

where ${{\boldsymbol{B}}}_{0}$ is the uniform background field, ${\boldsymbol{b}}$ is a zero-mean turbulent field being transverse to ${{\boldsymbol{B}}}_{0}$, and a Cartesian coordinate $({x}^{{\prime} },{y}^{{\prime} },{z}^{{\prime} })$ system is adopted with the ${z}^{{\prime} }$-axis to be parallel to the mean magnetic field ${{\boldsymbol{B}}}_{0}$. The turbulent magnetic fluctuations ${\boldsymbol{b}}$ are composed of a slab and 2D component (Matthaeus et al. 1990; Zank & Matthaeus 1992; Bieber et al. 1996; Gray et al. 1996; Zank et al. 2006),

$$\begin{eqnarray}&&{\boldsymbol{b}}(x^{\prime} ,y^{\prime} ,z^{\prime} )={{\boldsymbol{b}}}^{\mathrm{slab}}(z^{\prime} )+{{\boldsymbol{b}}}^{2D}(x^{\prime} ,y^{\prime} ),\end{eqnarray} \tag{ 6 }$$

where both of the two components are perpendicular to the large-scale background magnetic field, ${{\boldsymbol{B}}}_{0}$, and have Kolmogorov's spectra with a power index of $\nu =-5/3$ at the high wavenumber, k. The slab component is created with fast Fourier transform (FFT) in a box with a size of ${L}_{{z}^{{\prime} }}=1000{\lambda }_{\mathrm{slab}}$ along the mean magnetic field ${{\boldsymbol{B}}}_{0}$ and number of grids ${N}_{{z}^{{\prime} }}={2}^{21}\,=2097152$, where ${\lambda }_{\mathrm{slab}}$ is the correlation scale of the slab component and is set to be 0.02 au at 1 au. The 2D component is created with FFT in a 2D box with the correlation scale length of ${\lambda }_{2{\rm{D}}}=0.1{\lambda }_{\mathrm{slab}}$, the 2D box size of ${L}_{{x}^{{\prime} }}={L}_{{y}^{{\prime} }}=100{\lambda }_{\mathrm{slab}}$, and the number of grids ${N}_{{x}^{{\prime} }}={N}_{{y}^{{\prime} }}=4096$. We take the ratio of magnetic energy densities from different components in the turbulence to be ${E}_{\mathrm{slab}}\,:{E}_{2{\rm{D}}}=20\,:80$. For more details on the slab–2D magnetic field model, see Qin et al. (2002a, 2002b), Qin (2002), and Zhang et al. (2017). Note that the bulk velocity ${\boldsymbol{U}}$ and the mean magnetic field ${{\boldsymbol{B}}}_{0}$ satisfy the Rankine–Hugoniot (RH) conditions, but we do not consider RH for the magnetic turbulence for simplicity based on the small level and zero average of turbulence.

Zoom In Zoom Out Reset image size

In the diffusive shock acceleration model, energetic particles can gain energy by crossing the shock back and forth because of the magnetic field turbulence scatterings. In this work we adopt the time-backward method to calculate the trajectories of a large number (N₀) of test particles, which are released with the same target energy, E_t, measured in the spacecraft frame at the shock front. Note that N₀ varies from 6000 to 25000 in different events. Each particle orbit is followed with numerical calculations by a fourth-order Runge–Kutta method with an adjustable time step that maintains an accuracy of the order of ${10}^{-9}$ until the test particle's motion time is up to the preset time boundary t_b. At the time boundary, the particles have different energy E_b (with the corresponding momentum p_b), most of which are lower than the target energy E_t with the time-backward calculation because of shock acceleration. In addition, a few number (${N}^{{\prime} }$) of particles with E_b lower than a preset energy ${E}_{b0}$, which is larger than that associated with the background plasma speed, are abandoned. Since ${N}^{{\prime} }$ is much smaller than N₀, the fact that we abandon the ${N}^{{\prime} }$ particles does not affect the simulation results. Note that our simulation box is large enough to ensure that the particles do not escape the acceleration region, so the spatial boundary is not taken into account. We consider eight target energies according to the energy channels in the range of 47–4750 keV from the Low-Energy Magnetic Spectrometer 30 (LEMS30) and LEMS120 detectors in the Electron, Proton, and Alpha Monitor (EPAM) instruments on board the ${Advanced}\ {Composition}\ {Explorer}\ ({\text{}}{ACE})$ spacecraft. For each of the ${N}_{E}=8$ target energies, E_i, we consider the target pitch angle cosine with an uniform distribution in the spacecraft frame, so the distribution function ${f}_{a}({p}_{i})$ at the shock front is obtained from an initial upstream background distribution function ${f}_{0}(p)$ through

$$\begin{eqnarray}&&{f}_{a}({p}_{i})=\displaystyle \frac{1}{{N}_{i}}\displaystyle \sum _{j=1}^{{N}_{i}}{f}_{0}({p}_{{bij}}),\end{eqnarray} \tag{ 7 }$$

where $i=1,\ldots ,{N}_{E}$ identifies the channels of energy, p_bij indicates the boundary momentum of the jth test particle, and ${N}_{i}={N}_{0}-{N}_{i}^{{\prime} }$ denotes the number of injected particles that contribute to the statistics. The corresponding accelerated energy spectrum is given by

$$\begin{eqnarray}&&j({E}_{i})=f({p}_{i}){p}_{i}^{2}.\end{eqnarray} \tag{ 8 }$$

For the upstream background distribution of energetic particles we choose the kappa distribution as shown in Equation (2), and we also choose the Maxwellian distribution as shown in Equation (1) for the comparison.

In our work for each of the shock events with the preset energy, ${E}_{b0}$, we gradually increase the preset simulation time until the accelerated particle energy spectra do not change obviously, and thus we obtain steady-state spectra and the acceleration time, t_acc. We then obtain the injection energy, E_inj, from the best ${E}_{b0}$ to fit steady-state spectra from simulations with EPAM observations. The process for the fit is done by examining the dependence of the target energy spectrum on the injection energy, and in our simulations, the injection energy error is up to the adopted steps of 0.1 keV and 0.01 keV for kappa and Maxwellian distributions, respectively.

Our model adopts the test particle method, so energetic particles do not affect shock waves and magnetic fields. In addition, we adopt a particle-splitting algorithm to improve the statistics of the accelerated particles. This method has been used in previous studies by Giacalone et al. (1992). In our present time-backward simulations, when the energy of a particle drops below a preset threshold (0.75 of the initial energy), two new particles are created at the location of the mother particle with slight velocity vector deviation. The weight of the split particles is one half of the weight of the mother particle prior to the split.

Next, we describe the turbulence level settings. In Neergaard Parker et al. (2014) the steady-state diffusion transport equation was solved, and in Giacalone (2015) the time-dependent transport equation was solved using a diffusion coefficient corresponding to a turbulence level that is significantly larger than of that observed. There is a general idea that the turbulence would be enhanced when the position is close to the shock front. In this work, we employ a decaying turbulence level model from the shock front which has the form

$$\begin{eqnarray}b/{B}_{0}=\left\{\begin{array}{ll}c\left(\tfrac{1}{{b}_{1}-{az}}+{d}_{1}\right), & z\lt 0\\ c\left(\tfrac{1}{{b}_{2}+{az}}+{d}_{2}\right), & z\geqslant 0\end{array}\right.\end{eqnarray} \tag{ 9 }$$

where a, b₁, d₁, and b₂ are the parameters assumed according to magnetic turbulence variations, c is the fitting parameter varying with different events, and ${d}_{2}={d}_{1}+1/{b}_{1}-1/{b}_{2}$ to make $b/{B}_{0}$ continuous at z = 0. For simplicity, we assume for different events that parameters a, b₁, d₁, b₂, and d₂ do not change, with c being the only variable. We use 16 s averaged observation data of $b/{B}_{0}$ as a function of time, t, from the Magnetometer Instrument (MAG) on board ${\text{}}{ACE}$ preceding and following the shock arrival during different events, where positive and negative t indicate the time after and before the shock front crossing, respectively. Here we assume $b/{B}_{0}$ as a function of z, which is the distance to the shock front with $z={V}_{\mathrm{sh}}t$, where V_sh is the shock speed in the spacecraft frame.

Because the main ingredient in the accelerated ions from EPAM observations are protons, we only consider proton acceleration at the shocks. We select SEP events observed at 1 au during 1998–2005 at quasi-perpendicular shocks with ${\theta }_{\mathrm{Bn}}\geqslant 75^\circ$ that are identified from the ACE shock lists provided by the ACE Web site (http://www.ssg.sr.unh.edu/mag/ace/ACElists/obs_list.html). In this work we only include quasi-perpendicular forward shocks where there was no additional shock within the previous 20 hr by following Neergaard Parker et al. (2014).

First, similar to Neergaard Parker et al. (2014), we construct the upstream kappa particle distribution with $\kappa =4$ using density, velocity, and temperature values obtained by averaging the five data points from the 64 s averaged data from the Solar Wind Electron, Proton, Alpha Monitor (SWEPAM) instrument on board ${\text{}}{ACE}$ ahead of the shock arrival, as listed in Table 1. Second, we construct a Maxwellian distribution as shown in Equation (1) with density, velocity, and temperature derived from observed plasma parameters for modeling the upstream particle distribution.

Comparing the MAG data on board ${\text{}}{ACE}$ for different shock events and the decaying turbulence level model with Equation (9), we get the parameters of $a=0.0001$ km⁻¹, ${b}_{1}=0.2,{d}_{1}=1,{b}_{2}=0.25$, and ${d}_{2}=2$. In addition, we get c for each event. In Figure 2, we plot the ACE/MAG observations of the magnetic field and the ratio of the root mean square (rms) value to the total magnetic field during a representative shock event on 1998 August 26. It is shown that there is an abrupt jump of $b/{B}_{0}$ at the time of the shock passage. Besides the data from spacecraft, we also plot the results from our model with Equation (9). Table 2 summarizes the input parameters of the shock events selected, such as the shock obliquity ${\theta }_{\mathrm{Bn}}$, shock speed V_sh, upstream speed U₁, Mach number M_A, compression ratio s, upstream magnetic field strength B₀, the turbulence parameter c, and the upstream convective electric field. The upstream background magnetic field is from 320 s averaged data of MAG on board ACE.

Zoom In Zoom Out Reset image size

We first present the results for the event that occurred on 1998 August 26. Figure 3 shows the accelerated energy spectra from simulations with the increase of the acceleration time (colored solid circles) along with 5 minute sector averages (red diamonds) of ACE/EPAM LEMS30 and LEMS120 data taken immediately after the shock. The dashed curve shows the original upstream kappa distribution with $\kappa =4$. It is clearly seen that as the time increases, the accelerated spectrum is gradually hardened, and after about 11 minutes, the spectrum reaches the steady state which is in good agreement with the observation, thus we assume the acceleration time of this event t_acc to be 11 minutes. It is noted that the injection energy E_inj for this event is 1.6 keV as indicated by the vertical red line in Figure 3.

Zoom In Zoom Out Reset image size

Similarly, we can get the acceleration time, t_acc, and injection energy, E_inj, for all of the quasi-perpendicular shock events. Figure 4 shows the accelerated energy spectra at the shock front for all of the 11 shock events with ${t}_{\mathrm{acc}}\leqslant 1$ hr. The black solid circles are the accelerated energy spectra from simulations with EPAM observations (red diamonds) overlaid. For the four events occurring on 1998 August 26, 1998 November 8, 2000 June 23, and 2002 May 23, the accelerated spectra fit well with the observations. However, for the rest of the seven cases, the simulated energy spectra are harder than those of the observations. These results indicate that our model is very efficient in accelerating particles with energies from a few tens of keV to several MeV.

Zoom In Zoom Out Reset image size

Figure 5 is similar to Figure 4 except that the shock events are with the acceleration time t_acc ranging from 1 to 20 hr. The accelerated energy spectra from simulations for these events, in general, fit well with the observed spectra. As can be seen, the 1999 September 15 and the 2001 May 12 events show the best agreement between simulations and observations. It is worth noting that the obviously slower upstream speeds of these events, which are within 130–200 km s⁻¹, compared to those of most shocks in the first subcategory are responsible for the longer acceleration time.

Zoom In Zoom Out Reset image size

We also simulated a number of shock events with even slower upstream speeds, which are about tens of kilometers per second, and found that all of these shocks are very difficult to accelerate the particles to higher energies within 20 hr. In principle, simulations should be performed for the time long enough with an upper limit of several days considering the finite age of interplanetary shocks. However, we here set an upper limit of an acceleration time of 20 hr, i.e., a 20 hr cutoff exactly, in the simulations to ensure the single shocks in our model. We list three representative examples in Table 2 and plot the simulated results in Figure 6. From Figure 6 we can see that all of the three shocks do not have the ability to accelerate particles to energies of ∼3 MeV with the acceleration time of 20 hr; therefore, all of their acceleration time ${t}_{\mathrm{acc}}\gt 20$ hr.

Zoom In Zoom Out Reset image size

For all of the shock events, the acceleration time t_acc with the shock obliquities ${\theta }_{\mathrm{Bn}}$ and upstream speeds U₁ are shown in Table 3. We plot the relationship between t_acc and U₁ for the events with ${t}_{\mathrm{acc}}\lt 20$ hr in the top panel of Figure 7 as asterisks. The red curve in top panel of Figure 7 shows the fitting to an inverse proportional function in log-linear space with the form of

$$\begin{eqnarray}&&{t}_{\mathrm{acc}}={t}_{\mathrm{acc}0}\times {10}^{{U}_{10}/{U}_{1}},\end{eqnarray} \tag{ 10 }$$

where ${t}_{\mathrm{acc}0}=0.291$ minutes, U₁₀ = 413 km s⁻¹, and the correlation coefficient $r=-0.76$. The fact that $r\lt -0.7$ indicates a negative correlation between the acceleration time and the upstream speed. For the shock events with U₁ larger than 200 km s⁻¹, t_acc is less than 1 hr, while for most events with U₁ in the range of 130–200 km s⁻¹, t_acc is in the range of 1–20 hr. The shocks with even lower U₁, such as those during the 1999 March 10 and 1999 October 28 events, are found to be incapable of accelerating sufficient particles to high energies within 20 hr.

Zoom In Zoom Out Reset image size

Table 3.  Shock Obliquity, Upstream Speed, Acceleration Time, Injection Energy, and the Ratio of the Intensity for the Shock Events

Date of Shock	${\theta }_{\mathrm{Bn}}$	U1 (km s−1)	tacc (minutes)	Einj (keV)	Ra (47–65 keV)
1998 Aug 26	88	325	11	1.6	8
1998 Sep 24	87	297	3	2.0	12
1998 Nov 8	77	210	12	1.5	3
1999 Mar 10	78	127	>1200	2.7	1
1999 Sep 15	83	168	336	3.9	2
1999 Oct 28	83	64	>1200	1.9	1
2000 Feb 11	88	224	14	1.4	7
2000 Jun 23	85	201	31	1.5	2
2000 Jul 19	81	150	>1200	1.5	4
2000 Sep 6	85	131	396	1.2	5
2001 Apr 28	88	492	3	1.8	26
2001 May 12	84	141	486	3.0	1
2001 Aug 27	89	150	16	2.2	5
2002 Jan 10	75	206	35	3.6	8
2002 May 23	84	466	4	1.3	62
2002 Sep 7	89	231	8	1.0	16
2003 Nov 15	81	162	450	4.0	5
2004 Jan 22	89	254	7	1.7	19
2005 Jul 10	82	147	390	1.5	4
2005 Aug 1	87	142	186	1.9	3

Download table as:  ASCIITypeset image

There are some individual events deviating from the general trend discussed above, as the observed shock events are very complicated and the acceleration time can be affected by many shock parameters such as the shock-normal angle, upstream speed, magnetic field strength, compression ratio, and so on. The shock that occurred on 2001 August 27, which is classified into subcategory 1 as shown in Figure 4, has ${t}_{\mathrm{acc}}=16$ minutes but an upstream speed of U₁ = 150 km s⁻¹, which is slower than 200 km s⁻¹. The short acceleration time for this event may be related to other factors, such as the shock obliquity discussed below. This shock has a normal angle of 89°, which is very close to 90° and may compensate for the low U₁ to result in ${t}_{\mathrm{acc}}\lt 1\,\mathrm{hr}$. On the other hand, there are some events in subcategory 2 with U₁ around 150 km s⁻¹ and ${\theta }_{\mathrm{Bn}}$ 3°–9° deviated from 90° in which t_acc is about a couple of hours. Therefore, it is suggested that the shock geometry (${\theta }_{\mathrm{Bn}}$) also plays an important role in determining the acceleration time besides the upstream speed.

To better understand the role of upstream speed in the shock acceleration, we examine the energetic particles intensity during the shock crossing for all 20 of the shock events. Here we define a quantity R_a from spacecraft observations as the ratio of the flux of 47–65 keV at the time of the shock crossing to the background flux, which is determined from the averaged flux of 12 points around 2 hr preceding the shock arrival. The values of the flux ratio R_a are shown in column 7 of Table 3, and we plot the relationship between R_a and U₁ in the bottom panel of Figure 7 as asterisks. In addition, the red line in the bottom panel of Figure 7 shows a linear fitting in log-linear space with the form of

$$\begin{eqnarray}&&{R}_{a}={R}_{a1}\times {10}^{{U}_{1}/{U}_{11}},\end{eqnarray} \tag{ 11 }$$

where ${R}_{a1}=0.853$, ${U}_{11}=270$ km s⁻¹, and the correlation coefficient $r=0.83$. It is noted that the upstream speed is not the only parameter to control particle acceleration, and as the shock event conditions are complicated with many shock parameters varying, a relatively large upstream speed does not always correspond to a relatively high flux ratio. However, here we intend to focus on the general trend between U₁ and R_a. The fact that the correlation coefficient is large ($r=0.83$) indicates a positive correlation between R_a and U₁. Considering that the flux ratio R_a is actually a measure of the shock ability for accelerating particles to high energies, the observational results with Equation (11) are consistent with our simulation results, showing a strong acceleration ability for the shocks with large upstream speeds, which can be explained from the theoretical view: assuming a fixed background magnetic field, a large upstream speed will yield a high convective electric field, which contributes to a large flux at the shock crossing.

We also consider an upstream Maxwellian distribution in order for a comparison with the kappa distribution for each shock event. The injection energy ${E}_{\mathrm{inj}}^{K}$ and ${E}_{\mathrm{inj}}^{M}$ from kappa and Maxwellian distributions, respectively, are shown in Table 4. We find that in each shock event ${E}_{\mathrm{inj}}^{K}\gt {E}_{\mathrm{inj}}^{M}$. Furthermore, most of the values of ${E}_{\mathrm{inj}}^{K}$ and ${E}_{\mathrm{inj}}^{M}$ lie within the theoretical range of 1–10 keV, obtained in the condition of $\langle {b}^{2}\rangle /{B}^{2}=0.78$ in Zank et al. (2006). Therefore, the Maxwellian distribution can produce an accelerated distribution that matches the data. The higher injection energy for the kappa distribution is due to the fact that it is a more appropriate representation of the suprathermal population for quasi-perpendicular shocks.

Table 4.  Injection Energy for Kappa and Maxwellian Upstream Distributions for the Shock Events

Date of Shock	${E}_{\mathrm{inj}}^{K}$ (keV)	${E}_{\mathrm{inj}}^{M}$ (keV)
1998 Aug 26	1.6	1.43
1998 Sep 24	2.0	1.54
1998 Nov 8	1.5	1.42
1999 Mar 10	2.7	1.40
1999 Sep 15	3.9	2.65
1999 Oct 28	1.9	1.13
2000 Feb 11	1.4	1.30
2000 Jun 23	1.5	1.26
2000 Jul 19	1.5	1.45
2000 Sep 6	1.2	1.08
2001 Apr 28	1.8	1.60
2001 May 12	3.0	1.70
2001 Aug 27	2.2	1.66
2002 Jan 10	3.6	2.38
2002 May 23	1.3	1.20
2002 Sep 7	1.0	0.93
2003 Nov 15	4.0	3.05
2004 Jan 22	1.7	1.53
2005 Jul 10	1.5	0.99
2005 Aug 1	1.9	1.55

Download table as:  ASCIITypeset image

In this paper, we perform test particle simulations of particle acceleration associated with quasi-perpendicular interplanetary shocks with modeling magnetic turbulence by solving the Newton–Lorentz equation numerically with a backward-in-time method. We identify a set of quasi-perpendicular shocks from the ACE shock database at 1 au from 1998 to 2005, and use the observed solar wind parameters to construct kappa and Maxwellian functions as the background upstream particle distributions (Neergaard Parker et al. 2014). By comparing the accelerated energy spectra between simulations and observations, we find that the shocks are capable of accelerating thermal particles to high energies of the order of MeV with both kappa and Maxwellian upstream particle distributions. In addition, the injection energy and timescale of particle acceleration are obtained. Through examining the relationship between the acceleration time and the parameters such as upstream speed and shock-normal angle ${\theta }_{\mathrm{Bn}}$, we clarify the crucial shock features that are responsible for efficient particle acceleration.

It it noted that our simulation results indicate that we can get accelerated energy spectra with the Maxwellian upstream distribution from quasi-perpendicular shocks. Neergaard Parker et al. (2014) also calculated accelerated distribution from the Maxwellian using the final kappa injection energy. It was shown that in many cases, the obtained accelerated distribution from the Maxwellian is slightly less than or comparable to that produced by the kappa distribution. However, in our simulations for all cases the injection energy from the Maxwellian distribution is less than that produced by the kappa distribution. This difference may originate from the fact that we used test particle simulations by solving the Newton–Lorentz equation while Neergaard Parker et al. (2014) theoretically solved the steady-state transport equation, and we adopted a turbulence model in which the local background magnetic field has a component parallel to the shock normal, so the shock is not precisely quasi-perpendicular locally. Although in our simulations for the quasi-perpendicular shocks the shock normal is almost perpendicular to the large-scale background magnetic field, it is usually less perpendicular to the local background magnetic field, since in our magnetic turbulence model the local averaged magnetic field deviates from the large-scale averaged one. Therefore, the injection energy in theory of quasi-perpendicular shock acceleration should not be a problem in our model.

In the future, we plan to study the particle acceleration at quasi-parallel interplanetary shocks with similar models used in this work.

The authors thank Drs. Yang Wang, Gang Li, and Chao-Sheng Lian for useful discussions related to this study. We thank the ACE EPAM, SWEPAM, and MAG instrument teams, and the ACE Science Center for providing the ACE data. We also thank Qiang Hu and Vasiliy Vorotnikov for help with the ACE shock database. Our work are partly supported by grants NNSFC 41574172, NNSFC 41374177, and NNSFC 41125016.