A REVERSE SHOCK IN GRB 130427A

Although long-duration gamma-ray bursts (GRBs) have been associated with the deaths of massive stars (Woosley & Bloom 2006), the precise nature of their progenitors, the structure of their explosion environments, and the composition of their ejecta remain only partially explored. Studies of the afterglow (forward shock, FS) emission provide insight into the explosion energy, geometry, and the structure of the circumburst medium (e.g., Sari et al. 1998; Chevalier & Li 2000). On the other hand, the most useful probe of the initial bulk Lorentz factor and the ejecta composition is afforded by the reverse shock (RS), which is expected to produce a similar synchrotron spectrum as the FS, with well-defined properties relative to the FS (e.g., Sari & Piran 1999b; Kobayashi & Zhang 2003a; Zou et al. 2005).

The expected observational signature of the RS is early time flares in the optical and radio bands and several studies have found hints of excess early-time emission attributable to a RS like component (e.g., Akerlof et al. 1999; Sari & Piran 1999a; Kulkarni et al. 1999; Soderberg & Ramirez-Ruiz 2002, 2003; Kobayashi & Zhang 2003b; Berger et al. 2003; Chevalier et al. 2004). However, a detailed understanding of RS emission requires a careful decomposition of the afterglow spectral energy distribution (SED) into RS and FS components. Since the peak frequencies of the two components are related by a factor of Γ² ≳ 10⁴, such a decomposition requires multi-wavelength observations spanning several orders of magnitude in frequency. Similarly, the density profile of the environment affects the hydrodynamical evolution of both the FS and the RS, leading to discernible differences in the behavior of the light curves. Consequently, well-sampled light curves are also essential for any systematic study of RSs.

Such datasets have not been available to date primarily due to sensitivity limitations of radio and millimeter facilities. However, the recent upgrade of the Very Large Array to the Karl G. Jansky Very Large Array (VLA), with an order of magnitude improvement in sensitivity, has opened a new avenue for the study of RSs in GRBs. Here we present the first example of a RS detected in a multi-wavelength dataset spanning 1–100 GHz of the nearby energetic GRB 130427A. By combining our detailed radio and millimeter observations with X-ray data from Swift and UV/optical/near-infrared (NIR) observations from Swift and ground-based telescopes, we present the most comprehensive dataset for RS studies to date. We undertake a joint model fit to the entire dataset and extract parameters for the explosion and the environment and draw general conclusions on RS studies in the VLA and ALMA era.

2.1. GRB Discovery

2.2. Radio to X-Ray Observations

We interpret the observed behavior of the afterglow emission from radio to X-rays in the framework of the standard synchrotron model, described by three break frequencies (the self-absorption frequency: ν_a, the characteristic synchrotron frequency: ν_m, and the synchrotron cooling frequency: ν_c) and an overall flux normalization. In this model, there are well-defined relations between the temporal evolution and spectral indices that allow us to determine the location of the synchrotron frequencies as well as the density profile of the circumburst medium (interstellar medium (ISM) profile: ρ = const; wind: ρ∝r⁻²). The models are described in detail in Sari et al. (1998) and Chevalier & Li (2000).

A striking feature of the afterglow light curves is the unbroken power-law decline in the X-rays with α_X ≈ −1.35 at ≳ 200 s (F_ν∝t^α). Given the X-ray spectral index of β_X ≈ −0.76 during this time (F_ν∝ν^β), we expect an X-ray decline rate of α_X = −0.64 if ν_c < ν_X (independent of the circumburst density profile), while for ν_c > ν_X, we expect α_X ≈ −1.14 for an ISM profile or α_X ≈ −1.64 for a wind profile. Since the decline rate predicted in the ISM model is shallower than the observed value independent of the ordering of ν_c and ν_X, we conclude that a wind profile is required, with ν_c ∼ ν_X providing a reasonable match to the observed decline. This conclusion is further supported by the similar decline rate in the UV/optical/NIR bands at ≳ 0.3 days, α_opt ≈ −1.35, which indicates that the optical and X-ray bands are located on the same portion of the synchrotron spectrum (i.e., ν_m < ν_opt < ν_X ∼ ν_c). Indeed, as shown with the broadband SEDs in Figure 1, the NIR/X-ray spectral index at ≈2.0 days is β_{NIR, X} ≈ −0.70, in agreement with our conclusion that the optical and X-ray bands are located on the same synchrotron slope. On the other hand, the spectral index within the UV/optical/NIR bands, β_opt ≈ −0.85, is steeper, indicating that modest extinction is present. The observed optical and X-ray decline rates require an electron power-law index of p ≈ 2.2, consistent with the observed X-ray spectral index given the proximity of the cooling break. We note that in the wind model ν_c∝t^1/2 and hence once ν_c crosses the X-ray band, the X-ray flux is expected to decline as an unbroken power law, matching the observed behavior.

