The Electromagnetic Counterpart of the Binary Neutron Star Merger LIGO/Virgo GW170817. II. UV, Optical, and Near-infrared Light Curves and Comparison to Kilonova Models

The era of gravitational wave (GW) astronomy began on 2015 September 14 when the Advanced Laser Interferometer Gravitational-wave Observatory (LIGO) made the first direct detection of gravitational waves, resulting from the merger of a stellar mass binary black hole (BBH; GW150914; Abbott et al. 2016a). LIGO has since announced the detection of three additional BBH events (Abbott et al. 2016b, 2017a; LIGO Scientific Collaboration et al. 2017). There are currently no robust theoretical predictions for electromagnetic (EM) emission associated with such mergers.

By contrast, mergers involving at least one neutron star can produce a wide range of EM signals, spanning from gamma-rays to radio (e.g., Metzger & Berger 2012). In the optical/near-infrared (NIR) bands, the most promising counterpart is the kilonova (KN), a roughly isotropic thermal transient powered by the radioactive decay of rapid neutron capture (r-process) elements synthesized in the merger ejecta (Li & Paczyński 1998; Metzger et al. 2010; Roberts et al. 2011; Metzger & Berger 2012; Barnes & Kasen 2013; Tanaka & Hotokezaka 2013). The properties of the KN emission (luminosity, timescale, spectral peak) depend sensitively on the ejecta composition. For ejecta containing Fe-group or light r-process nuclei with atomic mass number $A\lesssim 140$, the KN emission is expected to peak at optical wavelengths at a luminosity ${L}_{{\rm{p}}}\sim {10}^{41}\mbox{--}{10}^{42}$ erg s⁻¹ on a short timescale of ${t}_{{\rm{p}}}\sim 1$ day (a so-called "blue" KN; Metzger et al. 2010; Roberts et al. 2011; Metzger & Fernández 2014). By contrast, for ejecta containing heavier lanthanide elements ($A\gtrsim 140$) the emission is predicted to peak at NIR wavelengths with ${L}_{{\rm{p}}}\sim {10}^{40}\mbox{--}{10}^{41}$ erg s⁻¹ over a longer timescale of ${t}_{{\rm{p}}}\sim 1$ week (a so-called "red" KN; Barnes & Kasen 2013; Kasen et al. 2013; Tanaka & Hotokezaka 2013).

The first direct detection of gravitational waves from the inspiral and merger of a binary neutron star (BNS) was made on 2017 August 17 (GW170817; LIGO Scientific Collaboration & Virgo Collaboration 2017a, 2017b, 2017c). This source was coincident with a short burst of gamma-rays detected by both Fermi/Gamma-ray Burst Monitor (GBM) and INTEGRAL (GRB 170817A; Blackburn et al. 2017; Goldstein et al. 2017a, 2017b; Savchenko et al. 2017a, 2017b; von Kienlin et al. 2017). Rapid optical follow-up by our Dark Energy Camera (DECam) program (Flaugher et al. 2015), starting just 11.4 hr after the GW trigger, led to the discovery of an associated optical counterpart in the nearby ($d\approx 39.5$ Mpc; Freedman et al. 2001) galaxy NGC 4993 (Allam et al. 2017; Soares-Santos et al. 2017). This optical source was independently discovered by several groups (Abbott et al. 2017b), and first announced as SSS17a by Coulter et al. (2017a, 2017b). The source has also been independently named DLT17ck (Valenti et al. 2017; Yang et al. 2017) and AT2017gfo.

Here we present rapid-cadence UV, optical, and NIR observations spanning from the time of discovery to 18.5 days post-merger. We construct well-sampled light curves and spectral energy distributions (SEDs) using data from DECam, along with Gemini-South/FLAMINGOS-2 (GS/F2) and Hubble Space Telescope (HST). We show that the data cannot be fit by a model with heating from ⁵⁶Ni radioactive decay and Fe-peak opacities (as in normal supernovae), but instead requires heating from r-process nuclei and at least two components consistent with lanthanide-poor and lanthanide-rich opacities. We further use the data to determine the ejecta masses and velocities for each component.

