Magnetic Field Amplification by Turbulence in A Relativistic Shock Propagating through An Inhomogeneous Medium
Yosuke Mizuno, Martin Pohl, Jacek Niemiec, Bing Zhang, Ken-Ichi Nishikawa, Philip E. Hardee
aa r X i v : . [ a s t r o - ph . H E ] N ov Magnetic Field Amplification by Turbulence in A Relativistic ShockPropagating through An Inhomogeneous Medium
Yosuke Mizuno , , Martin Pohl , , Jacek Niemiec , Bing Zhang , Ken-Ichi Nishikawa , , andPhilip E. Hardee ABSTRACT
We perform two-dimensional relativistic magnetohydrodynamic simulations of amildly relativistic shock propagating through an inhomogeneous medium. We showthat the postshock region becomes turbulent owing to preshock density inhomogeneity,and the magnetic field is strongly amplified due to the stretching and folding of fieldlines in the turbulent velocity field. The amplified magnetic field evolves into a fila-mentary structure in two-dimensional simulations. The magnetic energy spectrum isflatter than the Kolmogorov spectrum and indicates that a so-called small-scale dynamois occurring in the postshock region. We also find that the amount of magnetic-fieldamplification depends on the direction of the mean preshock magnetic field, and thetime scale of magnetic-field growth depends on the shock strength.
Subject headings: gamma-ray burst: general - MHD - methods: numerical - relativisticprocesses - shock waves - turbulence
1. Introduction
Relativistic shocks accelerate plasma to ultrarelativistic velocities. Emission from such shocksis observed in astrophysical sources on many length, time and energy scales, and is typically at-tributed to the synchrotron process. The composition of the plasma, the magnitude of the preshockmagnetic field, and the Lorentz factor of the shock are generally unknown. One of the best-studiedsystems is the radiative afterglow that follows the bright flash of γ -rays in gamma-ray bursts Center for Space Plasma and Aeronomic Research, University of Alabama in Huntsville, 320 Sparkman Drive,NSSTC, Huntsville, AL 35805, USA; [email protected] National Space Science and Technology Center, VP62, Huntsville, AL 35805, USA Institut fur Physik und Astronomie, Universit¨at Potsdam, 14476 Potsdam-Golm, Germany DESY, Platanenallee 6, 15738 Zeuthen, Germany Institute of Nuclear Physics PAN, ul. Radzikowskiego 152, 31-342 Krak´ow, Poland Department of Physics and Astronomy, University of Nevada, Las Vegas, NV 89154, USA Department of Physics and Astronomy, The University of Alabama, Tuscaloosa, AL 35487, USA ǫ B ∼ − − − of the internal energy density (e.g., Panaitescu & Kumar 2002, Yost etal. 2003; Panaitescu 2005). However, simple compressional amplification of the weak pre-existingmagnetic field of the circumburst medium can not account for such high magnetization (Gruzinov2001).The leading hypothesis for field amplification in GRB afterglows is the relativistic Weibelinstability that produces filamentary currents aligned with the shock normal; these currents are re-sponsible for the creation of transverse magnetic fields (e.g., Gruzinov & Waxman 1999; Medvedev& Loeb 1999). Recent plasma simulations using the particle-in-cell (PIC) method have demon-strated magnetic-field generation in relativistic collisionless shocks (e.g., Nishikawa et al. 2003,2009; Spitkovsky 2008). However, the size of the simulated regions is orders of magnitude smallerthan the GRB emission region. It remains unclear whether magnetic fields generated on scales oftens of plasma skin depths will persist at sufficient strength in the entire emission region, some10 plasma skin depths in size. On the other hand, in magnetohydrodynamic (MHD) processes, ifthe density of the preshock medium is strongly inhomogeneous, significant vorticity is produced inthe shock transition that stretches and deforms magnetic field lines leading to field amplification(Goodman & MacFadyen 2007; Sironi & Goodman 2007; Palma et al. 2008). For the long durationGRBs, i.e., associated with the iron core collapse of mass-losing very massive stars, density fluc-tuations in the interstellar medium may arise through several processes (e.g., Sironi & Goodman2007).There is also a direct observational motivation for relativistic turbulence in GRB outflows.Significant angular fluctuations by relativistic turbulence have been invoked to explain the largevariation of the γ -ray luminosity among bursts as well as the intraburst variability observed inmany afterglows. The afterglow polarization indicates a breaking of axial symmetry; the correlationof this polarization with large-amplitude variability suggests blastwave anisotropy as the source ofvariability. Narayan & Kumar (2009) and Lazar et al. (2009) have proposed a relativistic turbulencemodel instead of the well-known internal shock model as the production mechanism for variableGRB light curves and applied it to GRB 080319B (Kumar & Narayan 2009). Zhang & Yan (2010)have developed a new GRB prompt-emission model in the highly magnetized regime, namely, aninternal-collision-induced magnetic reconnection and turbulence (ICMART) model. Within thismodel the short-time variability “spikes” are related to turbulence in the magnetic dissipationregion while the broad variability component in the GRB lightcurves is related to central engineactivity.The innermost parts of the relativistic jets seen in Active Galactic Nuclei (AGN) are exploredby multi-wavelength observational campaigns that provide probes into the structure of the blazarzone (e.g., Marscher et al. 2008, 2010). Recently correlated X-ray/TeV gamma-ray flares with 3 –timescales from 15 minutes (Mrk 421) to a few hours (Mrk 501 and 1ES 1959+650) have beenobserved (Aharonian et al. 2003; Albert et al. 2007; Krawczynski et al. 2004). The fast variableflares may come from small regions, a few Schwarzschild radii in size. Several authors have proposedfast-moving needles within a slower jet or a jet-within-a-jet scenario to explain the fast variabilityof blazars (Levinson 2007; Begelman et al. 2008; Ghisellini & Tavecchio 2008; Giannios et al. 2008).On the other hand, Marscher & Jorstad (2010) have proposed that some short-term fluctuationscan be understood as the consequence of a turbulent ambient-jet-plasma that passes through shocksin the jet flow (Marscher et al. 1992).Supernova remnants (SNRs) are bound by an expanding nonrelativistic spherical blastwave.The synchrotron emission from SNRs is generally consistent with compression of the interstellarfield. However, recent discovery of the year-scale variability in the synchrotron X-ray emission ofSNRs suggests that the magnetic field needs to be amplified in the SNR up to the milli-Gauss level(Uchiyama et al. 2007, Uchiyama & Aharonian 2008), although the interpretation is disputed (e.g.,Bykov et al. 2008, Huang & Pohl 2008). The evidence for magnetic field amplification in SNRs hasbeen obtained also from the thickness of X-ray filaments that typically indicate ∼
100 micro-Gaussmagnetic fields (Bamba et al. 2003, 2005a, 2005b; Vink & Laming 2003, but see also Pohl et al.2005). Since the typical magnetic-field strength in the interstellar medium is on the order of a fewmicro-Gauss, amplification beyond simple shock compression is necessary to achieve a milli-Gausslevel magnetic field in SNRs.Recently, Giacalone & Jokipii (2007) have performed non-relativistic MHD shock simulationsincluding density fluctuations with a Kolmogorov power spectrum in the preshock medium. Theyobserved a strong magnetic-field amplification caused by turbulence in the postshock medium;the final rms magnetic-field strength is reportedly a hundred times larger than the preshock fieldstrength. Density fluctuations with a Kolmogorov power spectrum would be present in the preshockmedium from the transonic turbulence expected in the diffuse interstellar medium (ISM) (e.g., Hen-nebelle & Audit 2007, Inoue & Inutsuka 2008, Inoue et al. 2009). Inoue et al. (2009) have alsoobtained strong magnetic-field amplification by turbulence associated with a strong thermal in-stability driven shock wave propagating through an inhomogeneous interstellar two-phase medium.The interaction of SNRs with a turbulent, magnetized ISM has been studied by Balsara et al. (2001)via three-dimensional MHD simulations. They found that the interaction of strong shocks with theturbulent density fluctuations could be a strong source of helicity generation which provides a po-tential source of magnetic field amplification. Kim et al. (2001) studied the energetics, structure,and spectra of a supernova (SN) driven turbulent medium with higher resolution simulations. Bal-sara et al. (2004) have investigated in detail the role of magnetic field amplification in SN-driventurbulence and found that the helicity generation mechanism by SN-driven, turbulent, multiphaseISM processes was adequate to sustain the magnetic field growth (also Balsara & Kim 2005; Kim& Balsara 2006). Haugen et al. (2004) have performed simulations and found that magnetic fieldamplification can be sustained by shock-driven turbulence in an isothermal, compressible, shock-driven medium. All the above simulations are relevant to SNR shocks, which are non-relativistic. 4 –GRBs and AGN jets have relativistic shocks. Similar physical processes should apply, so that theseshocks would also experience strong magnetic-field amplification by turbulence.The magnetic-field amplification mechanism is by a so-called turbulent dynamo or small-scaledynamo, also known as a fast dynamo, and is different from the classical mean-field dynamo mech-anism. Mean-field dynamo theory (e.g., Steenbeck et al. 1966; Parker 1971; Ruzmaikin et al. 1988,Diamond et al. 2005) has been used to study the growth of the large-scale magnetic field in ourGalaxy. In general, mean-field dynamos are slow, involve almost incompressible motions, and thegrowth rate decreases with decreasing resistivity ( ν ). Kulsrud & Anderson (1992) pointed out thatthe magnetic Reynolds number (Re m ≡ ℓv/ν , where ℓ is typical length of flow) in our Galaxyis very large (Re m ∼ ), and such a high Re m leads to the development of turbulent velocityand magnetic-field fluctuations on small scales. While the growth of large-scale magnetic fieldsis probably driven by a mechanism similar to a mean-field dynamo, other mechanisms are likelyresponsible for the more rapid generation of field on smaller scales. The theory of fast dynamos hasbeen developed to investigate the rapid growth of magnetic field on smaller scales in the limit ofinfinite magnetic Reynolds number, and the growth rate does not depend on resistivity (e.g., Chil-dress & Gilbert 1995; Galloway 2003; Brandenburg & Subamanian 2005). Here the amplificationof the magnetic field arises from sequences of a “stretch, twist, and fold” nature. While the limit(Re m = ∞ ) is unrealistic in astrophysical simulations, published results of incompressible MHDsimulations (e.g., Schekochihin et al. 2004) and compressible MHD simulations (e.g., Balsara et al.2004; Balsara & Kim 2005; Kim & Balsara 2006) are nevertheless in agreement with the existenceof a small-scale turbulent driven dynamo.Here we report on magnetic-field amplification by turbulence in two-dimensional relativisticmagnetohydrodynamic (RMHD) simulations of a mildly relativistic shock wave propagating throughan inhomogeneous medium. In this paper we focus on the small-scale growth of the magnetic field.This paper is organized as follows: We describe the numerical method and setup used for oursimulations in §
2, present our results in §
3, and discuss the astrophysical implications in §
2. Numerical Method and Setup
To study the time evolution of a mildly relativistic shock propagating in an inhomogeneousmedium, we use the 3D GRMHD code “RAISHIN” in two-dimensional Cartesian geometry ( x − y plane). RAISHIN is based on a 3 + 1 formalism of the general relativistic conservation laws ofparticle number and energy-momentum, Maxwell’s equations, and Ohm’s law with no electricalresistance (ideal MHD condition) in a curved spacetime (Mizuno et al. 2006) . In the RAISHIN Constained transport schemes are used to maintain divergence-free magnetic field in the RAISHIN code. Thisscheme requires the magnetic field to be defined at the cell interfaces. On the other hand, conservative, high-resolutionshock capturing schemes (Godonov-type schemes) for convservation laws require the variables to be defined at thecell center. In order to combine variables defined at these different positions, the magnetic field calculated at the cell r into a weighted average of order 2 r − ρ and containing fluctuations δρ is established across the whole simulation region and uniformlyflows in the positive x -direction with speed v , where ρ is an arbitrary normalization constant(our simulations are scale free) and we choose ρ = 1 .
0. The density fluctuations are generatedaccording to a two-dimensional Kolmogorov-like power-law spectrum of the form P k ∝
11 + ( kL ) / , (1)where k is the wave number and L is the turbulence coherence length.A method used to generate fluctuations with a Kolmogorov-like turbulence spectrum wasproposed in Giacalone & Jokipii (2007) (also Giacalone & Jokipii 1999). Following their method,our turbulence is described by summing over a large number of discrete wave modes, δρ ( x, y ) = N m X n =1 A ( k n ) exp[ i ( k n cos θ n x + k n sin θ n y + φ n )] , (2)where A ( k n ) is the amplitude, θ n is the direction of the wavevector with magnitude k n , and φ n isthe phase. For each n a phase and a direction are randomly selected from the ranges 0 < φ n < π and 0 < θ n < π , respectively. The total number of modes is N m = 50. The amplitude in theequation above is given by A ( k n ) = a G ( k n ) " N m X n =1 G ( k n ) − , (3)and G ( k n ) = 2 πk n ∆ k n P ( k ) , (4) interfaces is interpolated to the cell center and as a result the scheme becomes non-conservative even though we solvethe conservation laws (Komissarov 1999). a is the wave variance. In the simulations we choose a = 0 .
02, which makes a fluctuationvariance of p < δρ > = 0 . ρ . We choose the maximum wavelength as λ max = 0 . L , one-half ofthe size of the simulation box in the y -direction, and the minimum wavelength as λ min = 0 . L .We consider a low gas-pressure medium with constant p = 0 . ρ c , where c is the speed oflight. The equation of state is that of an ideal gas with p = (Γ − ρe , where e is the specific internalenergy density and the adiabatic index Γ = 5 /
3. The specific enthalpy is h ≡ e/c + p/ρc .The pre-shock plasma carries a weak constant magnetic field. The magnetic field amplitude is B = 4 . × − p πρ c which leads to a high plasma- β ( β = p gas /p mag = 10 ). To investigatethe effect of the magnetic field direction with respect to the shock propagation direction, we choosetwo different directions, magnetic field parallel ( B x , θ Bn = 0 ◦ ) and perpendicular ( B y , θ Bn = 90 ◦ )to the shock normal. Using the mean plasma density ( ρ ), the sound speed c s /c = (Γ p/ρh ) / andthe Alfv´en speed v A /c = [ B / ( ρh + B )] / are c s = 0 . c and v A = 0 . c .The size of the computational domain is ( x, y ) = (2 L, L ) with an
N/L = 256 grid resolution.We impose periodic boundary conditions in the y -direction. In order to create a shock wave, arigid reflecting boundary is placed at x = x max . The fluid, which is moving initially in the positive x -direction, is stopped at this boundary, where the velocity v x is set to zero. As the density buildsup at x = x max , a shock forms and propagates in the − x -direction. In the frame of referenceof the simulation (contact-discontinuity frame), the plasma velocity downstream of the shock iszero on average. New fluid continuously flows in from the inner boundary ( x = 0) and densityfluctuations are advected with the flow speed. In the frame of reference of the inflowing fluid thedensity fluctuation pattern is maintained at each time step using eq.(2) but extended in the -xdirection at the flow speed, v , so that density fluctuations are inserted in the inflowing materialat the inner simulation boundary. We choose two different flow speeds, v x = v = 0 . c and 0 . c .The various different cases that we have considered are listed in Table 1.Table 1. Models and ParametersCase v θ Bn t max /c N/L A 0.4 0 10 256B 0.4 90 10 256C 0.9 0 5 256D 0.9 90 5 256 7 –
3. Results3.1. Shock structure
Figures 1 and 2 show 2-D images of the density and the total magnetic field strength at t s = 5 . t s = 10 .
0, where t s is in units of L/c with c = 1, for the parallel ( B x ) and per-pendicular ( B y ) magnetic field cases with slow flow velocity v = 0 . c When the preshock plasmaFig. 1.— Two-dimensional images of density ( upper ) and total magnetic field strength ( lower ) at t s = 5 . left ) and t s = 10 . right ) for case A (magnetic field parallel to the shock propagationdirection and v = 0 . c ). White arrows indicate the flow direction in the postshock region.with inhomogeneous density encounters the shock, the shock front is rippled, leading to significant,random transverse flow behind the shock. These ripples in the shock front introduce rotation andvorticity in the postshock region through a process similar to the Richtmyer-Meshkov instability(e.g., Brouillette 2002). This instability is excited when a sudden density jump encounters a shockand is a possible source of vorticity generation (e.g., Zabusky 1999) . Our simulation starts frompreexisting fluctuations having a broad range of scales in the preshock medium. The flow patternin the postshock region is initially highly nonlinear and comparison with linear instability analysisis not useful.Since the preexisting magnetic field is much weaker than the postshock turbulence, the tur-bulent velocity field can easily stretch and deform the frozen-in magnetic field, resulting in its Samtaney & Zabusky (1994) calculated scaling relations for vorticity generation via the Richtmyer-Meshkovinstability when shocks interact with density inhomogeneities in two dimensions and found that increasing shockMach number or density contrast enhances the generation of vorticity at the shocks. upper ) and total magnetic field strength ( lower ) at t s = 5 . left ) and t s = 10 . right ) for case B (magnetic field perpendicular to the shock propagationdirection and v = 0 . c ). White arrows indicate the flow direction in the postshock region.distortion and amplification. This creates regions with higher magnetic field intensity. Near theshock front, the vorticity scale size is small, but far away from the shock front the vorticity scalesize becomes larger through an inverse cascade of vortices, and the magnetic field is strongly am-plified with time. The amplified magnetic field evolves into a filamentary structure. We note thatthe magnetic-field amplification under study here does not increase the total magnetic flux, onlythe magnetic-field strength. The average turbulent velocity in the postshock region is ∼ . c at t s = 10 in both cases A and B. The average sound speed and Alfv´en velocity in the postshockregion at t s = 10 are ∼ . c and ∼ . c , respectively in case A and ∼ . c and ∼ . c ,respectively in case B. Thus the turbulent velocity is subsonic and super-Alfv´enic in most of thepostshock region in both cases A and B. This result is consistent with previous non-relativisticstudies (e.g., Balsara et al. 2001; Kim et al. 2001; Balsara et al. 2004; Kim & Balsara 2006; Hauganet al. 2004, Giacalone & Jokipii 2007; Inoue et al. 2009), although the turbulent velocity in thepostshock region is locally supersonic in these previous non-relativistic simulations.Figure 3 shows one-dimensional profiles along the x -axis at y/L = 0 . t s = 10 for casesA and B, involving slow flow velocity and parallel or perpendicular magnetic field. The density,normalized to ρ , the transverse velocity ( v y ), and the total magnetic field strength are plottedfrom top to bottom. The shock front is located at x = 0 . L . The left and right sides of theshock front are upstream and downstream, respectively. The shock propagation speed is around v sh ∼ . c , in good agreement with analytical calculations (see the Appendix). We parameterizethe shock strength by the relativistic Mach number M s ≡ γ ′ sh v ′ sh /γ s c s , where γ ′ sh and v ′ sh arethe shock propagation speed and the shock Lorentz factor measured in the inflow frame (see the 9 –Fig. 3.— One-dimensional profile of ( panels a, b ) the density normalized by the initial mean density( ρ ), ( panels c, d ) the transverse velocity ( v y ), and ( panels e, f ) the total magnetic field strength( B tot ) normalized by the initial magnetic field ( B ), all shown along the x -direction at y/L = 0 . t s = 10 . left , a parallel initial magnetic field and v = 0 . c ), and case B ( right , aperpendicular initial magnetic field and v = 0 . c ). The shock front is located at x/L = 0 .
31. Theleft and the right sides of shock front are upstream and downstream regions, respectively.Appendix) and γ s ≡ (1 − c s ) − / is the Lorentz factor associated with the sound speed. Theshock propagation speed obtained from the simulations and the sound speed in the preshock regionleads to M s ∼ . . c in both cases, which is subsonic in thepostshock region ( < c s > ≃ . c in slow flow cases). The total-magnetic-field profiles also showstrong fluctuations. In the perpendicular field case, the magnetic field is compressed by nearly afactor 3.5 behind the shock front, but further downstream the magnetic field locally reaches morethan 10 times the amplitude of the initial field. As a result of the density fluctuations the magneticfield compression behind the shock front is not identical to the density jump. In the parallel-fieldcase, magnetic compression at the shock does not occur. Instead, the magnetic field is amplifiedpurely by turbulent motion, locally reaching more than 5 times the initial field amplitude.Figure 4 shows a 2-D image of the density and the total magnetic-field strength at t s = 4 . B x ) and perpendicular ( B y ) magnetic field cases with fast flow velocity, v = 0 . c .When the flow speed increases, the shock becomes stronger. The shock propagation speed for v = 0 . c is more than two times faster than that in v = 0 . c cases. In the postshock region, 10 –Fig. 4.— Two-dimensional images of the density ( left ) and the total magnetic-field strength ( right )at t s = 4 . upper , parallel magnetic field, v = 0 . c ) and D ( lower , perpendicularmagnetic field, v = 0 . c ). White arrows indicate the flow direction in the postshock region.turbulent motion is created by ripples in the shock front, and the magnetic field is deformed andamplified by the turbulence. Filamentary magnetic fields are seen also in these cases. The averageturbulent velocity in the postshock region is ∼ . c at t s = 4 . t s = 4 . ∼ . c and ∼ . c incase C and ∼ . c and ∼ . c in case D. As in the slower velocity cases, the turbulent velocityis subsonic, but super-Alfv´enic in most of the postshock region.One-dimensional profiles along the x -axis for cases C and D (parallel or perpendicular initialfield with fast flow) are displayed in Figure 5. The shock front is located at x = 0 . L , and theshock speed is around v sh ∼ . c . This is in good agreement with analytical estimates and about2.5 times faster than that of v = 0 . c cases. The relativistic Mach number is M s ∼ . . c in both cases. This is also a subsonic speed in thepostshock region ( < c s > ≃ . c in fast flow cases). The total magnetic field also shows substantialfluctuations. In the case of an initially perpendicular magnetic field (D), the field is compressedby nearly a factor of 3 at the shock front. The field compression behind the shock front is notidentical to the density jump due to the density fluctuations. Further downstream the field isstrongly amplified, reaching more than 15 times the initial field strength. In case C, involving aninitially parallel magnetic field, the field is amplified purely by turbulence and locally reaches morethan 8 times the initial amplitude. For a fast flow, the field amplification is stronger than in thecase of a slow flow, because the stronger shock produces faster vortices in the postshock region. 11 –Fig. 5.— One-dimensional profiles along the x -direction at y/L = 0 . t s = 4 . panels a, b )the density normalized by ρ , ( panels c, d ) the transverse velocity ( v y ), and ( panels e, f ) the totalmagnetic-field strength ( B tot ) normalized by B . Results for case C, involving a parallel magneticfield and v = 0 . c , are shown on the left, whereas plots for case D, a perpendicular magnetic fieldand v = 0 . c , are displayed on the right.Because the turnover time of vortices is shorter, magnetic-field amplification is faster. Thus, themagnetic-field amplification time scale depends on the shock strength. In order to understand the statistical properties of the turbulent fluctuations in the postshockregion, it is helpful to observe their spectra. We define the kinetic energy spectrum E kin ( k ) as Z E kin ( k ) dk = h ργ ( γ − i , (5)where γ is the Lorentz factor and the electromagnetic energy spectrum E mag ( k ) as Z E mag ( k ) dk = h T EM i , (6)where T EM is the 00-component of the energy-momentum tensor for the electromagnetic field(Zhang et al. 2009). Figures 6 and 7 show spherically-integrated spectra of the kinetic and electro-magnetic energy for slow and fast flow parallel and perpendicular initial field cases, respectively. 12 –Fig. 6.— Spherically-integrated spectra of ( a, b ) the kinetic energy, ( c, d ) the electromagneticenergy, and ( e, f ) the ratio of electromagnetic and kinetic energy spectra for the case A ( left ,parallel magnetic field, v = 0 . c ) and B ( right , perpendicular magnetic field, v = 0 . c ). Differentlines are for different times: t s = 4 ( solid ), 6 ( dashed ), 8 ( dotted ), and 10 ( dash-dotted ). A dottedline, representing the power law E ( k ) ∼ k − / , is shown for comparison in the upper and middlepanels ( a, b, c, d ).The kinetic-energy spectra almost follow a Kolmogorov spectrum in all cases, E kin ( k ) ∝ k − (5 / − ( D − with D = 2 in two-dimensional systems. Initially the density in the preshock re- 13 –Fig. 7.— Spherically-integrated spectra of ( a, b ) the kinetic energy, ( c, d ) the electromagneticenergy, and ( e, f ) the ratio of electromagnetic and kinetic energy spectra for the case C ( left ,parallel magnetic field, v = 0 . c ) and D ( right , perpendicular magnetic field, v = 0 . c ). Differentlines are for different times: t s = 2 ( solid ), 3 ( dashed ), and 4 ( dotted ). Again, a dash-dotted linerepresenting the power law E ( k ) ∼ k − / is shown for comparison.gion is inhomogeneous with a Kolmogorov-like power spectrum, and this density spectrum stillexists in the postshock region. The slope and magnitude of the kinetic-energy spectra do notchange with time. 14 –The electromagnetic-energy spectrum increases over time in all cases, implying that magnetic-field amplification is ongoing and not yet saturated. However, the overall shape of the spectra arevery similar. The magnetic energy spectra are almost flat and strongly deviate from a Kolmogorovspectrum. Spectra flatter than a Kolmogorov spectrum are typical of a small-scale dynamo (e.g.,Childress & Gilbert 1995; Galloway 2003; Brandenburg & Subamanian 2005). In a small-scale dy-namo, a forward cascade of magnetic energy from large scales to intermediate scales, and an inversecascade from small scales to intermediate scales, are introduced. Therefore, the electromagnetic-field spectrum is flatter than Kolmogorov. A flat magnetic-energy spectrum is generally seen inturbulent-dynamo simulations for incompressible MHD (e.g., Schekochihin et al. 2004) and com-pressible MHD (e.g., Balsara et al. 2004; Balsara & Kim 2005; Kim & Balsara 2006). The sameproperties are also observed in relativistic MHD turbulence induced by the Kelvin-Helmholtz in-stability (Zhang et al. 2009).The ratio of the electromagnetic and kinetic energy spectra is shown in the lower panelsof Figures 6 and 7. This ratio increases with time in all cases. The electromagnetic energy islarge at small scales but the kinetic energy remains dominant at all scales, because magnetic-fieldamplification by turbulence has not yet saturated. The ratio of the electromagnetic and kineticenergy spectra will probably increase up to equipartition ( ∼ Figure 8 shows the time evolution of the magnetic field, both as a volume-averaged (mean)total magnetic field and as a peak total magnetic-field strength in the postshock region (whose sizeincreases with time as the shock propagates away from the wall), in all cases normalized by B .Note that the mean magnetic-field strength is still increasing when the simulation was terminated;apparently the magnetic-field strength takes some time to saturate. The mean postshock magneticfield is stronger for the perpendicular shock case ( B y ), for which we observe an increase by morethan a factor of 5 compared to B . For the parallel shock case ( B x ) an amplification by abouta factor of 2 is observed. In the perpendicular shock case, the perpendicular magnetic field iscompressed by a factor of 3 at the shock, and additional magnetic field amplification by turbulentmotion in the postshock region is almost the same as for the parallel magnetic field case (see inFig. 9). In both cases, the peak field strength is much larger than the mean magnetic field, about13 and 26 times for the slow flow velocity and 12 and 49 times for the fast velocity, where the twonumbers are for the parallel and perpendicular cases, respectively. The time scale for magnetic-fieldamplification in the fast-flow case is about half that in the slow-flow case, because vortices withfaster velocity can be produced.To investigate the amount of magnetic-field amplification via turbulent motion, the time evo- 15 –Fig. 8.— Time evolution of ( a, c ) the volume-averaged total magnetic field and ( b, d ) the maximumtotal magnetic-field strength in the postshock region normalized by B for the case of lower flowvelocity ( v = 0 . c : upper panels ) and higher velocity ( v = 0 . c : lower panels ). Different lines arefor different initial magnetic field orientation: parallel ( B x : solid ) and perpendicular ( B y : dashed )with respect to the shock propagation direction.lution of the mean fluctuation amplitudes of δB x and δB y in the postshock region normalizedby the mean magnetic-field strength ( B s ≃ . × − for the parallel-field case, ≃ .
014 for theperpendicular-field case) is shown in Figure 9. In all cases, the fluctuation amplitude increases withtime, showing no sign of saturation. Each fluctuating magnetic field component grows similarly,implying that the turbulence amplifies each magnetic field component equally. The growth rateis independent of the original magnetic-field direction. Therefore, the amplification of magneticfield by turbulent motion does not depend on the direction of a homogeneous magnetic field. Weshould note that this magnetic field fluctuation amplitude is normalized by the mean magneticfield in the postshock region, which is about 3 times larger in the perpendicular-field case than inthe parallel-field case on account of magnetic field compression by the shock. Therefore, the totalmagnetic-field amplification is larger for a perpendicular homogeneous field. 16 –Fig. 9.— Time evolution of the volume-averaged root mean square (rms) fluctuation amplitudes of δB x ( solid lines ) and δB y ( dashed lines ) in the postshock region normalized by the mean magneticfield strength in postshock region for the cases A, B, C, and D.
4. Summary and Discussion
We have performed two-dimensional relativistic MHD simulations of the propagation of amildly relativistic shock through an inhomogeneous medium. We have shown that the postshockregion becomes turbulent owing to the preshock density inhomogeneity, and the magnetic field isstrongly amplified by the turbulent motion in the postshock region. The amplified magnetic fieldevolves into filamentary structures in these two-dimensional simulations. This is consistent withearlier non-relativistic studies (e.g., Balsara et al. 2001, 2004; Giacalone & Jokipii 2007; Inoue etal. 2009). The magnetic energy spectrum is flatter than Kolmogorov, which is typical for a small-scale dynamo. We also found that the total magnetic-field amplification from the preshock valuedepends on the direction of the homogeneous magnetic field, and the time scale of magnetic fieldgrowth depends on the shock strength. The mean magnetic field strength in the postshock region isstill increasing when the simulations were terminated. Therefore, longer simulations with a longersimulation box are needed to follow the magnetic field amplification to saturation.In this paper we have performed simulations in two-dimensional geometry as a first step towardunderstanding magnetic-field amplification by turbulence in the relativistic MHD regime. However,in general the turbulence structure, for example the slope of the Kolmogorov-like power spectrum, isdifferent in two and three dimensions. To more realistically analyze three-dimensional phenomena,we will extend the current investigation to three-dimensional simulations in future work. 17 –Recently, the Fermi satellite observed GRB 080916C which showed a featureless smoothly-joint broken power-law spectra covering 6-7 decades in energy (Abdo et al. 2009). Zhang & Pe’er(2009) analyzed this burst in detail and argued that the non-detection of a thermal componentstrongly suggests that the outflow is Poynting-flux dominated with magnetization σ &
15 (thePoynting-to-kinetic energy ratio) at the central engine. Therefore, at least some GRB outflowsappear to be magnetically dominated. In our simulations, we assume the preshock medium isweakly magnetized. If we consider internal shocks in GRBs, which are internal collisions withinrelativistically propagating shells, the magnetization in the preshock medium can cover a wide rangeof values. In particular, it is expected that in the strongly-magnetized (high σ ) regime the MHDturbulence is anisotropic and has a different scaling in the directions along and perpendicular tothe homogeneous field, E ( k ⊥ ) ∝ k − / ⊥ (Kolmogorov-type) and E ( k k ) ∝ k − k (Goldreich & Sridhar1995; Cho & Lazarian 2003). While the kinetic power in the turbulence quickly decreases with k ,the power in the magnetic field does not significantly decrease with k . As a result, eddies of smallerscales are more stretched and appear elongated along the magnetic field. Current simulations donot show any anisotropy because the magnetic field in the turbulent postshock medium is still weak.MHD turbulence in the highly-magnetized (high σ ) regime is not well studied. A comprehensivestudy of relativistic turbulence with a wide range of magnetization is now underway and will bepresented in future publications.In the internal shock model, the GRB central engine ejects an unsteady outflow that canbe visualized as discrete magnetized shells with variable Lorentz factor and mass loading. Thesemagnetized shells collide with each other at the conventional internal shock radius, R IS ≃ γ cδt ≃ γ . δt cm, where γ . = γ/ . . Hereafter the convention Q s = Q/ s . The gamma-rayemission radius R can be in principle in the range of 10 -10 cm (between the photosphere radiusand the external shock deceleration radius). Various observational constraints suggest that thenon-thermal emission radius is at R ≥ cm (e.g., Kumar et al. 2007; Racusin et al. 2008; Shen& Zhang 2009; Zhang & Pe’er 2009). The thickness of the discrete shell is ∆ ′ ≃ R/γ ≃ R γ − . cm in the comoving frame. From a turbulence model of GRBs (Narayan & Kumar 2009; Lazar etal. 2009), the size of an eddy in the comoving frame is R e ≃ f /a R/γγ t , where f is a filling factorand γ t is a random Lorentz factor for an energy-bearing eddy. Lazar et al. (2009) suggest that a is between 2 and 3. Narayan & Kumar (2009) have pointed out that R e = R/γγ t if one requiresthat the eddies be of the maximal causally allowed size. Our simulation results show that eddieshave different sizes and randomly oriented flow directions. Clearly R e will be very small.There is a question as to whether a fluid treatment of the GRB plasma is justified. Thesame question arises in supernova remnants, to which the standard answer is that magnetizationof the plasma ensures a short effective mean free path. From the detailed discussion in Zhang& Yan (2010), the strong collision radius for Coulomb collisions may be defined by e /r col ∼ kT so that r col ∼ e /kT ∼ − /T . The comoving collision mean free path of electrons can beestimated as l ′ e,col = ( n ′ e πr col ) − ≃ cm, where n ′ e is the electron number density in the comovingframe. In order to have the plasma in the collisional regime, one needs to require l ′ e,col < ∆ ′ . 18 –Therefore GRB shocks must be collisionless. However, the GRB fluid stays together throughmagnetic interactions. A more relevant mean free path would be a magnetic gyro-radius andplasma skin depth. The gyro(cyclotron)-radii in the comoving frame are r ′ B,e = γ e m e c /eB ′ ∼ γ e is the electron Lorentz factor and B ′ is the comoving magnetic fieldstrength), and r ′ B,p = γ p m p c /eB ′ ∼ cm for protons (where γ p is the proton Lorentz factor).The comoving relativistic plasma skin depths are δ ′ e = c/ω ′ p,e = ¯ γ e m e c / πn ′ e e ∼ cm forelectrons, and δ ′ p = c/ω ′ p,p = ¯ γ p m p c / πn ′ p e ∼ cm for protons. Here ¯ γ e and ¯ γ p denote themean Lorentz factor of the relativistic electrons and protons, respectively. Physically, plasma skindepths are relevant to non-magnetized ejecta. For magnetized ejecta, the plasma skin depths arerelevant in the direction parallel to the magnetic field lines. The gyro-radii are more relevant inthe direction perpendicular to the magnetic field lines. Both magnetic gyro-radii and plasma skindepths are ≪ ∆ ′ . Therefore this justifies a fluid treatment of the GRB plasma.Our simulations show filamentary magnetic-field structure. A turbulent plasma with fast-moving magnetic filaments is likely a site for second order Fermi acceleration of charged particlesand a source of high energy cosmic rays. The spatial variability and large fluctuations of theturbulent magnetic field also have implications for nonthermal radiation with possibly observableconsequences. A simple light curve (the variability amplitude) can be obtained from integrating B ⊥ with respect to the light of sight measured in the laboratory frame over the postshock regionas a proxy of the electron emissivity. Within the GRB context, turbulence-powered lightcurveswere calculated by Narayan & Kumar (2009) and Lazar et al. (2009). In their calculations, therelativistic turbulent eddies have a characteristic Lorentz factor, and the electron emissivity can becharacterized by γ e B ′ along with a Lorentz boost due to the turbulent motion. As a result, thelightcurve shows a peak when an eddy enters into the line of sight (similar to synchrotron radiation).This can produce significant variability. In our case, the turbulence is mildly relativistic, so that atany epoch, one sees the contribution from many eddies since the Lorentz boost eddy effect is notsignificant. For our current simulations we would expect the eddy emission to average out, and asa result we would expect a roughly flat lightcurve without significant variability. A self-consistentcalculation of nonthermal radiation, which includes the acceleration of nonthermal electrons andthe amplification of magnetic fields would be highly valuable.Y.M. thanks UNLV for hospitality during the preparation of part of the work. Y.M., P.H.,and K.N. acknowledge partial support by NSF awards AST-0506719, AST-0506666, AST-0908010,and AST-0908040, and NASA awards NNG05GK73G, NNX07AJ88G, and NNX08AG83G. J.N. ispartially supported by MNiSW research project N N203 393034, and The Foundation for PolishScience through the HOMING program, which is supported by a grant from Iceland, Liechtenstein,and Norway through the EEA Financial Mechanism. BZ acknowledges NSF award AST-0908362and NASA awards NNX09AT66G and NNX10AD48G for support. The simulations were performedon the Columbia Supercomputer at the NAS Division of the NASA Ames Research Center, theSGI Altix (cobalt) at the National Center for Supercomputing Applications in the TeraGrid projectsupported by the NSF, and the Altix3700 BX2 at YITP in Kyoto University. 19 – A. Calculation of shock wave propagation speed
In our simulations, simulation box is the contact discontinuity frame and plasma moves into thecontact discontinuity (the boundary at x = x max ) from the left side with a flow speed v x = v = 0 . c or 0 . c (a Lorentz factor of γ = 1 . . − x direction. From the relativistic shock jump conditions (Blandford & McKee 1976) in the“flow” (primed) frame, i.e., the frame of the plasma injected at the boundary of x = x min (preshockframe), the shock Lorentz factor is given by γ ′ sh = ( γ + 1)[Γ( γ −
1) + 1] Γ(2 − Γ)( γ −
1) + 2 . (A1)The shock propagation speed in the contact discontinuity frame is obtained from the Lorentztransformation of the velocity as v sh = v ′ sh − v − v ′ sh v , (A2)where v ′ sh is the shock propagation speed in the flow frame. Using an adiabatic index Γ = 5 /
3, theshock propagation speed in the contact discontinuity frame is v sh = 0 . c for the v = 0 . c caseand v sh = 0 . c for the v = 0 . c case, respectively. REFERENCES
Abdo, A. A. et al. 2009, Science, 323, 1688.Aharonian, F. et al. 2003, A&A, 410, 813.Albert, J. et al. 2007, ApJ, 669, 862.Balsara, D. S., & Kim, J. S. 2005, ApJ, 634, 390Balsara, D. S., Benjamin, R. A., & Cox, D. P. 2001, ApJ, 563, 800.Balsara, D. S., Kim, J. S., Mac Low, M., & Mathews, G. J. 2004, ApJ, 617, 339Bamba, A., Yamazaki, R., Ueno, M., & Koyama, K. 2003, ApJ, 589, 827.Bamba, A., Yamazaki, R., Yoshida, T., Terasawa, T., & Koyama, K. 2005a, ApJ, 621, 793.Bamba, A., Yamazaki, R., & Hiraga, J. S. 2005b, ApJ, 632, 294.Begelman, M. C., Blandford, R. D., & Rees, M. J. 2008, MNRAS, 384, L19.Blandford, R. D., & McKee, C. F. 1976, Physics of Fluids, 19, 1130.Brandenburg, A. & Subramanian, K. 2005, Phys. Rep., 417, 1. 20 –Brouillette, M. 2002, Annu. Rev. Fluid Mech., 34, 445.Bykov, A.M., Uvarov, Y.A., Ellison, D.C. 2008, ApJL 689, 133Childress, S., & Gilbert, A. 1995, Stretch, Twist, Fold: The Fast Dynamo (Berlin: Springer).Cho, J., & Lazarian, A. 2003, MNRAS, 345, 325.Diamond, P. H., Hughes, D. W., & Kim, E. 2005, in Self-Consistent Mean Field Electrodynamics inTwo and Three Dimensions, ed. A. M. Soward, C. A. Jones, D. W. Hughes, & N. O. Weiss(Boca Raton: CRC Press), 145.Galloway, D. J. 2003, in Advances in Nonlinear Dynamos, ed. A. Ferriz-Mas & M. N´u˜nez (London:Taylor & Francis), 37.Ghisellini, G., & Tavecchio, F. 2008, MNRAS, 386, L28.Giannios, D., Uzdensky, D. A., & Begelman, M. C. 2009, ApJ, 395, L29.Goldreich, P. & Sridhar, H. 1995, ApJ, 438, 763.Goodman, J. & MacFadyen, A. 2007, J. Fluid. Mech., 604, 325.Gruzinov, A. 2001, ApJ, 563, L15.Gruzinov, A. & Waxman, E. 1999, ApJ, 364, 451.Haugan, N. E., Brandenburg, A., & Mee, A. 2004, MNRAS, 353, 947.Hennebelle, P. & Audit, E. 2007, A&A, 465, 431.Huang, C.-Y., & Pohl, M. 2008, Astropart. Phys. 29, 282Inoue, T. & Inutsuka, S. 2008, ApJ, 687, 303.Inoue, T., Yamazaki, R., & Inutsuka, S. 2009, ApJ, 695, 825.Jiang, G. S., & Shu, C. W. 1996, J. Comput. Phys. 126, 202.Kim, J. S., & Balsara, D. S. 2006, Astro. Nach., 327, 433Kim, J. S., Balsara, D. S., & Mac Low, M.-M. 2001, JKAS, 34, 333.Krawczynski, H. et al. 2004, ApJ, 601, 151.Kulsrud, R. M., & Anderson, S. W. 1992, ApJ, 396, 606.Kumar, P., & Narayan, R. 2009, MNRAS, 395, 472.Kumar, P., et al. 2007, MNRAS, 376, L57. 21 –Lazar, A., Nakar, E., & Piran, T. 2009, ApJ, 695, L10.Levinson, A. 2007, ApJ, 671, L29.Marscher A. P., & Jorstad, S. G. 2010, in Proc. Fermi meets Jansky - AGN in Radio and Gamma-Rays, ed. T. Savolainen, E. Ros, R. W. Porcas, & J. A. Zensus, in press.Marscher, A. P., Gear, W. K., & Travis, J. P. 1992, in Blazar Variability, ed. E. Valtaoja & M.Valtonen (Cambridge Univ. Press), 85.Marscher, A. P. et al. 2008, Nature, 452, 966.Marscher, A. P. et al. 2010, ApJ, 630, L5.Medvedev, M. V. & Loeb, A. 1999, ApJ, 526, 697.M´esz´aros, P. 2006, Rep. Prog. Phys. 69, 2259.Mizuno, Y., Nishikawa, K.-I., Koide, S., Hardee, P., & Fishman, G. J. 2006, ArXiv Astrophysicse-prints, 0609004.Mizuno, Y., Hardee, P. E., & Nishikawa, K.-I. 2010, ApJ, submitted.Narayan, R., & Kumar, P. 2009, MNRAS, 394, L117.Nishikawa, K.-I., Hardee, P., Richardson, G., Preece, R., Sol, H., & Fishman, G. J. 2005, ApJ, 622,927.Nishikawa, K.-I., Niemiec, J., Hardee, P. E., Medvedev, M., Sol, H., Mizuno, Y., Zhang, B., Pohl,M., Oka, M., & Hartmann, D. H. 2009, ApJ, 698, L10.Palma. G., Mignone, A., Vietri, M., & Del Zanna, L. 2008, ApJ, 686, 1103.Panaitescu, A. 2005, MNRAS, 363, 1409.Panaitescu, A., & Kumar, P. 2002, ApJ, 571, 779.Parker, E. N. 1971, ApJ, 163, 255Piran, T. 2005, Rev. Mod. Phys., 76, 114.Pohl, M., Yan, H., & Lazarian, A. 2005, ApJL 626, 101Racusin, J. L., et al. 2008, Nature, 455, 183.Ruzmaikin, A. A., Shukurov, A. M., & Sokolov, D. D. 1988, Magnetic Fields of Galaxies (Dordrecht:Kluwer).Samtaney, R., & Zabusky, N. J. 1994, J. Fluid Mech., 269, 45 22 –Schekochihin, A. A., Cowley, S. C., Taylor, S. F., Maron, J. L., & McWilliams, J. C. 2004, ApJ,612, 276.Shen, R.-F., & Zhang, B. 2009, MNRAS, 398, 1936.Sironi, L., & Goodman, J. 2007, ApJ, 671, 1858.Spitkovsky, A. 2008, ApJ, 682, L5.Steenbeck, M., Krause, F., & R¨adler, G. 1966, Z. Naturforsch, 21a, 369.Uchiyama, Y. et al. 2007, Nature, 449, 576.Uchiyama, Y., & Aharonian, F. 2008, ApJ, 677, L105.Vink, J., & Laming, J. M. 2003, ApJ, 584, 758.Yost, S. A., Harrison, F. A., Sari, R., & Frail, D. A. 2003, ApJ, 597, 459.Zhang, B., & Pe’er, A. 2009, ApJ, 700, L65.Zhang, B., & Yan, H. 2010, ApJ, submitted.Zhang, W., MacFadyen, A., & Wang, P. 2009, ApJ, 690, L40.Zabusky, N. J. 1999, Ann. Rev. of Fluid Mech., 31, 495.