EDP Sciences
Press Release
Free Access
Issue
A&A
Volume 617, September 2018
Article Number A104
Number of page(s) 17
Section Planets and planetary systems
DOI https://doi.org/10.1051/0004-6361/201732535
Published online 25 September 2018

© ESO 2018

1 Introduction

Extrasolar planetary systems have always shown a huge variety of orbital architectures, ever since the very first exoplanet orbiting a main sequence star was discovered (Mayor & Queloz 1995). In the following two decades the number of known exoplanets has increased to more than 35001, mainly discovered from transit (e.g. Kepler) and radial-velocity (RV; e.g., HARPS, HARPS-N) observations. The detected systems encompass a wide range of masses, radii, and orbital separations, including small mass rocky Earth-twins (e.g. Pepe et al. 2011).

Low-mass M dwarfs are the most common main sequence stars, comprising ~ 70% of the stars in the solar neighbourhood (Henry et al. 2006). Furthermore, M dwarfs have become the most promising ground for the hunt for low-mass, rocky planets (e.g. Dressing & Charbonneau 2013; Sozzetti et al. 2013; Astudillo-Defru et al. 2017b), due to their more advantageous mass and radius ratios compared to solar-type stars. There is a solid evidence, arising both from HARPS and Kepler observations, that super Earths and Neptunes are commonly found in multiple systems (e.g. Udry et al. 2017; Rowe et al. 2014, and references therein).

The nearby M1 dwarf Gl15A was studied by Howard et al. (2014), who find a short-period super-Earth orbiting the star: they measure a period of 11.44 d and an amplitude of 2.94 m s−1. They also study the activity signals of the host star, identifying the rotational period of the star to be 44 d, both from the SHK index analysis and also from the precise photometric light-curve, collected with the automatic photometric telescope (APT) at Fairborn Observatory (Eaton et al. 2003). Gl15A has also a known M3.5 binary companion, Gl15B, identified astrometrically by Lippincott (1972) from a small fragment of its orbit, with a measured orbital separation of 146 AU, and an orbital period of 2600 yr.

In this framework we present the clear detection of a long-period eccentric super-Neptune planet around Gl15A. We have used five years of high-precision Doppler monitoring with the High Accuracy Radial velocity Planet Searcher for the Northern emisphere (HARPS-N) high-resolution (resolving power R ~ 115 000) optical echelle spectrograph (Cosentino et al. 2012) at the Telescopio Nazionale Galileo (TNG), combined with 15 yr of archival RVs data from the Lick–Carnegie Exoplanet Survey (LCES) HIRES/Keck Precision Radial Velocity survey (Butler et al. 2017). The HARPS-N data were collected as part of the Harps-n red Dwarf Exoplanet Survey (HADES) programme, a collaboration between the Italian Global Architecture of Planetary Systems (GAPS; Covino et al. 2013; Desidera et al. 2013; Poretti et al. 2016) Consortium2, the Institut de Ciències de l’Espai de Catalunya (ICE) and the Instituto de Astrofísica de Canarias (IAC). We also confirm the presence and update the amplitude of the Gl15A b RV signal, significantly reducing its minimum mass, while confirming the orbital period measured by Howard et al. (2014). Our findings are also discussed in the light of the recent non-detection of Gl15A b in the CARMENES visual RV time series by Trifonov et al. (2018).

With tens of exoplanets orbiting binary stars observed to date (e.g. Eggenberger 2010), including a confirmed wide binary system with planetary companions orbiting both components (Desidera et al. 2014; Damasso et al. 2015), it is debated how the presence of stellar companions has influenced the distribution of planetary orbital parameters. For this reason, we derived new orbital parameters for the stellar companion Gl15B, using HARPS-N RV measurements along with astrometric measurements from the WDS (Washington Double Star; Mason et al. 2001) catalogue, and performed several numerical simulations to test the dynamical influence of Gl15B on the planetary system.

We describe the Doppler and photometric measurements collected for the analysis in Sect. 2, and then in Sect. 3 we describe the host star and binary companion updated properties. The complete analysis of the RVs and activity indices of the system is presented in Sect. 4. We then analyse the binary orbit and the perturbations produced on the two planets orbits in Sect. 5. We summarize and discuss our findings in Sect. 6.

2 Observations and HIRES catalogue data

As part of the HADES RV programme, Gl15A has been observed from BJD = 2456166.7 (27th August 2012) to BJD = 2457772.4 (18th January 2017). The total number of data points acquired was 115 over a time span of 1605 days. The HARPS-N spectra were obtained using an exposure time of 15 minutes, and achieving an average signal-to-noise ratio (S/N) of 150 at 5500 Å. Of the 115 epochs, 49 were obtained within the GAPS time and 67 within the Spanish time. Since the simultaneous Th-Ar calibration could contaminate the CaII H&K lines, which are crucial in the analysis of stellar activity of M dwarfs (Forveille et al. 2009; Lovis et al. 2011a), the observations were gathered without it. However, to correct for instrumental drift during the night, we used other GAPS targets spectra, gathered by the Italian team during the same nights as the Gl15A observations using the simultaneous Th-Ar calibration.

The data reduction and RV extraction were performed using the TERRA pipeline (Template-Enhanced Radial velocity Re-analysis Application; Anglada-Escudé & Butler 2012), which is considered to be more accurate when applied to M-dwarfs, with respectto the HARPS-N Data Reduction Software (DRS; Lovis & Pepe 2007). For a more thorough discussion of the DRS and TERRA performances on the HADES targets, see Perger et al. (2017). The rms of the TERRA RVs is 2.69 m s−1, while the mean internal error is 0.62 m s−1. The TERRA pipeline also corrected the RV data for the perspective acceleration of Gl15A, dvrt = 0.69 m s−1 yr−1.

We also used the HIRES-Keck binned data for Gl15A in our analysis, downloaded from the LCES HIRES/Keck Precision Radial Velocity Exoplanet Survey. The HIRES time series spans 6541 days, from BJD = 2450461.8 (13th January 1997) to BJD = 2457002.7 (11th December 2014). These data were newly reduced with respect of the original RVs used by Howard et al. (2014), and also new data have been observed with the HIRES spectrograph since the paper came out. For a full description of the catalogue and data reduction, see Butler et al. (2017). We discarded the last data point of the HIRES time series (BJD = 2457002.7), since it is almost ~10 m s−1 off with respect of the rest of the data. In our analysis we thus used 169 HIRES RVs, showing a variation of 3.26 m s−1 and an average internal error of 0.84 m s−1. The RV dataset from the Keck archive was already corrected for the perspective acceleration of the host star.

We combined the HARPS-N and HIRES RV datasets, obtaining a 290 points time series, spanning 7310 days. The complete RV time series has rms = 3.08 m s−1, mean error = 0.71 m s−1. Figure 1 shows the combined RV time series, with the respective mean RV subtracted from each dataset to visually compensate for the offset expected between the measurements of the two instruments.

As a part of the analysis of the binary orbit described in Sect. 5, we obtained 5 RV data points of the companion star Gl15B with HARPS-N, collected from BJD = 2457753.9 (31 December 2016) to BJD = 2457771.9 (18 January 2017). The rms of the time series is 1.26 m s−1, while the mean internal error is 2.42 m s−1.

