{'fraction_non_alphanumeric': 0.03573140209808829, 'fraction_numerical': 0.03626391181639065, 'mean_word_length': 4.624438454627134, 'pattern_counts': {'":': 0, '<': 5, '<?xml version=': 0, '>': 4, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
The chemical signatures of planetary engulfment events in binary systems

Planetary engulfment events involve the chemical assimilation of a planet into a star's external layer. This can cause a change in the chemical pattern of the stellar atmosphere in a way that mirrors the composition of the rocky object engulfed, with the refractory elements being more abundant than the volatiles. Due to these stellar chemical changes, planetary engulfment events can render the process of chemical tagging potentially inaccurate. A line-by-line differential analysis of twin stars in wide binary systems allows us to test the chemical homogeneity of these associations with typical individual stellar Fe I uncertainties of 0.01 dex and eventually unveil chemical anomalies that could be attributed to planetary engulfment events. Out of the 14 systems analysed here, we report the discovery of the most chemically inhomogeneous system to date (HIP34407/HIP34426). The median difference in abundances of refractory elements within the pair is 0.19 dex and the trend between the differential abundances and condensation temperature suggests that the anomaly is likely due to a planetary engulfment event. Within our sample, five other chemically anomalous systems are found.


INTRODUCTION
Within the last decade, radial velocity and transit surveys have discovered thousands of exoplanets around Sun-like stars. The global picture that has emerged exhibits a remarkable degree of diversity in terms of the architectures of these planetary systems (Winn & Fabrycky 2015). Presumably, this observed diversity has arisen as a result of dynamical processes acting since the first stage of planetesimal formation (Raymond et al. 2018;Chambers 2018). The fact that some systems have undergone complex phases of dynamical evolution is also attested to by the presence of planets on highly eccentric orbits (e.g., Kane et al. 2012) that are misaligned or even counter rotating with respect to the spin axes of their hosting star (e.g., Naoz et al. 2011), interstellar exoplanets (e.g., Mróz et al. 2018) or by the observation of dusty debris disks formed through exoplanet collisions (e.g., Kenyon & Bromley 2016). It is likely that in systems with evidence of a dynamical past, part of the planetary material has fallen into the hosting star (e.g., Martinez et al. 2019), polluting its atmosphere and producing a significant increase in the stellar metallicity, which can be reliably detected (e.g., Spina et al. 2015Spina et al. , 2018Oh et al. 2018;Tucci Maia et al. 2019). In fact, such a dilution will not yield an indiscriminate abundance rise of all the heavy elements, but instead will produce a characteristic chemical pattern that mirrors the composition observed in rocky materials (Chambers 2010;Yana Galarza et al. 2016;Kunitomo et al. 2018) with most refractory elements (e.g., those having condensation temperature T cond >1000 K) being over-abundant relative to the volatiles (T cond <1000 K).
Prior proposals to explain these signatures suggest that the anomalous volatile to refractory ratios may be due to the selective accretion of volatiles. This could occur after the formation of rocky planets around a star (Meléndez et al. 2009). However, as already discussed in Spina et al. (2018), this explanation seems unlikely, as the chemical signature would have been imprinted on the star when it was too young (age≤10 Myrs), diluting any accreted material due to the star's thick external layer.
Unveiling the chemical signatures of planetary engulfment events in stellar atmospheres is key to studying the frequency over which these catastrophic events occur, but it is also relevant to test the chemical homogeneity of stellar associations and probe the success of "chemical tagging". The concept of chemical tagging is to use the chemical information of stars to assign them to their progenitor cloud (Freeman & Bland-Hawthorn 2002). A critical assumption of this technique is that members of the same stellar associations, such as open clusters or binary systems, are chemically identical as they formed at the same time and from the same material. Therefore, if the abundance of any element is altered due to a planetary engulfment event, the results from chemical tagging will no longer be reflective of the star's progenitor cloud.
In this letter we examine the chemical homogeneity of 14 wide binary star systems, making use of highprecision abundance determinations. In Section 2 we describe the target selection, observations and method of analysis. The results of our study are presented in Section 3 and our conclusions in Section 4.

OBSERVATIONS AND ANALYSIS
The 28 target stars were selected from lists of wide binary systems (separation > 4") compiled by Martín et al. (2002), Desidera et al. (2004), and Fuhrmann & Chini (2015). Proper motions, radial velocities and parallaxes from Gaia DR2 (Gaia Collaboration 2018) are listed in Table 1 and confirm that the pairs are physical systems. For this project, we have chosen only the pairs formed by twin dwarf stars with spectral types similar to that of the Sun. Twin stars are objects with atmospheric parameters very similar to each other (e.g., ∆T eff 300 K and ∆log g 0.2 dex). It has been shown that a strictly differential line-by-line analysis of twin stars permits us to obtain differential abundances at the highest precision possible (e.g., Meléndez et al. 2009;Bedell et al. 2014Bedell et al. , 2018Liu et al. 2014;Spina et al. 2018). In fact, when studying samples of twin stars, most of the systematic uncertainties that plague chemical abundance analyses are so similar among stars in the same binary pair that a strict differential line-by-line analysis cancels them out. Also NLTE effects have been found negligible for all the elements considered by our analysis of solar twin stars (Meléndez et al. 2012;Melendez et al. 2016;Nissen 2015;Spina et al. 2016), and so no NLTE corrections were required. This leaves the observational noise as the main source of error. Thus, error bars can be made very small simply by acquiring very high-quality spectra (i.e. a resolving power R≥60,000 and signal-to-noise ratio S/N≥300).
Eight binary systems were observed by the Ultraviolet and Visual Echelle Spectrograph (UVES; Dekker et al. 2000) on the Very Large Telescope of the European Southern Observatory 1 . The observations have been performed with a resolving power R∼85,000 and a wavelength coverage between 330 -680 nm, though the DIC-1 (390+580). Another six pairs have been observed with the High Dispersion Spectrograph (HDS; Noguchi et al. 2002) on the Subaru telescope 2 . For these latter we have used a R∼80,000 and the Yc standard setup which covers the 439-705 nm spectral range. Thanks to the brightness of our targets (V≤9.3 mag), we achieved a SNR ranging from 300-400 pixel 1 at 600 nm, with a median of 350 pixel 1 . Solar spectra have been acquired both with UVES (S/N∼350) and HDS (S/N∼400) using the same instrument configurations described above.
All spectra have been normalized and Doppler-shifted using IRAF's continuum and dopcor tasks. Equivalent widths (EWs) of the atomic and molecular transitions reported in Meléndez et al. (2014) and listed in Table  2 3 have been measured with Stellar diff 4 . This code allows the user to select one or more spectral windows for the continuum setting around each line of interest. Ideally, these windows coincide with regions devoid of other absorption lines. We employed the same window settings to calculate continuum levels and fit the lines of interest with Gaussian profiles in every spectrum. Therefore, the same assumptions have to be taken in the choice of the local continuum around the lines of interest. This is expected to minimize the effects of an imperfect spectral normalization or unresolved features in the continuum that can lead to larger errors in the differential abundances (Bedell et al. 2014). Furthermore, Stellar diff is able to identify points affected by hotpixels or cosmic rays and remove them from the calculation of the continuum. The code delivers the EW of each line of interest along with its uncertainty.
The iron EW measurements are processed by the qoyllur-quipu (q2) code ) that performs a line-by-line differential analysis relative to the Solar spectrum and automatically estimates the stellar parameters (effective temperature T eff , surface gravity log g, metallicity [Fe/H], and microturbulence ξ) by iteratively searching for the three equilibria: ex-citation, ionisation, and the trend between the iron abundances and the reduced equivalent width. We assumed the nominal solar parameters, T eff =5777 K, log g=4.44 dex, [Fe/H]=0.00 dex and ξ =1.00 km s −1 (Cox 2000). The iterations are executed with a series of steps starting from a set of initial parameters and employing the Kurucz (ATLAS9) grid of model atmospheres (Castelli & Kurucz 2004). In each step the abundances are estimated using MOOG (version 2014(version , Sneden 1973. The errors associated with the stellar parameters are then evaluated by the code. This takes into account the dependence between the parameters in the fulfillment of the three equilibrium conditions (Epstein et al. 2010). We first run q2 adopting the Solar parameters as a first guess for each star. After q2 has converged to a set of stellar parameters, the differential abundances relative to the Sun for the following elements are calculated: C I, CH, NH, Na I, Mg I, Al I, Si I, S I, Ca I, Sc II, Ti I, Ti II, V I, Cr I, Cr II, Mn I, Fe I, Fe II, Co I, Ni I, Cu I, Zn I, Y II, Zr II, and Ba II. Through the blends driver in MOOG and adopting the line list from the Kurucz database, the q2 code corrected the abundances of V, Mn, Co, Cu, and Y for hyperfine splitting effects, by using the HFS components in the input line list. For each element, we performed a 3-sigma clipping on the abundances yielded by each EW measurement. This allowed us to remove the EW measurements affected by telluric lines or other unresolved blendings with adjacent lines. A second run of q2 using the restricted list of EW measurements yielded the stellar parameters listed in Table 3. Atmospheric T eff values range from 5515 -6295 K, while log g values range within 3.975 -4.523 dex. All pairs are composed of twin stars.
With these final parameters, we repeated the calculation of the differential abundances relative to the Sun. The resulting abundances are listed in Table 4, together with their uncertainties and the number of lines used for the abundance determinations. The error budget associated with each elemental abundance has been obtained by summing in quadrature the standard error of the mean among the lines, and the propagated effects of the uncertainties on the stellar parameters. The typical precision that we achieved in individual stellar Fe I abundances is 0.01 dex. We also determined the differential abundances within each pair which are listed in Table 5  tral noise or were contaminated by O 2 telluric lines. The measurement of this line requires particular care as it is small (typically 3 -6 mÅ) and it is blended by a Ni line with nearly the same wavelength (Allende Prieto et al. 2001). Thus, we have followed a procedure already tested by Nissen (2015). Namely, i) using Iraf's Splot task we measured the EW of [OI]+Ni line and its uncertainty by assuming different local continuum levels ; ii) through MOOG's ewfind task and [Ni/H] abundances in Table 4, we have calculated the EW contribution of the Ni line and its uncertainty; iii) we have subtracted the Ni contribution from the measured EW in order to estimate the EW of the [OI] line along with its uncertainty; iv) using the [OI] EW and the parameters listed in Table 3, we calculated the [O/H] and its uncertainty. The Oxygen differential abundances are also listed in Tables 4 and 5.
The differential abundances within each pair are plotted in Fig. 1 as a function of the condensation temperature.

RESULTS
In Table 6 we report the reduced chi-square value χ 2 red of the pairs, which is defined as where ∆[X i /H] is the differential abundance of the X th ielement within the pair, σ Xi its uncertainty and N is the number of elements detected in the components of the pair. In order to assess if the pair is chemically anomalous or not, we can compare the χ 2 red values with the value expected for the model where the pair is assumed to be chemically identical: χ 2 red =1. However, as noted by (Andrae et al. 2010), this value has an uncertainty σ χ =2/N which is given by the width of the χdistribution due to the random noise of the data (which we can assume to be Gaussian). The σ χ values are also listed in Table 6. Consequently, the six pairs of our sample with χ 2 red >1+3σ χ are all considered chemically anomalous.
The HIP34407/HIP34426 system has a χ 2 red value of 118.99, indicating that the two stars differ considerably in chemical composition. Five other systems are chemically anomalous, although to a lesser extent. They are HIP44858/HIP44862, HD98744/HD98745, HD103431/HD103432, HD105421/HD105422, and HIP70386A/HIP70386B, which have χ 2 red values equal to 2.47, 5.56, 4.92, 5.59, and, 1.95 respectively. Therefore, out of the 14 pairs considered in this letter, six could render chemical tagging potentially inaccurate at our level of precision. However, assuming uncertainties of 0.05 dex for each element, which is a typical value for large spectroscopic surveys (e.g., Smiljanic et al. 2014;Buder et al. 2018), we calculated the χ 2 Survey , also listed in Table 6. Based on this, only two pairs out of the 14 would result in being chemically inhomogeneous if targeted by a spectroscopic survey: HIP34407/HIP34426 and HD98744/HD98745. A rocky planet falling into a star and polluting its atmosphere would result in a selective enhancement of refractory elements (Chambers 2001). Therefore, in order to verify if the chemical anomalies identified above could be attributed to planetary engulfment events, we have performed a linear fit of the differential abundances as a function of the condensation temperature. The resulting slopes are also listed in Table 6. Among the initial six chemically anomalous pairs, three have a slope that is consistent with being zero, while another three pairs (i.e. HIP44858/HIP44862, HD105421/HD105422, and HIP34407/HIP34426) have a slope that is inconsistent with zero. However, it should be noted that most of the slopes listed in Table 6 are heavily driven by the differential abundances of the most volatile species, such as C, CH, O, and N (i.e. T Cond <200 K). In fact, the number of volatiles is much smaller than the number of refractory elements, as is visible in Fig. 1. In addition, the chemical abundances of volatile elements are often more uncertain than those of refractories.
In the case of HIP34407/HIP34426, the chemical anomaly is extremely large compared to all the other pairs. The elements with an intermediate condensation temperature, such as Na, S, and Zn (i.e. 500<T Cond <1000 K) are well aligned by the ∆[X/H]-T Cond relation plotted in Figure. 1. This neat linear relation between elemental abundances and condensation temperature for this system may suggest that the anomaly was caused by a planetary engulfment event. Similar results were independently obtained by Ramirez et al. (2019).
Interestingly, the HIP34407/HIP34426 system was determined to be the most chemically inhomogeneous system found to date among other pairs in binary systems or clusters, that show similar trends between abundances and condensation temperature (Spina et al. 2015(Spina et al. , 2018Tucci Maia et al. 2014Teske et al. 2015Oh et al. 2018;Biazzo et al. 2015;Saffe et al. 2016Saffe et al. , 2017. Its median differential abundance in the refractory elements is 0.19 dex. In addition, we note that the results from the HIP34407/HIP34426 system involve similar but not identical volatile elemental abundances. Out of the four most volatile species measured (i.e. T Cond <200 K), only one is consistent with zero, leading to the hy-pothesis that the slope is due to the engulfment of a gas giant. This is because such planets are capable of causing abundance changes in both volatiles (due to their gaseous outer layers) and refractories (due to their rocky cores and metallic inner layers). This hypothesis has been also proposed to explain a similar anomaly observed in the chemical composition of the 16 Cyg binary system (Tucci Maia et al. 2014.

CONCLUSIONS
The occurrence of chemically anomalous stars at the 0.01 dex precision level can be deduced from line-byline differential analysis of twin stars in binary systems. Out of the 14 systems measured using UVES and HDS, six were found to exhibit chemical inhomogeneities (i.e. χ 2 red >1+3σ χ ), rendering chemical tagging on these systems inaccurate. However, only two pairs were chemically inhomogeneous at the precision level typical of large spectroscopic surveys (i.e. 0.05 dex).
Of special note is that the HIP34407/HIP34426 system was determined to be the most chemically inhomogeneous system found to date. For this pair, the trend between differential abundances and condensation temperature is neat, and it may suggest that rocky material has polluted the atmosphere of HIP34407. We speculated that the planet engulfed by the star was a giant gaseous planet, as HIP34426 is also anomalously rich in the volatiles. The planet population around the stars of this pair is currently unknown. Follow-up observations are required to establish if the two stars also have two different architectures for their planetary systems. This is necessary to understand the origin of similar anomalies observed within other stellar associations (Spina et al. 2015(Spina et al. , 2018Tucci Maia et al. 2014Teske et al. 2015Oh et al. 2018;Biazzo et al. 2015;Saffe et al. 2016Saffe et al. , 2017, where the stars richer in refractory elements could also be the ones with a more chaotic architecture of their planetary systems. Finally, systems whose components have identical stellar parameters and chemical patterns were also observed, such as the HIP58298A/HIP58298B and HR4443/HR4444 systems. These systems are excellent laboratories to test if other quantities, such as the thickness of the convective zone, stellar rotation, and activity can differ (and to what extent) among stars with equal mass, age and chemical composition.
We thank the many scientists and engineers who made the UVES and HDS observations possible. It is a pleasure to acknowledge M. Asplund, A.R. Casey, J. Meléndez and D. Yong for helpful discussions. L. S. and A. I. K. acknowledge financial support from the Australian Research Council (Discovery Project 170100521).