HIGH PRECISION ABUNDANCES IN THE 16 Cyg BINARY SYSTEM: A SIGNATURE OF THE ROCKY CORE IN THE GIANT PLANET^*

It is common to assume that stars of multiple stellar systems have the same chemical composition, since they originated from the same natal cloud. However, some studies indicate that in binary systems, there may be small differences in the chemical composition of their components (Gratton et al. 2001; Laws & Gonzalez 2001; Desidera et al. 2004, 2006; Ramírez et al. 2011). One explanation for these anomalies is planet formation (e.g., Laws & Gonzalez 2001; Ramírez et al. 2011).

The binary system 16 Cygni is known for having a detected giant planet orbiting the B component, with a minimum mass of 1.5 M_Jup (Cochram & Hatzes 1997) and a probable true mass of about 2.4 M_Jup (Plávalová & Solovaya 2013). Even though the system has been monitored for small radial velocity variations for over two decades, so far no planets have been detected around the primary, which makes this system ideal to study the formation of giant planets. However, the chemical signatures of planet formation on the host star are expected to be very small, of only a few 0.01 dex (Meléndez et al. 2009; Ramírez et al. 2009; Chambers 2010), hence a high precision is needed to detect these effects.

Although earlier analyses of the 16 Cyg system suggested that 16 Cyg A is about 0.05 dex more metal-rich than 16 Cyg B (e.g., Gonzalez 1998), the difference is so small that it could be due to the relatively large abundance uncertainties of these earlier studies. In a pioneer precise line-by-line differential study of this binary, Laws & Gonzalez (2001) found a difference (A − B) of +0.025 ± 0.009 dex in the iron abundance of both components. Seeking for potential additional signatures of giant planet formation, Ramírez et al. (2011) performed a differential abundance determination of 25 elements and discovered significant differences among all chemical elements that were analyzed, with component A being more metal-rich by 0.04 ± 0.01 dex than B. In contrast, in a study published at about the same date, Schuler et al. (2011) found no difference in the chemical composition of these two stars. The intent of this work is to shed more light into this matter using better quality spectra and discuss the possible chemical signature caused by the formation of gas giant planets.

Spectra of 16 Cyg A and B were obtained with the Echelle SpectroPolarimetric Device for Observation of Stars (ESPaDOnS) on the 3.6 m Canada–France–Hawaii Telescope (CFHT) at Mauna Kea. The observations took place on 2013 June 6 on Queued Service Observing (QSO) mode. The observations were taken with the fiber only on the object (Spectroscopy, star o), that is the highest resolution (R = 81,000) on the instrument. Note that our resolving power is significantly higher than that used in the previous studies of Schuler et al. (2011) and Ramírez et al. (2011), R = 45,000 and R = 60,000, respectively.

The exposure times were 3 × 280 and 3 × 350 s on 16 Cyg A and B, respectively, with 16 Cyg B observed immediately after 16 Cyg A. We achieved a S/N ∼ 700 around 600 nm for each of the binary components. The asteroids Vesta and Ceres were also observed with the same spectrograph setup to acquire the solar spectrum that served as the reference in our differential analysis. A similar S/N (∼700) was achieved for both asteroids. Our S/N ratios are higher than those obtained by Ramírez et al. (2011), S/N ∼ 400, and about the same as that obtained by Schuler et al. (2011) for 16 Cyg A (S/N = 750).

We used the pipeline reduced spectra provided by CFHT, which passed through the usual reduction process including bias subtraction, flat fielding, spectral order extractions, and wavelength calibration. We performed the continuum normalization of the spectra using IRAF.

We used the line-by-line differential method to obtain stellar parameters and chemical abundances, as described in Meléndez et al. (2012) and Monroe et al. (2013). The 2002 version of the LTE code MOOG (Sneden 1973) was used with Kurucz ODFNEW model atmospheres (Castelli & Kurucz 2004).

The adopted line list is an updated version of that presented in Meléndez et al. (2012), with several dozen lines added. The equivalent width (EW) measurements were made by hand with the task splot in IRAF, using Gaussian profile fits. The local continuum was carefully selected by overplotting the spectra of both binary components and the solar spectrum for each line.

