A kilonova as the electromagnetic counterpart to a gravitational-wave source

S. J. Smartt, T. -W. Chen, A. Jerkstrand, M. Coughlin, E. Kankare, S. A. Sim, M. Fraser, C. Inserra, K. Maguire, K. C. Chambers, M. E. Huber, T. Kruhler, G. Leloudas, M. Magee, L. J. Shingles, K. W. Smith, D. R. Young, J. Tonry, R. Kotak, A. Gal-Yam, J. D. Lyman, D. S. Homan, C. Agliozzo, J. P. Anderson, C. R. Angus C. Ashall, C. Barbarino, F. E. Bauer, M. Berton, M. T. Botticella, M. Bulla, J. Bulger, G. Cannizzaro, Z. Cano, R. Cartier, A. Cikota, P. Clark, A. De Cia, M. Della Valle, L. Denneau, M. Dennefeld, L. Dessart, G. Dimitriadis, N. Elias-Rosa, R. E. Firth, H. Flewelling, A. Flors, A. Franckowiak, C. Frohmaier, L. Galbany, S. Gonzalez-Gaitan, J. Greiner, M. Gromadzki, A. Nicuesa Guelbenzu, C. P. Gutierrez, A. Hamanowicz, L. Hanlon, J. Harmanen, K. E. Heintz, A. Heinze, M. -S. Hernandez, S. T. Hodgkin, I. M. Hook, L. Izzo, P. A. James, P. G. Jonker, W. E. Kerzendorf, S. Klose, Z. Kostrzewa-Rutkowska, M. Kowalski, M. Kromer, H. Kuncarayakti, A. Lawrence, T. B. Lowe, E. A. Magnier, I. Manulis, A. Martin-Carrillo, S. Mattila, O. McBrien, A. Muller, J. Nordin, D. O'Neill, F. Onori, J. T. Palmerio, A. Pastorello, F. Patat, G. Pignata, Ph. Podsiadlowski, M. L. Pumo, S. J. Prentice, A. Rau, A. Razza, A. Rest, T. Reynolds, R. Roy, A. J. Ruiter, K. A. Rybicki, L. Salmon, P. Schady, A. S. B. Schultz, T. Schweyer, I. R. Seitenzahl, M. Smith, J. Sollerman, B. Stalder, C. W. Stubbs, M. Sullivan, H. Szegedi, F. Taddia, S. Taubenberger, G. Terreran, B. van Soelen, J. Vos, R. J. Wainscoat, N. A. Walton, C. Waters, H. Weiland, M. Willman, P. Wiseman, D. E. Wright, L. Wyrzykowski, O. Yaron

astro-ph.HE

References

1 Distance and Reddening

The host galaxy NGC4993 has been identified as a member of a group of 10 galaxies (LGG 332). The heliocentric recessional velocity of $2951\pm 26$ $\rm{\,km\,s^{-1}}$ , or $z=0.009843\pm 0.000087$ , is from optical data. The kinematic distance (correcting for various infall models and using $H_{0}=71\pm 2$ $\rm{\,km\,s^{-1}}$ Mpc-1) and the Tully Fisher distances to the group containing NGC4993 are in good agreement within the uncertainty of $d=$ $40\pm 4$ Mpc (distance modulus $\mu=$ $33.01\pm 0.20$ ), and we adopt this value. The foreground reddening values in the direction of NGC4993 and AT2017gfo (as reported in NED The NASA/IPAC Extragalactic Database (NED) is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration.) are adopted to be $A_{U}=0.54,A_{g}=0.39,A_{r}=0.28,A_{i}=0.21A_{z}=0.16,A_{y}=0.13,A_{J}=0.09,A_{H}=0.06,A_{K}=0.04$ (Landolt $U$ , Pan-STARRS1 $grizy_{\rm P1}$ and UKIRT $JHK$ ), or $E(B-V)=0.11$ mag. These reddening corrections were applied to the photometry to calculate absolute magnitudes and bolometric luminosities.

2 Hubble Space Telescope pre-discovery data

NGC4993 was observed by the Hubble Space Telescope using the Advanced Camera for Surveys (ACS) Wide Field Channel on 2017 April 28, less than four months prior to the discovery of AT2017gfo. 2 $\times$ 348 s exposures were taken with the F606W filter (comparable to Sloan $r^{\prime}$ ). As this is the deepest image of the site of AT2017gfo taken prior to discovery, we examined it for any possible pre-discovery counterpart.

