Microinstabilities in the Transition Region of Weak Quasi-Perpendicular Intracluster Shocks
DDraft version February 12, 2021
Typeset using L A TEX twocolumn style in AASTeX63
Microinstabilities in the Transition Region of Weak Quasi-Perpendicular Intracluster Shocks
Sunjung Kim, Ji-Hoon Ha, Dongsu Ryu, and Hyesung Kang Department of Physics, School of Natural Sciences UNIST, Ulsan 44919, Korea Department of Earth Sciences, Pusan National University, Busan 46241, Korea (Received; Revised; Accepted)
Submitted to The Astrophysical JournalABSTRACTMicroinstabilities play important roles in both entropy generation and particle acceleration in colli-sionless shocks. Recent studies have suggested that in the transition zone of quasi-perpendicular ( Q ⊥ )shocks in the high-beta ( β = P gas /P B ) intracluster medium (ICM), the ion temperature anisotropydue to the reflected-gyrating ions could trigger the Alfv´en ion cyclotron (AIC) instability and theion-mirror instability, while the electron temperature anisotropy induced by magnetic field compres-sion could excite the whistler instability and the electron-mirror instability. Adopting the numericalestimates for ion and electron temperature anisotropies found in the particle-in-cell (PIC) simulationsof Q ⊥ -shocks with sonic Mach numbers, M s = 2 −
3, we carry out a linear stability analysis for thesemicroinstabilities. The kinetic properties of the microinstabilities and the ensuing plasma waves onboth ion and electron scales are described for wide ranges of parameters, including β and the ion-to-electron mass ratio. In addition, the nonlinear evolution of the induced plasma waves are examined byperforming 2D PIC simulations with periodic boundary conditions. We find that for β ≈ − M ∗ AIC ≈ .
3. Also electron-scale waves are generated primarilyby the whistler instability in these high- β shocks. The resulting multi-scale waves from electron to ionscales are thought to be essential in the electron injection to diffusive shock acceleration in Q ⊥ -shocksin the ICM. Keywords: acceleration of particles – cosmic rays – galaxies: clusters: general – methods: numerical –shock waves INTRODUCTIONMajor mergers of galaxy clusters are known to driveweak shocks with sonic Mach numbers, M s (cid:46)
3, in thehot intracluster medium (ICM) of high β (e.g., Ryuet al. 2003; Skillman et al. 2008; Vazza et al. 2009; Honget al. 2014; Ha et al. 2018). Here, the plasma beta, β = P gas /P B , is the ratio of the gas pressure to the mag-netic pressure. The radiative signatures of such shockshave been detected in X-ray and radio observations (e.g.,Br¨uggen et al. 2012; Brunetti & Jones 2014). In the caseof so-called radio relics, the radio emission has been in-terpreted as the synchrotron radiation from relativis- Corresponding author: Sunjung [email protected] tic electrons accelerated via diffusive shock acceleration(DSA) in merger-driven shocks (see van Weeren et al.2019, for a review).To explain the origin of radio relics, this DSA modelrequires an electron preacceleration mechanism, becausepostshock thermal electrons do not have momenta largeenough to participate in the standard DSA process, inwhich cosmic ray (CR) electrons diffuse across the shock(Kang et al. 2012). Since the width of the shock transi-tion layer is comparable to the gyro-radius of postshockthermal ions, thermal electrons need to be energized tothe so-called injection momentum, p inj ∼ a few × p i , th .Here, p i , th = √ m i k B T i2 is the ion thermal momentumin the postshock gas of temperature T i2 , m i is the ionmass, and k B is the Boltzmann constant. For shocks inthe solar wind, the electron injection is observed prefer-entially at the quasi-perpendicular ( Q ⊥ ) configuration a r X i v : . [ a s t r o - ph . H E ] F e b Kim et al.
Table 1.
Linear Properties of the Instabilities driven by Perpendicular Temperature Anisotropiesinstability AIC whistler ion-mirror electron-mirrorfree energy source T i ⊥ > T i (cid:107) T e ⊥ > T e (cid:107) T i ⊥ > T i (cid:107) T e ⊥ > T e (cid:107) propagation angle with γ ma parallel parallel oblique obliquewavenumber ck/ω pi ≤ ck/ω pe ≤ ck/ω pi ≤ ck/ω pe ≤ < ω r < Ω ci Ω ci < ω r < Ω ce ω r = 0 ω r = 0wave polarization LHCP b RHCP b Non-propagating Non-propagating aγ m is the maximum growth rate. b LHCP (RHCP) stands for left-hand (right-hand) circular polarization. with θ Bn (cid:38) ◦ , where θ Bn is the shock obliquity anglebetween the shock normal and the upstream magneticfield direction (e.g., Gosling et al. 1989; Oka et al. 2006;Burgess 2006).The electron preacceleration has been a key outstand-ing problem in understanding the production of CR elec-trons in weak ICM shocks. Previous studies have shownthat, in low- M s , high- β , Q ⊥ -shocks, thermal electronscould be preaccelerated primarily through the Fermi-likeacceleration in the shock foot (Matsukiyo et al. 2011;Guo et al. 2014a,b; Kang et al. 2019) and the stochas-tic shock drift acceleration (SSDA) in the shock transi-tion region (Katou & Amano 2019; Niemiec et al. 2019;Trotta & Burgess 2019).These two preacceleration mechanisms rely on the var-ious microinstabilities triggered by the ion and elec-tron temperature anisotropies in the shock structure(Gary 1993). If T e (cid:107) > T e ⊥ , for example, the elec-tron firehose instability (EFI) can grow with the follow-ing two branches: the nonresonant, parallel-propagatingmode with left-hand circular polarization, and the reso-nant, nonpropagating, oblique mode (Gary & Nishimura2003). Hereafter, the subscripts (cid:107) and ⊥ denote theparallel and perpendicular directions to the backgroundmagnetic field, B , respectively. Under the condition of T e ⊥ > T e (cid:107) , by contrast, the whistler instability and theelectron-mirror (e-mirror) instability can be triggered(Scharer & Trivelpiece 1967; Gary 1992; Hellinger &ˇStver´ak 2018). The most unstable whistler mode prop-agates in the direction parallel to B with right-handcircular polarization, while the e-mirror mode is non-propagating and has the maximum growth rate at thewavevector direction oblique with respect to B . In thecase of T i ⊥ > T i (cid:107) , the Alfv´en ion cyclotron instabil-ity (AIC, or the proton cyclotron instability) and theion-mirror (i-mirror) instability may become unstable(Winske & Quest 1988; Gary 1993; Gary et al. 1997;Burgess 2006). The fastest-growing mode of the AIC in- stability propagates in the direction parallel to B withleft-hand circular polarization, while the i-mirror modeis non-propagating and has the maximum growth rate atthe wavevector direction oblique with respect to B . Ta-ble 1 summarizes these linear properties of the instabil-ities driven by perpendicular temperature anisotropies,which are relevant for the present study.In the foot of Q ⊥ -shocks, the shock-reflected electronsbackstream mainly along the upstream magnetic fieldand induce an electron parallel anisotropy ( T e (cid:107) > T e ⊥ ),which could trigger the EFI and facilitate the Fermi-likepreacceleration (Guo et al. 2014b; Kang et al. 2019; Kimet al. 2020). In the transition region behind the shockramp, on the other hand, the AIC and i-mirror instabili-ties can be triggered by the ion perpendicular anisotropy( T i ⊥ > T i (cid:107) ) mainly due to the shock-reflected ions ad-vected downstream, while the whistler and e-mirror in-stabilities can be excited by the electron perpendicularanisotropy ( T e ⊥ > T e (cid:107) ) mainly due to magnetic fieldcompression at the shock ramp (Guo et al. 2017; Katou& Amano 2019). Such multi-scale waves from electronto ion scales are essential in the electron preaccelerationvia the SSDA (Matsukiyo & Matsumoto 2015; Niemiecet al. 2019; Trotta & Burgess 2019).Using two-dimensional (2D) particle-in-cell (PIC) sim-ulations for β ≈ − Q ⊥ -shocks, Kang et al.(2019) showed that the Fermi-like preacceleration in-volving multiple cycles of shock drift acceleration (SDA)in the shock foot could be effective only in supercriticalshocks with M s greater than the EFI critical Mach num-ber, M ∗ ef ≈ .
3. However, they argued that the electronpreacceleration may not proceed all the way to p inj , be-cause the EFI-driven waves are limited to electron scales.Niemiec et al. (2019), on the other hand, performed aPIC simulation for M s = 3 shock with β = 5 in a 2D do-main large enough to include ion-scale fluctuations, andsuggested that electrons could be energized beyond p inj via the SSDA due to stochastic pitch-angle scattering icroinstabilities in the Shock Transition Region β ≈ M ∗ AIC ≈ .
5. In a sep-arate paper (Ha et al. 2021, HKRK2021, hereafter), wereport a similar study of β ≈ −
100 shocks, whichis design to explore through 2D PIC simulations howthe multi-scale waves excited mainly by the AIC andwhistler instabilities in the shock transition can assistthe electron injection to DSA in ICM shocks.In this paper, adopting the numerical estimates forthe temperature anisotropies in the transition region ofthe simulated shocks of HKRK2021, we first performa linear stability analysis for the microinstabilities forwide ranges of parameters such as M s = 2 − β =1 − m i /m e =50 − β = 50and m i /m e = 50, we carry out 2D PIC simulations withperiodic boundary conditions (periodic-box simulations,hereafter) to study the nonlinear evolution of the plasmawaves excited by such microinstabilities.Note that throughout the paper we refer two differentsets of PIC simulations: (1) The ‘periodic-box simula-tions’ are designed to study the nonlinear evolution ofthe excited plasma waves in the same set-up as in thelinear analysis, and will be presented in Section 3. (2)The ‘shock simulations’ reported in HKRK2021 providethe numerical estimates for the ion and electron temper-ature anisotropies in the shock transition zone.The paper is organized as follows. Section 2 describesthe linear analysis of the AIC, whistler, and mirror in-stabilities. In Section 3, we present the evolution of thewaves driven by these instabilities in 2D periodic-boxsimulations. In Section 4, the implication of our workon the shock criticality and shock surface ripples is dis-cussed. A brief summary is given in Section 5. LINEAR ANALYSIS2.1.
Basic Equations
We consider a homogeneous, collisionless, magnetizedplasma, which is specified by the density and tempera-ture of ions and electrons, n i , n e , T i , T e , and the back-ground magnetic field B . The linear dispersion relationof general electromagnetic (EM) modes is given asdet (cid:18) (cid:15) ij − c k ω (cid:0) δ ij − k i k j k (cid:1)(cid:19) = 0 , (1) Figure 1. (a) Coordinate system employed in the presentstudy. The background magnetic field, B = B ˆ z , is parallelto the +ˆ z direction, while the wavevector, k = k x ˆ x + k z ˆ z , liesin the x - z plane. The wave propagation angle, θ , is the anglebetween B and k . (b) Schematic configuration showing thevelocity ellipsoid of a bi-Maxwellian VDF with a temperatureanisotropy A a , where a denotes either ‘ion’ or ‘electron’. where the dielectric tensor, (cid:15) ij , is determined by theplasma parameters and the velocity distribution func-tions (VDFs) of particles. Here, k i and k j are the com-ponents of the wavevector k . Then, the complex fre-quency, ω = ω r + iγ , can be calculated as a functionof the wavenumber, k , and the propagating angle, θ ,between k and B . Without loss of generality, we set B = B ˆ z along the + z direction and k = k x ˆ x + k z ˆ z in the x - z plane, as schematically illustrated in Figure1(a).In order to compute (cid:15) ij , we adopt the VDFs with bi-Maxwellian distributions for ions and electrons: f a ( v ⊥ , v (cid:107) ) = n π / v T a ⊥ v T a (cid:107) exp (cid:32) − v ⊥ v T a ⊥ − v (cid:107) v T a (cid:107) (cid:33) , (2)where v ⊥ = (cid:113) v x + v y and v (cid:107) = v z . The subscript a denotes e or i defined as the electron or ion species, re-spectively. Here, n is the number density of electronsor ions, which satisfies the charge neutrality condition,i.e., n = n e = n i . The parallel and perpendicular (to B ) thermal velocities are v T a (cid:107) = (cid:112) k B T a (cid:107) /m a and v T a ⊥ = (cid:112) k B T a ⊥ /m a , respectively. Then, the perpen-dicular temperature anisotropy of each particle speciesis given as A a ≡ T a ⊥ /T a (cid:107) = v T a ⊥ /v T a (cid:107) . The schematicconfiguration of the thermal velocity ellipsoid for a bi-Maxwellian VDF with the temperature anisotropy A a is shown in Figure 1(b). As A a increases, the thermalvelocity surface in the velocity space deviates furtheraway from the spherical shape. Under these considera- The quantity i is the imaginary unit, not the coordinate compo-nent nor for ion species. Kim et al.
Table 2.
Model Parameters and Linear PredictionsModel Name M s β e a β i a A e a A i a m i /m e AIC b whistler c ion-mirror b electron-mirror c LM2.0 β
20 2.0 9.7 10.3 1.1 1.1 50 stable stable stable stableLM2.0 β
50 2.0 24 26 1.1 1.2 50 stable (0.013,0.20,0) quasi-stable (0.0029,0.14,69)LM2.0 β
100 2.0 48 52 1.1 1.2 50 quasi-stable (0.035,0.20,0) (0.015,0.14,62) (0.008,0.15,56)LM2.3 β
20 2.3 8.4 12 1.1 1.5 50 (0.041,0.21,0) (0.0016,0.26,0) (0.04,0.28,63) stableLM2.3 β
50 2.3 22 28 1.2 1.5 50 (0.048,0.15,0) (0.03,0.24,0) (0.054,0.21,61) (0.0074,0.18,67)LM2.3 β
100 2.3 44 56 1.2 1.5 50 (0.053,0.11,0) (0.056,0.24,0) (0.063,0.16,58) (0.016,0.17,61)LM3.0 β β β
20 3.0 7.5 13 1.2 2.0 50 (0.127,0.29,0) (0.0156,0.30,0) (0.094,0.32,56) (0.0016,0.17,74)LM3.0 β
50 3.0 19 31 1.2 2.0 50 (0.145,0.20,0) (0.059,0.29,0) (0.11,0.23,55) (0.015,0.22,64)LM3.0 β
100 3.0 38 62 1.2 2.0 50 (0.156,0.15,0) (0.10,0.29,0) (0.12,0.18,54) (0.03,0.21,56)LM3.0 β β a The quantities, β e , β i , A e , and A i , are obtained by averaging numerical values over the transition zone in the simulated Q ⊥ -shocks presented in HKRK2021. b Linear predictions for the fastest growing mode of the ion-driven instabilities, ( γ m / Ω ci , ck m /ω pi , θ m ), normalized with the iongyro and plasma frequencies. θ m is given in units of degree. c Linear predictions for the fastest growing mode of the electron-driven instabilities, ( γ m / Ω ce , ck m /ω pe , θ m ), normalized with theelectron gyro and plasma frequencies. θ m is given in units of degree. tions, (cid:15) ij is given as Equation (3) in Kim et al. (2020)without the bulk drift velocities.For n , B , T a (cid:107) , and T a ⊥ , we adopt the numerical val-ues estimated for the simulated shocks of HKRK2021,where the preshock conditions are specified with thetypical parameters of high- β ICM plasmas, n ICM =10 − cm − , k B T ICM = ( k B T i + k B T e ) / . β ICM = 20 − β a = 8 πn a k B T a /B , the plasma frequency, ω =4 πn a e /m a , and the gyro-frequency, Ω ca = eB /m a c ,for electrons and ions are used. Plasma waves are characterized with the growth rate, γ , and the real frequency, ω r , which are calculatedby solving the dispersion relation in Equation (1) forwavevector k . If the propagation angle of the wave withthe maximum growth rate, γ m , is θ m ≈ ◦ , the wavemode is called ‘parallel-propagating’. If θ m (cid:29) ◦ , it is‘oblique-propagating’. If the wave frequency, ω r ≈ P , can be estimated also using the solution of the dis- Note that in HKRK2021 the results are expressed in terms ofthe upstream parameters, n up0 ≈ n /r and B up0 ≈ B /r , where r is the density compression ratio across the shock ramp. So forexample ω uppa ≈ ω pa / √ r and Ω upca ≈ Ω ca /r . persion relation as follows: P ≡ sign( ω r ) | δE + | − | δE − || δE + | + | δE − | , (3)where δE ± ≡ δE x k ,ω ∓ iδE y k ,ω (Verscharen & Chandran2013). The left-hand circular polarization (LHCP) cor-responds to P = −
1, whereas the right-hand circularpolarization (RHCP) corresponds to P = +1. Wavesare in general elliptically polarized with P (cid:54) = ±
1. Inthe case of non-propagating mode ( ω r = 0), P = 0 (seeTable 1).2.2. Linear Properties of AIC, Whistler and MirrorInstabilities
In this section, we report the results of the linear sta-bility analysis for the microinstabilities triggered by theion and electron temperature anisotropies in the transi-tion region of high- β , Q ⊥ -shocks. The first column ofTable 2 lists the model name, which is assigned withthe two parameters, the shock Mach number, M s , and β ( ≈ β e + β i ) in the shock transition region. For ex-ample, the LM3.0 β
50 model has M s = 3 . β ≈ β e , β i , A e , and A i , arelisted in the 3 − β ≈ −
100 and the mass ration m i /m e = 50,they are obtained with n , B , T a (cid:107) , and T a ⊥ estimated icroinstabilities in the Shock Transition Region - 1 - 2 - 1 L M 3 . 0 b b e = 1 8 . 9 , b i = 3 1 . 1 , m i / m e = 1 8 3 6 , w p e / W c e = 2 6 ) g / W i A e = 1 . 2 3 , A i = 2 . 0 2 A e = 1 , A i = 2 . 0 2 A e = 1 . 2 3 , A i = 1 ( a ) A I C ( q m = 0 o ) ( b ) W I ( q m = 0 o ) e - m i r r o r ( c ) M i r r o r ( q m = 6 4 o ) i - m i r r o r ( d ) A I C ( q m = 0 o ) w r/ W ci c k / w p i ( e ) W I ( q m = 0 o ) c k / w p i ( f ) M i r r o r ( q m = 6 4 o ) c k / w p i Figure 2. (a)-(c): Linear growth rate, γ , at the propagation angle of the fastest growing mode, θ m , for the AIC, whistler,and mirror modes, respectively, as a function of the wavenumber k for the LM3.0 β A e = 1 . A i = 1 . A e = 1 . A i = 2 . A e = 1 . A i = 2 .
0. In panel (c) both γ and k areplotted in logarithmic scales. (d)-(f): Real frequency, ω r , for the same case as the black dashed lines in the upper panels. Notethat γ and ω r are normalized with Ω ci and k is normalized with ω pi /c , uniformly for both the ion and electron modes. by averaging the numerical values over the transitionregion in the simulated shocks with M s = 2 − β up = 20 −
100 of HKRK2021. Considering the un-certainties in averaging over nonlinear structures withovershoot/undershoot oscillations, they are given onlyup to two significant figures.For the models with higher mass ratios, LM3.0 β m i /m e = 100 and LM3.0 β m i /m e = 1836, the parameters for the LM3.0 β
50 model( β e = 19, β i = 31, A e = 1 .
2, and A i = 2 .
0) are usedonly for the linear analysis. Also we carried out twoadditional shock simulations for M3.0 β β = 1and M3.0 β β = 5, which were not considered inHKRK2021, in order to obtain the parameters to be usedfor LM3.0 β β
5. Our fiducial models have m i /m e = 50, which is adopted in order to ease the re-quirements of computational resources for the periodic- Note that β up for the shock models in HKRK2021 represents theplasma beta of the upstream, preshock plasmas, while β for thelinear analysis models in Table 2 is the plasma beta in the shocktransition zone. We found that β ≈ β up for the simulated shocks,although, in general, the plasma beta of the far downstream re-gion is higher than β up . box PIC simulations that will be described in Section3. The linear predictions for the AIC, whistler, i-mirror,and e-mirror instabilities are given in the 8 −
11 columnsof Table 2. The three numbers inside each parenthe-sis present the linear properties of the fastest growingmode: ( γ m / Ω ci , ck m /ω pi , θ m ) for the AIC and i-mirrorinstabilities, and ( γ m / Ω ce , ck m /ω pe , θ m ) for the whistlerand e-mirror instabilities. Here, k m is the wavenumberthat has the maximum growth rate γ m at θ m , and θ m is given in units of degree. For a clear distinction be-tween the ion and electron mirror modes, in the 10 − γ m of each mirror mode, obtained with eitherisotropic electrons (i.e., A e = 1, A i >
1) or isotropicions (i.e., A i = 1, A e > γ m <
0, and‘quasi-stable’ corresponds to γ m / Ω ci < − .Figure 2 shows the linear analysis results for theLM3.0 β ω pe / Ω ce = 26. Panels (a)-(c) display the growth (ordamping) rate at θ m of the AIC, whistler, and mirrorinstabilities, respectively, as a function of the wavenum-ber. To make a uniform comparison, γ and k are nor-malized with Ω ci and ω pi /c , respectively, for both theion-driven and electron-driven instabilities. Note that Kim et al. ( a ) A I C ( q m = 0 o ) g / W ci L M 3 . 0 b b b g / W ce ( f ) W I ( q m = 0 o ) g / W ci ( c ) i - m i r r o r ( q m = 5 5 o ) g / W ce ( d ) e - m i r r o r ( q m = 6 4 o )( e ) A I C ( q m = 0 o ) g / W ci c k / w p i L M 3 . 0 b
1 L M 3 . 0 b
5 L M 3 . 0 b b b g / W ce c k / w p e ( b ) W I ( q m = 0 o ) q m = 5 6 o g / W ci c k / w p i q m = 5 5 o ( g ) i - m i r r o r q m = 5 4 o q m = 6 3 o g / W ce c k / w p e ( h ) e - m i r r o r q m = 7 4 o q m = 6 4 o q m = 5 6 o Figure 3.
Dependence of the linear growth rate, γ , on m i /m e (top) and β (bottom); γ at the propagation angle of the fastestgrowing mode, θ m , is given as a function of the wavenumber k . The model parameters are listed in Table 2. Note that for themirror modes, θ m depends on β , but not on m i /m e . in panel (c) both γ and k are given in logarithmic scales,in order to show both the i-mirror and e-mirror modesin the same panel. To examine the effects of A i and A e separately and also their combination, we presentthe black dashed lines for the case with both the ionand electron anisotropies, the red solid lines with theion anisotropy only, and the blue solid lines with theelectron anisotropy only.The AIC instability induces quasi-parallel modes with θ m = 0 ◦ . Although A i > A e > A e > A i . The whistlermode is also quasi-parallel propagating with θ m = 0 ◦ .The mirror modes, on the other hand, are highly obliquewith θ m = 64 ◦ for LM3.0 β k ( ck/ω pi > .
3) grows much faster thanthe i-mirror mode (red) at low- k ( ck/ω pi < . A i > A e >
1, a mixture of the two mirrormodes appears in the intermediate- k regime ( ck/ω pi ∼ . β γ WI (cid:29) γ E − M (cid:29) γ AIC > γ I − M , (4)where γ WI , γ E − M , γ AIC and γ I − M are the maximumgrowth rates of the whistler, e-mirror, AIC and i-mirrorinstabilities, respectively. Note that in general γ WI > γ E − M (Gary & Karimabadi 2006), and γ AIC > γ I − M un-der space plasma conditions with low- β and large tem-perature anisotropies (Gary 1992, 1993).The real frequency, ω r / Ω ci , at θ m for the mixed case( A e = 1 . A i = 2 .
0) are shown in panels (d)-(f) ofFigure 2. The AIC-driven mode has ω r / Ω ci ∼ . − . ck/ω pi ∼ . − .
4, while the whistler mode has ω r / Ω ci ∼ −
350 for ck/ω pi ∼ −
20. The mir-ror modes are nonpropagating or purely growing with ω r = 0. Moreover, the polarization, calculated usingthe solutions of the dispersion relation, is P = −
1, +1,and 0 for the AIC, whistler, and mirror instabilities, re-spectively, as expected.2.3.
Parameter Dependence of Linear Properties
As listed in Table 2, we consider a number of modelsto explore the dependence on m i /m e and β . The upperpanels of Figure 3 show the linear predictions for themodels with M s = 3, β = 50, and m i /m e = 50 − M s = 3, m i /m e = 50, and β = 1 − − . Nevertheless, γ AIC and γ E − M are al-most independent of m i /m e . In the case of the whistlerand i-mirror instabilities, on the other hand, overall,the growth rates are slightly lower for smaller m i /m e .Also the damping rate for the whistler instability isslightly higher for smaller m i /m e in the small wavenum-ber regime ( ck/ω pe ∼ . icroinstabilities in the Shock Transition Region ( a ) A I C ( q m = 0 o ) g / W ci L M 2 . 0 b b b ( b ) W I ( q m = 0 o ) g / W ce ( c ) i - m i r r o r q m = 5 6 o q m = 6 3 o g / W ci q m = 7 7 o ( d ) e - m i r r o r q m = 8 1 o q m = 8 5 o q m = 7 4 o g / W ce ( e ) A I C ( q m = 0 o ) g / W ci L M 2 . 0 b b b ( f ) W I ( q m = 0 o ) g / W ce ( g ) i - m i r r o r g / W ci q m = 5 5 o q m = 6 1 o q m = 6 4 o ( h ) e - m i r r o r g / W ce q m = 6 7 o q m = 6 9 o q m = 6 4 o ( i ) A I C ( q m = 0 o ) g / W ci c k / w p i L M 2 . 0 b b b ( j ) W I ( q m = 0 o ) g / W ce c k / w p e ( k ) i - m i r r o r g / W ci c k / w p i q m = 5 4 o q m = 5 8 o q m = 6 2 o ( l ) e - m i r r o r g / W ce c k / w p e q m = 6 1 o q m = 6 3 o q m = 5 6 o Figure 4.
Dependence of the linear growth rate, γ , on M s and β ; γ at the propagation angle of the fastest growing mode, θ m ,is given as a function of the wavenumber k . In each panel, the black, red, and blue lines show the results for M s = 2 . , β = 20 (top), 50 (middle), and 100 (bottom). The model parameters are listed inTable 2. Note that for the mirror modes, θ m depends on β . suppressed in the shock simulations with reduced massratios. However, even in the case of m i /m e = 50, thiseffect is expected to be only minor, because the inequal-ity in Equation (4) is still valid and the changes of k m and θ m are negligible (see the table 2).The plasma beta is another important parameterthat affects the stability of the system. Note that theanisotropy parameters, A e and A i , are almost indepen-dent of β for β ≈ − β in the second digit to theright of the decimal point. In the low- β case (LM3.0 β A i = 1 . β models due to the strong magnetization of ions. Thisis because A i in the shock transition is closely relatedto the fraction of reflected ions. On the other hand, A e in the shock transition is not substantially affected by β , because it is mainly determined by the magnetic fieldcompression rather than the fraction of reflected elec-trons. Given the same temperature anisotropies, thegrowth of the instabilities tends to be suppressed by strong magnetic fields at low- β plasmas. As can be seenin the lower panels of Figure 3, the peak values of ei-ther γ m / Ω ci for the ion-driven modes or γ m / Ω ce for theelectron-driven modes increase with increasing β . Forthe AIC, whistler, and i-mirror instabilities, γ m / Ω ci or γ m / Ω ce occurs at smaller ck/ω pi or ck/ω pe , for higher β . But such a trend is not obvious in the case of thee-mirror mode.In the high- β cases ( β ≈ − M s = 3,all the AIC, whistler, i-mirror and e-mirror waves canbe triggered, as shown in the lower panels of Figure 3,leading to the generation of multi-scale waves from elec-tron to ion scales. On the other hand, in the LM3.0 β β M s , is the key parame-ter that determines the temperature anisotropies in thetransition of high- β ICM shocks ( β ≈ − Kim et al. x [ c / w p e ] z[c/ w pe] B B x [ c / w p e ] z[c/ w pe] z[c/ w pi] x [ c / w p i ] - 0 . 6 0- 0 . 4 5- 0 . 3 0- 0 . 1 50 . 00 . 1 50 . 3 00 . 4 50 . 6 0 d By/B0 x [ c / w p e ] z[c/ w pe] x [ c / w p e ] z[c/ w pe] ( f )( e )( d ) ( c )( b ) t ~ t A I C t ~ W I z[c/ w pi] x [ c / w p i ] - 0 . 1 0- 0 . 0 7 5- 0 . 0 5 0- 0 . 0 2 50 . 00 . 0 2 50 . 0 5 00 . 0 7 50 . 1 0 d Bz/B0 ( a ) t ~ t W I B k Figure 5.
Magnetic field fluctuations, δB y (top) and δB z (bottom), in the periodic-box simulation for the LM3.0 β
50 model,plotted in the x - z plane. At early times, t ∼ − τ WI , shown in panels (a)-(b) and (d)-(e), electron-scale waves are excitedby the whistler and e-mirror instabilities, while ion-scale waves are generated by the AIC and i-mirror instabilities at t ∼ τ AIC shown in panels (c) and (f). Note that the 2D domain with [84 . × . c/ω pe ) is shown in panels (a)-(b) and (d)-(e), while the2D domain with [84 . × . c/ω pi ) is shown in panels (c) and (f). The black arrows indicate the direction of the backgroundmagnetic field, B , while the blue arrow in panel (f) shows the direction of the wave propagation, k , for the i-mirror mode withthe maximum growth rate. pression are closely related to M s . Figure 4 shows thegrowth rates of the instabilities for M s = 2 . β = 20 (top),50 (middle) and 100 (bottom). As M s increases, both A e and A i increase, so all the modes grow faster and k m shifts towards larger k , regardless of β .Note that the AIC and whistler modes have γ m at θ m = 0 independent of M s , whereas θ m decreases withincreasing M s for the i-mirror and e-mirror modes (seealso Table 2). In LM2.0 β
50 and LM2.0 β β
20, allthe instabilities are stable (see black lines in top panels).In the models with M s = 2 . − NONLINEAR EVOLUTION OF INDUCEDWAVES IN PERIODIC-BOX SIMULATIONS3.1.
Numerical Setup
To investigate the development and nonlinear evolu-tion of the instabilities, we performed 2D PIC simula-tions with periodic boundary conditions for the three fiducial models, LM2.0 β
50, LM2.3 β
50 and LM3.0 β β e , β i , A e , A i given in Table 2. As noted before,here m i /m e = 50 is employed due to the computationallimitations, but at least the early, linear-stage develop-ment of the plasma instabilities under consideration isexpected to depend rather weakly on the mass ratio.The simulations were carried out using a paral-lelized EM PIC code, TRISTAN-MP (Buneman 1993;Spitkovsky 2005). The simulation domain is a squareof box size L x = L z = 84 . c/ω pi = 600 c/ω pe in the x - z plane, which consists of the grid cells of ∆ x = ∆ z =0 . c/ω pe . In each cell, 32 particles (16 for ions and 16for electrons) are placed. The time step of the simula-tions is ∆ t = 0 . ω − , and the simulations were runup to t end = 130Ω − .3.2. Results of Periodic-Box Simulations
With the inequality in Equation (4), we expect thatthe whistler mode grows much faster than the othermodes, resulting in the relaxation of A e during the earlystage. As the whistler and e-mirror modes grow andthen decay on the time scale of τ WI ≡ /γ WI , the AIC icroinstabilities in the Shock Transition Region A I Ct ~ t A I C t ~ 3 t W I c k l l / w p i c k l l / w p i ck ^ / w pick ^ / w pick ^ / w pi c k l l / w p e c k l l / w p e c k l l / w p e c k l l / w p e ck ^ / w peck ^ / w peck ^ / w peck ^ / w pe c k l l / w p e ck ^ / w pe c k l l / w p e ck ^ / w pe t ~ t W I i - m i r r o r
W I W I i - m i r r o r A I C( i ) M 2 . 0( h ) M 2 . 0( g ) M 2 . 0 ( f ) M 2 . 3( e ) M 2 . 3( d ) M 2 . 3 ( c ) M 3 . 0( b ) M 3 . 0 c k l l / w p i - 5 - 4 - 3 - 2 d B(k)2/B02 ( a ) M 3 . 0
Figure 6.
Power spectra of the magnetic field fluctuations, δB y ( k ), in the period-box simulations for LM3.0 β
50 (top), LM2.3 β β
50 (bottom), plotted in the k (cid:107) - k ⊥ (that is, k z − k x ) plane. The results are shown at t ∼ τ WI (left), t ∼ τ WI (middle), and t ∼ τ AIC (right). See the text for the remarks on τ AIC for LM2.0 β
50. The gray star symbol marks the locationof the maximum linear growth rate, γ m , estimated from the linear analysis. In the models with M s ≥ .
3, AIC, whistler andi-mirror waves appear, while those waves do not grow substantially in the model with M s = 2. and i-mirror modes become dominant later on the timescale of τ AIC ≡ /γ AIC .Figure 5 shows the magnetic field fluctuations, δB y (upper panels) and δB z (lower panels), in the x - z plane (simulation plane) at three different times in theLM3.0 β
50 model. Here, the growth time scales, τ WI and τ AIC , are estimated by γ m of each mode in Ta-ble 2. At t ∼ τ WI , the transverse component, δB y ,appears on electron scales and the waves containing itpropagate parallel to B in panel (a), but the longitu-dinal component, δB z , does not grow significantly inpanel (d). In this early stage, the dominant mode is thewhistler mode, while the e-mirror mode is much weak tobe clearly manifested. As A e decreases in time due to the electron scattering off the excited waves, the whistlerwaves decay as shown in panel (b). On the time scale of τ AIC , both the AIC and i-mirror instabilities grow andbecome dominant. It is clear that the AIC-driven waves,shown in panel (c), are parallel-propagating, while thei-mirror-driven waves, shown in panel (f), are oblique-propagating; the blue arrow in the bottom-left corner ofpanel (f) denotes the wavevector of the i-mirror-drivenmode with the maximum growth rate.Figure 6 shows the time evolution of the power spec-trum for the magnetic field fluctuations, δB y ( k ), forLM2.0 β
50, LM2.3 β
50, and LM3.0 β
50 at t ∼ τ WI , t ∼ τ WI , and t ∼ τ AIC . Again, the growth time scale ofeach mode is estimated with γ m listed in Table 2, except0 Kim et al.
Table 3.
Shock Criticality of the Simulated Shock Models and Stability of the Linear Analysis ModelsSimulated Shock Model Shock Criticality Linear Analysis Model AIC WI ion-mirror electron-mirrorM2.0 β
20 sub LM2.0 β
20 stable stable stable stableM2.0 β
50 sub LM2.0 β
50 stable unstable quasi-stable unstableM2.0 β
100 sub LM2.0 β
100 quasi-stable unstable unstable unstableM2.3 β
20 super LM2.3 β
20 unstable unstable unstable stableM2.3 β
50 super LM2.3 β
50 unstable unstable unstable unstableM2.3 β
100 super LM2.3 β
100 unstable unstable unstable unstableM3.0 β β β β β
20 super LM3.0 β
20 unstable unstable unstable unstableM3.0 β
50 super LM3.0 β
50 unstable unstable unstable unstableM3.0 β
100 super LM3.0 β
100 unstable unstable unstable unstable for the LM2.0 β
50 model, in which the AIC instabilityis stable, and so the output time of panel (i) is chosenat the evolutionary stage similar to that of LM2.3 β M s = 2 . t ∼ τ WI . After the initial linear stage, the energy ofthe whistler waves is transferred to smaller wavenum-bers and the waves gradually decay, as shown in panels(b) and (e). On the time scale of ∼ τ AIC , AIC wavesand i-mirror waves appear dominantly at quasi-paralleland highly oblique angles, respectively, as shown in pan-els (c) and (f). This is consistent with the evolutionarybehavior which we have described with Figure 5. Forthe AIC and whistler instabilities, the linear predictionsfor k m with the maximum growth rate (gray star sym-bols) agree reasonably well with the peak locations ofthe magnetic power spectrum realized in the PIC simu-lations. But the linear estimates for the i-mirror modeare slightly off, because γ m is obtained without the elec-tron anisotropy, as stated in Section 2.2. In summary,the results of the periodic-box simulations are quite con-sistent with the linear predictions described earlier. Alsowe note that the results of our PIC simulations are ingood agreement with those of Ahmadi et al. (2016),in which PIC simulations were carried out to explorethe evolution of the instabilities due to the temperatureanisotropies in space plasmas with β ∼
1. The bottompanels of Figure 6 confirm that waves do not grow no-ticeably in the LM2.0 β
50 model.In these periodic-box simulations, the electron-scalewaves develop first and then decay as A e is relaxed inthe early stage, followed by the growth of the ion-scalewaves due to A i . In the shock transition region, by con-trast, temperature anisotropies are to be supplied con-tinuously by newly reflected-gyrating ions and magnetic field compression, hence multi-scale plasma waves fromelectron to ion scales are expected to be simultaneouslypresent. IMPLICATIONS FOR SHOCK SIMULATIONS4.1.
Shock Criticality
As mentioned in the introduction, the Fermi-like ac-celeration, which relies on the upstream waves excitedby the EFI, is effective only for supercritical shocks with M s ≥ M ∗ EFI ≈ . β ≈ −
100 plasma (Guo et al.2014b; Kang et al. 2019). The SSDA, which dependson the multi-scale waves excited mainly by the AIC andwhistler instabilities, is thought to occur in supercriticalshocks with M s ≥ M ∗ AIC ≈ . β ≈ M s ≥ M ∗ AIC ≈ . β ≈ − M ∗ EFI and M ∗ AIC are related to the sonic critical Mach number, M ∗ s , for ion reflection, since the structure of collisionlessshocks is governed primarily by the dynamics of shock-reflected ions.Table 3 summarizes the shock criticality of the simu-lated shock models and the stability of the linear anal-ysis models. The first column lists the name of the sim-ulated shock models considered in HKRK2021, and thetwo additional models for low- β shocks performed forthis study. The shock criticality of each model is givenin the second column. The name of the correspondinglinear analysis models is given in the third column, whilethe last four columns show the stability for the four in-stabilities (see also Table 2). We note that the nameof the shock models includes β up in the preshock, up-stream plasmas, while that of the linear analysis modelsincludes β in the shock transition zone given in Table 2.As discussed in Sections 2 and 3, in β ≈ −
100 plas-mas, the AIC instability operates for M s (cid:38) .
3, while icroinstabilities in the Shock Transition Region M s . In theM3.0 β M s = 5 and β = 5shock reported earlier by Niemiec et al. (2019). In theM3.0 β A i is smaller than that ofhigh- β models, and all the four instability modes aresuppressed by strong magnetization. This is consistentwith the results of M ∗ AIC ≈ . β ≈ Shock Surface Rippling
Another important feature of the supercritical shocksabove M ∗ AIC is the shock surface rippling. According toprevious shock simulations (e.g., Lowe & Burgess 2003;Matsukiyo & Matsumoto 2015; Niemiec et al. 2019;Trotta & Burgess 2019), the rippling has the charactersof AIC waves with the fastest growing mode at θ m ∼ ∼ λ AIC ( ≈ c/ω pi ).The parallel-propagating AIC and whistler waves inhomogeneous plasmas are purely electromagnetic andincompressible with both the electric and magnetic wavevectors pointing normal to B . The fluctuating mag-netic fields of oblique mirror modes, on the other hand,have a substantial longitudinal component, that is, δ B has a significant component parallel to B (Gary 1993).Since the density fluctuations are proportional to theparallel electric and magnetic field fluctuations (Hojoet al. 1993), we expect to see only weak ion density fluc-tuations due to the i-mirror mode in our 2D periodic-boxsimulations.Panel (b) of Figure 7 displays the variations of theion density, [ (cid:104) n i − n (cid:105) x , avg /n ] ≈ ± .
01, averaged overthe x -direction in the periodic-box simulation for theLM3.0 β
50 model. Panel (a) shows the fluctuations ofthe transverse component of B , [ (cid:104) B y − B (cid:105) x , avg /B ] ≈± .
4, which have a relatively large amplitude due to theAIC-driven waves. It shows that even after the AIC-driven waves have fully grown, they have little effectson the ion density fluctuations.However, the ion density fluctuations of the ripplingwaves propagating along the shock surface behind theshock ramp are rather significant in the shock simulationfor the M3.0 β
50 model in HKRK2021. Panel (c) showsthat both the variations of [ (cid:104) n i − n (cid:105) x , avg /n ] ≈ ± . (cid:104) B y − B (cid:105) x , avg /B ] ≈ ± . n i are much larger than those of the lin-ear prediction expected for the parallel-propagating AICmode. Note that here the quantities are averaged alongthe x -direction over the shock transition zone includingthe first and second overshoot oscillations behind the t ~ 4 6 W c i - 1 t ~ 7 9 W c i - 1 t ~ 1 3 0 W c i - 1 z [ c / w p i ]z [ c / w p i ] t ~ 4 6 W c i - 1 t ~ 7 9 W c i - 1 t ~ 1 3 0 W c i - 1 s h o c k s i m u l a t i o n ( M 3 . 0 b W c i - 1 ) y [ c / w p i ] ( a ) < B y - B > x , a v g / B ( b ) < n i - n > x , a v g / n ( c ) < n i - n > x , a v g / n , < B y - B > x , a v g / B Figure 7. (a)-(b): Variations in the transverse compo-nent of magnetic field, (cid:104) B y − B (cid:105) x , avg /B , and the ion-density, (cid:104) n i − n (cid:105) x , avg /n , averaged over the x -domain inthe 2D periodic-box simulation for the LM3.0 β
50 model,plotted along B ( z -direction) at three different times. (c):Variations in the longitudinal component of magnetic field, (cid:104) B y − B (cid:105) x , avg /B (red), and the ion-density, (cid:104) n i − n (cid:105) x , avg /n (black), averaged along the x -direction over the shock transi-tion zone in the 2D shock simulation for the M3.0 β
50 modelin HKRK2021, plotted along the y -direction. Note that thepreshock magnetic field, B up0 , lies in the x - y plane, and theobliquity angle between B up0 and the y -axis is θ Bn = 63 ◦ inthe shock simulation. ramp. Hence, the basic assumptions of the linear the-ory, such as the homogeneous background, charge neu-trality, and zero net-current, are likely to be violated inthis region.We point that such large-amplitude fluctuations of n i ,comparable to the fluctuations of B y , were previouslyrecognized in the 2D hybrid simulations of supercriti-cal, perpendicular shocks presented in Winske & Quest(1988). The authors suggested that large compressivewaves might result from nonlinear effects in addition to2 Kim et al. oblique i-mirror modes. The effects due to nonlinearcouplings between various wave modes could be signifi-cant as well (e.g. Shukla & Stenflo 1985; Verscharen &Marsch 2011; Marsch & Verscharen 2011). Therefore,the pure AIC-driven waves in the shock transition couldhave been modified by such possible nonlinearities, lead-ing to the enhancement of ion density fluctuations. SUMMARYIn supercritical Q ⊥ -shocks, a substantial fraction ofincoming ions and electrons are reflected, and the trans-verse components of magnetic fields are amplified atthe shock ramp. The reflected-gyrating ions and theamplified magnetic fields induce the ion and electronperpendicular anisotropies, A i and A e , respectively, inthe shock transition region (Guo et al. 2017). Theyin turn trigger various microinstabilities including theAIC, whistler, i-mirror, and e-mirror instabilities (e.g.,Winske & Quest 1988; Lowe & Burgess 2003; Guo et al.2017). The kinetic properties of these four instabili-ties are summarized in Table 1. The multi-scale plasmawaves generated by these microinstabilities are thoughtto be crucial for the electron preacceleration via theSSDA (e.g., Katou & Amano 2019; Niemiec et al. 2019;Trotta & Burgess 2019).In this work, adopting the numerical estimates for theion and electron temperature anisotropies found in the2D PIC simulations of Q ⊥ -shocks with M s = 2 − β = 1 −
100 and m i /m e = 50 − γ m , k m , θ m , of eachinstability are given in Table 2. In addition, in orderto investigate the development and nonlinear evolutionof the waves induced by the microinstabilities, we haveperformed the 2D PIC simulations with periodic bound-ary conditions for the three fiducial models, LM2.0 β β
50, and LM3.0 β
50. Finally, the results were alsocompared with the 2D PIC simulations for ICM shocksreported in HKRK2021.The main results can be summarized as follows:1. In the LM3.0 β γ WI (cid:29) γ E − M (cid:29) γ AIC > γ I − M (Fig. 2). Hence, the parallel-propagating AIC andwhistler waves are expected to be more dominant thanthe oblique-propagating mirror waves.2. In the LM2.0 β
50 model, which represents a subcrit-ical ICM shock, by contrast, the AIC mode is stable(Table 2), so mainly the electron-scale whistler wavesare generated. 3. The maximum growth rates for the AIC and e-mirrorinstabilities, γ AIC and γ E − M , are almost independent of m i /m e , while γ WI for the whistler instability and γ I − M for the i-mirror instability are slightly lower for smaller m i /m e (Fig. 3).4. For all the four instabilities, the maximum growthrates increase with increasing β (Fig. 3).5. As the sonic Mach number M s increases, both A e and A i increase, all the modes grow faster, and k m ofeach mode shifts towards larger k , regardless of β in therange of β ≈ −