We obtained the abundance of 18 elements: C, O, Na, Mg, Al, Si, S, Ca, Sc, Ti, V, Cr, Mn, Fe, Co, Ni, Cu, and Zn. All abundances were differentially determined line-by-line using the Sun as standard in a first approach and then using 16 Cyg B as reference to obtain the 16 Cyg A − B ratios. The differential method minimizes errors due to uncertainties in the line transition probabilities and shortcomings of model atmospheres, allowing thus an improved determination of stellar parameters and chemical abundances. The elements V, Mn, Co, and Cu had their abundances corrected for hyperfine structure (HFS). For this calculation the blends driver in MOOG was used adopting the HFS data from Meléndez et al. (2012).

The atmospheric parameters for 16 Cyg A and B were obtained by differential excitation equilibrium (for T_eff) and differential ionization equilibrium (for log g), using as reference solar abundances for Fe i and Fe ii lines. First, we determined absolute abundances for the Sun using the solar atmospheric parameters of 5777 K for T_eff and 4.44 for log g, and adopting an initial microturbulence velocity of v_t = 0.9 km s⁻¹. Then, we estimated v_t by the usual method of requiring zero slope in the absolute abundances of Fe i lines versus reduced EW. We obtained a final v_t = 0.86 km s⁻¹ for the Sun, and computed our reference solar abundances for each line.

The next step was the determination of stellar parameters for the 16 Cygni stars. Initially, we used model atmospheres with the parameters published in Ramírez et al. (2011): T_eff = 5813 K, log g = 4.28 and [Fe/H] = 0.10 for 16 Cyg A, and T_eff = 5749 K, log g = 4.33 and [Fe/H] = 0.06 for 16 Cyg B. Then, we iteratively changed the stellar parameters of 16 Cyg A and B until we achieved the differential excitation and ionization equilibrium, and also no trend in the differential Fe i abundances with reduced EW (to obtain v_t), changing the metallicity of the models at each iteration until reaching convergence.

Our derived stellar parameters using the Sun as a standard are T_eff = 5830 ± 11 K, log g = 4.30 ± 0.02, v_t = 0.98 ± 0.02 km s⁻¹ and [Fe/H] = 0.101 ± 0.008 dex for 16 Cyg A, and T_eff = 5751 ± 11 K, log g = 4.35 ± 0.02, v_t = 0.90 ± 0.02 km s⁻¹ and [Fe/H] = 0.054 ± 0.008 dex for 16 Cyg B. These errors take into account the errors in the measurements and the degeneracy of stellar parameters. A similar procedure was repeated but using 16 Cyg B as the reference star instead of the Sun to perform the differential spectroscopic equilibrium (Figure 1), and fixing the stellar parameters of the B component to our results from the differential analysis relative to the Sun. The resulting atmospheric parameters for the A component are the same as when the Sun is used as a reference, but with smaller errors for T_eff (±7 K) and v_t (±0.01 km s⁻¹). The final Δ(Fe) difference for 16 Cyg A minus 16 Cyg B is 0.047 ± 0.005 dex, confirming that there is indeed a difference in the metallicity between the two stars of this binary system.

Zoom In Zoom Out Reset image size

Our stellar parameters are in very good agreement with the ones determined by Ramírez et al. (2011). We obtain somewhat higher both effective temperatures and log g by 17 K and 0.02 dex in the case of 16 Cyg A, and by 2 K and 0.03 dex for 16 Cyg B. The relative difference between the components A and B is in even better agreement, with our results showing a ΔT_eff higher by 15 K than Ramírez et al. (2011), and the difference in the Δlog g is 0.00 dex. Using the photometric IRFM calibrations of Ramírez & Meléndez (2005) for 13 different optical and infrared colors (Table 1) from the Jhonson, Cousins, Vilnius, Geneva, DDO, Tycho and 2MASS systems (Taylor 1986; Mermilliod et al. 1997; Høg et al. 2000; Cutri et al. 2003), with the corresponding [Fe/H] for each binary component, we determined average, median and trimean³ effective temperatures (Table 1) for the binary pair, resulting in a temperature difference of $\Delta T_{\rm eff}^{\rm phot}$ (A−B) = 58 ± 10, 78 ± 10, 73 ± 10 K, for the difference of average, median and trimean temperatures. The two robust indicators, median and trimean, are in excellent agreement with our spectroscopic $\Delta T_{\rm eff}^{\rm spec}$ (A−B) = 79 ± 7 K, and also in agreement with the results from Ramírez et al. (2011), who found $\Delta T_{\rm eff}^{\rm spec}$ (A−B) = 64 ± 25 K. Compared to Schuler et al. (2011), our ΔT_eff and Δlog g are higher by +36 K and +0.03 dex, respectively. Note that according to the trigonometric log g (Ramírez et al. 2011), Δlog g should be 0.05 dex between the components, that is the value found in our work and by Ramírez et al. (2011), but Schuler et al. (2011) found a lower Δlog g = 0.02 dex, although our results are in agreement with Schuler et al. (2011) within their error bars.