As for all targets in the HADES survey, Gl15A was also monitored photometrically by the EXORAP (EXOplanetary systems Robotic APT2 Photometry) programme, to estimate the stellar rotation from the periodic modulation in the differential light curve. The observations were performed by INAF-Catania Astrophysical Observatory at Serra la Nave on Mt. Etna (+14.973°E, +37.692°N, 1725 m a.s.l.) with the APT2 telescope, which is an 80 cm f/8 Ritchey–Chretien robotic telescope. The detector is a 2 k × 2 k e2v CCD 230–42, operated with standard Johnson-Cousins UBV RI filters. The IDL routine aper.pro was used to implement the aperture photometry. The data were collected from the 13th of August 2013 to the 15th of June 2017, for a total 242 and 233 data points, for B and V photometry, respectively.

thumbnail Fig. 1

Combined HARPS-N and HIRES RV time series. As described in Sect. 4.3, we considered the HIRES data as splitinto two separate datasets for the purpose of the RV analysis.

Open with DEXTER

3 Stellar properties of Gl15A and Gl15B

The star Gl15A is a high proper motion nearby (π = 280.74 ± 0.31 mas) early-M dwarf, of type M1. We used the stellar parameters published by Maldonado et al. (2017), which were calculated applying the empirical relations by Maldonado et al. (2015) on the same HARPS-N spectra from which we derived the RV time series. This technique calculates stellar temperatures from ratios of pseudo-equivalent widths of spectral features, and calibrate the metallicity on combinations and ratios of different features. Although such techniques are mainly used for solar-type stars, Maldonado et al. (2015) proved them to be equally effective for low-mass stars.

Howard et al. (2014) derived a rotational period of 43.82 ± 0.56 d, both from their Keck-HIRES measurements of the SHK index and their APT photometric observations at Fairborn Observatory. Recently, Suárez Mascareño et al. (2018) analysed the potential signatures of magnetic activity in the CaII H&K and Hα activity indicators of the HADES M-dwarfs sample. For Gl15A they computed a mean level of chromospheric emission , and a rotation period Prot = 45.0 ± 4.4 d, fully consistent with the value found by Howard et al. (2014). The definition of can be extended for application on M-dwarfs spectra (Suárez Mascareño et al. 2018, and references therein).

The stellar parameters for Gl15A used in this work are summarized in the left column of Table 1. We can see that most are fully consistent with the values used by Howard et al. (2014).

The stellar companion of Gl15A, Gl15B, is a type M3.5 dwarf whose orbit was measured astrometrically by Lippincott (1972). For the purpose of our orbital analysis in Sect. 5, we took five HARPS-N spectra of Gl15B during the last observing season at TNG (see next section) and we applied on them the Maldonado et al. (2015) techniques to calculate updated stellar properties. The derived values are listed in the right column of Table 1.

Table 1

Stellar parameters for the two stars Gl15A and Gl15B.

Photometric analysis

In order to identify the potential modulation in the stellar photometry due to the presence of active regions, we used the GLS periodogram (Zechmeister & Kürster 2009) to analyse the B and V time series collected within the EXORAP survey, composed of 242 and 233 points, respectively, taken over five seasons from 13 August 2013 to 15 June 2017. No evident periodicity is found in either time series, in contrast with the findings of Howard et al. (2014) who found a clear signal at 43.82 days (see their Fig. 4), identified as the rotation period of the star, with a corresponding signal seen in the CaII activity index. The discrepancy between the two analysis may be due to the different photometricprecision of the observations: the amplitude of the signal found by Howard et al. (2014) was only of 6 mmag, which is below the sensitivity of the EXORAP survey for targets in the magnitude range of Gl15A.

4 Spectroscopic data analysis

In Fig. 1 an evident long-period signal can be seen. Howard et al. (2014) identified in their HIRES data a hint of a long period decreasing trend, which they included in their model as a constant negative acceleration of − 0.26 ± 0.09 m s−1 yr−1, and concluded it to be consistent with the orbit of the Gl15AB system as calculated by Lippincott (1972). But considering also the new HARPS-N data, it becomes clear how the long-term RV variation cannot be modelled as a linear trend.

Even if it could no longer be treated as a constant acceleration, the long-period signal could still be due to the binary reflex motion of Gl15A, so we investigate this possibility. To model the possible RV signal due to the stellar companion, we follow the procedure used by Kipping et al. (2011) to study the presence of a long-period companion in the HAT-P-31 system: they modelled the long-period signal as a quadratic trend, and then derived a range of orbital parameters from the quadratic coefficients. The ratio between the semi-amplitude, KB, and period, PB, of the companion signal is given by the second-order term of the trend, , as (their Eq. (4)) (1)

and, since KB depends on the orbital period and on the mass of the two stars, the orbital period can be derived as a function of the second-order term and the masses .

From the fit of the complete RV time series, we obtain m s−1 day−2. Using the stellar masses listed in Table 1, the resulting period is PB ≃ 680 yr, which would correspond to a semi-major axis aB ≃ 63 AU. This solution is clearly unphysical, since the presumed semi-major axis is less than half the observed orbital separation of the two objects.

Moreover, several high-contrast imaging surveys ruled out the presence of additional co-moving stellar objects close to Gl15A: van Buren et al. (1998) excluded the possibility of objects with M > 0.084M at separations of 9−36 AU, while Tanner et al. (2010) ruled out objects up to a magnitude contrast of ΔKs ≃ 6.95 mag within 1 (3.57 AU) and ΔKs ≃ 10.24 mag within 5 (17.8 AU).

The long-period signal observed in the RV time series is therefore very unlikely to be caused either by Gl15B or by additional stellar companions. Instead, a long-period planetary-mass companion orbiting around Gl15A could be a possible explanation of the observed signal. We now investigate this hypothesis, by analysing both the potential presence of a Keplerian signal in the RVs and the stellar activity signals in the chromospheric indicators, which could also cause long-period variations due to the star magnetic cycle.

4.1 The MCMC model

Bearing in mind that a signal tightly linked to the stellar rotation period is clearly present in the RV data, as shown by Howard et al. (2014), we have selected the Gaussian process (GP) regression as a useful, and physically robust, tool both to investigate the presence of periodicities in the chromospheric activity indicators and to mitigate the stellar activity contribution to the RV variability. The GP regression is becoming a commonly used method to suppress the stellar activity correlated “noise” in RV time series (e.g. Dumusque et al. 2017, and references therein), especially when adopting a quasi-periodic covariance function. This function is described by four parameters, called hyperparameters, and it can model some of the physical phenomena underlying the stellar noise through a simple but efficient analytical representation. The quasi-periodic kernel is described by the covariance matrix: (2)

where t and t′ indicate twodifferent epochs.

This kernel is composed of a periodic and an exponential decay term, so that it can model a recurrent signal linked to stellar rotation, also taking into account the finite-lifetime of the active regions. Such approach is therefore particularly suitable when modelling signals of short-term timescales, as those modulated by the stellar rotation period.

About the covariance function hyperparameters, h is the amplitude of the correlations; θ represents the rotation period of the star; w is the length scale of the periodic component, linked to the size evolution of the active regions; and λ is the correlation decay timescale, that can be related to the active regions lifetime. In Eq. (2), σinstr, data(t) is the data internal error at time t for each instrument; σinstr, jit is the additional uncorrelated “jitter” term, one for each instrument, that we added in quadrature to the internal errors in the analysis of the RV datasets to take into account additional instrumental effects and noise sources neither included in σinstr, data(t) nor modelled by the quasi-periodic kernel; is the Kronecker delta function.