We localised the position of AT2017gfo on the ACS image by aligning this to the GROND $i^{\prime}$ images taken on each of the nights from 2017 Aug 18 to 21. Nine point sources common to both the GROND and ACS images were matched, and the final position on the ACS image has an uncertainty of 28, 50 mas in $x$ and $y$ respectively, determined from the scatter among the positions as measured on different GROND images.

No sources were detected by the DOLPHOT photometry package at a significance of 3 $\sigma$ or higher, within a radius $>3\times$ the positional uncertainty. We determined the limiting magnitude at the position of AT2017gfo to be F606W $>27.5$ (VEGAMAG), based on the average magnitude of sources detected at $3\sigma$ within a 100 $\times$ 100 pixel region centred on the position of AT2017gfo. For our adopted distance modulus and foreground reddening, this implies that any source at the position of AT2017gfo must have an absolute magnitude F606W $>-5.8$ .

3 ATLAS system and observational data and upper limit to the rate of kilonova events

The Asteroid Terrestrial-impact Last Alert System (ATLAS) , is a full-time near Earth asteroid survey. It is currently running two 0.5 m f/2 wide-field telescopes on Haleakala and Mauna Loa. The ATLAS sensor is a single thermoelectrically-cooled STA1600 detector with 1.86 arcsecond per pixel platescale (10.56k $\times$ 10.56k pixels) giving a 29.2 square degree field of view. The two units work in tandem to survey the entire visible sky from $-40^{\circ}<\delta<80^{\circ}$ with a cadence of two to four days, depending on weather. The ATLAS unit on Haleakala has been working in scientific survey mode since April 2016 and was joined by the Mauna Loa unit in March 2017.

ATLAS observes in two wide-band filters, called “cyan” or “ $c$ ”, which roughly covers the SDSS/Pan-STARRS $g$ and $r$ filters, and “orange” or “ $o$ ”, which roughly covers the SDSS/Pan-STARRS $r$ and $i$ . The observing cadence for identifying moving asteroids is typically to observe each footprint 4-5 times (30 s exposures, slightly dithered) within about an hour of the first observation of each field. All data immediately go through an automatic data processing pipeline. This produces de-trended, sky-flattened images which are astrometrically corrected to the Gaia stellar reference frame and photometrically corrected using Pan-STARRS1 reference stars. Difference images are produced using a static-sky template and source extraction is carried out on both the target and difference images using DOPHOT on the target frames and a custom written package for PSF-fitting photometry, which we call TPHOT (on the difference frames). Sources found on the difference images are then cataloged in a MySQL database and merged into astrophysical objects if there are at least three detections from the five (or more) images. These objects are subject to a set of quality filters, a machine-learning algorithm and human scanning.

Our database did not contain any astrophysical object at the position of AT2017gfo between MJD 57380.64463 and 57966.26370. The position was observed 414 times and on each of these we forced flux measurements at the astrometric position of the transient on the difference image. We measured 5 $\sigma$ flux limits and any epochs with greater than 5 $\sigma$ detections. The 5 $\sigma$ flux limits were in the range $o>18.6\pm 0.5$ (AB mag, median and standard deviation) and $c>19.3\pm 0.4$ (see Extended Data Figure A kilonova as the electromagnetic counterpart to a gravitational-wave source). We found 44 images which formally had flux detections greater than 5 $\sigma$ , but on visual inspection we rule out these being real flux variability at the transient position. They all appear to be residuals from the host galaxy subtraction. With ATLAS, we rule out any variability down to 18.6 to 19.3 (filter dependent) during a period 601 to 16 days before discovery of AT2017gfo.

We can estimate an approximate upper limit to the rates of these kilonovae, without a GW trigger from the ATLAS survey. Extended Data Figure A kilonova as the electromagnetic counterpart to a gravitational-wave source implies that we would be sensitive to objects like AT2017gfo to 60 Mpc. ATLAS typically surveys 5000 sq deg per night, 4-5 times, which provides a sampled volume of $10^{-4}$ Gpc3 within 60 Mpc. If we assume that a kilonova lightcurve is visible for 4 days and we have observations every 2 - 4 days, and observe 60% of clear time then the control time is 0.9 yr. We have no candidates, therefore the simple Poisson probabilities of obtaining a null result are 50%, 16% and 5% when the expected values are 0.7, 1.8 and 3.0 $\times 10^{4}$ Gpc ${}^{-3}\,$ yr-1. Therefore the 95% confidence upper limit to the rate of kilonovae is $<3.0\times 10^{4}$ Gpc-3 yr-1. This simple approach is in broad agreement with the upper limit from the Dark Energy Survey and the LIGO Scientific collaboration for NS-NS mergers and a more sophisticated calculation is warranted for the ATLAS data.

