Predictions for radio afterglows of binary neutron star mergers: a population study for O3 and beyond

Predictions for radio afterglows of binary neutron star mergers: a population study for O3 and beyond

R. Duque Sorbonne Université, CNRS, UMR 7095, Institut d’Astrophysique de Paris, 98 bis boulevard Arago, 75014 Paris, France    F. Daigne & R. Mochkovitch Sorbonne Université, CNRS, UMR 7095, Institut d’Astrophysique de Paris, 98 bis boulevard Arago, 75014 Paris, France
Key Words.:
Methods: statistical – Stars: neutron – Gravitational waves – Gamma-ray burst: general – Gamma-ray burst: individual: GRB170817A

Context:Following the historical observations of GW170817 and its multi-wavelength afterglow, more radio afterglows from neutron star mergers are expected in the future as counterparts to gravitational wave inspiral signals.

Aims:We wish to describe these events using our current knowledge of the intrinsic population of neutron star mergers coming from gamma-ray burst science, and taking into account the sensitivities of current and future gravitational wave and radio detectors.

Methods:This is done by coupling a semi-analytical model for the jet-dominated radio afterglow from a neutron star merger, an analytical simplified model for the gravitational wave signal from such events, and a population model prescribing the energetics, external density and other relevant parameters of the mergers.

Results:We report the expected distributions of observables (distance, orientation, afterglow peak time/flux, etc.) from future events and study how these can be used to further probe the population of binary neutron stars, their mergers and related outflows during future observing campaigns. In the case of the O3 run of the LIGO-Virgo Collaboration, the radio afterglow of one third of gravitational-wave-detected mergers should be detectable–and detected if the source is localized thanks to the kilonova counterpart–by the Very Large Array, and these events should have viewing angles similar to that of GW170817.

Conclusions:These findings confirm the radio afterglow as a powerful insight on these events, though some key afterglow-related techniques–such as Very Long Base Interferometry imaging of the merger remnant–may no longer be possible as the gravitational wave horizon increases.

1 Introduction

The first detection of the gravitational waves (GW) from the inspiral phase of a binary neutron star merger (Abbott et al., 2017, 2019) was followed by all three electromagnetic counterparts expected after the coalescence: a short gamma-ray burst (GRB) (Goldstein et al., 2017; Savchenko et al., 2017), its multi-wavelength afterglow (AG) in the X-ray (Haggard et al., 2017; Margutti et al., 2017; Troja et al., 2017), radio (Hallinan et al., 2017) and optical bands (Lyman et al., 2018) and an optical-IR thermal transient source (Evans et al., 2017; Arcavi et al., 2017; Lipunov et al., 2017; Soares-Santos et al., 2017; Valenti et al., 2017; Tanvir et al., 2017; Cowperthwaite et al., 2017; Drout et al., 2017; Kasliwal et al., 2017; Coulter et al., 2018), which allowed to localize the event with sub-arcsecond precision in the S0-type galaxy NGC4993 at a distance of  Mpc (Palmese et al., 2017; Cantiello et al., 2018). This thermal transient showed evidence for heating from -process nucleosynthesis (Pian et al., 2017; Smartt et al., 2017), classifying it as a “kilonova” (KN), the first with such detailed observations.

According to estimates on joint GW-GRB events by Beniamini et al. (2019), such a combined detection of GW+GRB+AG+KN should remain rare. Indeed, GRB170817A would not have been detected at a distance larger than 50 Mpc (Goldstein et al., 2017) or from a viewing angle larger than .

What can we expect regarding the afterglows of future events during the present O3, and future observing runs of the LIGO-Virgo Collaboration (LVC)? What will be the future rates of GW and joint GW-AG observations from binary neutron star mergers? How will the future events distribute themselves in terms of observables such as the viewing angle , distance , afterglow peak time and peak flux and and proper motion ?

In this paper111Preliminary results were presented at the 8th Fermi Symposium in Baltimore (October 2018); see presentation here., we consider these questions from the angle of a population study. Precisely, we generate a sample of mergers with an intrinsic variability in parameters such as kinetic energy, external medium density and shock microphysical conditions, etc. These reflect prior knowledge from GRB science and observations of the 170817 event. We then calculate the afterglows arising from these mergers in the context of a jet-dominated afterglow emission, and finally apply criteria of current and future GW and EM detectors to determine those events which will be detected jointly in GW and EM domains. These make up a population of events likely to be observed, for which we can study the distributions of observables.

Such a population model for short GRB afterglows was already considered by Saleem et al. (2018a, b), with detailed numerical models for the multi-band afterglow and GW signals, but without clear astrophysical motivation for the event parameter distributions. Moreover, the discussion was centered on the parameter-space constraints one can derive from detections or non-detections of the afterglow in various electromagnetic bands, and on the conditions to observe afterglows from off-axis lines-of-sight. Here, we will rather discuss the effect of the population parameters and detector configurations on the expected population of afterglows.

More recently, this approach was also followed by Gottlieb et al. (2019), and again our study differs in that it accounts for some prior knowledge on short GRBs (such as their luminosity function) and observations from GW170817 to sharpen our predictions on the events to come. This allows a detailed study of the impact of the population parameters on the predicted observations.

Our study relies on the following assumptions:

  • In most cases the merger produces a successful central jet, whose contribution dominates at the peak of the afterglow. This is clearly a strong assumption that will have to be validated by future observations. Recent studies based on the population of short GRBs (Beniamini et al., 2019) and the dynamics of jets interacting with merger ejecta (Duffell et al., 2018) suggest that successful jets could be common in these phenomena.

  • In this case, the afterglow peak flux depends on the jet’s parameters (isotropic kinetic energy , opening angle of the central core jet, shock microphysics parameters , and ), on its environment (external medium density ), and its viewing conditions (the distance and viewing angle). The distributions of most of these quantities remain quite uncertain. They rely on afterglow fitting of a limited sample of short GRB afterglows with known distance (e.g. Berger, 2014; Fong et al., 2015), on the expectation that short GRBs generally occur in a low density environment and on the results obtained on the (up to now) single event GRB170817A (e.g. van Eerten et al., 2018).