In the GP framework, the log-likelihood function to be maximized by the Markov chain Monte Carlo (MCMC) procedure is (3)

where n is the number of the data points, K is the covariance matrix built from the covariance function in Eq. (2), and is the data vector.

We use the publicly available emcee algorithm (Foreman-Mackey et al. 2013) to perform the MCMC analysis, and the publicly available GEORGE Python library to perform the GP fitting within the MCMC framework (Ambikasaran et al. 2015). We used 150 random walkers to sample the parameter space. The posterior distributions have been derived after applying a burn-in as explained in Eastman et al. (2013, and references therein). To evaluate the convergence of the different MCMC analyses, we calculated the integrated correlation time for each of the parameters, and ran the code a number of steps around 100–1000 times the autocorrelation times of all the parameters, depending on the complexity of each analysis. This ensured that the emcee walkers thoroughly sample the parameter space (Foreman-Mackey et al. 2013).

4.2 Analysis of the activity indexes

We first investigated both the HIRES and HARPS-N CaII H&K and Hα index time series, in order to test the potential stellar origin of the long period variation seen in the combined RV time series shown in Fig. 1. No long-period trend is found in either the HIRES or HARPS-N datasets. Suárez Mascareño et al. (2018) found a magnetic-cycle type periodicity at Pcycle = 2.8 ± 0.5 yr in the HARPS-N CaII H&K and Hα indexes, even if we find no similar signal in the respective HIRES time series. Nevertheless, the period of this cycle is far from the time span of the long-period signal we investigate, so it could not be the origin of it. Another clue of the stellar origin of thelong-period modulation of the RVs could be the correlation between the activity indexes and RV time series, which we computed via the Spearman’s rank correlation coefficients. No significant correlation was identified ( for all the indexes). These facts suggest that the long-period signal is not due to the star magnetic cycle or chromospheric activity, and reinforces the hypothesis of it be due to a wide-orbit planetary companion.

To further test the effect of the activity of Gl15A, we then investigated the stellar activity behaviour over the four seasons covered by HARPS-N observations by analysing chromospheric activity indicators based on the CaII H&K, and Hα spectroscopic lines. They were extracted using the definition of the line cores and of the reference intervals given in Gomes da Silva et al. (2011). We analysed measurements obtained only from HARPS-N spectra because they represent an homogeneous dataset. By spanning more than 1600 days, the HARPS-N data alone can provide robust insights into the long- and short-term stellar activity variability. The average S/N was 18 for the CaII H&K index and 213 for the Hα.

We performed an analysis of the activity indicators based on a GP regression, as detailed in the previous section. By adopting the same covariance function (Eq. (2)), we used to model the stellar contribution present in the RV variations in the following section, our primary goal is to investigate some properties of the stellar activity and use them to constrain some parameter priors in the analysis of the radial velocities. This represents a reasonable expectation, because neither the activity indicators nor the RVs have been pre-whitened, and the variability patterns in the former could be present with similar properties in the latter, as was noticed by Affer et al. (2016) during the analysis of another HADES target. In this sense, results from the analysis of the activity indicators can be used to train the GP regression of the RVs, by keeping unchanged the way the stellar activity contribution is modelled. Here we have used Eq. (2) to describe the variability correlated with the stellar rotation period Prot, by adopting a uniform prior for the corresponding hyperparameter θ which is constrained between 40 and 60 days (the list of priors is shown in Table 2).

The MCMC analysis was stopped after running for around 1000 times the highest autocorrelation time. The posterior distributions of the model parameters are shown in Fig. 2, and the best-fit estimates are listed in Table 2. We note that the estimates of the stellar rotation periods are well constrained but slightly different for the two activity indicators. The value found from the CaII H&K index is very similar to that found by Howard et al. (2014) for the CaII H&K S-index derived from HIRES spectra. Thehyper-parameter θ found in theanalyses of the two time series is in very good agreement with the rotation period found by Suárez Mascareño et al. (2018) with a GLS analysis of the activity indexes time series of Gl15A.

Table 2

Priors and best-fit results for the Gaussian process regression analysis of the chromospheric activity indicators extracted from the HARPS-N spectra of Gl15A.

thumbnail Fig. 2

Posterior distributions of the fitted (hyper)parameters of the GP quasi-periodic model applied to the time series of CaII H&K (upper panel) and Hα (lower panel) activity indexes. The vertical lines denote the median (solid) and the 16th–84th percentiles (dashed).

Open with DEXTER

4.3 Analysis of the combined HIRES and HARPS-N RV time series

First, we studied the HARPS-N RV time series for confirmation of the presence of the Gl15A b signal found by Howard et al. (2014). We performed a GLS analysis to identify any periodicity in our data in P ∈[2, 100] d. In the periodogram, shown in Fig. 3, we can see that the P = 11.44 d period of planet Gl15A b is clearly recovered. As a confirmation of the coherence of the signal throughout the entire HARPS-N time series, Fig. 4a shows the evolution of the period, amplitude and phase recovered by the GLS periodogram as a function of the number of observations: we can clearly see how from ~ 70 forward the period and phase of the signal remain stable, with small oscillations of the order of the final error on the parameters, even if the remainder of the observations covers almost one and a half years. We also studied the periodogram of the CaII H&K and Hα HARPS-N time series, and found no significant signal at periods close to the inner planet period Pb, thus confirming its planetary nature as announced by Howard et al. (2014).

We note that, as can be seen in Fig. 3, the final amplitude recovered by the GLS periodogram on the HARPS-N dataset, Kb = 1.82 ± 0.31 m s−1, is smaller than the one published by Howard et al. (2014), Kb,How = 2.94 ± 0.28 m s−1. The decreasing behaviour of the amplitude recovered by GLS with increasing number of observation is not unexpected, as the sampling of the signal strongly influence the periodogram structure and fit, and we tested that a similar behaviour can be observed also in simulated datasets with the same epochs as our HARPS-N time series but containing only the planetary signal and white noise, as shown in Fig. 4b. The lower amplitude value than that found by Howard et al. (2014) can be explained similarly, since the HARPS-N time series is composed of five years of continuous high-cadence observations of the target, while the dataset used by Howard et al. (2014) consisted in mostly sparse measurements with an intensive high-cadence monitoring only in the last seasons, which could affect the signal recovery. Moreover, it is worth noticing that the HIRES data from the now public LCES HIRES/Keck archive have been reduced with a different and more effective technique (Butler et al. 2017) than the data used by Howard et al. (2014). Performing a GLS analysis of the two time series we obtain the same peak periodicity but an amplitude Kb ~ 23% smaller in the new archive data. Since no other short-period signal emerges from the GLS analysis of the RVs, we focussed our attention on the study of the long-period signal, which we assume to be due to a planetary-mass wide-orbit companion.

To estimate the orbital and physical parameters of the known planet Gl15A b and the candidate companion, hereafter Gl15A c, we have performed an MCMC analysis of the RVs. Following the prescription of Howard et al. (2014), we modelled the RVs dividing the HIRES dataset in a pre-upgrade and a post-upgrade sublist, due to the HIRES CCD upgrade occurred in August 2004. Each subsample was then treated as an independent dataset, with its own zero-point offset and additional jitter term. In this case, the in Eq. (3) represents the RV residuals, obtained by subtracting the Keplerian signal(s) from the original RV dataset.

