No Surviving Companion in Kepler's Supernova

The supernova of 1604, observed by Johannes Kepler and other European, Korean, and Chinese astronomers, is one of the five "historical" supernovae (SNe) that have been classified as belonging to Type Ia (thermonuclear), the other four being SN 1572 (Tycho Brahe's SN), SN 1006, SN 185 (supposed to have created the remnant RCW86), and the recently discovered youngest SN Ia G1.9+03 that occurred in our Galaxy as recently as around 1900 but was not discovered due to dust extinction and being observable only from the Southern Hemisphere.

As is well known, Type Ia SNe (hereafter SNe Ia) are well explained by the thermonuclear explosion of a C+O white dwarf (WD) star in a close binary system. There are three proposed channels to bring the WD to the point of explosion, depending on the nature of the companion star. In the single-degenerate (SD) channel, the companion is a star that is still evolving thermonuclearly (Whelan & Iben 1973; Nomoto 1982); in the double-degenerate (DD) channel it is another WD, either a C+O WD or an He WD (Iben & Tutukov 1984; Webbink 1984). Another possible channel, known as the core-degenerate (CD) scenario, involves a C+O WD and its merging with the core of an asymptotic giant branch (AGB) star, following a common-envelope episode (Livio & Riess 2003; Soker et al. 2014; Aznar-Siguán et al. 2015).

Significant progress has been made in the identification of progenitors of SNe Ia. For instance, there has been proof of a specific scenario that works to give rise to SNe Ia. This is the double detonation scenario studied theoretically by Fink et al. (2010), Sim et al. (2012), and others. In this scenario, the CO WD accumulates a helium-rich layer on in its surface. The detonation of the helium-rich layer ignites the CO WD. This seems to be the explosion mechanism involved in MUSSES1604D (Jiang et al. 2017) and similar events. An He WD companion seems to be favored. The donor He-rich WD might survive in particular cases studied by Shen & Schwab (2017).

However, the observed double detonation scenario, as seen from the effect in the very early light curves of the SNe Ia, cannot account for more than a small percentage of the SNe Ia observed (Jiang et al. 2017).

Another possible path to SNe Ia is a WD–WD collision. In this new DD path, no surviving companion is expected. A study of this mechanism to give rise to SNe Ia shows that it would account for <1% of the observed events (see Soker 2018 for an overview).

So, mainly in the single-degenerate scenario, a surviving companion should remain after the explosion. There is no surviving companion, but merging of the two components of the system in the DD scenario and in the CD scenario. A surviving companion might be identified from its kinematics (large radial velocity and/or proper motion, fast rotation), anomalous luminosity, or contamination of its surface layers by the SN ejecta (Wang & Han 2012 and Ruiz-Lapuente 2014 review those effects). Detailed simulations can be found in Marietta et al. (2000), Podsiadlowski (2003), Pakmor et al. (2008), Pan et al. (2012, 2013, 2014), Liu et al. (2012, 2013), Shappee et al. (2013), and Shen & Schwab (2017). We will compare observations with their predictions.

The central regions of the supernova remnant (SNR) of Tycho's SN (Ruiz-Lapuente et al. 2004; González Hernández et al. 2009; Kerzendorf et al. 2009, 2013; Bedin et al. 2014) and of SN 1006 (González Hernández et al. 2012; Kerzendorf et al. 2012, 2018) have already been explored, to search for a possible surviving companion of the SN, as have extragalactic remnants such as SNR 0509-67.5 (Schaefer & Pagnotta 2012), SNR 0509-68.7 (Edwards et al. 2012), and N103B (Pagnotta & Schaefer 2015; Li et al. 2017). Studies of other SNRs are in progress or have been proposed. Through comparison of the work done by various authors on those SNRs, the double-degenerate scenario seems favored in several SNe Ia.

The classification of Kepler's SN, SN 1604, as an SN Ia was a matter of debate for a long time; some authors classified it as a core-collapse SN, in spite of its position, quite above the Galactic plane. The question has been settled by X-ray observations of the remnant (Cassam-Chenaï et al. 2004), showing an O/Fe ratio characteristic of SNe Ia (Reynolds et al. 2007).

There are indications (Vink 2008) that one component of the binary system giving rise to the SN might have created a detached circumstellar shell with a mass ∼1 M_⊙, expanding into the interstellar medium. More recently, Katsuda et al. (2015) have deduced that the shell should have lost contact with the binary years before the explosion. It has been suggested (Chiotellis et al. 2012; Vink 2017) that the companion star was an AGB star that lost its envelope at the time of the explosion.

The distance to the remnant has also been the subject of discussion, the estimates ranging between 3 and 7 kpc. Thus, Reynoso & Goss (1999), based on the H i absorption toward the remnant, estimated 4.8 < d < 6.2 kpc. Later, Sankrit et al. (2005) found d = ${3.9}_{-0.9}^{+1.4}$ kpc from the proper motion of the optical filaments. But very recently, Sankrit et al. (2016) have revisited their method and give d = ${5.1}_{-0.7}^{+0.8}$ kpc. Even more recently, Ruiz-Lapuente (2017), from the reconstruction of the optical light curve of the SN based on the historical records, also infers a distance d = 5.0 ± 0.61 kpc, in agreement with Sankrit et al. (2016). We thus adopt here a distance d ∼ 5.0 ± 0.7 kpc to Kepler's SN. At that distance, given the Galactic latitude of the SNR, b = 68, it lies ≃590 pc above the Galactic plane.

The aim of this paper is to address the progenitor system that led to Kepler's supernova, SN 1604. A first paper on the possible progenitor of the Kepler supernova suggested a marginal possibility that there was a donor, but only tentatively (Kerzendorf et al. 2014, K14 hereafter). At that time, the Hubble Space Telescope (HST) proper motions were not analyzed and the stellar parameters of the stars were unknown. Here we provide a complete analysis of a survey using the FLAMES instrument at the SN ESO VLT-UT2 and we add all the proper-motion information from HST with a baseline of 10 years.

The present paper is organized as follows. Section 2 describes the search in Kepler's SNR and what can be obtained from it. Section 3 describes the observations made with the Very Large Telescope (VLT) using the FLAMES instrumentation and the reduction of those observations. It also presents the proper motions obtained from data from the HST archive, from programs GO-9731 and GO-12885 (PI: Sankrit). Section 4 presents the method used to derive the stellar parameters and the results. Section 5 presents the estimated distances to the stars, and discusses the radial velocities obtained, comparing them to those in previous studies. Section 6 compares the candidate stars with a kinematical model of the Galaxy. Section 7 discusses the results and Section 8 provides a summary of the conclusions.

Our survey has a limiting apparent magnitude m_R = +19 mag. The visual extinction, A_V, in the direction to the remnant of Kepler's SN is A_V = 2.7 ± 0.1 mag, and A_R/A_V = 0.748. Thus, we have reached down to an absolute magnitude M_R = +3.4 mag. That corresponds to a luminosity L = 2.6 L_⊙. For the spectroscopic observations we used FLAMES (Pasquini et al. 2002) mounted to the UT2 of the VLT. For the measurement of the proper motions, we used archival data from the HST.

As a comparison, Kerzendorf et al. (2014) performed a shallower survey of possible survivors down to L > 10 L_⊙, according to them (L > 6 L_⊙, if we take our reliable, newly determined distance to SN 1604).

The remnant of SN 1604 has an average angular diameter of 225 arcsec. Our survey is complete down to m_R ≤ 19 within 24 arcsec of the center of the SNR (≃20% of its radius: blue circumference in Figure 1(a)) at α_J2000 = 17^h30^m4125, δ_J2000 =−21°29'3295 (Vink 2008). Additional fibers were used to extend the search beyond 20% of its radius (the green circumference encompasses 38 arcsec of the radius), although the supplementary stars are not very relevant, due to their distances from the center of the SNR. The radius of the search area of 24 arcsec, at a distance of 5 kpc, corresponds to a transverse displacement from the center of the SNR by 0.58 pc. That is the path that a possible companion star would have travelled in 400 years, moving at v = 1460 km s⁻¹ perpendicularly to the line of sight. A total of 32 stars were observed. They are listed, with their coordinates, in Table 1.

Zoom In Zoom Out Reset image size

While preparing the final version of this work, a new analysis of the X-ray knots of Kepler's SN by Sato & Hughes (2017) has provided new estimates of the expansion center. Both are very close to the center that we used, by Vink (2008), thus it does not impact the results of the stars included in our search. There is an estimate that does not take into account a possible deceleration of the knots. This places the center at α_J2000 = 17^h30^m41189 ± 36 and δ_J2000 = −21°29'2463 ±35. The center that takes into account a deceleration coincides practically with that of Vink (2008). The newly determined center taking into account a model for the deceleration of the knots is α_J2000 = 17^h30^m41321 ± 44 and δ_J2000 = −21°29'3051 ± 43. We include these new two centers in our Figure 1(a). We include our search area in relation to the whole SNR in Figure 1(b).

K14 explored the central region of Kepler's SNR. The search has been photometric and spectroscopic, covering a square field of 38'' × 38'' around the center of the SNR determined by Katsuda et al. (2008), down to m_V ≃ 18 mag. They have used, for their spectroscopy, the 2.3 m telescope of the Australian National University (ANU), and archival HST images for the photometry. The WiFeS-spectrograph is an image slicer with 25 slitlets of size 38'' × 1'' and 05 sampling in the spatial direction on the detector. They chose this instrument for its large field of view. However, that instrumentation did not allow them to determine the stellar parameters. Apparently they noted that, due to sky subtraction errors, the continuum placement in their data was unreliable. Without determined stellar parameters, it is not possible to estimate distances, because the absolute magnitudes of the stars then remain unknown. They also had problems in estimating rotational velocities to better than 200 km s⁻¹, due to the resolution and quality of the spectra (K14).

Here, we present a study that includes the stellar parameters (T_eff, log g, and [Fe/H]), and we also have rotational velocities and radial velocities, apart from the proper motions from HST. The conclusion on the supernova companion is thoroughly tested.

3.1. Spectral Observations and Reductions

3.2. Proper Motions

To derive proper motions we used data from two HST programs collected at two epochs separated by almost 10 years.

The first epoch is the data set from GO-9731 (PI: Sankrit), and it was collected in 2003 August 28–29 with the Wide Field Channel (WFC) of the Advanced Camera for Surveys (ACS). The images are in three narrow-band filters F502N, F660N, F658N, and in one medium-band filter, F550M.

The second epoch is from GO-12885 (also PI: Sankrit), and it was collected in 2013 July 1–3 with the UV-VISual Channel (UVIS) of the Wide Field Camera 3 (WFC3). The images are in F336W, F438W, F547M, F656N, F658N, and F814W. We refer here to Table 2 for more details.

Table 2.  HST Data Used in This Work

ACS/WFC	Epoch 1	2003.65842–2003.66141
F502N	4 × ∼1100 s	GO-9731
F502N	4 × ∼500 s	GO-9731
F660N	8 × ∼500 s	GO-9731
F658N	4 × ∼1250 s	GO-9731
F658N	4 × ∼600 s	GO-9731
F550M	4 ×  210 s	GO-9731
WFC3/UVIS	Epoch 2	2013.50023–2013.50666
F336W	4 × ∼750 s	GO-12885
F438W	4 × ∼270 s	GO-12885
F547M	4 × ∼470 s	GO-12885
F656N	8 × ∼1220 s	GO-12885
F658N	4 × ∼970 s	GO-12885
F814W	4 × ∼240 s	GO-12885

Download table as:  ASCIITypeset image

To derive proper motions we used 28 images in each of the two epochs 10 years apart. In Table 2 we give a log of the images used.

As the charge transfer efficiency (CTE) losses have a major impact on astrometric projects (Anderson & Bedin 2010), in this work every single ACS/WFC and WFC3/UVIS image employed was treated with the pixel-based correction for imperfect CTE developed by Anderson & Bedin (2010). The improved corrections are directly included in the MAST data products.¹³

Raw positions and fluxes were extracted in every image using the software and the spatially variable effective point-spread functions (PSFs) produced by Anderson & King (2006). All libraries of PSFs for both ACS/WFC and WFC3/UVIS are publicly available.¹⁴

The raw positions were then corrected for the average geometric distortion of these instruments using the prescriptions by Anderson (2002, 2007) for ACS/WFC and by Bellini & Bedin (2009) and Bellini et al. (2011) for WFC3/UVIS.

The methodology and the procedures followed to derive proper motions were extensively described in a previous work with similar goals (Bedin et al. 2014). In the following we give a brief description.

We first defined a reference frame with respect to which we will measure all the relative positions. This was built using the four images in filter WFC3/UVIS/F814W from the second epoch, because they have the best image quality and highest number of stars with high signal. In the reference frame we used only relatively bright, unsaturated, isolated stars, with a stellar profile (the quality fit q described in Anderson et al. 2008).

The positions from all these suitable stars measured in the four F814W images were then transformed into a common distortion-corrected reference system and their clipped means taken as the final positions of the reference frame. The field is a relatively dense one, and the final reference frame contains over 20,000 stars, giving a typical distance of ∼30 pixels to one such reference star. The consistency in the positions of the four F814W images gives a positional accuracy of ∼0.01 pixels for the brightest unsaturated sources, perfectly consistent with the accuracy given in the geometric distortion of Bellini et al. (2011), and as achieved in other work (for example, Bedin et al. 2013).

To calibrate our reference frame we need to derive the astrometric zero-points, plate scale, orientation, and skew terms, which are needed to bring its coordinate systems into an absolute astrometric system. To achieve this we used 266 sources in common between our reference frame and the Two Micron All Sky Survey (2MASS) catalog (Skrutskie et al. 2006).

Hereafter our absolute coordinates will be given in equatorial coordinates at equinox J2000, with positions given at the reference epoch of the reference frame, 2013.50. Finally, we note that our absolute accuracy will be the same as that of the 2MASS catalog, about 0.2 arcsec, while our relative precisions should be much better than that, down to a few times 0.1 mas.

As in Bedin et al. (2014), we have also produced image stacks, which often provide a useful representation of the astronomical scene. There are 4 + 6 such stacks, one per filter/epoch, and they give a critical inspection of the region surrounding each star at any given epoch.

We also produced trichromatic exposures for both epochs.

As extensively explained in Bedin et al. (2014), even the best geometric distortion available is always an average solution. There are always deviations from that, typically as large as several milliarcseconds. To remove these residual systematic errors in our geometric distortion, we used the "bore-sight" correction described in detail in Bedin et al. (2014), which is essentially a local approach (Bedin et al. 2003). We used a network of at least 15 stars at no more than 500 pixels from target stars, used only stars with consistent positions between the two epochs better than 0.75 pixels, and used stars brighter than at least 250 photoelectrons above the local sky.

Note that although most of the field objects are background Galactic bulge stars at ∼8 kpc, with an intrinsic dispersion of about ∼3 mas yr⁻¹ (e.g., Bedin et al. 2003), our proper-motion precision in the catalog is significantly superior to this limit. Indeed, motions of the reference local network are averaged over N reference stars, implying systematics lowered by the factor ∼1/$\sqrt{(N-1)}$, with N typically a few hundred and hard-coded to be at least 15.

To compute proper motions, we divide the measured displacements of sources between epoch 2 and epoch 1 by the time baseline after having transformed the displacements in the astrometric reference frame provided by the sources in our reference frame and 2MASS.

The proper motions are shown in Figure 2 and in Table 3. The uncertainties were computed as the sum in quadrature of the rms of positions observed within each of the two epochs.

Zoom In Zoom Out Reset image size

Table 3.  Proper Motions of the Targeted Stars

Name	dRAa (mas yr−1)	dDCb (mas yr−1)	sRAc (mas yr−1)	sDCd (mas yr−1)
T01	0.96333105	3.42910032	0.12081037	0.07071602
T02	−2.33240031	3.86571195	45.22749329	54.34668660
T03	1.65005841	−0.00972573	0.03339045	0.08114443
T04	−4.85102345	−1.31896956	0.10271324	0.21623776
T05	2.78755994	0.98653555	0.06256493	0.08762565
T06	−0.14035101	0.80693946	0.06580537	0.08934278
T07	2.22968225	0.63820899	0.08311504	0.09635703
T08	−0.64593882	0.48014508	45.22705027	54.34598682
T09	1.15779958	0.38160776	0.11226873	0.09974361
T10	−1.17098183	−1.93729322	0.06072457	0.10494832
T11	1.07273546	1.85830604	0.07653632	0.12997579
T12	−1.84466939	−2.21716030	0.05498542	0.08621992
T13	0.82006479	3.46512561	0.05292586	0.06884338
T14	−3.85388310	1.59710464	0.10770873	0.10911011
T14b	2.19854365	−1.00379486	0.03274381	0.09097165
T15	−5.93215420	1.39310415	0.03482210	0.08710643
T16	0.29501700	5.03375451	0.02070014	0.05802598
T17	0.31106430	0.16393080	0.06300106	0.07123875
T18	25.06353304	−43.80675171	0.02033413	0.08582691
T19	−3.14334844	−4.36036079	0.05860900	0.08929045
T20	−0.81856961	1.70037126	0.06403092	0.04892372
T21	−3.56172476	−2.59434617	0.05555466	0.04844967
T22	1.63830144	−2.74355827	0.05628767	0.07366438
T23	7.54610667	−1.28547796	0.04192214	0.05333274
T24	−0.54450471	3.25364568	0.08436434	0.06721759
T25	0.16070871	−0.27626484	45.22705027	54.34598682
T26	−3.39922823	0.75818276	0.08496424	0.07949148
T27	2.52236058	4.58473159	0.04709384	0.08656006
T28	2.09224795	2.28248031	0.03391564	0.04773037
T29	2.81992271	−3.87317107	0.08391351	0.07822889
T30	2.37310197	7.78087136	45.22739769	54.34606942
T31	−2.61147871	−2.23889319	0.01148096	0.08836027

Notes.

^adRA stands for relative proper motion in R.A. with respect to the reference frame (see the text). ^bdDC stands for relative proper motion in decl. ^csRA is the standard error (1σ) in relative proper motion in R.A. ^dsDC is the standard error (1σ) in relative proper motion in decl.

Download table as:  ASCIITypeset image

To assess the completeness of the measured sample in our field we calibrated our instrumental magnitudes and then performed artificial star tests.

The photometric calibration was performed only in filters ACS/WFC/F550M (narrow V) and WFC3/UVIS/F814W (wide I) and was obtained following the detailed procedures described in Bedin et al. (2005). We used the most updated zero-points, which were adopted following STScI instructions.¹⁵ The calibrated color–magnitude diagram of the detected sources is shown in left panel of Figure 3.

Zoom In Zoom Out Reset image size

Artificial star tests were performed as described in great detail in Anderson et al. (2008) for ACS/WFC, and more recently applied to WFC3/UVIS detectors (e.g., Bedin et al. 2015). Briefly, we built a fiducial line of the "main" main sequence (MS) of the Galactic field stars as representative of the sources in the field (middle panel of Figure 3), and added those sources in the individual images, which were reduced in an identical manner to the real sources. The ratio between the number of recovered stars and the number of inserted artificial stars at various instrumental magnitude levels then provides the completeness curve, which is shown in the right panel of Figure 3.

We can be reasonably confident that our sample is virtually complete down to m_F814W ≃ 22.5, and 50% complete down to ∼23.4. However, completeness always has a statistical value, as individual stars can be missed for several reasons.

The stellar parameters T_eff, log g, and [Fe/H] have been determined from the spectra obtained with the Giraffe spectrograph in the HR15n setup and from the UVES spectra in the 580 nm setup, through the spectroscopic indices defined by Damiani et al. (2014). The method is based on a set of narrow-band spectral indices, and calibrated stellar parameters are derived from suitable combinations of these indices. Those indices sample the amplitude of the TiO bands, the Hα core and wings, and temperature- and gravity-sensitive sets of lines at several wavelength intervals. The latter are the close group of lines between 6490 and 6500 Å, which are sensitive to gravity, and Fe i lines falling within the range covered by the HR15 setup, which are sensitive to temperature. Further gravity-sensitive features are found in the region 6750–6780 Å.

Two global indices, τ (sensitive to temperature) and γ (sensitive to gravity), are computed from the former ones. A further composite index, ζ, measures stellar metallicity.

Tests and calibrations of those indices have been performed (Damiani et al. 2014), based on photometry and reference spectra from the UVES Paranal Observatory Project (Bagnulo et al. 2003) and the ELODIE 3.1 Library (Prugniel & Soubiran 2001).

The method works well for stars in the approximate temperature range 3000 K ≲ T_eff ≲ 9000 K. The values of T_eff, log g, and [Fe/H], with their errors, for our targeted stars are given in columns 4–9 of Table 4. There were four stars (T17, T19, T25, and T26) whose surface gravities could not be determined in this way, due to the low S/N of their spectra. Their distances are thus left largely undetermined. We know their effective temperatures, however, and so we can estimate the distance depending on the luminosity classes to which they can belong.

The resulting set of stellar parameters suggests that the stars come from a rather ordinary mixture of field stars (mostly giants). A few of the stars seem to have low [Fe/H] (<−1), although with large errors; they are all consistent with being metal-poor giants. The radial velocities and the rotational v sin i values are very well determined. Radial velocities were measured from several lines in the spectra of the targeted stars. All spectra gave a clear and narrow peak in the cross-correlation function (CCF) with template spectra. Even when spectral lines are poorly defined, the CCF may be sharply peaked, since all lines add up. Radial and rotational velocities were measured from the CCF, as above, using a set of template spectra covering the relevant T_eff range, taken from the Gaia–ESO sample studied in Damiani et al. (2014); the template giving the highest CCF peak for each program star was used for the velocity determinations, an approach well-tested within the Gaia–ESO Survey. The radial velocities and v sin i values are very well determined and their errors (i.e., uncertainties on CCF peak center and width) are lower in percentage than the errors in the stellar parameters. Uncertainties on radial and rotational velocities based on Giraffe data were carefully studied by Jackson et al. (2015), and found to be dependent on S/N, T_eff, and v sin i (besides, obviously, spectral resolution). For the ranges of these parameters relevant to this work, this implies typical uncertainties of 1–2 km s⁻¹ on radial velocities, and 10–15 km s⁻¹ on v sin i.

T18 is a clear outlier, having very fast motion across the line of sight. However, this star is just a cold M-dwarf, located only about half a kiloparsec away.

The distances to the targeted stars are determined by comparison of the stellar parameters with the observed apparent magnitudes in different bands. We use, for that, the isochrones of Marigo et al. (2017) (see Figure 4 for an example), which give, for each combination of the parameters T_eff, log g, and [Fe/H], the absolute magnitudes M_V, M_R, M_J, M_H, and M_K. Apparent magnitudes of the stars in our sample are taken from the NOMAD catalog, and compared with the absolute magnitudes in the different photometric bands, taking into account the corresponding extinctions. The results are given in Table 5. The distances indicated, with their errors, are weighted averages over the different bands.

Zoom In Zoom Out Reset image size

Table 5.  Absolute Magnitudes of the Targeted Stars, Observed Apparent Magnitudes, and Inferred Distances

Star	MV	mV	MR	mR	MJ	mJ	MH	mH	MK	mK	d (kpc)
T01	${1.1}_{-2.2}^{+1.5}$	⋯	${0.5}_{-2.2}^{+1.4}$	17.9	$-{0.5}_{-2.1}^{+1.3}$	15.4	$-{1.0}_{-2.1}^{+1.2}$	14.8	$-{1.1}_{-2.0}^{+1.3}$	14.5	${11.5}_{-5.1}^{+17}$
T02	$-{4.3}_{-0.7}^{+6.0}$	⋯	$-{5.2}_{-0.7}^{+6.4}$	16.0	$-{7.1}_{-0.6}^{+7.1}$	⋯	$-{7.9}_{-0.4}^{+7.3}$	⋯	$-{8.1}_{-0.4}^{+7.4}$	⋯	>20
T03	${1.2}_{-2.8}^{+2.2}$	⋯	${0.5}_{-2.8}^{+2.2}$	16.2	$-{0.8}_{-2.7}^{+2.1}$	14.3	$-{1.5}_{-2.6}^{+2.0}$	13.6	$-{1.6}_{-2.6}^{+2.0}$	13.3	${8.2}_{-2.3}^{+20}$
T04	${3.3}_{-0.7}^{+1.9}$	⋯	${2.9}_{-0.9}^{+1.9}$	15.9	${2.2}_{-1.2}^{+1.9}$	15.1	${1.7}_{-1.3}^{+1.9}$	14.5	${1.6}_{-1.3}^{+1.9}$	14.5	${3.3}_{-1.9}^{+2.7}$
T05	$-{4.6}_{-1.8}^{+4.4}$	⋯	$-{5.4}_{-1.8}^{+4.3}$	17.5	$-{7.0}_{-1.6}^{+4.1}$	14.9	$-{7.7}_{-1.5}^{+4.0}$	14.1	$-{7.9}_{-1.5}^{+4.0}$	14.0	>20
T06	${0.8}_{-2.6}^{+2.4}$	16.3	${0.1}_{-2.5}^{+2.3}$	⋯	$-{1.3}_{-2.4}^{+2.3}$	13.7	$-{2.0}_{-2.5}^{+2.5}$	13.0	$-{2.1}_{-2.5}^{+2.3}$	12.7	${8.0}_{-5.3}^{+17}$
T07	${1.5}_{-1.1}^{+1.9}$	⋯	${1.2}_{-1.1}^{+1.8}$	19.3	${0.6}_{-1.0}^{+1.8}$	⋯	${0.3}_{-1.0}^{+1.8}$	⋯	${0.3}_{-1.0}^{+1.8}$	⋯	16 ${}_{-9}^{+11}$
T08	${0.6}_{-1.0}^{+1.1}$	15.9	$-{0.3}_{-0.7}^{+0.7}$	12.8	$-{1.1}_{-0.9}^{+0.9}$	13.1	$-{1.9}_{-0.7}^{+0.7}$	12.3	$-{2.0}_{-0.7}^{+0.7}$	12.0	${5.5}_{-1.5}^{+2.0}$
T09	${4.7}_{-1.1}^{+4.5}$	⋯	${4.3}_{-1.1}^{+4.0}$	⋯	${3.5}_{-1.2}^{+2.7}$	15.4	${3.2}_{-1.3}^{+2.4}$	14.6	${3.1}_{-1.3}^{+2.3}$	14.3	${1.5}_{-1.0}^{+2.0}$
T10	${0.7}_{-2.9}^{+2.3}$	⋯	${0.0}_{-2.9}^{+2.2}$	15.9	$-{1.4}_{-2.8}^{+2.1}$	14.9	$-{2.1}_{-2.8}^{+2.2}$	13.7	$-{2.2}_{-2.8}^{+2.1}$	13.4	12 ${}_{-4.0}^{+32}$
T11	$-{2.7}_{-3.0}^{+2.5}$	16.2	$-{3.5}_{-2.9}^{+2.5}$	13.0	$-{5.3}_{-2.8}^{+2.5}$	14.0	$-{6.0}_{-2.7}^{+2.4}$	13.2	$-{6.2}_{-2.7}^{+2.4}$	13.0	>20
T12	$-{4.0}_{-2.8}^{+3.1}$	16.4	$-{4.8}_{-2.4}^{+3.3}$	⋯	$-{6.4}_{-2.6}^{+2.7}$	14.8	$-{7.1}_{-2.8}^{+2.6}$	14.2	$-{7.3}_{-2.6}^{+2.7}$	14.1	>20
T13	$-{2.7}_{-2.6}^{+2.7}$	⋯	$-{3.3}_{-2.4}^{+2.6}$	17.7	$-{4.5}_{-2.3}^{+2.5}$	15.0	$-{5.1}_{-2.2}^{+2.5}$	14.2	$-{5.2}_{-2.2}^{+2.5}$	14.0	>20
T14	$-{1.1}_{-1.6}^{+1.6}$	16.7	$-{1.8}_{-1.6}^{+1.6}$	13.9	$-{3.3}_{-1.6}^{+1.6}$	14.0	$-{3.9}_{-1.6}^{+1.5}$	13.2	$-{4.1}_{-1.6}^{+1.6}$	12.9	${21}_{-10}^{+25}$
T14b	${2.8}_{-1.0}^{+2.3}$	⋯	${2.3}_{-1.1}^{+2.4}$	17.4	${1.1}_{-1.1}^{+2.1}$	15.4	${0.6}_{-1.2}^{+2.1}$	14.8	${0.5}_{-1.2}^{+2.1}$	14.5	${5.5}_{-3.4}^{+4.0}$
T15	${5.7}_{-3.2}^{+2.8}$	16.7	${5.3}_{-3.3}^{+2.3}$	13.9	${4.3}_{-3.3}^{+1.6}$	⋯	${3.9}_{-3.5}^{+1.3}$	13.2	${3.8}_{-3.5}^{+1.2}$	12.9	${0.6}_{-0.3}^{+2.2}$
T16	$-{2.0}_{-1.8}^{+2.5}$	16.1	$-{2.6}_{-1.8}^{+2.5}$	15.3	$-{5.8}_{-1.8}^{+2.4}$	13.6	$-{4.3}_{-1.8}^{+2.3}$	12.9	$-{4.4}_{-1.8}^{+2.4}$	12.6	22 ${}_{-15}^{+27}$
T17	⋯	⋯	⋯	17.1	⋯	15.3	⋯	14.4	⋯	14.2	⋯
T18	${8.6}_{-0.4}^{+0.1}$	⋯	${7.6}_{-0.5}^{+0.2}$	17.6	${5.6}_{-0.5}^{+0.4}$	15.1	${5.0}_{-0.5}^{+0.4}$	14.5	${4.8}_{-0.5}^{+0.4}$	14.2	${0.6}_{-0.1}^{+0.2}$
T19	⋯	⋯	⋯	17.9	⋯	15.6	⋯	14.9		14.6	⋯
T20	${0.0}_{-5.7}^{+3.8}$	⋯	$-{0.8}_{-5.7}^{+3.7}$	17.2	$-{2.5}_{-5.6}^{+3.6}$	15.3	$-{3.3}_{-5.5}^{+3.6}$	14.5	$-{3.4}_{-5.5}^{+3.5}$	14.3	>20
T21	$-{1.6}_{-1.3}^{+2.6}$	17.7	$-{2.2}_{-1.3}^{+2.6}$	16.7	$-{3.4}_{-1.2}^{+2.6}$	14.7	$-{3.9}_{-1.2}^{+2.6}$	13.9	$-{4.0}_{-1.1}^{+2.6}$	13.6	${29}_{-27}^{+20}$
T22	$-{2.5}_{-0.9}^{+1.7}$	16.3	$-{2.9}_{-0.9}^{+1.6}$	17.0	$-{3.9}_{-0.8}^{+1.5}$	14.1	$-{4.3}_{-0.8}^{+1.5}$	13.4	$-{4.4}_{-0.7}^{+1.5}$	13.2	${28}_{-14}^{+11}$
T23	${9.9}_{-7.8}^{+4.3}$	16.3	${8.9}_{-7.1}^{+4.3}$	15.4	${6.6}_{-5.5}^{+4.1}$	14.6	${6.0}_{-5.1}^{4.0+}$	14.0	${5.8}_{-4.9}^{+4.0}$	13.8	${0.3}_{-0.2}^{+3.0}$
T24	${4.3}_{-2.0}^{+2.1}$	⋯	${3.2}_{-1.9}^{+2.0}$	17.8	${2.0}_{-1.8}^{+2.0}$	15.4	${1.4}_{-1.8}^{+2.0}$	14.8	${1.3}_{-2.2}^{+1.9}$	14.4	${3.8}_{-2.2}^{+6.8}$
T25	⋯	⋯	⋯	⋯	⋯	12.7	⋯	12.0	⋯	11.8	⋯
T26	⋯	⋯	⋯	⋯	⋯	15.6	⋯	15.1	⋯	14.8	⋯
T27	$-{5.6}_{-3.3}^{+5.4}$	16.2	$-{6.1}_{-3.2}^{+5.2}$	⋯	$-{7.0}_{-3.1}^{+5.0}$	14.2	$-{7.4}_{-3.6}^{+4.8}$	13.7	$-{7.5}_{-3.1}^{+4.8}$	13.5	>20
T28	${0.0}_{-2.2}^{+1.1}$	⋯	$-{0.7}_{-1.9}^{+1.1}$	⋯	$-{2.0}_{-2.0}^{+1.0}$	13.9	$-{2.7}_{-2.0}^{+1.1}$	13.3	$-{2.8}_{-2.0}^{+1.1}$	13.1	${13}_{-5.0}^{+20}$
T29	${4.1}_{-1.6}^{+1.7}$	⋯	${3.0}_{-1.6}^{+1.5}$	⋯	${1.7}_{-1.5}^{+1.3}$	13.8	${1.1}_{-1.5}^{+1.2}$	13.0	${1.0}_{-1.5}^{+1.1}$	12.8	${2.1}_{-0.8}^{+2.1}$
T30	$-{3.0}_{-3.0}^{+3.4}$	15.4	$-{3.7}_{-3.0}^{+3.3}$	14.9	$-{5.2}_{-2.9}^{+3.3}$	14.0	$-{5.9}_{-2.8}^{+3.3}$	13.7	$-{6.0}_{-2.9}^{+3.2}$	13.4	>20
T31	${2.1}_{-3.0}^{+0.7}$	⋯	${1.4}_{-2.9}^{+0.8}$	17.8	${0.3}_{-2.5}^{+0.7}$	15.8	$-{0.3}_{-2.2}^{+0.8}$	14.9	$-{0.3}_{-2.2}^{+0.7}$	14.9	${9.6}_{-2.7}^{+17}$

Download table as:  ASCIITypeset image

Peculiar radial velocities, as referred to the average velocities of the stars at the same position in the Galaxy, can be one characteristic of a surviving companion of the SN, the excess velocity coming from the orbital motion of the star before the binary system is disrupted by the explosion, plus the kick imparted by the collision with the SN ejecta.

The results are given in the first and second columns in Table 6, and they are compared with those of K14 for stars common to the two surveys in the third and fourth columns. The velocities are in the heliocentric system. There is good agreement for some stars but not for all of them. In our case, we have two or three spectra of the same star, so we can be sure about the radial velocities and exclude possible binarity of those stars. As an example, we can quote the velocity of our star T18, named A in K14. We measure its radial velocity with an uncertainty of 1–2 km s⁻¹. The star moves at 42.5 km s⁻¹, while K14 give −69.07 km s⁻¹.

We will address the point of radial velocities and proper motions as compared to the kinematics of the Galaxy in Section 6.

Being the surviving companion star of an SN Ia means having a peculiar velocity, referred to the average velocity of the stars at the same position in the Galaxy, due to the orbital motion in the binary progenitor of the SN, plus the kick velocity caused by the impact of the SN ejecta. An estimate of the expected velocities, depending on the type of companion (MS, subgiant, red giant, supergiant) was made by Canal et al. (2001), and more recently by Han (2008). The highest peculiar velocities (∼450 km s⁻¹) would correspond to MS companions and the smallest ones (∼100 km s⁻¹) to red giants.

As the reference for the average velocities of the stars, depending on the location within the Galaxy and on the stellar population considered, we adopt the Besançon model of the Galaxy (Robin et al. 2003). We have run the model to find the distributions of both radial velocities and proper motions in the direction of the center of Kepler's SNR and within the solid angle subtended by our search, and including all stellar populations. The same model has been taken as the reference in K14.

In Figure 5, the 1σ, 2σ, and 3σ regions of the radial-velocity distribution (in the heliocentric reference system) are shown. We see that the average velocities first steadily increase, with positive values. That corresponds to the differential rotation of the Galactic disk. The dispersion also increases as we move from the thin to the thick disk, given the direction of the line of sight. Then, at a distance of ∼7 kpc, both the slope and the dispersion increase when reaching the Galactic bulge, and start decreasing beyond ∼9 kpc.

Zoom In Zoom Out Reset image size

In the same figure we compare the measured radial velocities with the distribution predicted by the Besançon model. We see that no star deviates significantly from the model distribution.

In Figure 6, the same is done for the proper motions perpendicular to the Galactic plane. In Figure 7, the same is done for the proper motions in Galactic latitude.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We see the same increase in the dispersion as in Figure 5, when reaching the Galactic bulge. Again, no star is a significant outlier with respect to the theoretical distribution. Without knowing the distances to the targets in Table 5, T18 appears outstanding from its proper motion (total proper motion μ = 50.5 mas yr⁻¹), but that corresponds in fact to the short distance to the star, of only ${0.6}_{-0.2}^{+0.1}$ kpc. For that distance, the velocity perpendicular to the line of sight is v_perp ≃ ${144}_{-48}^{+44}$ km s⁻¹, which falls within the range of the model predictions. The same applies to T23, with a proper motion μ = 7.6 mas yr⁻¹, at a distance of only ${0.3}_{-0.2}^{+3.0}$ kpc. Table 7 summarizes the main conclusions about the stars within 10 kpc and shows the parameters in relation to the model of the Galaxy.

Table 7.  Radial Velocities, Proper Motions, and Distances to the 12 Stars Located at Less than 10 kpc

Star	vr (km s−1)a	μl (mas yr−1)	μb (mas yr−1)	d (kpc)
T03	140.8	1.40 ± 0.00	−0.38 ± 0.08	${8.2}_{-2.3}^{+20}$
T04	8.2	−3.75 ± 0.03	2.41 ± 0.24	${3.3}_{-1.9}^{+2.7}$
T06	94.4	−0.34 ± 0.03	−0.75 ± 0.10	${8.0}_{-5.3}^{+17}$
T08	−84.6	−0.68 ± 23.70	−0.31 ± 63.01	${5.5}_{-1.5}^{+2.0}$
T09	−15.8	0.88 ± 0.07	−0.64 ± 0.11	${1.5}_{-1.0}^{+2.0}$
T14b	142.3	2.13 ± 0.00	0.45 ± 0.09	${5.5}_{-3.4}^{+4.0}$
T15	−23.2	−5.39 ± 0.00	0.06 ± 0.09	${0.6}_{-0.3}^{+2.2}$
T18	42.4	32.94 ± 0.01	36.31 ± 0.09	${0.6}_{-0.1}^{+0.2}$
T23	50.8	6.73 ± 0.02	−0.53 ± 0.06	${0.3}_{-0.2}^{+3.0}$
T24	38.2	−1.33 ± 0.05	−3.00 ± 0.09	${3.8}_{-2.2}^{+6.8}$
T29	45.4	3.42 ± 0.05	3.06 ± 0.10	${2.1}_{-0.8}^{+2.1}$
T31	−67.5	−1.61 ± 0.02	2.77 ± 0.09	${9.6}_{-2.7}^{+17}$

Note.

^aThe errors in the radial velocities are of 1–2 km s⁻¹.

Download table as:  ASCIITypeset image

K14 note that, according to Blair et al. (1991) and Sollerman et al. (2003), Kepler's SNR has a systemic radial velocity of −180 km s⁻¹. As can be seen from Figure 5, that lies between 2σ and 3σ from the average radial velocity of the stars at ∼5 kpc from us, in the direction of the SNR. K14 suggest that a possible surviving companion of the SN should have −180 km s⁻¹ added to the radial velocity component of its orbital motion at the time of the explosion, which would make it more easily identifiable. The systemic velocity of the SNR, however, is that of the exploding WD, which is the sum of the velocity of the center of mass of the system plus the orbital motion of the WD at the time of the explosion. We do not know what the velocity of the center of mass was, and thus we cannot just add −180 km s⁻¹.

Concerning theoretical predictions, we can compare our results with the observational features expected from numerical simulations. All groups that have simulated the impact of the ejecta of a supernova on the companion star (Marietta et al. 2000; Pakmor et al. 2008; Liu et al. 2012, 2013; Pan et al. 2012) find that the companion star should have survived the explosion and gained momentum from the disruption of the binary system. Their predictions vary on how luminous the untied companion would be, on how much mass it would have lost due to the impact of the supernova ejecta, and on how fast it would rotate and move away from the center of mass of the original system.

Let us start with the rotational velocities of the post-explosion companions. Liu et al. (2013) did binary population synthesis after performing 3D hydrodynamic simulations of the impact of the ejecta on a main-sequence star with different orbital periods and separations from the exploding WD. They obtained the expected distribution of rotational velocities for the surviving companion. This leaves room for a wide range in this parameter, unlike previous assumptions that the post-impact star should have a very high velocity. Pan et al. (2012) also found that angular momentum of the companion would have been lost with the stripped material. In the case of the stars studied in the survey for the companion of SN 1604, all of them have rotational velocities lower than 20 km s⁻¹. This is not uniquely interpreted as a sign of the absence of companions, though.

Concerning the luminosity discussion, there are some differences in the surviving companions from the supernova explosion. Podsiadlowski (2003) found that, for a subgiant companion, the object ∼400 years after the explosion might be significantly overluminous or underluminous relative to its pre-SN luminosity, depending on the amount of heating and the amount of mass stripped by the impact of the SN ejecta. More recently, Shappee et al. (2013) have also followed the evolution of luminosity for years after the impact of the ejecta on the companion. The models in which there is mass loss rise in temperature and luminosity, peaking at 10⁴ L_⊙, then starting to cool and dim down to 10 L_⊙ some 10⁴ yr after the explosion. Around 500 days after explosion the companion luminosity is 10³ L_⊙. Pan et al. (2012, 2013, 2014) found lower luminosities for the companions than these previous authors. They found luminosities of the order of 10 L_⊙ for the companions, several hundred days after the explosion. It is interesting to see that they predict, for the surviving companions, effective temperatures T_eff in the range 5000–9000 K. This allows us to discard possible candidates below 5000 K in our sample. Only four stars in our sample have T_eff higher than 5000 K (see Tables 4 and 7). T04 has 5545 ± 251 K and log g = 3.8 ± 0.6. The distance is uncertain but consistent with that of Kepler's SN: a distance of ${3.3}_{-1.9}^{+2.7}$ kpc. The heliocentric velocity, however, is v = 8.2 km s⁻¹ and the proper motion moderate, μ_l = 3.75 ± 0.03 mas yr⁻¹ and μ_b = 7.41 ± 0.24 mas yr⁻¹. So, it is within the expectations of the Besançon model. There is only one more target at a distance compatible with Kepler's SN explosion and with T_eff higher than 5000 K. This target is T23, with T_eff = 5436 ± 336 K and log g = 5.3 ± 1.6, at a distance of ${0.3}_{-0.2}^{+3.0}$ kpc. The radial velocity is 50.8 km s⁻¹ and the proper motion is μ_l = 6.73 ± 0.02 mas yr⁻¹ and μ_b = −0.53 ± 0.06 mas yr⁻¹.

At a distance compatible with Kepler's SN there are no main-sequence stars in our sample (let us recall that we go down to 2.6 L_⊙). There are four subgiants and two giants at distances compatible with the Kepler distance. Their stellar parameters, radial velocities, rotational velocities, and proper motions are within what is expected of a sample field at the Kepler position in our Galaxy.

Overall, there are no stars showing any peculiarity. All of them have rotational velocities around 10–20 km s⁻¹ or less, since Giraffe HR15n is not able to measure values lower than that. Their radial velocities are within those expected for field stars.

The predictions by Pan et al. (2012, 2013) were that 400 yr after the SN Ia explosion, the luminosities of the companion stars would still be 10 times higher than those before receiving the impact of the ejecta of the SN. They have extended those predictions to MS companion masses down to 0.656 M_⊙ and He WDs down to 0.697 M_⊙ (Pan et al. 2014). We have gone below the luminosities predicted for surviving companions of the kind examined by these authors and the predicted T_eff values are higher than those found in our sample. They have calculated the post-impact evolution of MS companions and of He–WD companions of very low mass at the time of the explosion.

The He WDs at the time of the explosion (Table 1 in Pan et al. 2014) have runaway velocities within the range 490–730 km s⁻¹, which would correspond, for purely transverse motions at a distance of 5 kpc, to proper motions μ between 21 and 31 mas yr⁻¹ or, if assumed to make a 45° angle with the line of sight, to radial velocities between 350 and 516 km s⁻¹ and proper motions between 15 and 22 mas yr⁻¹. Those proper motions would have been detected by the HST astrometry, even for objects fainter than our targets down to ${m}_{{\rm{F}}814{\rm{W}}}$ ∼22.5 mag.

There are several channels through which WDs could be surviving companions of SN Ia explosions, apart from the He WD companions above mentioned. One is dynamically stable accretion on a CO WD from an He WD or from a lower-mass CO WD (Shen & Schwab 2017). In that case, an He-shell detonation could induce a core explosion (Shen & Bildsten 2014). The mass-donor WD might survive. One salient characteristic of those companions is that, due to their extreme closeness to the exploding WD and to their strong gravitational fields, they should capture part of the radioactive material (⁵⁶Ni) produced by the SN.

Shen & Schwab (2017) study the effects of the decays of ⁵⁶Ni to ⁵⁶Co and of ⁵⁶Co to ⁵⁶Fe, for different masses of material captured by WDs of masses between 0.3 M_⊙ and 0.9 M_⊙. The decays, in the physical conditions prevailing at the surfaces of those WDs, drive persistent winds and produce residual luminosities that, 400 yr after the explosion, are higher than ∼10 L_⊙ in all cases (see Figure 4 in Shen & Schwab 2017). Furthermore, the surviving WDs should be running away from the site of the explosion at velocities of ∼1500–2000 km s⁻¹. A search for such WD companions has recently been made by Kerzendorf et al. (2018), in the central region of the remnant of SN 1006, with a negative result. We have not detected faint hot surviving WDs moving at high speed. We are at a greater distance than SN 1006 and the exploration does not go so deep, though (see our discussion on completeness).

Another possible channel producing a surviving WD companion is the spin-up, spin-down model (Di Stefano et al. 2011; Justham 2011): the WD, spun up by mass accretion from the companion star, can grow beyond the Chandrasekhar mass; then, when the accretion ceases, it has to lose angular momentum before reaching the point of explosion. During this last time interval, the companion might have evolved past the AGB stage and become a cool WD. The timescale for spin-down is hard to determine theoretically, but Meng & Podsiadlowski (2013) empirically obtain an upper limit of a few times 10⁷ yr for progenitor systems that contain a red giant donor and for which circumstellar material has been detected. We must note, however, that the spin-up, spin-down model should mostly produce super-Chandrasekhar explosions, since there is nothing there to tell the system to stop mass transfer just when the WD has reached that mass. In the case of Kepler's SN, reconstruction of its light curve (Ruiz-Lapuente 2017) clearly indicates that the SN was in no way overluminous.

From all the preceding evidence, we can exclude MS, subgiants, giants, and to a certain extent stars below the solar luminosity.

As an interesting point, no one has yet attempted to calculate how much and for how long the impact of the SN Ia ejecta would affect the luminosity of a WD companion in the spin-up, spin-down case. One cannot just assume that the WD would be cold and dim and remain so after the explosion. This has not been proved by any hydrodynamic simulation. It has only been done for closer pairs of WDs, as mentioned above. The typical separation between the two WDs, at the time of the explosion, should be larger that in the cases considered by Pan et al. (2014) and Shen & Schwab (2017), so less radioactive material would be captured and the runaway velocities would also be smaller, but the narrowing of the orbit by the emission of gravitational waves during the cooling stage of the companion WD might not be negligible, and the loss of angular momentum by the system during a likely common-envelope episode preceding the formation of the detached WD–WD system would also have considerably narrowed, previously, the separation between the two objects. We encourage others to perform these hydrodynamical calculations.

We present a study that includes the first detailed stellar parameters—T_eff, log g, v_rot—as well as accurate radial velocities of the stars, and proper motions, using the HST, of possible companions of Kepler's SN within 20% of the remnant center. This last part of the research is very important, since one does not know whether the peculiar velocities expected for the surviving companion will mostly be along the line of sight or perpendicular to it. No attempt to measure the proper motions of the stars in the core of Kepler's SNR had ever been made before.

We have determined luminosities and distances to the candidate companions of Kepler's SN. Any companion would have a luminosity more than twice the solar luminosity, which is the lowest luminosity of our sample. The radius of our search is 24 arcsec, which is 20% of the average radius of the SNR. Our stars correspond to stellar parameters and velocities consistent with being from a mixture of stellar populations in the direction of Kepler's SNR.

From our study, we conclude that the SD scenario is disfavored in the case of Kepler's supernova. The idea that Kepler's SN could come from the merging of two stars within a common envelope seems plausible. It would explain why the SN is surrounded by a large circumstellar medium. The idea of the core-degenerate scenario (Kashi & Soker 2011), whereby an already existing WD and a degenerate red giant core merge inside an AGB envelope, appears very likely in this case.

This analysis makes relevant intensive studies to detect surviving companions in very nearby SNe Ia remnants. There are many good cases for study in our Galaxy and in nearby ones. There are cases in our Galaxy far away enough that Gaia cannot estimate their proper motion, the stars being too dim. The HST plays a key role here. In addition, telescope time in 10 m class telescopes and in the coming generation of large telescopes with high-resolution spectrographs is the key to determining the nature of the surviving companions of SNe Ia.

As mentioned in the discussion, more hydrodynamical simulations are needed to compare predictions with observational results.

Based on observations made with ESO Telescopes at the La Silla Paranal Observatory under program ID 093.D-0384(A). Based on archival images from the HST programs GO-9731 and GO-12885. The scientific results reported in this article make use of observations by Chandra X-ray Observatory and published previously in cited articles. This research was supported by the Munich Institute for Astro- and Particle Physics (MIAPP) of the DFG cluster of excellence "Origin and Structure of the Universe." P.R.-L. is supported by AYA2015-67854-P from the Ministry of Industry, Science and Innovation of Spain and the FEDER funds. J.I.G.H. acknowledges financial support from the Spanish MINECO under the 2013 Ramon y Cajal program MINECO RyC-2013-14875, and also from the Spanish Ministry Project MINECO AYA2014-56359-P. L.G. was supported in part by the US National Science Foundation under grant AST-1311862. We thank the referee, Wolfgang Kerzendorf, for his very useful report.