Table 1. Colors of 16 Cyg A and B on Different Photometric Systems, and their Corresponding Mean Effective Temperatures using the Calibrations of Ramírez & Meléndez (2005)

Color	16 Cyg A	16 Cyg B
(B − V)	0.644	0.663
(b − y)	0.410	0.416
(Y − V)	0.569	0.575
(V − S)	0.557	0.569
(B2 − V1)	0.398	0.402
(B2 − G)	0.109	0.117
(V − R)C	0.357	0.363
(V − I)C	0.698	0.706
(R − I)C	0.341	0.343
C(42-45)	0.648	0.669
C(42-48)	1.671	1.698
(BT − VT)	0.722	0.732
(V − K2)	1.533	1.577
$T_{\rm eff}^{\rm average}$ (K)	5726	5668
$T_{\rm eff}^{\rm median}$ (K)	5737	5659
$T_{\rm eff}^{\rm trimean}$ (K)	5734	5661
σ (K)	29	23
s.e. (K)	8	6

Download table as:  ASCIITypeset image

Once the stellar parameters of the 16 Cygni stars were set using iron lines, we computed abundances for all remaining elements. In Table 2 we present the final differential abundances of 16 Cyg A relative to 16 Cyg B, and their respective errors, while in Table 3 we present the abundances and errors for 16 Cyg A and B using the Sun as standard. We present both the observational errors and systematic errors due to uncertainties in the stellar parameters, as well as the total error obtained by adding quadratically both errors.

The differential abundances of the 16 Cyg pair relative to the Sun are shown in Figure 2. Both 16 Cyg A and B show abundances that have a clear trend with condensation temperature, as already shown by Ramírez et al. (2011) and Schuler et al. (2011). There is a reasonable agreement with the mean trend of 11 solar twins relative to the Sun by Meléndez et al. (2009), shown by solid lines in Figure 2, after a vertical shift is applied to match the refractory elements. Interestingly, the same qualitative pattern as in Ramírez et al. (2011) is found for individual volatile elements in both components, with O somewhat higher than C, and Zn somewhat higher than S. Thus, the variations among the volatile elements are likely real.

Zoom In Zoom Out Reset image size

From Figure 2, it is already noticeable that there are abundance differences between the two 16 Cygni components, with 16 Cyg A being more metal-rich. The differential abundances of 16 Cyg A relative to 16 Cyg B, plotted in Figure 3, shows this more clearly. As already found by Ramírez et al. (2011), all elements seem enhanced in 16 Cyg A, but now this is more evident due to our higher precision. This is contrary to the results obtained by Schuler et al. (2011), who found no chemical difference in the binary pair. The differential analysis of Takeda (2005) also showed both components to have the same iron abundance, but the S/N of his spectra (S/N ∼ 100) is too low for a precise analysis.

Zoom In Zoom Out Reset image size