The general form for the models that we tested in this work is given by the following equation: (4)

where nplanet = 1, 2; ν is a function of time t, time of the inferior conjuntion T0,j, orbital period Pj, eccentricity e and argument of periastron ωj; γinstr is the RV offset, one for each instrument; GP is the stellar noise modelled with the Gaussian Process. The term is the residual acceleration of the system, which can include also the RV acceleration due to the orbital motion within the Gl15AB binary system, which is difficult to predict due to the large uncertainties in the orbital parameters (as we will discuss in Sect. 5.1). Instead of fitting separately ej and ωj, we use the auxiliary parameters and to reduce the covariance between ej and ωj. Howard et al. (2014) find the eccentricity of planet Gl15A b to be consistent with zero, and concluded a circular orbit to be the best fit for the data. To confirm this conclusion or to evaluate the real value of the eccentricity eb, we assumed in our analysis a Keplerian orbit for both planets. The eventual eccentricity orbit for Gl15A b would also be of interest as a potential consequence of the dynamical interaction of planet b with Gl15B, as we will discuss in Sect. 5.

With the exception of very few parameters, we assumed uniform priors in our analysis. Our choice for the range of λ is justified by the results obtained from the GP analysis of the CaII and Hα spectroscopic activity indexes (see Table 2), and from a preliminary, quick MCMC test on the data, which showed that the chains were converging towards values not far from the expected stellar rotation period. For the semi-amplitudes K of the Doppler signals, we used the modified invariant scale prior: (5)

following the prescription in Gregory (2010). For the orbital period of Gl15A b there was a strong constraint from the results of Howard et al. (2014), confirmed by our GLS analysis previously discussed. Nevertheless, since the same dataset from Howard et al. (2014) is used in our analysis, it would be improper to use this results as an informative prior. We instead adopted an uninformative prior over a small period interval centred on the value found by Howard et al. (2014), e.g. Pb ∈ [10.9, 12.0]. For the orbital period of the candidate planet Gl15A c, we used a logarithmic prior to assume all the orders of magnitude equally likely.

The MCMC analysis was stopped after running for around 300 times the highest autocorrelation time. Table 3summarizesthe best-fit values for each of the 19 free parameters of our model, together with the details of the prior distributions used to draw the samples. The posterior distributions of the model parameters are shown in Fig. 5.

About the best-fit values for the free parameters, we first note that the semi-amplitude of the Doppler signal induced by Gl15A b is even lower than the estimate obtained from the GLS periodogram, differing from the estimate of Howard et al. (2014) (2.94 ± 0.28 m s−1) for more than 3.5-σ, resulting in a lower minimum mass. In addition to the previously discussed effect, this is a consequence of our global model, where planetary signals are fitted jointly with a term describing the stellar correlated noise, and also explains the lower values of σHIRES, jit required by our fit.

Our dataset does not cover a complete orbit of the outer planet Gl15A c, then we cannot reliably constrain the eccentricity, which appears to be significant within less than 1.5-σ (e < 0.40 at a 68% level of confidence). Our estimate for the minimum mass of Gl15A c places this companion in the group of the so-called super-Neptunes (as planets in the [20, 80]M mass range are generally referred to).

The eccentricity of Gl15A b also resulted consistent with zero within 1.5-σ, with a much lower value (e < 0.13 at a 68% level of confidence). Our analysis thus agrees with the results by Howard et al. (2014), who concluded that there was no evidence for an eccentric orbit to be preferable to a circular one.

Figure 6 shows the RV curves folded at the best-fit orbital periods for the known planet Gl15A b and the detected outer companion Gl15A c. In Fig. 7 we show the RV residuals, after the two Keplerians have been removed, with superposed the best-fit stellar correlated, quasi-periodic signal.

As for the GP hyper-parameters, the stellar rotation period θ assumes a value which is in agreement with that expected, and the higher uncertainties than those for the activity indicators (Table 2) should be due to the longer time span covered by the RVs, which include HIRES data. The longer timespan should also explain why the active-region evolutionary timescale λ sets on a value less than those found for the activity indicators. Without any ancillary data available as photometry, covering the same time span of the RVs, we cannot conclude if the value for λ is actually physically robust, but a value not far from the rotation period seems not unreliable, as already observed in studies of different systems (e.g. Affer et al. 2016). Also, the h hyper-parameter, corresponding to the mean amplitude of the stellar signal, is fully compatible with the expected value derived by Suárez Mascareño et al. (2018) from the mean activity level of the star, , that is Kexp = 1.9 ± 0.4 m s−1.

thumbnail Fig. 3

GLS periodogram of the complete HARPS-N RV dataset.

Open with DEXTER
thumbnail Fig. 4

Panel a: evolution of the GLS orbital parameters as a function of the number of HARPS-N RV measurements; panel b: evolution of the GLS orbital parameters as a function of the number of simulated data. The error bars on the final points are shown as references.

Open with DEXTER
Table 3

Planetary parameter best-fit values obtained through a joint modelling of Keplerian signals and correlated stellar noise, using Gaussian process regression.

thumbnail Fig. 5

Posterior distributions of the fitted (hyper)parameters of the two-planet model, where the stellar correlated noise has been modelled with a GP regression using a quasi-periodic kernel. On the y-axis is shown the logarithm of the product between the likelihood and the prior. The vertical lines denote the median (solid) and the 16th–84th percentiles (dashed).

Open with DEXTER

4.4 Comparison with CARMENES results

In a very recent work, Trifonov et al. (2018) publish an RV analysis at optical wavelengths of the Gl15A system as part of the CARMENES survey. The RVs derived from the visible arm of CARMENES show no evidence of the presence of Gl15A b, thus the authors conclude it to be an artefact of the stellar noise. The authors also identify a long-term trend in the combined Keck+CARMENES dataset, which they propose to be due to a distant low-mass companion. As discussed in Sect. 4.3, our HARPS-N data alone clearly confirm the presence of the 11.44 d period due to Gl15A b, albeit with a reduced amplitude and mass, and the analysis of the combined HARPS-N+Keck datasets makes a decisive case for the existence of the long-period low-mass giant Gl15A c based on a five-year time span of HARPS-N observations, partly overlapping with the Keck time-series. The orbital elements for Gl15A c derived in this work have larger formal uncertainties than those reported in Trifonov et al. (2018) due to the likely much longer orbital period inferred for the planet, but the inferred companion mass is fully in agreement with the preliminary CARMENES estimate.

In Fig. 8a we can see the last season of HARPS-N observations, that overlaps with most of the CARMENES data, compared toour two-Keplerian plus correlated stellar noise model. This season is clearly dominated by activity-related variations in the optical HARPS-N spectra, while the internal errors for HARPS-N are typically three times smaller than those of CARMENES (0.6 m s−1 vs. 1.8 m s−1). We nonetheless tested the effect of combining the HARPS-N and CARMENES datasets: the GLS periodogram, shown in Fig. 8b, clearly peaks on the 11.44 d period of Gl15A b, and remarkably resembles Fig. 3, recovering the same amplitude as that recovered on the HARPS-N dataset alone.