Zoom In Zoom Out Reset image size

Finally, a 24.5 GHz radio observation at Δt ≈ 29 days with a flux density of 0.16 ± 0.03 mJy also strongly rules out the ISM environment, since ν_m < ν_NIR at 0.3 days and the NIR K-band flux density is ≈1.6 mJy at this time, which together predict a flux density ≳ 0.9 mJy at 24.5 GHz at ≈29 days in the ISM model. We defer a complete analysis of the observations after ≈15 days to a future work.

Whereas the X-ray light curve follows a single power-law decline, the UV/optical/NIR light curves display a clear change in slope at ≈0.1 days (see Figure 2). The initial slope is shallower, with α_opt ≈ −0.8, and the time of the break is chromatic, occurring later in the redder filters. In the wind model, the passage of ν_m through an observing band results in a transition from α = 0 to α = (1–3p)/4 ≈ −1.4 (for p = 2.2), in clear contrast with the observed evolution. This result indicates that a different emission component dominates the UV/optical/NIR light curves at ≲ 0.1 days.

Zoom In Zoom Out Reset image size

Even stronger evidence for a distinct emission component is provided by the radio and millimeter data. The relatively flat spectral index between the radio/millimeter and optical bands at all times, β_{radio, opt} ≈ −0.25, is inconsistent with a single power-law extrapolation from the optical. This shallow slope cannot be caused by the location of ν_m between the radio and optical bands because all light curves below ν_m should be flat in the wind model, while the observed radio and millimeter light curves clearly decline at all frequencies spanning 6.8–90 GHz. Similarly, the expected spectral index below ν_m is β = 1/3 or 2 (the latter if ν_a is located above an observing band), while we observe instead β_radio ≈ −0.2 at Δt ≈ 2–5 days (Figure 1).

One way to flatten the spectral slope between the radio and optical bands is to introduce a break in the electron energy distribution at γ_b such that $N(\gamma)\propto \gamma ^{-p_{2}}$ (with p₂ < 2) for γ_m < γ < γ_b and such that ν_syn(γ_m) ≲ 1 GHz at Δt ≲ 9.7 days, where ν_syn(γ) is the characteristic synchrotron frequency corresponding to electrons with a Lorentz factor γ. For the observed value of β_{radio, opt}, p₂ ≈ 1.4 is required. However, in this scenario, the radio spectral index can only remain constant (before the passage of ν_syn(γ_b) through the radio band) or steepen with time (after the passage of ν_syn(γ_b)). Instead, we observe that the radio spectral index becomes shallower from β_radio = −0.25 ± 0.02 at 2.0 days to β_radio = −0.09 ± 0.07 at 9.7 days. This fact argues against a break in the electron spectrum causing the shallow radio/optical spectral slope.

Instead, the distinct spectral and temporal behavior in the radio and millimeter clearly requires a different emission component, which is also required to explain the UV/optical/NIR data at ≲ 0.1 days. We associate this component with the RS and show in detail in the next section that it can explain the radio and millimeter data as well as the early optical data. To gain insight into the RS spectrum, we note that the SED at 0.67 days (Figure 1) requires ν_m to be located between 10 and 90 GHz, with an optically thin spectrum extending to the optical. In addition, the spectral index between 5.1 and 6.8 GHz, β ≈ 1.8, is indicative of self-absorption, with ν_a ≈ 7 GHz. The optically thin spectrum observed in the radio at about 2.0 days indicates that ν_a ≲ 5 GHz at that time. Additionally, since the emission from a RS cuts off exponentially above ν_c, whereas we invoke some contribution from this component to the UV flux at ≲ 0.1 days, we require ν_c ≳ 10¹⁵ Hz for the RS at ≈0.1 days.

To summarize, the optical data at ≳ 0.1 days, when combined with the full X-ray light curve at ≳ 200 s, require a wind medium with p ≈ 2.2, ν_c ∼ ν_X, and ν_m ∼ ν_opt at ∼0.1 days. In addition, the radio and millimeter emission and the early optical emission at ≲ 0.1 days cannot be accommodated by same synchrotron spectrum and therefore require a separate emission component. This emission component has ν_a ≈ 7 GHz and ν_m ≈ 10–90 GHz at 0.67 days, as well as ν_c ≈ 10¹⁵ Hz at 0.1 days. In the next section, we expand on these results with full broadband modeling of the FS and RS emission.