4 The Pan-STARRS1 system and observational data

The Pan-STARRS1 system comprises a 1.8 m telescope with a 1.4 Gigapixel camera (called GPC1) mounted at the Cassegrain $f/4.4$ focus. This wide-field system is located on the summit of Haleakala on the Hawaiian island of Maui. The GPC1 is composed of sixty Orthogonal Transfer Array devices (OTAs), each of which has a detector area of 48460 $\times$ 48680 pixels. The pixels are 10 microns in size (0.26 arcsec) giving a focal plane of 418.88 mm in diameter or 3.0 degrees. This corresponds to field-of-view area of 7.06 square degrees, and an active region of about 5 square degrees. The filter system (which we denote $grizy_{\rm P1}$ ) is similar to the SDSS and is described in detail in two papers. Images from Pan-STARRS1 are processed immediately with the Image Processing Pipeline. The existence of the Pan-STARRS1 3 $\pi$ Survey data provides a ready made template image of the whole sky north of $\delta=-30^{\circ}$ , and we furthermore have proprietary $i_{\rm P1}$ data in a band between $-40^{\circ}<\delta<-30^{\circ}$ , giving a reference sky in the $i_{\rm P1}$ band down to this lower declination limit. Images in $i_{\rm P1}<math xmlns="http://www.w3.org/1998/Math/MathML" display="block"><semantics><mrow><msub><mi>z</mi><mrow><mi mathvariant="normal">P</mi><mn>1</mn></mrow></msub></mrow><annotation encoding="application/x-tex">z_{\rm P1}</annotation></semantics></math>zP1y_{\rm P1}$ were taken on 7 nights, at high airmass due to the position of AT2017gfo.

A series of dithered exposures were taken in the three filters during the first available night (starting 2017 Aug 18 05:33:01 UT), and we placed the target on a clean detector cell. We repeated the $i_{\rm P1}<math xmlns="http://www.w3.org/1998/Math/MathML" display="block"><semantics><mrow><msub><mi>z</mi><mrow><mi mathvariant="normal">P</mi><mn>1</mn></mrow></msub></mrow><annotation encoding="application/x-tex">z_{\rm P1}</annotation></semantics></math>zP1y_{\rm P1}$ for two subsequent nights until the object became too low in twilight and we switched to $z_{\rm P1}$ and $y_{\rm P1}$ and then only $y_{\rm P1}$ . Frames were astrometrically and photometrically calibrated with standard Image Processing Pipeline steps. The Pan-STARRS1 3 $\pi$ reference sky images were subtracted from these frames and photometry carried out on the resulting difference image.

5 ePESSTO and Xshooter observational data

EFOSC2 consists of a combined 2048 $\times$ 2048 pixel CCD imaging camera and low-dispersion spectrograph, mounted at the Nasmyth focus of the 3.58 m New Technology Telescope (NTT) at La Silla, Chile. The SOFI instrument has a 1024 $\times$ 1024 pixel near infra-red array for long-slit spectroscopy and imaging, and is also mounted at the NTT on the other Nasmyth focus. All EFOSC2 spectra were taken at the parallactic angle using the configurations listed in Table A kilonova as the electromagnetic counterpart to a gravitational-wave source, and reduced using the PESSTO pipeline. Spectroscopic frames were trimmed, overscan and bias subtracted, and divided by a normalised flat field. In the case of the Gr#16 spectra, a flat field was obtained immediately after each spectrum to enable fringing in the red to be corrected. Spectra were wavelength calibrated using arc lamps, and the wavelength solution checked against strong sky emission lines. Cosmic rays were masked in the two-dimensional spectra using the LACosmic algorithm, before one-dimensional spectra were optimally extracted from each frame. Flux calibration of the spectra was done using an average sensitivity curve derived from observations of several spectrophotometric standard stars during each night, while the telluric features visible in the red were corrected using a synthetic model of the absorption.