We also tested our two-planet + stellar noise model described in Sect. 4.3 on the complete dataset obtained combining the HIRES, HARPS-N, and CARMENES time series (running the MCMC code for around 100 times the maximum autocorrelation time), and we obtained values of the system parameters entirely in line with those presented in Table 3. In particular, the recovered Doppler semi-amplitude and minimum mass for the inner planet are m s−1 and , respectively. These values are well within 1-σ of the results from the analysis of only the HIRES and HARPS-N datasets: no significant drop in the amplitude of Gl15A b signal is observed. We also note that the residual jitter found in the analysis is small ( m s−1) signifying that again the GP is able to model the stellar noise for this third dataset as well. The GP parameters are almost unchanged with respect to the values presented in Table 3 ( m s−1, d, and d), meaning that the activity signal description derived from the first two datasets is robust also in modelling the CARMENES data. Given the intrinsically higher quality of the HARPS-N data that drive the GP regression modelling and in which the coherence of the signal from the inner planet Gl15A b is clearly present, we decide to stick to the results in Table 3 for the purpose of the dynamical interaction analysis described in Sect. 5.2.

Trifonov et al. (2018) also stated, as an argument in favour of the non-Keplerian interpretation of Gl15A b, that the 11.44 d signal disappeared when analysing the last two years of HIRES observations, subsequent to the time series used by Howard et al. (2014). Studying the same time span of data, we observed that, when ignoring the outlier described in Sect. 2, the 11.44 days signal is still clearly visible in the GLS periodogram.

thumbnail Fig. 6

Phase-folded RV curves for panel a Gl15A b and panel b Gl15A c. Each curve shows the residuals after the subtraction of the other planet and the stellar correlated signal. The red curve represents the best-fit Keplerian orbit, while the red dots and error bars represent the binned averages and standard deviations of the RVs.

Open with DEXTER
thumbnail Fig. 7

Upper panel: best fit stellar quasi-periodic signal (blue line) compared to the RV residuals. Lower panels: magnification of the high-cadence HIRES/KECK observations (left panel) and of the last HARPS-N observing season (right panel).

Open with DEXTER
thumbnail Fig. 8

Panel a: comparison between the last season HARPS-N data (black triangles) and the overlapping CARMENES data (grey dots). The red solid line represents our two-planet plus stellar noise model, with the pink shaded area representing the 1-σ uncertainties on the GP hyper-parameters; panel b: GLS periodogram of the combined HARPS-N and CARMENES RV datasets.

Open with DEXTER

5 Binary orbital interaction

The orbital eccentricity of the outer planet, Gl15A c, is quite uncertain: as shown in Fig. 5, the posterior distribution has no clear peak, only constraining the eccentricity towards small values, with a 68% probability to be <0.40 and 95% probability to be <0.72. This would normally point towards the adoption of a circular orbit as best fit solution for the system, lacking a significant evidence of eccentricity, but this would be to reckon without the stellar companion. The presence of Gl15B cannot be ignored, especially when studying a wide orbit such as that of Gl15A c.

We thus investigated how the dynamical interaction with the companion star could influence the Gl15A system. One of the main mechanisms to excite exoplanets eccentricities is the Lidov–Kozai effect, in which the presence of an external perturber causes oscillations of the eccentricity, e, and the inclination, i, with the same period but opposite phase. It was originally studied by Lidov (1962) and Kozai (1962) to compute the orbits of high inclination small solar system bodies, like asteroids and artificial satellites, and it is strongly dependent on the eccentricities and mutual inclination of the involved objects. Thus, in order to better estimate the strength of Lidov–Kozai oscillations in the Gl15A planetary system, we need to understand as precisely as possible the orbit of the stellar companion Gl15B.

5.1 Orbital modelling from astrometry and RV data

The first obstacle in the dynamical analysis of the system was the poorly constrained orbital parameters of the companion. Lippincott (1972) estimated a period of 2600 yr from ~ 100 yr of astrometric measurements, which coversless than 4% of the orbit.

Dealing with a similar case of an eccentric planet hosted by a wide binary system, Hauser & Marcy (1999) developed a technique for constraining long-period binary orbital parameters, combining astrometric and radial velocity measurements. The method is based on the fact that, being Newtonian-mechanics deterministic, knowing the instantaneous full position and velocity vector, [x, y, z, Vx, Vy, Vz], of one mass with respect to the other, you can compute the exact orbit of the system. Therefore, even by observing a small fragment of the orbit, we should be able to gather all the orbital parameters of the stellar companion. Of course astrometry, which is restrained in the plane of the sky, cannot provide the complete 3D information needed for this analysis, so additional data from radial velocity, to compute the third component of the velocity vector, and parallax distance, to convert astrometric positions in Cartesian coordinates, is necessary.

Following the procedure of Hauser & Marcy (1999), we downloaded 122 astrometric observations of the Gl15 system from the WDS catalogue, spanning from 1860 to 2015. The variations of the position angle θ of Gl15B relative to Gl15A and the angular separation ρ are shown in Fig. 9. To derive θ and ρ at a specific time, along with their derivatives dt and dt, which we need to calculate the velocity component in the plane of the sky, we fitted the data with a second-order polynomial:

From these we can easily derive fitdt and fitdt as

Adopting the parallax value from Table 1, πP = 280.74 ± 0.31 mas, we can derive the Cartesian position and velocity using their Eqs. (1)–(4):

Thus, we gained 4 of the desired physical components, [x, y, Vx, Vy], as a function of time t, through θfit and ρfit (Eqs. (6) and (7)).

To obtain the third component of the velocity vector, Vz, we used the combined Doppler information of the two stars: from the same HARPS-N spectra we used to obtain the RV time series for the planet detection, collected as illustrated in Sect. 2, we extracted the absolute radial velocities with the DRS pipeline. We used the DRS pipeline instead of TERRA since the latter only produce relative RV, which cannot be used in the comparison of two different objects. For Gl15A we take all the 115 HARPS-N epochs, and subtract from the absolute RV the planetary and stellar signals, as derived in Sect. 4.3. For Gl15B we used the five spectra taken on January 2017, as described in Sect. 2. The two datasets are illustrated in Fig. 10. To derive the binary orbit, we need to know all the position and velocity components at the same instant. Therefore, we fitted the two RV time series with first-order polynomials. From the difference of the two linear fits, we obtained the relative line-of-sight velocity Vz of Gl15B with respect to Gl15A. We were then able to then select an epoch, and compute the values of [x, y, Vx, Vy, Vz] for that time. The selected epoch to derive the binary orbit is BJD = 2457754.5, that is 1st January 2017, which is well in-between all the datasets, and close to the Gl15B RV time series, which is the shortest and most uncertain.

The last missing piece of the puzzle is, of course, the line-of-sight separation z, which cannot be measured, but can be constrained imposing the condition of a bound orbit for the binary system. This is expected due to the similar spectral type and proximity of the two stars. The condition of bound orbit translate into a condition on the total energy of the system: (14)

where μ = G(MA + MB), and, of course, .

From this we get a range of acceptable z values − 400 ≲ z ≲ 400 AU. From every value of z is possible to compute all the orbital parameters [P, e, i, ω, Ω, TP] for the binary system (see Hauser & Marcy 1999, Eqs. (7)–(15)). The results are shown in Fig. 11. As we can see, there is a wide variety of possible orbits, with completely different eccentricities and orientations, even with the bound orbit constraint, and this is still insufficient for any meaningful dynamical analysis.