Ramírez et al. (2011) found a roughly constant difference of about 0.04 dex in the differential abundances (A − B) of volatiles and refractories. However, while in our study the volatile elements show a difference of about 0.03 dex, the refractories show larger differences and a trend with condensation temperature (Figure 3). A similar trend has been reported in a short note added in proof by Laws & Gonzalez (2001), where based on the analysis of 13 elements, a correlation with condensation temperature is found, with a slope of 1.4 ± 0.5 × 10⁻⁵ dex K⁻¹, however, no further details are given. Interestingly, the same slope of A minus B (1.4 ± 2.8 × 10⁻⁵) is found by Schuler et al. (2011). In this work, we obtain a slope for the refractories of 1.88 ± 0.79 × 10⁻⁵ dex K⁻¹, in reasonable agreement with the results by Laws & Gonzalez (2001) and Schuler et al. (2011). Note that the abundance difference that we find here for 16 Cyg A − B, is very distinct from the mean trend for the 11 solar twins of Meléndez et al. (2009), shown by a dot-dashed line in Figure 3 after a shift has been applied to fit highly refractory elements.

The overall deficiency in the abundances of 16 Cyg B (compared to 16 Cyg A), could be attributed to the formation of its giant planet, as the metals missing in 16 Cyg B could have been taken from the proto-planet disk to form its gaseous giant. More interestingly, for the refractories we seem to detect a trend with condensation temperature, which may represent the detection, for the first time, of a signature of the rocky accretion core of a giant planet. In the core accretion model (e.g., Papaloizou & Terquem 2006), first a rocky core forms through coagulation of planetesimals until it is massive enough for accretion of a gas envelope, to form a giant planet. The higher deficiency of refractories in 16 Cyg B, means that the giant planet 16 CygBb may have an excess of refractories, which could be due to its rocky accretion core. Another important hint is that the break in condensation temperature between volatiles and refractories seem to occur around T_cond ∼ 500–700 K. This means that most likely the rocky core was not formed in the inner disk regions (equivalent to the Earth–Sun distance), but it was formed at a larger distance, where giant planets are more likely to form.

Following Chambers (2010), we estimate the mass of the rocky core by adding a mixture of the composition of the Earth and CM chondrites to the convective zone of 16 Cyg B. Assuming a similar convection zone as in the Sun (0.023 M_☉), we can reproduce the trend of the refractories (Figure 3) by adding about 1.5 M_⊕ of rocky material with the above composition. Notice that this is the minimum mass of the rocky core because we do not know the exact size of the convection zone at the time the giant planet was formed (see discussion in Ramírez et al. 2011). Assuming a convection zone of about 0.1 M_☉ at the time of the formation of the giant planet's core, the mass of the rocky core would be higher (6 M_⊕). Our estimates (1.5–6 M_⊕) are consistent with Jupiter's core mass of ∼5 ± 5 M_⊕ (Guillot 2005).

We found significant differences in the chemical abundances of 16 Cyg A relative to 16 Cyg B, for all analyzed elements. The abundance differences range from 0.03 dex for the volatiles up to 0.06 dex for the refractories.

The 16 Cygni system is so far a unique case where high-precision abundance analyses show a sharp distinction in the chemical composition of the binary components. A recent study of the binary pair HAT-P-1, where the secondary hosts a giant planet of 0.53 M_Jup (Bakos et al. 2007) but no planets have been detected so far around the primary, reveals no abundance contrast (Liu et al. 2014). One explanation for the lack of abundance differences between the binary components of HAT-P-1 could be that the mass of the planet is much smaller than the planet in the 16 Cygni system, which has about 2.4 M_Jup (Plávalová & Solovaya 2013). Scaling by the mass of the planets, the observed difference of about 0.04 dex in the chemical abundances of the 16 Cygni pair, would imply in a dissimilarity of only 0.009 dex for the HAT-P-1 binary, which would be challenging to detect. Another recent study of a binary pair, HD 20781/HD 20782, where HD 20782 has a Jupiter-mass planet and HD 20781 hosts two Neptune-mass planets, show zero abundance differences (0.04 ± 0.07 dex) within the error bars (Mack et al. 2014).

In any case, our findings could be interpreted as due to the formation of the giant planet around 16 Cyg B. Within that scenario, we have detected, for the first time, the signature of the rocky accretion core of the giant planet 16 Cyg Bb, with a mass of ∼1.5–6 M_⊕. Our study opens new windows on the study of the planet–star connection.

M.T.M. acknowledges support by CAPES. J.M. acknowledges support by FAPESP (2012/24392-2) and CNPq (Bolsa de produtividade).

Facility: CFHT (ESPaDOnS) - Canada-France-Hawaii Telescope