The Xshooter instrument on the ESO Very Large Telescope was used for two epochs of spectra. The observational setup and spectral reductions were similar to those previously employed in and detailed in several publications, with the custom-built T. Krühler reduction pipeline used for the reduction and flux calibration and $molecfit$ package used for telluric correction. All spectra were scaled to contemporaneous photometric flux calibrations. Images with the NACO and VISIR instruments on the ESO Very Large Telescope were taken in the $L-$ band (NACO) and $N-$ band (VISIR) in the mid-infrared. These were kindly made public by ESO to all collaborating groups working with the LIGO-Virgo follow-up programmes and are publicly available through the ESO archive. We found no detection of the transient in either instrument. The host galaxy NGC4993 was faint, but visible in the $L-$ band NACO images. With only one standard star, at a vastly different airmass from the target we could not reliably determine an upper limit. Similarly, no flux was visible in the VISIR $N-$ band data.

The EFOSC2 and SOFI images were reduced using the PESSTO pipeline. All EFOSC2 images were overscan and bias subtracted, and divided by a flat-field frame created from images of the twilight sky. Individual images taken at each epoch were then aligned and stacked. The SOFI images were cross-talk and flat-field corrected, sky subtracted, aligned, and merged. The transient had faded below the detection limit in the $g^{\prime}r^{\prime}i^{\prime}z^{\prime}$ GROND images obtained on 2017 August 26.97 UT, and the $U$ EFOSC2 image observed on 2017 August 21.05 UT. The VISTA Hemisphere survey $JK_{\rm s}$ images observed on 2014 April 10 were used as references for the SOFI $JK_{\rm s}$ images. No VISTA archive images were available in $H-$ band, therefore we used the GROND $H-$ band on 2017 Aug 29.99 UT as the reference. Template image subtraction to remove the contribution from the host galaxy was carried out based on the ISIS2.2 package, and the subtractions were of good quality. Point-spread function (PSF) fitting photometry was carried out on each stacked and template-subtracted image. An empirical model of the PSF was made for each image from sources in the field, and fitted to the transient to determine its instrumental magnitude. In the case where the transient was not detected, artificial star tests were used to set a limiting magnitude. The photometric zeropoint for each image was determined through aperture photometry of Pan-STARRS1 or 2MASS sources in the field of the EFOSC2 and SOFI images, respectively, and used to calibrate the instrumental magnitudes onto a standard system. Three further epochs were taken with the Boyden 1.52-m telescope in South Africa, giving extra time resolution coverage over the first 72 hrs. The Boyden 1.52-m telescope, is a 1.52 m Cassegrain reflector combined with an Apogee 1152 $\times$ 770 pixel CCD imaging camera, providing a field of view of 3.7 arcmin $\times$ 2.5 arcmin. Observations were carried out during twilight and the early hours of the night at low altitude using 30 sec exposures. Observations were reduced and analysed using a custom pipeline for this telescope. All photometric observations were taken using a clear filter and then converted to SDSS r using four Pan-STARRS1 reference stars.

6 GROND system and observational data

Observations with GROND at the 2.2 m Max-Planck telescope at La Silla ESO started on 2017 Aug 18 23:15 UT. Simultaneous imaging in $g^{\prime}r^{\prime}i^{\prime}z^{\prime}JHK_{\rm s}$ continued daily, weather allowing until 2017 Sept 4 (see Tables A kilonova as the electromagnetic counterpart to a gravitational-wave source and A kilonova as the electromagnetic counterpart to a gravitational-wave source). GROND data were reduced in the standard manner using pyraf/IRAF. PSF photometry of field stars was calibrated against catalogued magnitudes from Pan-STARRS1 for $g^{\prime}r^{\prime}i^{\prime}z^{\prime}$ images and 2MASS for $JHK_{\rm s}$ images. The images were template subtracted using the ISIS2.2 package. GROND $g^{\prime}r^{\prime}i^{\prime}z^{\prime}$ images from 2017 Aug 26.97 UT and $JHK_{\rm s}$ images from 2017 Aug 29.99 UT were used as reference images. These were the best quality images we had with no detection of the source. The photometry results in typical absolute accuracies of $\pm$ 0.03 mag in $g^{\prime}r^{\prime}i^{\prime}z^{\prime}$ and $\pm$ 0.05 mag in $JHK_{\rm s}$ .

7 Spectral and lightcurve comparisons