This publication is organized as follows. In Sec. 2 we describe our models to calculate the afterglow emission and detail the input parameters of our population model. In Sec. 3 we give the criteria for GW and radio afterglow detection we apply to obtain the observed population from the intrinsic population. In Sec. 4 we report general results on the observed population, describe its characteristics (viewing angle, peak time and peak flux, proper motion distributions) and study their sensitivity to the detector configurations and the population model parameters. In Sec. 5, we discuss our results in the context of future multi-messenger campaigns and how they may help to interpret observations therein.

2 A population model for binary neutron star merger afterglows

2.1 Peak flux and peak time of the radio afterglow

The multi-wavelength afterglow of GW170817 is associated with the deceleration of a structured relativistic jet emitted by the central source formed after the merger, and more precisely to the synchrotron emission of electrons accelerated by the forward shock propagating in the external medium. It was observed for more than 300 days. The time evolution is similar at all wavelengths, indicating that the emission is produced in the same spectral regime of the slow cooling synchrotron process, i.e. , and that the non-thermal distribution of shock accelerated electrons has a non-evolving slope (Nakar et al., 2018; Mooley et al., 2018b; van Eerten et al., 2018). The observed slow rise of the afterglow is a clear evidence of the lateral structure of the outflow, which is observed off-axis with a viewing angle and the observed flux is initially dominated by the relativistically beamed emission from mildly-relativistic/mildly-energetic material coming towards the observer (Mooley et al., 2018c). Due to the deceleration, the relativistic beaming becomes less and less efficient and regions closer to the jet axis start to contribute to the flux. The peak is reached when the core jet, with an opening angle of a few degrees, becomes visible, when its Lorentz factor has decelerated down to . The alternative possibility, a quasi-spherical mildly relativistic ejecta with a steep radial structure (Kasliwal et al., 2017; D’Avanzo et al., 2018; Gottlieb et al., 2018; Hotokezaka et al., 2018; Mooley et al., 2018c; Nakar & Piran, 2018; Gill & Granot, 2018; Troja et al., 2018) is strongly disfavored by radio VLBI observations, which confirm the emergence of a relativistic core jet at the peak of the afterglow, showing an apparent superluminal motion with (Mooley et al., 2018a) and a compact angular size of the radio source below (Ghirlanda et al., 2019). The lateral structure of the jet is constrained by observations: best fits are obtained for a sharp decay of the kinetic energy and the Lorentz factor with the angle from the jet axis (Lamb & Kobayashi, 2017; Lazzati et al., 2017a, b; Troja et al., 2017, 2018; Resmi et al., 2018; Gill & Granot, 2018; D’Avanzo et al., 2018; Margutti et al., 2018). Typically, for a top-hat core jet with an opening angle , an isotropic equivalent kinetic energy and a Lorentz factor surrounded by a sheath with a power-law structure for the kinetic energy per solid angle and the Lorentz factor , the required initial distributions are




with steep slopes (Gill & Granot, 2018; Ghirlanda et al., 2019).

Figure 1: Left: Radio light curve at of the afterglow of a structured jet with a sharp power-law structure at a distance , with a viewing angle , an external density , a core jet with , , , surrounded by a sheath with a power-law lateral structure with and (see Eq. 1-2), and with microphysics parameters , and . The radio observations of GW170817 are also plotted and are compiled from Hallinan et al. (2017); Alexander et al. (2018); Margutti et al. (2018); Mooley et al. (2018b, c); Dobie et al. (2018). The respective contributions of the core jet (), the sheath and the total are plotted in dashed, dotted and solid red lines. The flux of the core jet computed with the simplified treatment described by Eq. 3 is also plotted in dashed blue line for comparison. Bottom right: Peak flux of the same structured jet as a function of the viewing angle (solid red line), peak flux of the light curve from the core jet only (dashed red line) or from the sheath only (dotted red line), and peak flux of the core jet as computed using the scaling law in Eq. 5 (thin black line). Upper right: Ratio of the peak flux from the whole outflow (core jet+sheath) to that of the core jet only.

An example of the radio lightcurves obtained for such a structured jet is plotted in Fig. 1 (left) using parameters typical of the various fits to the data that have been published (see e.g. Gill & Granot, 2018; van Eerten et al., 2018; Ghirlanda et al., 2019). To compute this lightcurve, the dynamics of the material at a latitude are computed independently of those at other latitudes (neglecting any lateral motion), with a deceleration radius which is constant in the core and slowly increases with outside the core ().

Then, the synchrotron emission of shock accelerated electrons in the shock comoving frame is assumed to follow the standard synchrotron slow-cooling spectrum (Sari et al., 1998) including self-absorption. Finally, the observed lightcurve is computed by summing the contributions of all latitudes on equal-arrival time surfaces, taking into account the relativistic beaming and Doppler boosting.

The separate contributions of the core jet and the sheath are plotted in Fig. 1 (left) and the emergence of the core at the peak is clearly visible. The evolution of the peak flux of the same structured jet is plotted in Fig. 1 (right) as a function of the viewing angle, as well as the ratio of the peak flux to the peak flux from the core jet only. Interestingly, this ratio is almost constant (), except for when the core jet is seen on-axis and is more dominant.

As long as the kinetic energy and the Lorentz factor decay steeply with , the core jet is expected to dominate the flux at the peak whatever the viewing angle, even if the precise value of the total-to-core ratio at the peak may slightly vary depending on the details of the assumed lateral structure. Therefore, as our population model considers only the properties of the afterglow at the peak, i.e. when it can more easily be detected, it is enough in the following to compute only the contribution of the core jet, recalling that it may slightly underestimate the peak flux by a factor . For a top-hat core jet as assumed above, this can be done very efficiently by computing the dynamics with the same simplification as above, and by computing the flux using the approximation suggested by Granot et al. (2002) to avoid the full integration over equal-arrival time surfaces.

Precisely, we compute the flux of a spherical ejecta with initial Lorentz factor and kinetic energy , and we correct it by


with and where is the source frame time. The latter correction, corresponding to the ratio of on-axis/off-axis arrival times is not included in Granot et al. (2002) and is found to improve the quality of the approximation. The corresponding light curve is plotted in Fig. 1 (left) and agrees well with the exact calculation.