Based on the basic considerations described in the previous section, we construct a model SED for the radio to X-ray emission at 0.67 days (Figure 1). The model is composed of two emission components: (1) a FS, which peaks between the millimeter and optical bands, fits the NIR to X-ray SED, and provides negligible contribution in the radio/millimeter and (2) a RS, which fits the radio to millimeter SED and provides negligible contribution at higher frequencies. The synchrotron parameters of the RS are ν_{a, RS} ≈ 7 GHz, ν_{m, RS} ≈ 20 GHz, ν_{c, RS} ≈ 10¹³ Hz, and F_{νm, RS} ≈ 10 mJy. The synchrotron parameters describing the FS are ν_{a, FS} < 5 GHz, ν_{m, FS} ≈ 4 × 10¹³ Hz, ν_{c, FS} ≈ 2 × 10¹⁷ Hz, and F_{ν, m, FS} ≈ 3 mJy. Both SEDs are in the slow cooling regime with the standard ordering of the synchrotron break frequencies, ν_a < ν_m < ν_c. We find that this combined RS plus FS model with a common value of p ≈ 2.2 completely describes the observed SED at 0.67 days (Figure 1).

We evolve both emission components in time to the other three epochs where we have extracted multi-wavelength SEDs. The evolution of the RS spectrum depends upon whether the shock is relativistic or Newtonian in the frame of the unshocked ejecta. Zou et al. (2005) derive the temporal evolution of the RS break frequencies, both before and after the ejecta crossing time (the so-called deceleration time, t_dec) in the two asymptotic regimes of relativistic and Newtonian evolution. An important constraint on the shock evolution is provided by an R-band flux density of about 77 mJy measured by the Faulkes Telescope North at ≈4.3 minutes (Melandri et al. 2013). We find that a relativistic RS overpredicts this early optical emission by about a factor of five. In addition, in the relativistic case, both ν_{m, RS} and ν_{c, RS}⁸ evolve as t^−15/8, leading to a predicted decline rate of t^{−3(5p + 1)/16} ∼ t^−2.25 above ν_{m, RS} (for p = 2.2), which is significantly steeper than the observed decline rate in the millimeter and radio bands, α_mm ≈ −1.5 and α_radio ≈ −1.0. Similarly, the decline rate below ν_{m, RS} is predicted to be t^−1/2 in the relativistic RS model, which is much shallower than the observed decline rates. The Newtonian RS model, on the other hand, offers an additional degree of freedom via the profile of the shocked ejecta, γ∝r^−g, where we expect 1/2 < g < 3/2 from theoretical considerations (Zou et al. 2005). We treat g as a free parameter and find a good match to the SED evolution with g ≈ 5. In the following, we focus on the Newtonian RS case.

Our implementation of the FS follows the smoothly connected power-law synchrotron spectrum for the wind environment of Granot & Sari (2002), where we compute the break frequencies and normalizations using the standard microphysical parameters (_e and _B), the explosion energy (E_{K, iso}), and the circumburst density in the wind profile (A_*). We also use the Small Magellanic Cloud extinction curve (Pei 1992) to model the extinction in the host galaxy (A_V). We use the flux density of the host galaxy in the griz filters as measured in the SDSS as an additive component to the relevant filters. Having determined a set of values for the break frequencies of the RS spectrum at 0.67 days, we evolve these parameters based on Newtonian evolution and fit separately for the FS contribution to the total flux at all frequencies. We determine the FS parameters by fitting all of the available photometry simultaneously with a combination of RS and FS spectra. To efficiently and rapidly sample the available parameter space, we carry out a Markov Chain Monte Carlo (MCMC) analysis using a Python implementation of the ensemble MCMC sampler emcee (Foreman-Mackey et al. 2013). For a detailed discussion of our modeling scheme, see Laskar et al. (2013).

We find that the self-absorption frequency of the FS (ν_{a, FS}) declines below about 300 MHz 0.5 hr after the burst and is not directly observed thereafter. Consequently, the derived blastwave parameters are degenerate with respect to ν_{a, FS}; we quote all results in terms of this frequency, scaled to 14 MHz at 1 day (Table 2). Imposing the theoretical restrictions _e, _B < 1/3, we further constrain the blastwave parameters to the following narrow ranges: 12 MHz < ν_{a, FS} < 16 MHz, 0.25 < _e < 0.33, 0.14 < _B < 0.33, 2.1 × 10⁻³ < A_* < 3.7 × 10⁻³, and 4.4 < E_{K, iso, 52} < 5.8. Marginalizing the posterior density functions, we find the following median values of the FS parameters: p = 2.23, _e = 0.30, _B = 0.20, A_* = 2.9 × 10⁻³, E_{K, iso, 52} = 4.9, and A_V = 0.18 mag. The FS spectrum transitions from fast cooling to slow cooling around Δt ≈ 1200 s. The best-fit combined RS and FS model is shown for the multi-epoch SEDs in Figure 1 and for all available radio to X-ray light curves in Figure 2. The model fully captures the observed evolution across nine orders of magnitude in frequency and over three orders of magnitude in time.