A comparison of our spectra with a sample of SNe is shown in Extended Data Fig. A kilonova as the electromagnetic counterpart to a gravitational-wave source. Both the spectral shape and features present differ significantly, with AT2017gfo showing a significantly redder SED than those of either Type Ia or Type II-P SNe within a few days of explosion. The spectra of AT2017gfo also lack the typical absorption features of intermediate-mass elements that are normally seen in early-time SN spectra. Fig. A kilonova as the electromagnetic counterpart to a gravitational-wave source also shows a comparison with optical spectra from a sample of some of the faintest and fastest evolving Type I SN discovered to date.

PESSTO has spectroscopically classified 1160 transients, and monitored 264, and none are similar to AT2017gfo. Volume-limited samples of supernovae (having samples of around 100-200 SNe within 30-60 Mpc) have never uncovered a similar transient. In the ATLAS survey, during the period up to Aug 2017, we have found 75 transients (all supernovae) in galaxies within 60 Mpc and no objects like AT2017gfo. This implies that objects like AT2017gfo have a rate of around 1% or less of the local supernova rate, justifying our probability calculation in the main text.

8 Bolometric light curve calculation

Firstly, the broad-band magnitudes in the available bands ( $U$ , $g$ , $r$ , $i$ , $z$ , $y$ , $J$ , $H$ , $K_{s}$ ) were converted into fluxes at the effective filter wavelengths, and then corrected for the adopted extinctions (see Methods - Distance and reddening). For completeness at early phases, we ensured consistency with the values for ultra-violet flux reported from the Swift public data in the bands $uvw2$ , $uvm2$ , $uvw1$ , and $U$ . An SED was then computed over the wavelengths covered. Fluxes were converted to luminosities using the distance previously adopted. We determined the points on the bolometric light curve at epochs when $K$ -band or ultra-violet observations were available. Magnitudes from the missing bands were generally estimated by interpolating the light curves using low-order polynomials ( $n\leq 2$ ) between the nearest points in time. We also checked that the interpolated/extrapolated magnitudes were consistent with the available limits. Finally, we fitted the available SED with a black-body function and integrated the flux from 1000 Å to 25000 Å. This provides a reasonable approximation to the full bolometric light curve but we caution that flux beyond 25000Å may contribute. It is not clear that the spectral energy distribution at this phase is physically well represented by a black body, and therefore we chose not to integrate fully under such a spectrum. Therefore the bolometric flux that we estimate at 8 days and beyond could be higher. For reference we report the temperature and radius evolution, together with uncertainties, from the SED fitting in Table A kilonova as the electromagnetic counterpart to a gravitational-wave source, although we again note that a black body assumption may not be valid at later times.

9 Light curve modelling - parameter range estimation

We compare the light curve data with the models by Arnett and Metzger using a Bayesian framework . The likelihood in our case is defined as $\mathcal{L}=e^{-\chi^{2}/2}$ . The time of the kilonova (used on both models) is defined to be that of the gravitational-wave trigger time. For both the Metzger and Arnett models considered in this analysis, we choose a log uniform prior of $-5\leq\log_{10}(M_{\rm ej})\leq 0$ for the ejecta mass, a uniform prior of $0\leq v_{\rm ej}\leq 0.3$ c for the ejecta velocity, and a uniform prior of $-1\leq\log_{10}(\kappa)\leq 2$ $\textrm{cm}^{2}\textrm{g}^{-1}$ for the opacity. Specifically for the Metzger model, we choose a uniform prior of $0\leq\alpha\leq 10$ for the slope of the ejecta velocity distribution. The power-law slope for radioactive powering given in the Arnett model is given a prior of $-5\leq\beta\leq 5$ .

We sample this given posterior using a nested sampling approach using the MultiNest implementation through a Python wrapper. Figure A kilonova as the electromagnetic counterpart to a gravitational-wave source shows the posterior of the Arnett model. Figure A kilonova as the electromagnetic counterpart to a gravitational-wave source shows the posterior of the Metzger model.

Systematic error for mass is dominated by uncertainty in the heating rate per mass of the ejecta. This consists of the product of intrinsic decay power, and thermalization efficiency. For the intrinsic decay power, we find values of $(1-3)\times 10^{10}$ erg g-1 s-1 in the literature , $1.9\times 10^{10}$ erg g-1 s-1 is our default value. There are only small uncertainties associated with nuclear mass models during the first few days, but this grows to a factor of $\sim$ 2-3 at later times.