This emission of the core jet corresponds to the standard calculation of the afterglow of a top-hat jet. Therefore analytical expressions for the properties at the peak are available. They depend slightly on the assumptions for a possible late jet lateral expansion. If lateral expansion is taken into account like in standard GRB afterglow theory, the peak flux of the radio light curve scales as (Nakar et al., 2002):


as long as the spectral regime remains . This leads to the following expression222We neglect the effect of redshift, as the events at play here are within the GW horizon distance of first generation interferometers, i.e. with . of the peak flux at :


where erg, rad, , , , and , and where we assume . This scaling law is plotted in Fig. 1 (right) and agrees well with the detailed calculation. The expression of the corresponding peak time is also given by Nakar et al. (2002), assuming lateral expansion of the jet. However, even if late observations of GRB 170817A may give some evidence for lateral expansion of the core jet (with a temporal decay index , as suggested by Mooley et al. 2018b), it is not clear if this expansion should be as strong for a core jet embedded in a sheath as for a “naked” top-hat jet. Therefore we also consider a limit case without lateral expansion. This only slightly affects the peak flux and Eq. 5 above remains a good approximation in both cases. On the other hand, it modifies the peak time, leading to:




where is the standard jet break time,


The ratio of the peak times without and with lateral expansion is simply


In the following, except if mentioned otherwise, we use Eqs. 5-7 to estimate the peak flux and peak time of the radio afterglow from a BNS merger.

Parameter Symbol Probability Distribution Function
Jet isotropic equivalent kinetic energy Broken power law (Eq. 10), with ,
Jet half-opening angle
External medium density Log-normal distribution, central value , standard deviation
Electron redistribution parameter Fixed at 0.1
Magnetic field redistribution parameter Same Log-normal as , restricted to
Electron population spectral index Fixed at 2.2
Table 1: Fiducial model parameter distributions.

2.2 Physical parameters of the jet and their intrinsic distribution

To generate a population of jet afterglows we have to fix the various parameters appearing in Eqs. 57. We first define a set of “fiducial distributions” for these parameters. We will present in Sec. 4 the corresponding predicted distribution of observables, and discuss the impact of some possible variations around this reference case. The parameters of our fiducial model are summarized in Tab. 1 and have been chosen as follows:

  • Jet isotropic equivalent kinetic energy : we adopt a broken power law distribution, which we directly deduce from the gamma-ray luminosity function of cosmological short GRBs, assuming a standard rest frame duration s and a standard efficiency (Beniamini et al., 2016) so that , leading to a density of probability


    where is normalized to unity. We fix and . Luminosity functions for short GRBs have been deduced from observations by several groups (e.g. D’Avanzo et al., 2014; Wanderman & Piran, 2015; Ghirlanda et al., 2016). For our fiducial model, we adopt that of Ghirlanda et al. (2016) with , and a break luminosity leading to .

  • Jet opening angle : for simplicity we assume a single value , which appears to be representative of the typical value found by fitting the afterglow light curve of GW170817 (e.g. Gill & Granot, 2018; van Eerten et al., 2018), by the VLBI observations of the remnant (Mooley et al., 2018a) and also consistent with the results of prior short GRB afterglow fitting (Fong et al., 2015)and short GRB/BNS merger rate comparisons (Beniamini et al., 2019).

  • External density : most BNS mergers are expected to occur in low density environments due to a long merger time, as assessed by the significant offsets of short GRB sources from their host galaxies (e.g. Nakar, 2007; Troja et al., 2008; Fong et al., 2010; Church et al., 2011; Berger, 2014). In the case of GW170817, afterglow fitting leads typically to , in agreement with the study of the HI content estimation of the host galaxy NGC 4996 (Hallinan et al., 2017). The external densities observed in short GRBs cover the interval (Berger, 2014), but this is probably biased towards high densities, which favor afterglow detection. In this work we consider  cm as typical and therefore take for the fiducial density distribution a log-normal of mean with standard deviation of .

  • Microphysics parameters: we adopt the common value generally used in GRB afterglow models (Wijers & Galama, 1999; Panaitescu & Kumar, 2001; Santana et al., 2014; Beniamini & van der Horst, 2017; D’Avanzo et al., 2018) and also adopted in most fits of the afterglow of GW170817 (e.g. van Eerten et al., 2018; Margutti et al., 2018; Nakar et al., 2018; Xie et al., 2018). We take , which is also a common value used in GRB afterglow models and is in agreement with shock acceleration theory in the relativistic regime (e.g. Sironi et al., 2015). In the case of GW170817 this value can be directly measured from multi-wavelength observations (Mooley et al., 2018b, c; van Eerten et al., 2018). Finally, more diversity is generally expected for . We assume a log-normal distribution similar to the one adopted for with the additional constraint that is restricted to the interval .

In Sec. 4.2 we consider possible alternatives to these fiducial distributions. In particular, we discuss the impact of the large uncertainties on the distribution of kinetic energy, , either due to the difficulty to measure the luminosity function of short GRBs, or to our simplifying assumptions regarding the duration and the efficiency of the prompt GRB emission. We also explore different values of the jet opening angle . For the external density, , we study the effect of an increase of the mean value of the log-normal distribution. Finally, we discuss the impact of our assumptions for the microphysics parameters.

3 From the intrinsic to the observed population: GW+radio joint detection

3.1 GW detection criterion

The signal-to-noise ratio (SNR, denoted by ) of an inspiral signal in a single LIGO-Virgo-type interferometer can be written as (Finn & Chernoff, 1993), where is the luminosity distance to the binary, is the chirp mass of the binary, is a quantity depending only the sensitivity profile of the interferometer, and represents the dependence of the signal to the binary sky position and orientation with respect to the plane of the interferometer. Again, is the angle of the line-of-sight to the binary polar direction, which we will also suppose is the jet-axis direction. admits a global maximum of corresponding to an optimally positioned and oriented source (binary at the zenith and polar axis orthogonal to the instrument’s arms). This optimal binary can be detected out to a distance known as the horizon of the instrument, which depends on the SNR threshold for detection, taken to be for the LVC network. We may thus rewrite the SNR in the following manner:


In a full population study, one would have to draw all four angles and evaluate this criterion in every case. We choose to reduce the number of parameters to two: and , as these are the ones relevant to the afterglow. We must thus review this criterion by averaging on sky-position and polarization angle . This is readily done from the expression given in Finn & Chernoff (1993) (Eq. 3.31) and it is found analytically that:


Hence, we use the following criterion to determine those binaries of inclination and distance which are detectable in GW on average in the sky:


which is:


where we have denoted the sky-position-averaged horizon.

Averaging on inclination angle and imposing a SNR threshold as in Eq. 13 would lead to the simple detection criterion of , where is the range of the instrument, i.e. the maximum distance to which a binary can be detected on average in sky-position and in orientation. The range is linked to the horizon by .

The criterion of Eq. 14 is valid for the detection using a single instrument. Instead, GW detection by the interferometer network is based on multi-instrument analysis, and is thus more complicated than that described by our criterion. Furthermore, as it was illustrated in the case of GW170817, true joint GW+AG detection requires the pin-pointing of the source, which in turn involves a small enough localization map from the GW data333and the detection of the kilonova counterpart, which proved impossible in the recent merger candidates S190425z and S190426c., and our criterion should incorporate this. We choose a simple localization criterion by supposing that a source is localized if it is detectable (once again on average) by the two most sensitive interferometers of the network. Thus our detection-localization criterion is that of Eq. 14, with the horizon of the second most sensitive instrument of the network, i.e. the LIGO-Hanford instrument. Tab. 2 reports the corresponding values taken for our study for the LVC O2, O3 runs and for the design interferometers.

In this framework, the mean viewing angle of GW-detected events is regardless of the horizon value. Also, the fraction of GW-detected events among all mergers is , independently of the horizon. This is thus an absolute maximum to the fraction of mergers to be jointly detected in both radio and GW channels.

LVC Run (Mpc) (Mpc)
O2 136 86
O3 226 143
Design 429 272
Table 2: Horizons used for the GW detection-localization criterion in our study. They are deduced from the published BNS ranges of the instruments (Abbott et al., 2018) by and considering the value for LIGO-Hanford (see text for details of this choice).
Figure 2: Left: Detected fraction of radio afterglows among gravitational wave events as a function of the horizon distance . Right: Expected number of detections normalized to the case of a horizon distance of Mpc (O3) and a radio detection limit of 10 Jy (VLA). In both panels the full (resp. dashed) lines correspond to a (resp. ) flux limit. The three dots represent the results for the O2, O3 ( limit) and design ( limit) configurations.

3.2 Radio detection criterion

The detection criterion in the radio band is simply that the 3 GHz peak flux as determined by our models be larger than the sensitivity of the radio array available for follow-up at the time considered, which we will denote by . We note that this is a criterion for detectability rather than detection. True detection only occurs if observations are pursued for long enough after the GW merger signal, as the flux rises to its maximum. We also note that events which are marginally detectable will provide only poor astrophysical output because the inference of jet parameters requires the fitting of an extended portion of the radio afterglow light-curve.

We will take two typical radio sensitivities reflecting the present and future capabilities of arrays: the Karl Jansky Very Large Array (VLA), at , and the Square Kilometer Array (SKA) or the Next Generation VLA (ngVLA) at . The combination of this radio detection criterion with that in the GW domain (Eq. 14) constitute a criterion for joint detection of mergers.

4 Results: detection rates and properties of the detected population

4.1 Fiducial model

In this section, we describe the population of events detected jointly in the GW and radio domains. For this purpose, we use a Monte Carlo approach where we simulate binary systems within the sky-averaged horizon using the parameter distributions of the fiducial model described in Sec. 2.2. Applying the detection criteria described in Sec. 3, We obtain systems detected by the GW interferometer network and systems which are also detectable by radio-telescopes through their afterglows. Then we study the properties of the sample of joint GW+radio detectable events.

4.1.1 Rate of joint GW+Radio detectable events

Starting from the fraction of binary systems detected by the LVC (Sec. 3.1), we will now estimate which fraction of those events also produce a detectable radio afterglow.

The results are shown in Fig. 2 (left) for O2 and O3 assuming the present VLA sensitivity, and for the design configuration of the LVC network assuming the ngVLA limits. In the design configuration of the LVC network, ngVLA can detect times more events than VLA.

The number of detections (normalized to the value for Mpc and a threshold Jy, i.e. the O3+VLA configuration) is represented in Fig. 2 (right). The number of joint detections behaves approximately as , with because of the reduction of the radio detection efficiency when the distance increases (Fig. 2, left). We find (resp. ) for a radio sensitivity (resp. ).

As illustrated in Fig. 2 (left) we find that the fraction of detected events decreases from 47% (O2) to 34% (O3) with a limiting sensitivity of 10 Jy. However the absolute number of detections will increase (Fig. 2, right). With the design+SKA/ng-VLA configuration, we recover a fraction of about half of the cases that can lead to joint detections, like for O2+VLA. Indeed the product is almost the same in the two configurations.

For completeness, the number of GW and jointly detected events per continuous year of GW network operation for future detector configurations can be found in Tab. 3, where we have used the horizons of Tab. 2 and our fiducial model.

4.1.2 Distance and viewing angle

LVC Run Radio Sensitivity GW Events Joint Events
O3 10
Design 10
Design 1
Table 3: Number of detectable events per continuous year of GW network operation. The expected rate of GW-only events is taken from Abbott et al. (2018) and values for joint events are derived using our fiducial population model. The uncertainties here are due to the uncertainty on the merger rate in the local Universe derived from GW observations of the LVC O1 and O2 runs. The additional uncertainties related to the population model will be discussed in Sec. 4.2.

The distribution in distance is shown in Fig. 3 for the O2 and O3 configurations. It can be seen that as a result of the GW and radio detection thresholds, sources are progressively lost when the distance increases so that this distribution exhibits a near linear increase (as opposed to of the intrinsic homogeneous population). Up to , all events are detected in GW, regardless of their orientation444It is only after this distance that events are detected in GW only if with . The differential distribution of distances thus transitions from before to a non-quadratic form afterwards, producing the small peak seen in Fig. 3. This of course is only the consequence of our simplified GW sky-averaged detection criterion (Eq. 14). and the only selection is due to the radio sensitivity. Near the horizon, the maximum viewing angle allowing GW-detection is , which produces a strict decrease of event density.