Table 2. Forward Shock Parameters for the Best Fit Wind Model

Parameter	Value
p	2.23
e	0.30(νa/14 MHz)5/6
b	0.20(νa/14 MHz)−5/2
A*	3.0 × 10−3(νa/14 MHz)5/3
EK, iso, 52	4.9(νa/14 MHz)−5/6
AV	0.18 mag
Γ	130
tjet	≳ 15 days
θjet	>25
νa, FS	(12–16)× 106 Hz
νm, FS	2.2 × 1013 Hz
νc, FS	2.8 × 1017 Hz
νa, RS	5.7 × 109 GHz
νm, RS	1.2 × 1010 Hz
νc, RS	5.9 × 1012 Hz

Notes. All frequencies in this table are calculated at Δt = 1 day. The self-absorption frequency of the forward shock at 1 day, ν_a, is constrained to 12 MHz <ν_a < 16 MHz by the requirement _e, _b  <  1/3.

Download table as:  ASCIITypeset image

The bulk Lorentz factor of the ejecta, Γ, can be calculated from a knowledge of the deceleration time and the parameters of the explosion. In the internal–external shock model, we expect the deceleration time to roughly match the duration of the burst. Since the Swift/BAT T₉₀ ≈ 163 s and the optical flux is already declining at ≈258 s (Melandri et al. 2013), we take t_dec ≈ 200 s. Using the relation $t_{\rm dec}\approx 29(1+z)E_{\rm K,iso,52}\Gamma _{1.5}^{-4} A_{*}^{-1}$ s (Zou et al. 2005), we find Γ ≈ 120–150 at t_dec, where the range corresponds to the uncertainty in E_{K, iso, 52} and A_* due to the uncertainty in ν_{a, FS}. For the median values of E_{K, iso, 52} and A_* reported above, we obtain Γ ≈ 130.

In the preceding analysis, we independently determined the RS break frequencies and fit for the parameters of the FS. However, the two synchrotron spectra are related since the two shocks propagate in opposite directions from the contact discontinuity. In particular, we expect that at the deceleration time ν_{c, RS} ∼ ν_{c, FS}, ν_{m, RS} ∼ ν_{m, FS}/Γ² and F_{ν, m, RS} ∼ ΓF_{ν, m, FS}. For the parameters given above, we find that the FS parameters at the deceleration time are ν_{c, FS} ≈ 7.8 × 10¹⁵(t_dec/200 s)^1/2 Hz,⁹ν_{m, FS} ≈ 2.2 × 10¹⁷(t_dec/200 s)^−3/2 Hz, and F_{ν, m, FS} ≈ 13(t_dec/200 s)^1/2 mJy. For the RS, we have ν_{c, RS} ≈ 1.5 × 10¹⁶(t_dec/200 s)^−1.3 Hz, Γ²ν_{m, RS} ≈ 5.3 × 10¹⁷(t_dec/200 s)^−1.3(Γ/130)² Hz, and F_{ν, m, RS}/Γ ≈ 13(t_dec/200 s)^−0.9(Γ/130)⁻¹ mJy at the deceleration time, where the power-law indices for ν_{c, RS}, ν_{m, RS}, and F_{ν, m, RS} are derived from the g-dependent expressions in Zou et al. (2005). We thus confirm that the expected relations between the RS and FS parameters are satisfied to within a factor of two, confirming our basic assumption that the two required emission components indeed correspond to the RS and FS. Whereas our results are consistent with the shock microphysical parameters being the same in the FS and RS regimes, we note that these observations can accommodate magnetization parameters (Zhang et al. 2003; Harrison & Kobayashi 2013) in the range $1\le \sqrt{\epsilon _{\rm B, RS}/\epsilon _{\rm B,FS}}\lesssim 5$ for our fiducial values of Γ and t_dec. These limits are unfortunately degenerate with Γ and t_dec and we are unable to constrain these more strongly at this time.

To summarize, a model that includes emission from the RS and FS consistently explains all the available data from the radio to X-rays over a timescale of ∼200 s to ∼10 days. The resulting bulk Lorentz factor of the outflow at the deceleration time is Γ ≈ 130.