The procedure by Hauser & Marcy (1999) provides the orbital solution for every single value of z, but does not in any way distinguish between the more probable configurations. But those orbital configurations are not all equally likely: from theory and observations of binary systems we know the expected distributions for different orbital parameters. From this information we can extract some priors to help us identify the most probable orbital configuration for the system, that is the best fit value of the line-of-sight separation z.

To do this we performed a Monte Carlo simulation in which the priors on the orbital parameters are injected via rejection sampling. Not all the a priori distributions of the orbital parameters are known, so we restricted the prior selection to the parameters that have a strong impact on the outcome and for which a good information is available. Due to the central role played by the eccentricity of the perturber in the Lidov–Kozai perturbation, we applied a prior on e, to have the best possible constraint on it. We selected the eccentricity distribution from Tokovinin & Kiyaeva (2016), who studieda sample of 477 wide binaries within 67 pc from the Sun with median projected separation of ~ 120 AU, very close to that of the Gl15AB system. The aforementioned prior is (15)

The second choice is a prior on the orbital period, in order to penalize long period poorly bound orbits. The chosen prior is the one suggested by Duchêne & Kraus (2013) for low-mass binary stars, that is a log-normal distribution, with ā ≈ 5.2 AU and σlogP ≈ 1.3.

The results of the Monte Carlo simulation are illustrated in Fig. 12. As we can see, there are three main peaks in the distribution. The third one, on the far right, represents orbits that are almost unbound, so it can be ignored, both because the two stars are expected to be in a binary system, and because even if bound, such wide orbits would have no influence whatsoever on the dynamics of the planetary system we intendedto study. The latter can also be said of the peak on the left, which corresponds to a period of yr, and a semi-major axis of AU, again too distant from the planetary system to have a significant influence.

We thus used the central peak of Fig. 12 to derive the orbital parameters best solutions and error bars. To do this we fitted a truncated normal distribution in the range z ∈ [−200, 200] AU, and use the mean value μz to calculate the corresponding best solution orbital parameters as described before. The upper and lower errors on the orbital parameters were calculated by taking μz + σz and μzσz and deriving the corresponding values of [P, e, i, ω, Ω, TP]. The orbital parameters solutions and errors are listed in Table 4.

Assuming the best fit orbital parameters listed in Table 4, we were able to calculate the RV signal caused by Gl15B on Gl15A, and compare it with the value of the residual acceleration found by our MCMC analysis. Gl15B RV signal is approximately linear, as is to be expected due to the small fraction of the orbit covered by the time series, and correspond to an acceleration m s−1 d−1, which is fully compatible with the result from our MCMC listed in Table 3.

thumbnail Fig. 9

Position angle (panel a) and angular separation (panel b) of Gl15B with respect to Gl15A. Position angle is measured from north towards east. The red lines in both panels represent the second order polynomial fits.

Open with DEXTER
thumbnail Fig. 10

Radial velocities of (panel a) Gl15A (after subtracting the planetary and stellar signals) and (panel b) Gl15B. The solid lines show the respective first order polynomial fits.

Open with DEXTER
thumbnail Fig. 11

Orbital parameters of Gl15B as a function of the line-of-sight separation z with Gl15A.

Open with DEXTER
thumbnail Fig. 12

Distribution of z resulting from the Monte Carlo simulation with the e and P priors injected.

Open with DEXTER
Table 4

Best orbital parameters for the Gl15 binary system from the MC simulation with priors on period and eccentricity.

5.2 Lidov–Kozai Interaction modelling

The results of the previous section show that the most likely orbit of the stellar companion Gl15B has a high eccentricity, e = 0.53. This is a further clue that strong orbital perturbation could affect the planetary system.

The interaction mechanism studied is the Eccentric Lidov–Kozai effect (commonly reffered to with the literature-coined acronym EKL). This mechanism applies to hierarchical triple-body systems, and consists in eccentricity and inclination oscillations on timescales much larger than the orbital period of the influenced body. We can safely ignore the mutual interaction between the two planets of the Gl15A system, except in the event of close encounters, and thus treat their interaction with the binary separately, as three-body systems.

The EKL mechanism is very sensitive to the mutual orientations of the orbits of the perturber, (the companion star), and of the influenced body (the planet). We have derived i, ω, Ω of Gl15B, but we do not knoweither i or Ω of the two planets Gl15A b and Gl15A c. Thus, some assumptions are to be made about their orbital orientation. To compare the results of the dynamical interaction model with the posterior distribution found by the MCMC analysis in Sect. 4.3, we calculated the fraction of time spent by Gl15A c below e = 0.40 and e = 0.72 (f(e < 0.40) and f(e < 0.72)), which can be considered a proxy of the probabilities (68% and 95%, respectively) to observe the system with eccentricities under those thresholds (e.g. Anderson & Lai 2017). The same can be done to study the interaction of Gl15B with Gl15A b, in which case the 68% and 95% probability levels correspond to e = 0.13 and e = 0.25.

We performed some preliminary tests to verify the various mechanisms which could be involved in the dynamical evolution of the system: we proved that the inner planet, Gl15A b, is too distant from Gl15B for any significant interaction. This is in good agreement with the extremely low level of eccentricity found in our MCMC analysis. We thus focussed our efforts on studying the orbit of the newly discovered Gl15A c; we also checked for the influence of dissipative tides, which could lessen the orbital eccentricity, and find that they act on much longer timescales than the EKL mechanism, so we neglected them in our analysis.

We denote the orbital parameters of Gl15B with the subscriptB and the ones of Gl15A c withc. All the EKL integrations were performed with a timescale of ~10 Myr. Since we cannot rule out planet–planet scattering events to have occurred in the early phases of the system’s evolution, we consider different values of the initial eccentricity ec,0 = 0.0, 0.1, 0.2, 0.3, 0.4. The other unknown is the longitude of the ascending node of the outer planet Ωc. However, the longitude of the ascending node influences the EKL interaction mainly by changing the relative inclination the two orbits, θ, which derives from the parameters of the two orbits as

and thus can be ignored so long as θ is conveniently sampled, which can be obtained by changing ic. For our analysis we fixed Ωc at 0° and varied the planet orbital inclination in the interval ic ∈ [−5°, 90°].

The results are shown in Fig. 13. In some cases the EKL oscillations are so extreme that the numerical integration has to be stopped due to the planet passing too close to the host star; the corresping systems are clearly unstable, and so we considered f(e < 0.40) = f(e < 0.72) = 0 to represent the incompatibility with the observed case. We also considered as unstable the cases in which the outer planet’s orbit becomes too close to the inner planet’s orbit, possibly producing planet–planet scattering events. In other terms, when (17)

where ec,max is the maximum eccentricity reached during the EKL oscillations. In these cases we also considered f(e < 0.40) = f(e < 0.72) = 0.

As we can see in both the panels of Fig. 13, there are regions in the ic space that are unstable regardless of the initial eccentricity of the planet ec,0: between 15° and 30° and around 0°. Only the solid black line, corresponding to ec,0 = 0, behaves differently, showing the system to be stable between 15° and 30° and unstable at larger inclinations. We can also see that the Lidov–Kozai interaction is weak for ic ~ 90° and strengthen as ic decreases. The top panel of Fig. 13 shows that for ic ∈ [75°, 90°] the resulting eccentricity ranges are compatible with the observed values. The constraints for the e < 0.72 threshold are somewhat looser, but pointing in the same direction.