All magnitudes presented in this work are given in the AB system and corrected for Galactic reddening⁷⁷ with $E(B-V)=0.105$, applying the calibration of Schlafly & Finkbeiner (2011). We assume a negligible reddening contribution from the host (Blanchard et al. 2017).

A summary of the observations and photometry described in this section is available in the Appendix.

2.1. DECam

2.2. HST

2.3. GS/F2

2.4. Swift/UVOT

3.1. Light Curves

Our UV/optical/NIR light curves are shown in Figure 1. The data span from 0.47 to 18.5 days post-merger, with bluer bands fading below the detection limits at earlier times. The light curve coverage was truncated by the proximity of the source to the Sun. We first note that the light curves are not well described by a power law, indicating minimal contribution from a GRB optical afterglow over the timescale of our observations. This is consistent with modeling of the afterglow based on X-ray and radio observations (Margutti et al. 2017; Alexander et al. 2017).

Zoom In Zoom Out Reset image size

The light curves exhibit a rapid decline in the bluest bands (ug), an intermediate decline rate in the red optical bands (rizY), and a shallow decline in the NIR (HK_s). However, while the u- and g-band light curves decline by ≈2 mag day⁻¹ starting with the earliest observations, the redder optical bands exhibit a more complex behavior: they exhibit a comparatively slow decline (≈0.3 mag day⁻¹) over the first 1.5 days, develop a shoulder at about 4 days, and subsequently begin to decline at about 8 days.

We find a similar rapid evolution in the colors of the transient (Figure 2). In particular, the u − g and g − r colors become redder by about 1 mag between about 1.5 and 3.5 days. The colors in the redder optical bands exhibit slower evolution, with $r-i\approx 0.5\mbox{--}1$ mag, $i-z\approx 0\mbox{--}0.5$ mag, and $z-Y\approx 0.3$ mag. These colors are significantly redder than those of known supernovae near explosion (e.g., Folatelli et al. 2010; Bianco et al. 2014; Galbany et al. 2016).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

3.2. Spectral Energy Distribution

3.3. Bolometric Light Curve

3.4. Qualitative Comparisons to Kilonova Emission

We test the conjecture that the UV/optical/NIR transient is an r-process KN by fitting several isotropic, one-zone, gray opacity models to the light curves. For each model, we assume a blackbody SED which evolves assuming a constant ejecta velocity until it has reached a minimum temperature, at which point the photosphere has receded into the ejecta and the temperature no longer evolves. A similar temperature "floor" is predicted in Barnes & Kasen (2013), and we include this minimum temperature as a fitted parameter. We additionally fit for a "scatter" term, added in quadrature to all photometric errors, which roughly accounts for additional systematic uncertainties that are not included in our model.

We use MOSFiT⁸² (Guillochon et al. 2017; Nicholl et al. 2017), an open-source light curve fitting tool that utilizes a Markov Chain Monte Carlo (MCMC) to sample the model posterior. For each model, we ensure convergence by enforcing a Gelman-Rubin statistic $\lt 1.1$ (Gelman & Rubin 1992). We compare models using the Watanabe–Akaike Information Criteria (WAIC; Watanabe 2010; Gelman et al. 2014), which accounts for both the likelihood score and the number of fitted parameters. The best-fit parameters, uncertainties, and WAIC scores for each model are provided in Table 1.

Table 1.  Kilonova Model Fits

Model	${M}_{\mathrm{ej}}^{\mathrm{blue}}$	${v}_{\mathrm{ej}}^{\mathrm{blue}}$	${\kappa }^{\mathrm{blue}}$	${M}_{\mathrm{ej}}^{\mathrm{purple}}$	${v}_{\mathrm{ej}}^{\mathrm{purple}}$	${\kappa }^{\mathrm{purple}}$	${M}_{\mathrm{ej}}^{\mathrm{red}}$	${v}_{\mathrm{ej}}^{\mathrm{red}}$	${\kappa }^{\mathrm{red}}$	${f}^{\mathrm{Ni}}$	WAIC
	(${M}_{\odot }$)	(c)	(${\mathrm{cm}}^{2}{{\rm{g}}}^{-1}$)	(${M}_{\odot }$)	(c)	(${\mathrm{cm}}^{2}\,{{\rm{g}}}^{-1}$)	(${M}_{\odot }$)	(c)	(${\mathrm{cm}}^{2}\,{{\rm{g}}}^{-1}$)