Due to the dominance of the post-diffusion tail in the fits, the mass scales roughly inversely with the powering level. Thus, if this is a factor of two higher than assumed our mass range declines by a factor of 2. However, the vast majority of decay models are close to our value, so we favour the $\sim$ 0.04 $M_{\odot}$ solutions over the $\sim$ 0.02 $M_{\odot}$ ones. We note also that even the high-opacity models fitting the later data points have $M\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{$ \sim $}}\hss}\hbox{$ > $}}}0.02\,M_{\odot}$ , so this should be a robust lower limit to the ejecta mass.

10 TARDIS Modelling details

For the temperature implied by the black-body like SED, Fe would be expected to be primarily in its neutral or singly-ionized state: in either case, detectable features would be expected. In particular, the lack of evidence for Fe ii features (e.g. the Fe ii $\lambda$ 5064 multiplet) in the blue part of our spectrum places a strong limit on the presence of this ion (simple TARDIS modelling suggests $<10^{-3}M_{\odot}$ of Fe ii can be present in the spectral forming region). This lack of Fe partly argues against ejecta compositions dominated by Fe-peak elements. Equivalent constraints on Ni, however, are weaker.

As noted in the main text, the combination of limited atomic data and simplistic modelling means that we cannot derive reliable elemental masses from the analysis carried out so far. However, we note that our model for the $+$ 1.4 d spectrum invokes ion masses of only $\sim 10^{-9}$ M⊙ and a few times $10^{-3}$ M⊙ for Cs i and Te i, respectively, at ejecta velocities above the adopted photosphere (i.e. $v>0.2$ c). In both cases, these are only lower limits on elemental masses, since the ions in question are expected to be sub-dominant at the conditions present in the ejecta (this is a particularly important consideration for Cs i, owing to its low ionization potential of only 3.9 eV). Nethertheless, these mass limits are consistent with the ejecta masses suggested in our light curve model.

11 Kilonova simulations

Kilonova simulations predict two distinct ejecta components: dynamic ejecta and disk winds. The dynamic ejecta is expelled directly in the merger. Starting from neutron star material with $Y_{e}\sim 0.03$ , it experiences some moderated de-neutronization by positron captures, but likely ends with $Y_{e}\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{$ \sim $}}\hss}\hbox{$ < $}}}0.2$ (as described in simulations). Such composition is predicted to produce all heavy r-process elements, including lanthanides and actinides. It is thus expected to lead to a high-opacity red component peaking on time scales of days/weeks. The disk wind has two components, a radiation driven wind and a dynamic torus ejection. These are exposed to neutrino irradiation, which can produce a larger variation in $Y_{e}$ . This component can thus be largely lanthanide and actinide free, and have low opacity, in particular for $0.2<Y_{e}<0.4$ . Dynamic and wind ejecta have similar heating rates. Thus, their contribution to the bolometric light curve is largely proportional to their masses. The compilation by Wu et al. shows that current simulations predict similar masses of the two components, but uncertainty of a factor few for their mass ratio.

The data suggest that we have detected the lower-opacity disk wind component, and that this has a $Y_{e}$ in the range giving low opacity (giving constraints on the poorly understood $Y_{e}$ setting processes). Whether a dynamic component is present as well is harder to ascertain. The whole light curve is reasonably well fit by a single disk wind component. Our models are too simplistic to warrant exploration of two-component scenarios. Assuming we have detected a disk wind of several times 0.01 $M_{\odot}$ , it is not easy to make this component drop away enough at late times to leave much flux for a dynamic ejecta component. Perhaps the opacity in the dynamic ejecta is as high ( $\kappa\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{$ \sim $}}\hss}\hbox{$ > $}}}100$ cm2 g-1) as speculated , and it then remains too dim to be seen compared to the wind for at least the first 20 days. Alternatively, this kilonova may simply have $M_{wind}\gg M_{ejecta}$ . The only circumstance which could substantially change these conclusions is if the first 2–3 data points are caused by a GRB afterglow. Then, a dynamic component with $\kappa\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{$ \sim $}}\hss}\hbox{$ > $}}}10$ cm2 g -1, can reasonably well fit the later data points. However, as we discuss in the main paper we find several arguments against this scenario, such as the chromatic light curve evolution and the absence of a strong X-ray afterglow, and assess that the early light is caused by a blue kilonova.