As previously mentioned, we do not know the initial eccentricity of the system, so a more robust way to consider the dynamical evolution is to consider the average of the results of the single integrations. This can be seen in Fig. 14, which confirms the trends discussed above, with initial inclinations around 0° and between 15° and 30° leading to an unstable system, and inclinations higher than 75° producing a good agreement with the observed orbital configuration.

thumbnail Fig. 13

Fraction of time spent below the e = 0.40 (panel a) and e = 0.72 (panel b) thresholds. The different linestyles and symbols correspond to different values of ec,0: solid and plus signs – ec,0 = 0; dotted and asterisks – ec,0 = 0.1; dashed and diamonds – ec,0 = 0.2; dash dot andtriangles – ec,0 = 0.3; dash dot dot and squares ec,0 = 0.4. The horizontal dashed grey lines indicate the corresponding probability from the MCMC – panel a: 68%; panel b: 95%).

Open with DEXTER
thumbnail Fig. 14

Fraction of time spent below the e = 0.40 (solid grey) and e = 0.72 (dotted black) thresholds, averaged on the initial eccentricity ec,0. The errorbars show the standard deviation.

Open with DEXTER

6 Discussion and conclusions

We present in this paper the fifth planet detected by the HADES programme conducted with HARPS-N at TNG. This long-period planet was found orbiting the planet-host M1 star Gl15A, from the analysis of high precision, high-resolution RV measurements collected as part of the survey in conjunction with archive RV data from the HIRES/Keck spectrograph.

The different trends observed in the two datasets suggest the presence of a long-period companion, which is confirmed by the homogeneous Bayesian analysis of the combined RV time series. The known inner planet Gl15A b is also recovered. The minimum masses derived from our analysis are and for the inner and outer planets, respectively. The mass we find for Gl15A b is much smaller than the value found by Howard et al. (2014), which was almost double due to the higher signals amplitude. The smaller value that we find can be easily explained by the additional information brought by the high-precision HARPS-N RVs, along with the new calibration of the archival HIRES data published by Butler et al. (2017). The combined dataset is almost twice as large that the one analysed by Howard et al. (2014), stretched on a significantly longer time span, with better sampling and precision. This, together with the simultaneous modelling of the stellar activity signal, can explain the much smaller uncertainty on the minimum mass of Gl15A b. It also highlights the importance of taking into account the chromospheric stellar activity for the correct identification of planetary signals. It is worth noticing that instead the orbital period Pb is almost unchanged from the previous estimate. Our fit also places an even lower upper limit to the eccentricity of Gl15A b than that found by Howard et al. (2014;e < 0.13 at a 68% level of confidence), thus confirming their conclusion that the planet’s motion is best described by a circular orbit.

With its period of ≃21 yr, Gl15A c is the longest-period sub-Jovian planet detected up to date with the RV method3, the second being HD 10180 h (Lovis et al. 2011b) with a period of approximately six years, and a minimum mass of 65.74M (Kane & Gelino 2014). With the confirmed presence of two widely spaced planetary mass companions, Gl15A is now the multi-planet system closest to our Sun, at a distance of only 3.5620 ± 0.0039 pc, slightly smaller than the distance to YZ Cet (Astudillo-Defru et al. 2017a), 3.69 ± 0.11 pc.

We also compared the results of our analysis of the HIRES and HARPS-N RV data with the recent results by Trifonov et al. (2018) based on CARMENES high-cadence monitoring of the target. Trifonov et al. (2018) found no evidence for the presence of planet b, while we showed in Sect. 4.4 that the signal can be clearly detected when combining the HIRES, HARPS-N and CARMENES data, without any loss in significance. The reason why CARMENES does not see the signal of Gl15A b is not fully clear, but the higher quality, and much longer time span, of the HARPS-N data, combined with our modelling of the stellar activity quasi-periodic signal could be a possible explanations for the non-detection based on the CARMENES data alone. Another factor contributing to this non-detection could be the sampling of the CARMENES RV data, which appears not to be homogeneously distributed over the 11.44 d period (see Fig. 10 from Trifonov et al. 2018).

The CARMENES visual arm contains a spectral region that extends all the way to 0.95μm; in other words, a significantly redder spectral range than that covered by HIRES and HARPS-N. As the amplitude of activity induced RV variations is known to be chromatic (Reiners et al. 2010), the non-detection of the 11.44 d signal in the CARMENES time series could be an indication of a wavelength dependent-amplitude of the signal, which would clearly indicate its stellar origin. However, the CARMENES visual arm spectral range still significantly overlaps with that of HIRES and HARPS-N and thus it must be affected by stellar activity in a similar way. Based on the higher RV precision of HARPS-N, allowing a detailed modelling of quasi-periodic stellar signal (as shown in Fig. 7), and on the coherence of the period and phase of the signal over the 20 yr time span covered by the combined HIRES and HARPS-N time series, the Keplerian origin of the signal still seems the most straightforward explanation for the observed data. It would, however, be interesting to carry out a systematic study on both CARMENES and HARPS-N time series of this target, adopting the same strategy outlined by Feng et al. (2017) for Tau Ceti, that is to study separately the RVs derived from different regions of the spectra, for a sistematic investigation of potential differences in the amplitude of the 11.44 d signal, as a function of the wavelength, but this lies well beyond the scope of this paper.

Dwarf stars are known to turn up much more frequently in multiple systems than they do in isolation, with a binary fraction as high as ≃57% for Sun-like stars (Duquennoy & Mayor 1991) and somewhat lower for M dwarfs (Bergfors et al. 2012). Many young binaries possess either circumstellar or circumbinary disks (e.g. Monin et al. 2007), and the existence of stable planetary orbits in binary systems was postulated well in advance of the first exoplanets discoveries (Dvorak 1982).

Early studies proposed different mass-period relations for planets around binaries and single stars (Zucker & Mazeh 2002) but in the following years the evidence of such diversity decreased (e.g. Desidera & Barbieri 2007; Eggenberger 2010), until most recently Ngo et al. (2017) claimed there not to be any difference of planetary properties between the two kind of systems. On the other hand, recent works such as Moutou et al. (2017) find statistical evidence for a much higher binary fraction in extrasolar systems hosting eccentric exoplanets than in the ones hosting only circular planets: this points towards the confirmation of the role of stellar multiplicity in orbital excitation of planetary systems, as predicted by theoretical studies which suggested a strong orbital influence of stellar companions on planetary systems, via mechanisms such as the eccentric Lidov–Kozai (EKL) oscillations.

Our numerical analysis of the EKL effect proved the strong influence of the Gl15B on the planetary system. We show that for a narrow range of initial inclination, 75°−90°, the outer planet maintains a low eccentricity orbit, regardless of the initial status in which the system was due to possible planet–planet scattering events. We also pointed out the presence of a forbidden ranges of inclination, 15°− 30° and ~ 0°, in which the Lidov–Kozai interaction become so strong that no stable orbit can be achieved, regardless of the initial eccentricity of Gl15A c.

The orbital parameters of Gl15A c have still large uncertainties due to the observation time span shorter than the orbital periods, and the semi-amplitude Kc is significant only at a 3-σ level, although the strong combined observational evidence from RV and imaging leaves no doubt as to the presence of a long-period planetary-mass companion. Additional RV observations in the years to come will, however, be very helpful in better constrainingthe orbit, and thus the mass of Gl15A c.