2-Comp	${0.014}_{-0.001}^{+0.002}$	${0.266}_{-0.002}^{+0.007}$	(0.5)	⋯	⋯	⋯	${0.036}_{-0.002}^{+0.001}$	${0.123}_{-0.014}^{+0.012}$	${3.349}_{-0.337}^{+0.364}$	⋯	−102
3-Comp	${0.014}_{-0.001}^{+0.002}$	${0.267}_{-0.011}^{+0.006}$	(0.5)	${0.034}_{-0.002}^{+0.002}$	${0.110}_{-0.010}^{+0.011}$	(3.0)	${0.010}_{-0.001}^{+0.002}$	${0.160}_{-0.025}^{+0.030}$	(10.0)	⋯	−106
56Ni	${0.008}_{-0.001}^{+0.007}$	${0.260}_{-0.031}^{+0.034}$	(0.1)	⋯	⋯	⋯	⋯	⋯	⋯	${0.749}_{-0.203}^{+0.214}$	17
Blue	${0.032}_{-0.004}^{+0.002}$	${0.180}_{-0.002}^{+0.002}$	(0.1)	⋯	⋯	⋯	⋯	⋯	⋯	⋯	17
Red	⋯	⋯	⋯	⋯	⋯	⋯	${0.026}_{-0.008}^{+0.010}$	${0.271}_{-0.002}^{+0.008}$	(10)	⋯	153
1-Comp	⋯	⋯	⋯	${0.040}_{-0.007}^{+0.002}$	${0.274}_{-0.093}^{+0.007}$	${0.817}_{-0.135}^{+0.146}$	⋯	⋯	⋯	⋯	11

Note. Model parameters and WAIC scores. Numbers in parentheses indicate fixed parameters of the model. The errors represent the 1σ confidence interval. Both the 2-component ("2-Comp") and 3-component ("3-Comp") models have significantly smaller WAIC scores (indicating better fits) compared to the four single-component models.

Download table as:  ASCIITypeset image

We first attempt a simple supernova model, namely heating by the radioactive decay of ⁵⁶Ni and Fe-peak opacity of $\kappa =0.1$ cm² g⁻¹ (see Villar et al. 2017). The model parameters are the ejecta mass and velocity, and the ⁵⁶Ni mass fraction in the ejecta (as well as the temperature floor and scatter). The best-fit model has ${M}_{\mathrm{ej}}\approx 0.01\,{M}_{\odot }$, ${v}_{\mathrm{ej}}\approx 0.26\,{\rm{c}}$, and ${f}^{\mathrm{Ni}}\approx 0.75$. The parameters are comparable to those we inferred from blackbody fits to the flux and SEDs in the previous section, but the overall fit is poor. In particular, this model severely underestimates the NIR light curves, while the high ⁵⁶Ni fraction is inconsistent with the optical spectra (Nicholl et al. 2017). We therefore conclude that the transient is not powered by the radioactive decay of ⁵⁶Ni.

We next turn to r-process heating, using the model outlined in Metzger (2017) and implemented in Villar et al. (2017). This model includes the ejecta mass, ejecta velocity, and opacity as fitted parameters (as well as the temperature floor and scatter). Within this context we first assume an Fe-peak opacity of κ = 0.1 cm² g⁻¹ (our "blue" model; e.g., as assumed historically in Li & Paczyński 1998) and fit for the ejecta mass and velocity. This model, with ${M}_{\mathrm{ej}}\approx 0.03\,{M}_{\odot }$ and ${v}_{\mathrm{ej}}\approx 0.18\,{\rm{c}}$, adequately describes the early light curves (≲3 days), but again is a poor fit to the NIR light curves. More recent calculations indicate that lanthanide-rich ejecta are expected to have a much higher opacity of $\kappa =10$ cm² g⁻¹, leading to a "red" KN (e.g., Barnes & Kasen 2013). However, such a model (our "red" model), with best-fit values ${M}_{\mathrm{ej}}\approx 0.03\,{M}_{\odot }$ and ${v}_{\mathrm{ej}}\approx 0.27\,{\rm{c}}$, produces a poor fit to the data as well. In particular, the model light curves exhibit an initial rise for $\approx 4$ days, in contrast to the observed rapid decline at early times, especially in the UV and blue optical bands. Finally, we allow the opacity to vary as a free parameter, finding a best-fit value of $\kappa \approx 0.82$ cm² g⁻¹, and an associated ${M}_{\mathrm{ej}}\approx 0.04\,{M}_{\odot }$ and ${v}_{\mathrm{ej}}\approx 0.27\,{\rm{c}}$. However, this model again fails to reproduce the initial rapid decline in the UV, as well as the NIR light curves. We therefore conclude that r-process heating with a single value for the opacity cannot explain the observed light-curve evolution and colors. The final light curves for these models can be seen in Figure 3.