The reduced, calibrated spectral data presented in this paper are openly available on the Weizmann Interactive Supernova data REPository (https://wiserep.weizmann.ac.il) and at the ePESSTO project website http://www.pessto.org. The raw data from the VLT, NTT and GROND (for spectra and imaging) are available from the ESO Science Archive facility http://archive.eso.org. The raw pixel data from Pan-STARRS1 and the 1.5m Boyden telescope are available from the authors on request.

The lightcurve fitting code described here is publicly available at the following website: https://star.pst.qub.ac.uk/wiki/doku.php/users/ajerkstrand/start. A code to produce the posteriors in this paper is available at: https://github.com/mcoughlin/gwemlightcurves. TARDIS is an open-source Monte Carlo radiative-transfer spectral synthesis code for 1D models of supernova ejecta and is publicly available here https://tardis.readthedocs.io/en/latest/. Standard software within the IRAF environment was used to carry out the spectral, and imaging reductions and photometry.

This work is based on observations collected at the European Organisation for Astronomical Research in the Southern Hemisphere, Chile as part of ePESSTO, (the extended Public ESO Spectroscopic Survey for Transient Objects Survey) ESO program 199.D-0143 and 099.D-0376. We thank ESO staff for their excellent support at La Silla and Paranal and for making the NACO and VISIR data public to LIGO-Virgo collaborating scientists. We thank Jacob Ward for permitting a time switch on the NTT. PS1 and ATLAS are supported by NASA Grants NNX08AR22G, NNX12AR65G, NNX14AM74G and NNX12AR55G. Part of the funding for GROND was generously granted from the Leibniz-Prize to Prof. G. Hasinger (DFG grant HA 1850/28-1). Pan-STARRS1 and ATLAS are supported by NASA Grants NNX08AR22G, NNX12AR65G, NNX14AM74G and NNX12AR55G issued through the SSO Near Earth Object Observations Program. We acknowledge the excellent help in obtaining GROND data from Angela Hempel, Markus Rabus and Régis Lachaume on La Silla. The Pan-STARRS1 Surveys (PS1) were made possible by the IfA, University of Hawaii, the Pan-STARRS Project Office, the Max-Planck Society, MPIA Heidelberg and MPE Garching, Johns Hopkins University, Durham University, University of Edinburgh, Queen’s University Belfast, Harvard-Smithsonian Center for Astrophysics, Las Cumbres Observatory Global Telescope Network Incorporated, National Central University of Taiwan, Space Telescope Science Institute, National Science Foundation under Grant No. AST-1238877, University of Maryland, and Eotvos Lorand University (ELTE) and the Los Alamos National Laboratory. We acknowledge EU/FP7-ERC Grants , and STFC funding through grant ST/P000312/1 and ERF ST/M005348/1. AJ acknowledges Marie Sklodowska-Curie grant No 702538. MG, AH, KAR and ŁW thank Polish NCN grant OPUS 2015/17/B/ST9/03167, JS is funded by Knut and Alice Wallenberg Foundation. CB, MDV., NE-R., AP and GT are supported by the PRIN-INAF 2014. MC is supported by the David and Ellen Lee Prize Postdoctoral Fellowship at the California Institute of Technology. MF is supported by a Royal Society - Science Foundation Ireland University Research Fellowship. MS and CI acknowledge support from EU/FP7-ERC grant no . PGJ acknowledges the ERC consolidator grant number . GREAT is funded by VR. JDL gratefully acknowledges STFC grant ST/P000495/1. TWC, PS and PW acknowledge the support through the Alexander von Humboldt Sofja Kovalevskaja Award. JH acknowledges financial support from the Vilho, Yrjö and Kalle Väisälä Foundation. JV acknowledges FONDECYT grant number 3160504. LG was supported in part by the US National Science Foundation under Grant AST-1311862. MB acknowledges support from the Swedish Research Council and the Swedish Space Board. AG-Y is supported by the EU via ERC grant No. 725161, the Quantum Universe I-Core program, the ISF, BSF and by a Kimmel award. LS acknowledges IRC grant GOIPG/2017/1525. AJR is supported by the Australian Research Council Centre of Excellence for All-sky Astrophysics (CAASTRO) through project number CE110001020. IRS was supported by the Australian Research Council Grant FT160100028. We acknowledge Millennium Science Initiative grant IC120009. This paper uses observations obtained at the Boyden Observatory, University of the Free State, South Africa.