We present a detailed multi-wavelength study of the bright afterglow of GRB 130427A spanning radio to X-rays. From the optical and X-ray data, we conclude that the progenitor exploded into a wind environment, pointing to a massive star. The radio and millimeter observations present a spectrum and temporal evolution that cannot be explained by emission from the FS alone. We show that this emission is consistent with synchrotron radiation from a Newtonian RS. With the available multi-band data, GRB 130427A is by far the best example to date of RS emission in the radio/millimeter, particularly when compared with previous detections that were based only on only one to two epochs at single frequencies (Kulkarni et al. 1999; Berger et al. 2003; Soderberg & Ramirez-Ruiz 2003). Our derived value of the index (g ≈ 5) of the Lorentz factor profile of the shocked ejecta formally lies outside the expected range. The exquisite dataset available for this GRB provides a unique opportunity for a detailed study of RSs in wind environments, requiring the combined power of numerical simulations and multi-wavelength afterglow modeling. While beyond the scope of this work, such a study might potentially reconcile the large value of g with theoretical expectations, as well as shed light on the question of ejecta magnetization.

Using multi-wavelength model fitting of the rich afterglow dataset, with the well-sampled light curves spanning over three orders of magnitude in time and nine orders of magnitude in frequency, we determine the properties of the explosion and the circumburst medium. In particular, we find a low circumburst density, A_* ≈ 3 × 10⁻³, corresponding to a mass loss rate of $\dot{M}=3\times 10^{-8}$ M_☉ yr⁻¹ (for a wind velocity of 1000 km s⁻¹). The low density causes the RS to be in the slow cooling regime (ν_{m, RS} < ν_{c, RS}), resulting in long-lived radio and millimeter emission. We note that a low density was also inferred for previous GRBs with likely radio RS emission (990123, 020405, and 021211; Chevalier et al. 2004), suggesting that low density is a requisite criterion for observable emission from a RS; in a high-density environment, the RS emission will decline rapidly due to efficient cooling of the radiating electrons (Kobayashi & Zhang 2003a).

From the derived properties of the explosion and environment, we obtain the Lorentz factor of the outflow at the deceleration time, Γ(200 s) ≈ 130, and show that the spectra of the FS and RS at the deceleration time are consistent with theoretical expectations. We infer E_{K, iso} ≈ 7 × 10⁵² erg. However, we note that the FS spectrum transitions from fast to slow cooling at ≈1200 s and radiative corrections between this time and t_dec ≈ 200 s are a factor of a few (Sari 1997; Dai & Lu 1998). This fact suggests that E_{K, iso} ≈ 2 × 10⁵³ erg, which is ≈20% of the isotropic-equivalent gamma-ray energy, implying a large radiative efficiency. The lack of an obvious break in the X-ray light curves to ≈15 days implies a lower limit of only θ_jet ≳ 25 on the opening angle of the jet, and hence E_K ≳ 2 × 10⁵⁰ erg and E_γ, ≳ 10⁵¹ erg. Thus, due primarily to the combination of a large isotropic energy and a low circumburst density, the lack of a jet break at ≲ 15 days is not surprising.

We conclude by noting that GRB 130427A is likely to become the benchmark for RS studies in the VLA and ALMA era. Our study demonstrates that a complete analysis of the explosion and ejecta properties requires detailed multi-wavelength modeling to leverage the anticipated exquisite datasets.

We thank Shiho Kobayashi, Yuanchuan Zou, Xuefeng Wu, Zigao Dai, Richard Harrison, and the anonymous referee for their helpful comments on the manuscript. The Berger GRB group is supported by the National Science Foundation under Grant AST-1107973. The MMT Observatory is a joint facility of the Smithsonian Institution and the University of Arizona. The National Radio Astronomy Observatory is a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc. This work made use of data supplied by the UK Swift Science Data Centre at the University of Leicester. CARMA observations were taken as part of the CARMA Key Project, c0999, "A Millimeter View of the Transient Universe" (PI: B. A. Zauderer) and VLA observations were taken as part of programs 13A-046 (PI: E. Berger), 13A-411 (PI: A. Corsi), SE0851 (PI: A. Fruchter), and S50386 (PI: S. B. Cenko). Support for CARMA construction was derived from the states of California, Illinois, and Maryland, the James S. McDonnell Foundation, the Gordon and Betty Moore Foundation, the Kenneth T. and Eileen L. Norris Foundation, the University of Chicago, the Associates of the California Institute of Technology, and the National Science Foundation. Ongoing CARMA development and operations are supported by the National Science Foundation under a cooperative agreement, and by the CARMA partner universities. The GMRT is run by the National Centre for Radio Astrophysics of the Tata Institute of Fundamental Research.