Figure 3: Left: Differential distribution of the distances (normalized to the sky-position-averaged horizon) of events detected through GW only (black), and jointly in the O2+VLA (blue) and O3+VLA (red) configurations. Right: Cumulative distribution of the same events.
Viewing angle.

Fig. 4 shows the distribution of the viewing angles. As a result of the rapid decline of the peak flux with angle () the peak of the distribution takes place at small viewing angles and about 40% (resp. 50%) of the events are seen with a viewing angle for O2 (resp. O3), the angle with which the GW 170817 event was seen. As shown in Fig. 5 (left), the mean viewing angle of the jointly observed events strongly depends on the GW horizon and the radio sensitivity. It appears that even if the radio sensitivity is not largely improved while the GW interferometers reach their design horizons, there will remain a significant number of events seen off-axis. The increasing and decreasing phases of these events’ afterglows will allow to better study the jet structure which is better revealed by off-axis events. Fig. 5 (right) shows the fraction of on-axis events () within the joint detections. It increases from 5.6% (O2) to 8.6% (O3).

Figure 4: Left: Differential distribution of the viewing angles of events detected through GW only (black), and jointly in the O2+VLA (blue) and O3+VLA (red) configurations. Right: Cumulative distribution of the same events.
Figure 5: Left: Mean viewing angle of the jointly detected events, as a function of the horizon and the radio sensitivity . The values of for O2, O3 and design instruments as well as for VLA and SKA/ng-VLA are indicated in dashed line. Right: Fraction of on-axis events (defined by ) among jointly detected events in the same diagram.

4.1.3 Peak times, peak fluxes and synchrotron spectral regime

The distributions of peak times and peak fluxes are represented in Fig. 6. The fraction of events with a peak time smaller than 150 days (as observed in GRB170817A) rises from about 54% without jet expansion to 79% with lateral expansion. The distribution in peak flux is shown for all sources which are detected in gravitational waves within the sky-position-averaged horizon of O2, O3 and design instruments. It appears clearly that for the present VLA sensitivity, most radio afterglows cannot be detected. This explains why improving the sensitivity of the future ngVLA has such an impact on the joint detection rates, as discussed in Fig. 2. However, even in the O3 configuration, the ngVLA sensitivity is still above the predicted averaged value of the radio peak flux.

The scaling law used to compute the radio emission at the peak assume that the observing frequency remains in the same spectral regime of the slow cooling synchrotron spectrum, i.e. , and above the absorption frequency. Using our more detailed calculation of the radio afterglow from the core jet (Eq. 3 and text thereabove), we could check that this condition was fulfilled for the bulk of the population. However, for mergers in high density environments (larger than ), this condition is not longer met as the absorption frequency and the injection frequency may be larger than the radio frequency at early times. For an event with and (the GW mean angle), the radio frequency typically meets the injection (resp. absorption) break around 30 days (resp. 100 days). In this case, chromatic lightcurves are expected.

Figure 6: Left: Cumulative distribution of the peak times of the radio afterglow for jointly detected events in the O3+VLA configuration, assuming a jet with (solid line) or without (dashed line) lateral expansion. Middle: Differential distribution of the peak fluxes of GW-detected events assuming the horizon for the O2 (blue), O3 (red) and design (black) configurations of the GW interometers. The radio sensitivities for VLA and SKA/ngVLA are indicated for comparison. Right: Differential distribution of the maximum proper motion for jointly detected events assuming the O2+VLA (blue), O3+VLA (red), design+SKA (black). The dashed vertical line shows the lower limit for detection of the proper motion.

4.1.4 Proper motion

Another independent observable is the proper motion of the remnant, which can be readily accessible to VLBI measurements as illustrated in the case of GRB170817A (Mooley et al., 2018a; Ghirlanda et al., 2019). An upper limit of the proper motion is obtained assuming that the core jet contribution dominates at the peak of the radio light curve, and that the shocked region at the front of this core jet constitutes the visible remnant of the merger. In first approximation, afterglows from jets peak when the observer enters the focalization cone of the radiation, i.e. . For an off-axis observer, we therefore estimate the maximum proper motion to be:




The distribution of is shown in Fig. 6.

For a given angular resolution, the VLBI network can only measure proper motions above a limit . Thus, for a given value of the viewing angle, there is a maximum distance beyond which the measurement of the proper motion is not possible. It is given by


where we have adopted as/day, inspired by the uncertainties quoted in Mooley et al. (2018a). For O2 (resp. O3) the measurement of the proper motion is possible in only 72% (resp. 54%) of the cases.

We finally represent in Fig. 7 various plots connecting two observable quantities: , , and . The grey dots correspond to events detectable in GW and the red ones to those detectable in both GW and radio. In the diagram, the limiting distance for the measurement of the proper motion (Eq. 17) is shown. This figure illustrates the various observational biases discussed in this section, and especially the bias towards small viewing angles, mainly due to the limiting sensitivity of radio-telescopes.

Figure 7: Predicted populations of GW-detected (grey) and jointly detected (red) mergers in the planes of various pairs of observables, for the O3+VLA configuration. Upper left: Distance and viewing angle, as well as the maximum distance to which the proper motion can be measured (, see Eq. 17, blue line). Upper right: Peak flux and viewing angle. Lower left: Time of radio afterglow peak (assuming lateral expansion of the jet) and viewing angle. Lower right: Maximum remnant proper motion and distance, as well as the limiting proper motion (blue line).

4.1.5 Kinetic energy, external density and microphysics parameters

The distributions of the kinetic energy and the external density , which are not directly observable but may result from fits of the afterglow, are shown if Fig. 8. In each case, the distributions of jointly detected events are compared to the GW-only detections. As the GW selection is independent of the kinetic energy and density, these distributions are the intrinsic distributions of our population model (Tab. 1). As could be expected, the detection of mergers leading to a relativistic core jet with a large is favored. In detected mergers, the energy distribution approximately remains a broken power-law but the respective indices below and above the break decrease from to and to respectively. The external density distribution is simply shifted to higher densities, as the mean value is increased by a factor of . The distribution of the microphysics parameter is very similar to that of the density, because both parameters appear with the same power in the peak flux (Eq. 5).