Inspired by the multi-component observed SED (Figure 2) and by the failure of single-component models to capture both the early rapid decline and the late-time red colors, we explore two multi-component models: (i) a two-component "blue" ($\kappa =0.5$ cm² g⁻¹) plus "red" (κ as a free parameter) model; and (ii) a three-component "blue" ($\kappa =0.5$ cm² g⁻¹) plus "purple" ($\kappa =3$ cm² g⁻¹) plus "red" ($\kappa =10$ cm² g⁻¹) model. These values were recently shown by Tanaka et al. (2017) to roughly capture the detailed opacity from radiative transfer simulations. For each component, we leave ${M}_{\mathrm{ej}}$ and ${v}_{\mathrm{ej}}$ as free parameters.

First, we explore the two-component model (with eight free parameters); we vary the ejecta masses, ejecta velocities, and temperature floors, the red component opacity, and a single scatter term. We find that the "blue" component has ${M}_{\mathrm{ej}}^{\mathrm{blue}}\approx 0.01\,{M}_{\odot }$ and ${v}_{\mathrm{ej}}^{\mathrm{blue}}\approx 0.27{\rm{c}}$ (with errors of roughly 10%), in good agreement with our inference from the SED at early times (Section 3.2). The "red" component has a much larger mass of ${M}_{\mathrm{ej}}^{\mathrm{red}}\approx 0.04\,{M}_{\odot }$ but a slower velocity of ${v}_{\mathrm{ej}}^{\mathrm{red}}\approx 0.12\,{\rm{c}}$. The best-fit opacity of this component is κ ≈ 3.3 cm² g⁻¹, lower than expected for lanthanide-rich ejecta. We find that most of the parameters are uncorrelated, with the exception of the red component's opacity and ejecta velocity, which have a Pearson correlation coefficient of ∼0.67. The resulting parameters and uncertainties from the MCMC fitting are summarized in Table 1.

For the three-component model (with 10 free parameters) we find similar values for the "blue" component (${M}_{\mathrm{ej}}^{\mathrm{blue}}\approx 0.01\,{M}_{\odot }$ and ${v}_{\mathrm{ej}}^{\mathrm{blue}}\approx 0.27\,{\rm{c}}$) and the "purple" component (${M}_{\mathrm{ej}}^{\mathrm{purple}}\,\approx 0.03\,{M}_{\odot }$ and ${v}_{\mathrm{ej}}^{\mathrm{purple}}\approx 0.11\,{\rm{c}}$). The "red" component is sub-dominant with ${M}_{\mathrm{ej}}^{\mathrm{red}}\approx 0.01\,{M}_{\odot }$ and ${v}_{\mathrm{ej}}^{\mathrm{red}}\approx 0.16\,{\rm{c}}$); see Table 1. These ejecta parameters are consistent with those determined from independent modeling of the optical and NIR spectra (Chornock et al. 2017; Nicholl et al. 2017).

Both sets of models are shown in Figure 1 and are essentially indistinguishable. Both provide a much better fit to the data than the single-component models described above, capturing both the initial blue colors and rapid decline, as well as the later redder colors and NIR light curves. Their similar WAIC scores suggest that neither model is statistically preferred. The two models differ most drastically at $\lesssim 5$ days in the K_s-band, where the two-component model is double-peaked, while the three-component model is single peaked. While neither model fully captures every feature of the light curves, it is remarkable that these simplified semi-analytic models produce such high-quality fits over a wide range of wavelength and time.