Our knowledge of this system will be also greatly improved by the results published in future Gaia data releases. For a circular orbit and assuming the minimum mass value for Gl15Ac, the expected astrometric signature on the primary is 570 μas. Gaia astrometry will only cover ~ 20%−25% of the full orbit. However, based on the Torres (1999) formalism curvature effects in the stellar motion should typically amount to 20− 30μas yr−2, thus they should be easily revealed by Gaia, that for such a bright star as Gl15A will be able to deliver end-of-mission proper motion accuracies ≲10μas yr−1 (e.g. Gaia Collaboration 2016b). Even with the orbital solutions now available, our analysis shows how interesting dynamical studies can be performed on the system, which, due to the presence of the eccentric binary companion, is an excellent playground to test orbital interaction mechanisms and their influence on the evolution of planetary systems.

Acknowledgements

GAPS acknowledges support from INAF through the “Progetti Premiali” funding scheme of the Italian Ministry of Education, University, and Research. The research leading to these results has received funding from the European Union Seventh Framework Programme (FP7/2007-2013) under Grant Agreement No. 313014 (ETAEARTH). The HARPS-N Project is a collaboration between the Astronomical Observatory of the Geneva University (lead), the CfA in Cambridge, the Universities of St. Andrews and Edinburgh, the Queen’s University of Belfast, and the TNG-INAF Observatory. M.P. thanks E.T. Russo for the helpful discussions in the final stages of the analysis. J.I.G.H., R.R.L., A.S.M., and B.T.P. acknowledge financial support from the Spanish Ministry project MINECO AYA2014-56359-P. J.I.G.H. also acknowledges financial support from the Spanish MINECO under the 2013 Ramón Cajal programme MINECO RYC-2013-14875. A.S.M also acknowledges financial support from the Swiss National Science Foundation (SNSF). We acknowledge the computing centres of INAF – Osservatorio Astronomico di Trieste/Osservatorio Astrofisico di Catania, under the coordination of the CHIPP project, for the availability of computing resources and support. This research has made use of the Washington Double Star Catalog maintained at the U.S. Naval Observatory. We gratefully acknowledge an anonymous referee for her/his insightful comments that materially improved an earlier version of this manuscript.

Appendix A Observation log for Gl15A

In this section we report the observational data collected with the HARPS-N spectrograph as part the HADES project and used in the present study. We list the observation dates (barycentric Julian date or BJD), the radial velocities (RVs), and the Hα and SHK indices, calculated by the TERRA pipeline. The RV values are corrected for perspective acceleration. The RV errors reported are the formal ones, not including the jitter term, while the Hα and SHK errors are due to photon noise through error propagation.

Table A.1

Data of the 115 observed HARPS-N spectra of Gl15A.

References


All Tables

Table 1

Stellar parameters for the two stars Gl15A and Gl15B.

Table 2

Priors and best-fit results for the Gaussian process regression analysis of the chromospheric activity indicators extracted from the HARPS-N spectra of Gl15A.

Table 3

Planetary parameter best-fit values obtained through a joint modelling of Keplerian signals and correlated stellar noise, using Gaussian process regression.

Table 4

Best orbital parameters for the Gl15 binary system from the MC simulation with priors on period and eccentricity.

Table A.1

Data of the 115 observed HARPS-N spectra of Gl15A.

All Figures

thumbnail Fig. 1

Combined HARPS-N and HIRES RV time series. As described in Sect. 4.3, we considered the HIRES data as splitinto two separate datasets for the purpose of the RV analysis.

Open with DEXTER
In the text
thumbnail Fig. 2

Posterior distributions of the fitted (hyper)parameters of the GP quasi-periodic model applied to the time series of CaII H&K (upper panel) and Hα (lower panel) activity indexes. The vertical lines denote the median (solid) and the 16th–84th percentiles (dashed).

Open with DEXTER
In the text
thumbnail Fig. 3

GLS periodogram of the complete HARPS-N RV dataset.

Open with DEXTER
In the text
thumbnail Fig. 4

Panel a: evolution of the GLS orbital parameters as a function of the number of HARPS-N RV measurements; panel b: evolution of the GLS orbital parameters as a function of the number of simulated data. The error bars on the final points are shown as references.

Open with DEXTER
In the text
thumbnail Fig. 5

Posterior distributions of the fitted (hyper)parameters of the two-planet model, where the stellar correlated noise has been modelled with a GP regression using a quasi-periodic kernel. On the y-axis is shown the logarithm of the product between the likelihood and the prior. The vertical lines denote the median (solid) and the 16th–84th percentiles (dashed).

Open with DEXTER
In the text
thumbnail Fig. 6

Phase-folded RV curves for panel a Gl15A b and panel b Gl15A c. Each curve shows the residuals after the subtraction of the other planet and the stellar correlated signal. The red curve represents the best-fit Keplerian orbit, while the red dots and error bars represent the binned averages and standard deviations of the RVs.

Open with DEXTER
In the text
thumbnail Fig. 7

Upper panel: best fit stellar quasi-periodic signal (blue line) compared to the RV residuals. Lower panels: magnification of the high-cadence HIRES/KECK observations (left panel) and of the last HARPS-N observing season (right panel).

Open with DEXTER
In the text
thumbnail Fig. 8

Panel a: comparison between the last season HARPS-N data (black triangles) and the overlapping CARMENES data (grey dots). The red solid line represents our two-planet plus stellar noise model, with the pink shaded area representing the 1-σ uncertainties on the GP hyper-parameters; panel b: GLS periodogram of the combined HARPS-N and CARMENES RV datasets.

Open with DEXTER
In the text
thumbnail Fig. 9

Position angle (panel a) and angular separation (panel b) of Gl15B with respect to Gl15A. Position angle is measured from north towards east. The red lines in both panels represent the second order polynomial fits.

Open with DEXTER
In the text
thumbnail Fig. 10

Radial velocities of (panel a) Gl15A (after subtracting the planetary and stellar signals) and (panel b) Gl15B. The solid lines show the respective first order polynomial fits.

Open with DEXTER
In the text
thumbnail Fig. 11

Orbital parameters of Gl15B as a function of the line-of-sight separation z with Gl15A.

Open with DEXTER
In the text
thumbnail Fig. 12

Distribution of z resulting from the Monte Carlo simulation with the e and P priors injected.

Open with DEXTER
In the text
thumbnail Fig. 13

Fraction of time spent below the e = 0.40 (panel a) and e = 0.72 (panel b) thresholds. The different linestyles and symbols correspond to different values of ec,0: solid and plus signs – ec,0 = 0; dotted and asterisks – ec,0 = 0.1; dashed and diamonds – ec,0 = 0.2; dash dot andtriangles – ec,0 = 0.3; dash dot dot and squares ec,0 = 0.4. The horizontal dashed grey lines indicate the corresponding probability from the MCMC – panel a: 68%; panel b: 95%).

Open with DEXTER
In the text
thumbnail Fig. 14

Fraction of time spent below the e = 0.40 (solid grey) and e = 0.72 (dotted black) thresholds, averaged on the initial eccentricity ec,0. The errorbars show the standard deviation.

Open with DEXTER
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.