We find that the Lorentz factor of the jet at the time of peak flux has a median value of (i.e. for a jet with seen with a mean viewing angle found for the O3+VLA configuration). This justifies to keep microphysics parameters representative of shock acceleration in the ultra-relativistic regime for the whole population.

Figure 8: Differential distribution of the core jet kinetic energy (left) and external medium density (right) of the GW-only (blue) and the jointly detected events (red), in the O3+VLA configuration.
Figure 9: Left: Fraction of joint events among GW events under the assumption of the WP15 short GRB luminosity function for different values of the minimum kinetic energy , and for our fiducial model (dashed lines). The results are plotted for different values of the jet opening angle , and rad (red, black and blue lines respectively). Middle: Mean viewing angle (same color coding). Right: Fraction of joint events with peaks earlier than 150 days post-merger, the time of peak being calculated with the hypothesis of jet lateral expansion.

4.2 Impact of the population model uncertainties on the predicted population

4.2.1 Jet energy distribution and opening angle

We now consider alternatives to our fiducial population model. With a distribution of jet energy that increases below the break (i.e. taking in Eq. 10) as would be obtained from the short GRB luminosity function of Wanderman & Piran (2015) () and still assuming the same standard values for the duration and efficiency, the fraction of events detected in radio would typically be three times smaller: 15% for O2 and 10% for O3 compared to 47% and 34% with . The viewing angle and peak times would also be strongly affected, showing peaks for at only (resp. ) and 96% (resp. 98%) of afterglows peaking before 150 days.

These results were obtained assuming a typical value for the gamma-ray efficiency but one cannot exclude that may vary from from less than 1% to more than 20%. Estimates of from afterglow fitting have been indeed found to cover a large interval (e.g. Lloyd-Ronning & Zhang, 2004; Zhang et al., 2007). Thus, the three orders of magnitude in isotropic gamma-ray luminosity could partially result from differences in efficiency, the jet isotropic kinetic energy being restricted to a smaller interval. Fig. 9 shows the effect of increasing the minimum kinetic energy from to  erg in the energy function deduced from the short GRB luminosity function of Wanderman & Piran (2015) on the detected fraction, the mean viewing angle and the fraction of afterglows peaking before 150 days. As could be expected, when the minimum kinetic energy is increased above erg, the results become closer to those of the fiducial model where the energy distribution function decreases below the peak.

We also show in Fig. 9 the effect of changing the typical value of the jet opening angle from 0.1 to 0.05 or 0.2 rad. When (resp. ), the detected fraction decreases (resp. increases) by about 50%, the mean viewing angle decreases (resp. increases) by , and the peak of the light curve takes place later (resp. earlier).

4.2.2 External density distribution and microphysics

The external medium acts as a revelator for the jet afterglow, and the peak flux is greatly dependent on its density, as can be seen from Eq. 5. By changing the central value of the distribution of external media, the fraction of joint events changes. This is illustrated in Tab. 4, where we represent this fraction for different external density distribution central values.

If there is indeed an excess of fast-merging NS binaries, as indicated by e.g. Matteucci et al. (2014); Hotokezaka et al. (2015); Vangioni et al. (2016); Beniamini & Piran (2016), these should merge in higher density environments, producing brighter afterglows and can notably contribute to the observed population even if there are intrinsically less numerous. With a density distribution centered at cm and the other parameters taking their fiducial values, 97% of the GW events produce detectable radio afterglows.

The microphysics parameter enters in Eq. 5 with the same power as the external density. Moreover our fiducial choices for the distribution of these two quantities are nearly the same (with simply a limit in range for ) so that changing its central value affects the detected fraction in the same way.

5 Discussion and conclusion

We studied the population of binary neutron star mergers to be observed jointly through GW and radio afterglow detection in future multi-messenger observing campaigns. For this we have assumed a likely population of mergers inspired by prior short GRB science, and simulated the GW and jet-dominated radio afterglow emissions from these. We have made predictions on the rates of such future events, and on how the observables from these events are expected to distribute themselves.

In the case of the ongoing O3 run of the LVC and assuming our fiducial set of parameters (Tab. 1), we predict that of GW events should have a detectable radio afterglow. These joint events should have mean viewing angles of , and should peak earlier than 150 days post-merger in of cases.

The fraction of these detectable events that can eventually be detected depends strongly on the accuracy of the localization, which is related to the size of the error box provided by the GW network, the possible detection of an associated GRB and the efficiency of the search for a kilonova, in terms of covered sky area and limiting magnitude.

-4 22 %
40 %
63 %
83 %
94 %
97 %
Table 4: Fraction of joint events among GW events for different intrinsic density distribution central values, for the O3+VLA detector configuration, as calculated with the semi-analytical model (Eq. 3).

Indeed, we have applied a threshold on the GW SNR and radio afterglow peak flux to label an event as jointly detected. Thus our study concerns the class of detectable events. As we know, a proper joint detection requires localization of the event, through the UV-IR detection of the kilonova counterpart. Also, once this localization is acquired, a continuous monitoring of the remnant up to may be necessary to detect the peak of the afterglow in the case of marginally detectable events. Even then, only events with peaks somewhat larger than the radio threshold should yield observations of astrophysical interest because extended observations of the rise, peak and decay of the afterglow are necessary to resolve the structure and dynamics (e.g. expansion) of the ultrarelativistic jet, as illustrated by the case of GRB170817A.

Furthermore, as also illustrated in this event, the lateral material may also contribute to the peak flux of the radio afterglow. Fig. 1 (right) shows that this may increase the peak flux by a factor of up to 1.5. In this case, our O3 predictions are revised to a fraction of 39% of events with detectable afterglow (instead of 34%), and a similar mean viewing angle of .

Also, VLBI measurements were instrumental to resolving the structured jet riddle in the case of GRB170817A. Thus an even more restrictive criterion for full event characterization could be the detectability of the remnant proper motion. We have shown in Sec. 4.1.4 that this would largely decrease the fraction of events, especially in the context of design GW detectors.