In the multi-component models, we can interpret each component as arising from distinct physical regions within the merger ejecta. In both models, the high velocity of the blue KN ejecta suggests that it originates from the shock-heated polar region created when the neutron stars collide (e.g., Oechslin et al. 2007; Bauswein et al. 2013; Sekiguchi et al. 2016). This dominant blue component is also seen in early-time optical spectra (Nicholl et al. 2017). By contrast, the low velocity red KN component in our three-component model could originate from the dynamically ejected tidal tails in the equatorial plane of the binary (e.g., Rosswog et al. 1999; Hotokezaka et al. 2013), in which case the relatively high ejecta mass $\approx 0.01\,{M}_{\odot }$ suggests an asymmetric mass ratio of the merging binary ($q\lesssim 0.8;$ Hotokezaka et al. 2013).

In both multi-component models we find that the κ ≈ 3 cm² g⁻¹ ejecta dominates by mass. The lower velocity of this component suggests an origin in the post-merger accretion disk outflow. Our inferred ejecta mass is consistent with that expected for a massive ∼0.1 ${M}_{\odot }$ torus (e.g., Just et al. 2015; Siegel & Metzger 2017). Similarly, the disk outflow composition is predicted to be dominated by ${Y}_{e}\sim 0.3$ matter that produces the $\kappa \approx 3$ cm² g⁻¹ component of the KN emission (Tanaka et al. 2017) as we observe. The fitted opacity indicates that the hyper-massive neutron star remnant is relatively short-lived (∼30 ms; Fernández & Metzger 2013; Just et al. 2015; Kasen et al. 2015). We additionally find that in both models the total kinetic energy is roughly $(1\mbox{--}2)\times {10}^{51}$ erg.

The fact that our multi-component models fit the data well provides strong evidence for the production of both light and heavy r-process nuclei, addressing one of the long-standing mysteries in astrophysics (Burbidge et al. 1957; Cameron 1957). We quantify this statement by comparing our blue and red ejecta masses to those necessary to reproduce the Milky Way (MW) r-process production rate. For heavy r-process elements (red KN), the MW inferred production rate is ${\dot{M}}_{\mathrm{rp},{\rm{A}}\gtrsim 140}\approx {10}^{-7}\,{M}_{\odot }$ yr⁻¹ (Bauswein et al. 2014). For light r-process elements (blue KN), the production rate is ${\dot{M}}_{\mathrm{rp},{\rm{A}}\gtrsim 100}\approx 7\times {10}^{-7}\,{M}_{\odot }$ yr⁻¹ (Qian 2000). Using a conservative estimate on the local BNS merger rate estimated by Abbott et al. (2016), ${R}_{0}\approx 1000$ Gpc⁻³ yr⁻¹, and a volume density of MW-like galaxies of $\approx 0.01$ Mpc⁻³, we estimate the MW rate of KN as ${R}_{\mathrm{MW}}\approx 100$ Myr⁻¹. Using this MW rate, we find that the average ejecta mass for a red KN is ${\dot{M}}_{\mathrm{rp},{\rm{A}}\gtrsim 140}/{R}_{\mathrm{MW}}\approx 0.001\,{M}_{\odot }$ and for a blue KN it is ${\dot{M}}_{\mathrm{rp},{\rm{A}}\gtrsim 100}/{R}_{\mathrm{MW}}\approx 0.007\,{M}_{\odot }$. These order-of-magnitude estimates are smaller than our inferred ejecta masses for this event, although the discrepancy can potentially be mitigated when properly taking into account the fraction of r-process materials that remains in a gas phase in the ISM and galactic halo. Nevertheless, this exercise suggests that BNS mergers can reproduce the r-process yields found in the MW and may be a dominant source of cosmic r-process nucleosynthesis.

We have presented a comprehensive UV, optical, and NIR data set for the first electromagnetic counterpart to be associated with a gravitational wave event. Analysis of these data reveals that the emission is due to an r-process-powered KN consisting of both "blue" and "purple/red" ejecta components. Models with ⁵⁶Ni heating, Fe-peak opacities, or a single component of r-process opacity fail to match the observations.

Our models indicate that the total ejecta mass is $\approx 0.05\,{M}_{\odot }$, with a high velocity ($v\approx 0.3c$) blue component and a slower ($v\approx 0.1\mbox{--}0.2c$) purple/red component. The presence of both components and the relatively large ejecta mass suggests that binary neutron star mergers (like GW170817) dominate the cosmic r-process nucleosynthesis.