We have studied the sensitivity of the expected observables’ distributions to the population parameters as well as the radio and GW detector configurations. In particular, we have shown the uncertainties on the rates and typical viewing angles of these events inherited from the uncertainties on the luminosity function of short GRBs and of the density of the media hosting the merger.

In Sec. 4.1.1 and 4.1.2 we have discussed the combined influences of increasing the radio and GW detector sensitivity on the rate and viewing angle of the events. At this stage, both improvements would be beneficial as the predicted rate in the current O3+VLA is well below the critical 50% of joint events among GW events. Also, as the GW horizon increases, more events will be seen on-axis. In the case of an opening angle of 0.1 rad and in the likely short-term evolution of the detector configuration from O3+VLA (current) to design+VLA (early 2020s), we predict that 6% to 10% of events should be seen on-axis, increasing the probability of a GRB counterpart in addition to the GW and afterglow. This figure is fully consistent with that announced in Beniamini et al. (2019) for events with a GW+GRB association. In fact, as a GW+GRB association is equivalent to an on-axis event, one may measure the opening angle of typical GRB jets by considering the ratio of GW events with afterglow counterpart only to events with both afterglow and GRB, thus exploiting the new possibility of observing afterglows without detecting the GRB.

More generally, a population of observed events corresponds to an intrinsic population of mergers. Thus a comparison of a statistical number of joint event observations to our predicted distributions is a means of measuring fundamental parameters of the population of mergers, and of constraining GRB quantities such as their energy function.

Finally, future facilities such as the SKA may bring detections of orphan radio afterglows through deep radio surveys. These would not be subject to the GW criterion and may probe another sub-population of mergers, bringing yet another class of constraints on these phenomena.

In conclusion, regardless of the evolution of GW and radio detectors, mutli-wavelength afterglows will remain instrumental in the study of binary neutron star mergers as revelators of both their environment and their role as progenitors of short GRBs. They will bring precious insight at the level of the population of mergers and on an event-to-event basis, provided the sources may be accurately localized in the sky.


  • Abbott et al. (2018) Abbott, B. P., Abbott, R., Abbott, T. D., et al. 2018, Living Reviews in Relativity, 21, 3
  • Abbott et al. (2019) Abbott, B. P., Abbott, R., Abbott, T. D., et al. 2019, Physical Review X, 9, 011001
  • Abbott et al. (2017) Abbott, B. P., Abbott, R., Abbott, T. D., et al. 2017, Physical Review Letters, 119, 161101
  • Alexander et al. (2018) Alexander, K. D., Margutti, R., Blanchard, P. K., et al. 2018, ArXiv e-prints, arXiv:1805.02870
  • Arcavi et al. (2017) Arcavi, I., Hosseinzadeh, G., Howell, D. A., et al. 2017, Nature, 551, 64
  • Beniamini et al. (2016) Beniamini, P., Nava, L., & Piran, T. 2016, MNRAS, 461, 51
  • Beniamini et al. (2019) Beniamini, P., Petropoulou, M., Barniol Duran, R., & Giannios, D. 2019, MNRAS, 483, 840
  • Beniamini & Piran (2016) Beniamini, P. & Piran, T. 2016, MNRAS, 456, 4089
  • Beniamini & van der Horst (2017) Beniamini, P. & van der Horst, A. J. 2017, MNRAS, 472, 3161
  • Berger (2014) Berger, E. 2014, Annual Reviews in Astronomy and Astrophysics, 52, 43
  • Cantiello et al. (2018) Cantiello, M., Jensen, J. B., Blakeslee, J. P., et al. 2018, ApJ, 854, L31
  • Church et al. (2011) Church, R. P., Levan, A. J., Davies, M. B., & Tanvir, N. 2011, MNRAS, 413, 2004
  • Coulter et al. (2018) Coulter, D., Foley, R., Kilpatrick, C., et al. 2018, in APS April Meeting Abstracts, Vol. 2018, D16.002
  • Cowperthwaite et al. (2017) Cowperthwaite, P. S., Berger, E., Villar, V. A., et al. 2017, ApJ, 848, L17
  • D’Avanzo et al. (2018) D’Avanzo, P., Campana, S., Salafia, O. S., et al. 2018, A&A, 613, L1
  • D’Avanzo et al. (2014) D’Avanzo, P., Salvaterra, R., Bernardini, M. G., et al. 2014, MNRAS, 442, 2342
  • Dobie et al. (2018) Dobie, D., Kaplan, D. L., Murphy, T., et al. 2018, ApJ, 858, L15
  • Drout et al. (2017) Drout, M. R., Piro, A. L., Shappee, B. J., et al. 2017, Science, 358, 1570
  • Duffell et al. (2018) Duffell, P. C., Quataert, E., Kasen, D., & Klion, H. 2018, ApJ, 866, 3
  • Evans et al. (2017) Evans, P. A., Cenko, S. B., Kennea, J. A., et al. 2017, Science, 358, 1565
  • Finn & Chernoff (1993) Finn, L. S. & Chernoff, D. F. 1993, Phys. Rev. D, 47, 2198
  • Fong et al. (2010) Fong, W., Berger, E., & Fox, D. B. 2010, ApJ, 708, 9
  • Fong et al. (2015) Fong, W., Berger, E., Margutti, R., & Zauderer, B. A. 2015, ApJ, 815, 102
  • Ghirlanda et al. (2019) Ghirlanda, G., Salafia, O. S., Paragi, Z., et al. 2019, Science, 363, 968
  • Ghirlanda et al. (2016) Ghirlanda, G., Salafia, O. S., Pescalli, A., et al. 2016, A&A, 594, A84
  • Gill & Granot (2018) Gill, R. & Granot, J. 2018, MNRAS
  • Goldstein et al. (2017) Goldstein, A., Veres, P., Burns, E., et al. 2017, ApJ, 848, L14
  • Gottlieb et al. (2018) Gottlieb, O., Nakar, E., & Piran, T. 2018, MNRAS, 473, 576
  • Gottlieb et al. (2019) Gottlieb, O., Nakar, E., & Piran, T. 2019, arXiv e-prints, arXiv:1903.08173
  • Granot et al. (2002) Granot, J., Panaitescu, A., Kumar, P., & Woosley, S. E. 2002, ApJ, 570, L61
  • Haggard et al. (2017) Haggard, D., Nynka, M., Ruan, J. J., et al. 2017, ApJ, 848, L25
  • Hallinan et al. (2017) Hallinan, G., Corsi, A., Mooley, K. P., et al. 2017, Science, 358, 1579
  • Hotokezaka et al. (2018) Hotokezaka, K., Kiuchi, K., Shibata, M., Nakar, E., & Piran, T. 2018, ApJ, 867, 95
  • Hotokezaka et al. (2015) Hotokezaka, K., Piran, T., & Paul, M. 2015, Nature Physics, 11, 1042
  • Kasliwal et al. (2017) Kasliwal, M. M., Nakar, E., Singer, L. P., et al. 2017, Science, 358, 1559
  • Lamb & Kobayashi (2017) Lamb, G. P. & Kobayashi, S. 2017, MNRAS, 472, 4953
  • Lazzati et al. (2017a) Lazzati, D., Deich, A., Morsony, B. J., & Workman, J. C. 2017a, MNRAS, 471, 1652
  • Lazzati et al. (2017b) Lazzati, D., López-Cámara, D., Cantiello, M., et al. 2017b, ApJ, 848, L6
  • Lipunov et al. (2017) Lipunov, V. M., Gorbovskoy, E., Kornilov, V. G., et al. 2017, ApJ, 850, L1
  • Lloyd-Ronning & Zhang (2004) Lloyd-Ronning, N. M. & Zhang, B. 2004, ApJ, 613, 477
  • Lyman et al. (2018) Lyman, J. D., Lamb, G. P., Levan, A. J., et al. 2018, Nature Astronomy, 2, 751
  • Margutti et al. (2018) Margutti, R., Alexander, K. D., Xie, X., et al. 2018, ApJ, 856, L18
  • Margutti et al. (2017) Margutti, R., Berger, E., Fong, W., et al. 2017, ApJ, 848, L20
  • Matteucci et al. (2014) Matteucci, F., Romano, D., Arcones, A., Korobkin, O., & Rosswog, S. 2014, MNRAS, 438, 2177
  • Mooley et al. (2018a) Mooley, K. P., Deller, A. T., Gottlieb, O., et al. 2018a, Nature, 561, 355
  • Mooley et al. (2018b) Mooley, K. P., Frail, D. A., Dobie, D., et al. 2018b, ApJ, 868, L11
  • Mooley et al. (2018c) Mooley, K. P., Nakar, E., Hotokezaka, K., et al. 2018c, Nature, 554, 207
  • Nakar (2007) Nakar, E. 2007, Phys. Rep, 442, 166
  • Nakar et al. (2018) Nakar, E., Gottlieb, O., Piran, T., Kasliwal, M. M., & Hallinan, G. 2018, ApJ, 867, 18
  • Nakar & Piran (2018) Nakar, E. & Piran, T. 2018, MNRAS, 478, 407
  • Nakar et al. (2002) Nakar, E., Piran, T., & Granot, J. 2002, ApJ, 579, 699
  • Palmese et al. (2017) Palmese, A., Hartley, W., Tarsitano, F., et al. 2017, ApJ, 849, L34
  • Panaitescu & Kumar (2001) Panaitescu, A. & Kumar, P. 2001, ApJ, 560, L49
  • Pian et al. (2017) Pian, E., D’Avanzo, P., Benetti, S., et al. 2017, Nature, 551, 67
  • Resmi et al. (2018) Resmi, L., Schulze, S., Ishwara-Chandra, C. H., et al. 2018, ApJ, 867, 57
  • Saleem et al. (2018a) Saleem, M., Pai, A., Misra, K., Resmi, L., & Arun, K. G. 2018a, MNRAS, 475, 699
  • Saleem et al. (2018b) Saleem, M., Resmi, L., Misra, K., Pai, A., & Arun, K. G. 2018b, MNRAS, 474, 5340
  • Santana et al. (2014) Santana, R., Barniol Duran, R., & Kumar, P. 2014, ApJ, 785, 29
  • Sari et al. (1998) Sari, R., Piran, T., & Narayan, R. 1998, ApJ, 497, L17
  • Savchenko et al. (2017) Savchenko, V., Ferrigno, C., Kuulkers, E., et al. 2017, ApJ, 848, L15
  • Sironi et al. (2015) Sironi, L., Keshet, U., & Lemoine, M. 2015, Space Sci. Rev., 191, 519
  • Smartt et al. (2017) Smartt, S. J., Chen, T.-W., Jerkstrand, A., et al. 2017, Nature, 551, 75
  • Soares-Santos et al. (2017) Soares-Santos, M., Holz, D. E., Annis, J., et al. 2017, ApJ, 848, L16
  • Tanvir et al. (2017) Tanvir, N. R., Levan, A. J., González-Fernández, C., et al. 2017, ApJ, 848, L27
  • Troja et al. (2008) Troja, E., King, A. R., O’Brien, P. T., Lyons, N., & Cusumano, G. 2008, MNRAS, 385, L10
  • Troja et al. (2018) Troja, E., Piro, L., Ryan, G., et al. 2018, MNRAS, 478, L18
  • Troja et al. (2017) Troja, E., Piro, L., van Eerten, H., et al. 2017, Nature, 551, 71
  • Valenti et al. (2017) Valenti, S., David, Sand, J., et al. 2017, ApJ, 848, L24
  • van Eerten et al. (2018) van Eerten, E. T. H., Ryan, G., Ricci, R., et al. 2018, arXiv e-prints, arXiv:1808.06617
  • Vangioni et al. (2016) Vangioni, E., Goriely, S., Daigne, F., François, P., & Belczynski, K. 2016, MNRAS, 455, 17
  • Wanderman & Piran (2015) Wanderman, D. & Piran, T. 2015, MNRAS, 448, 3026
  • Wijers & Galama (1999) Wijers, R. A. M. J. & Galama, T. J. 1999, ApJ, 523, 177
  • Xie et al. (2018) Xie, X., Zrake, J., & MacFadyen, A. 2018, ApJ, 863, 58
  • Zhang et al. (2007) Zhang, B., Liang, E., Page, K. L., et al. 2007, ApJ, 655, 989
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description