The data presented in this paper (and others in this series) represent by far the best observations of an r-process powered KN, and it is remarkable how well the observations match theoretical models. This event also marks the true beginning of joint GW--EM multi-messenger astronomy. We expect that this event will serve as a benchmark for future efforts to model and understand the behavior of these transients, and for the first time allow the development of data-driven KN models. The next Advanced LIGO/Virgo observing run (starting in Fall 2018) is expected to detect many more BNS events (Abbott et al. 2016c). The follow-up of these events will provide further understanding of the ubiquity of the features seen in this event, the relationship between event and host properties, and place even stronger constraints on r-process enrichment from BNS mergers.

The Berger Time-Domain Group at Harvard is supported in part by the NSF through grants AST-1411763 and AST-1714498, and by NASA through grants NNX15AE50G and NNX16AC22G. P.S.C. is grateful for support provided by the NSF through the Graduate Research Fellowship Program, grant DGE 1144152. V.A.V. acknowledges support by the National Science Foundation through a Graduate Research Fellowship. D.A.B. is supposed by NSF award PHY-1707954. We thank the University of Copenhagen, DARK Cosmology Centre, and the Niels Bohr International Academy for hosting R.J.F. during the discovery of GW170817/SSS17a, where he was participating in the Kavli Summer Program in Astrophysics, "Astrophysics with gravitational wave detections." This program was supported by the Kavli Foundation, Danish National Research Foundation, the Niels Bohr International Academy, and the DARK Cosmology Centre. The UCSC group is supported in part by NSF grant AST–1518052, the Gordon & Betty Moore Foundation, the Heising-Simons Foundation, generous donations from many individuals through a UCSC Giving Day grant, and from fellowships from the Alfred P. Sloan Foundation and the David and Lucile Packard Foundation to R.J.F. E.B. acknowledges financial support from the European Research Council (ERC-StG-335936, CLUSTERS). D.K. and E.Q. were funded in part by the Gordon and Betty Moore Foundation through Grant GBMF5076.

We are grateful for the heroic efforts of the entire staff at Gemini-South in order to continue obtaining observations of this target in evening twilight at high airmass as the object was setting.

This project used data obtained with the Dark Energy Camera (DECam), which was constructed by the Dark Energy Survey (DES) collaboration. Funding for the DES Projects has been provided by the DOE and NSF (USA), MISE (Spain), STFC (UK), HEFCE (UK), NCSA (UIUC), KICP (U. Chicago), CCAPP (Ohio State), MIFPA (Texas A&M), CNPQ, FAPERJ, FINEP (Brazil), MINECO (Spain), DFG (Germany), and the collaborating institutions in the Dark Energy Survey, which are Argonne Lab, UC Santa Cruz, University of Cambridge, CIEMAT-Madrid, University of Chicago, University College London, DES-Brazil Consortium, University of Edinburgh, ETH Zürich, Fermilab, University of Illinois, ICE (IEEC-CSIC), IFAE Barcelona, Lawrence Berkeley Lab, LMU München, and the associated Excellence Cluster Universe, University of Michigan, NOAO, University of Nottingham, Ohio State University, University of Pennsylvania, University of Portsmouth, SLAC National Lab, Stanford University, University of Sussex, and Texas A&M University.

Based in part on observations at Cerro Tololo Inter-American Observatory, National Optical Astronomy Observatory (NOAO Prop. ID: 2017B-0110, PI: E. Berger), which is operated by the Association of Universities for Research in Astronomy (AURA) under a cooperative agreement with the National Science Foundation.

Based in part on observations obtained at the Gemini Observatory (Program IDs GS-2017B-Q-8 and GS-2017B-DD-4; PI: Chornock), which is operated by the Association of Universities for Research in Astronomy, Inc., under a cooperative agreement with the NSF on behalf of the Gemini partnership: the National Science Foundation (United States), the National Research Council (Canada), CONICYT (Chile), Ministerio de Ciencia, Tecnología e Innovación Productiva (Argentina), and Ministério da Ciência, Tecnologia e Inovação (Brazil).