EDP Sciences
Free Access
Issue
A&A
Volume 599, March 2017
Article Number A90
Number of page(s) 10
Section Planets and planetary systems
DOI https://doi.org/10.1051/0004-6361/201629484
Published online 06 March 2017

© ESO, 2017

1. Introduction

One of the emerging scenarios of the exoplanets population, as revealed from the observations, is that planets are preferably organized in multiple systems (see e.g. Mayor et al. 2011; Latham et al. 2011). Nevertheless, no evidence of a strict analog of our solar system has been found so far, but hundreds of systems with a large variety of architectures (Fabrycky et al. 2014). The study of the planet distribution in multiple systems allows us to put constraints on theories for their complex formation, dynamics and evolution. Radial velocity (RV) surveys can play a crucial role in this framework for several reasons. The availability of new generation instruments, able to push down the minimum-mass detection limit, allows nowadays to search for additional low-mass companions in systems already known to host giant planets. This provides a more complete view of their architecture as well as indications on the frequency of systems similar to the solar one. A full characterization of multiplanet systems requires great observational effort, as the number of components increase. In addition, several sources of astrophysical noise affect the RV (stellar oscillations, granulation and activity) and could impact on the measurements, either mimicking the presence of a planetary companion or tangling the interpretation of the signal (Dumusque et al. 2011). Since August 2012 the GAPS (Global Architecture of Planetary Systems, see e.g. Covino et al. 2013) observing programme started its operations thanks to the high performances of the HARPS-N spectrograph (Cosentino et al. 2012), mounted at the Italian telescope TNG in La Palma, Canary Islands. Within GAPS we studied the G9V star HD 108874 in the framework of a dedicated sub-programme focused on the characterization of systems with already known planets. Butler et al. (2003) claimed the presence of a giant planet (Msini = 1.71 MJ) with a period of 397.5 ± 4.0 days (d), while Vogt et al. (2005) found an additional companion on a wider orbit (Msini = 1.02 ± 0.3MJ, Porb = 1605.8 ± 88.0 d) through the analysis of ~70 spectra from HIRES (Keck telescope). The presence of these planets was subsequently confirmed by other authors (Wright et al. 2009, hereafter Wr09, and Wittenmyer et al. 2009, hereafter Wi09), while further studies (Goździewski et al. 2006; Libert & Henrard 2007; Veras & Ford 2010) analysed the dynamics of this system in the presence of a 4:1 mean motion resonance between them. The planet HD 108874 b is probably located in the habitable zone of its host star and this stimulated Schwarz et al. (2007) to investigate possible stable regions for Earth-like Trojan planets around it. Here we present an updated analysis of the HD 108874 system based on a three-year intensive monitoring with HARPS-N: we describe the observations, the data reduction and the derivation of stellar parameters in Sects. 2 and 3; we present our RV and stellar activity analysis in Sects. 4 and 5; a new orbital solution is shown in Sect. 6; the dynamical analysis is described in Sect. 7; a discussion on the detection limits, the system architecture and its habitability is provided in Sect. 8.

2. Observations and data reduction

The GAPS observations of HD 108874 with HARPS-N at TNG lasted three seasons, from December 2012 to July 2015. During the first half of the observations the spectrograph was affected by a small defocus, evident in the time series of some parameters (FWHM and contrast of the Cross-Correlation Function, CCF) but negligible for others, like the bisector or the RVs, only causing a small increase in the estimated errors. The problem was successfully fixed in March 2014. The total number of collected spectra is 94, spread over 82 nights. The simultaneous Th-Ar calibration was used to obtain the required RV precision, and the median of the instrumental drift, monitored through the second fiber of HARPS-N, was 0.39 m s-1, with a rms of the total time series of 0.76 m s-1. The spectra were obtained with 900s integration time and the average value of the signal-to-noise ratio (S/N) is 110 per pixel on the extracted spectrum at 5500 Å. The data reduction and the RV measurements were performed by means of the data reduction software of HARPS-N (Pepe et al. 2002). Through the method of the CCF of the acquired spectrum by using a mask that depicts spectral features of a G2 star, we obtained the RV measurements listed in Table A.1. The median of the internal errors is σRV = 0.7 m s-1.

3. Stellar parameters

The atmospheric parameters of HD 108874 are measured as in Biazzo et al. (2011), based on the line equivalent widths measurements, by using the 2013 version of the MOOG code (Sneden 1973) and the line list in Biazzo et al. (2012). The analysis is performed on a merged spectrum, obtained by coadding the available spectra of the target after the correction of the corresponding radial velocity shift, showing a S/N of ~1200 at 5500 Å. A summary of the extracted parameters is presented in Table 1. We used the web interface PARAM1 (da Silva et al. 2006) which is based on isochrones by Bressan et al. (2012) for the estimation of the stellar mass, radius, and age. Besides the effective temperature and the metallicity, we also included the parallax (15.97 ± 1.07 mas, van Leeuwen 2007) and the V mag of the star (8.76 ± 0.02, Yoss & Griffin 1997) as input quantities. A general agreement with the literature is found for all of our estimates (e.g. Valenti & Fischer 2005, Wi09, Torres et al. 2010). The spectral analysis yields a value of vsini = 1.6 ± 0.5 km s-1 (see D’Orazi et al. 2011). A calibration of the FWHM of the CCF, using stars with transiting planets with known photometric rotation periods observed with HARPS-N, provides a value of vsini = 1.36 ± 0.26 km s-1, which we adopt in the following. The average value of the activity index , provided by the HARPS-N pipeline (see Sect. 5.1), is equal to –5.05 and indicates that HD 108874 is less active than the Sun. The rotation period measured in this work is also reported (see the following sections).

Table 1

Stellar parameters of HD 108874.

4. Radial velocity analysis

We extended our RV dataset by considering the data available in the literature from HIRES at Keck and HRS at HET (see Table 2), obtaining a total time span of 16 yr.

Table 2

Summary of the datasets.

We performed a two-planets fit, including a zero-point correction of the RVs (see Sect. 6), obtaining the model overplotted to the three sets of data in panel a of Fig. 1. Panel b shows only HARPS-N data.

thumbnail Fig. 1

Panel a): RV dataset for HD 108874. Green dots: HIRES at Keck; orange dots: HRS at HET; blue dots: HARPS-N at TNG; black solid line: global fit of the three sets. Panel b): RV time series for HARPS-N data only. Panel c): residuals after the subtraction of the two planets fit.

Open with DEXTER

After the subtraction of this fit from the original time series, the resulting residuals show an rms of 4.1 m s-1 for HIRES data, 5.9 m s-1 for HRS (smaller than the typical errors) and 2.8 m s-1 for HARPS-N. The rms of the residuals for the whole dataset is 4.03 m s-1. The Lomb-Scargle periodogram of the RV residuals of HARPS-N (Fig. 2) shows a clear periodicity at ~0.0248 d-1 with an uncertainity of 1.0 × 10-4 d-1 evaluated with the relation in Montgomery & O’Donoghue (1999), corresponding to 40.20 ± 0.15 d in the domain of periods, with a normalized power equal to 16.3 (Horne & Baliunas 1986) and a confidence level higher than 99.99%, obtained after 100 000 bootstrap random permutations. The spectral window of the HARPS-N dataset has been computed (inset of Fig. 2) as in Deeming (1975): the first relevant peak in our region of interest is at 0.0025 d-1, corresponding to one cycle per year and is due to the visibility of the star. Many other small peaks are also present between 0.045 d-1 and 0.065 d-1 (i.e., 15–22 d), probably due to the scheduling of the GAPS observing runs. Therefore, it does not seem that the periodicity of 40.2 d is due to an aliasing effect with the orbital periods of the two planets.

thumbnail Fig. 2

Periodogram of RV residuals of HARPS-N dataset: both frequency and period domains are shown. Horizontal dashed lines represent the confidence level of the peaks. The inset shows the window function due to the temporal sampling of the data.

Open with DEXTER

Being this periodicity very robust, we performed a three-planets fit to the data, with the addition of a further Keplerian function aiming to model the signal in the RV residuals. The resulting r.m.s of the time series decreased from 4.03 to 2.2 m s-1 and the periodogram does not show any residual power around 40.2 d. The fit implies a minimum mass of about 17 M and a circular orbit for the planet candidate. A few non-significant peaks are found around 40 d in the periodogram of residuals for HIRES and HRS, explainable by sub-optimal temporal sampling, but there is also the possibility that this signal is not driven by Keplerian motion. For this reason we analysed the single observing seasons of HARPS-N in order to verify the presence of this periodicity from one year to the next (Fig. 3). In the figure, the location of P = 40.2 d is indicated as a reference, and the second and the third harmonics of the main peak are also shown (20.1 d, 13.4 d and 10.05 d). The signal at 40.2 d is present only in the first season (panel a, with a confidence level of 99.97%), and marginally in the third one (panel c), together with its harmonics. In the second season (panel b) the 40.2 d peak is suppressed: we only observe two large features at 25 and 73 d probably associated to each other.

thumbnail Fig. 3

Periodogram of HARPS-N RV residuals subdivided into each observing season. The red solid line marks the location of the 40.2 d periodicity, while the red dotted lines represent the 1st, 2nd and the 3rd harmonics of the main peak.

Open with DEXTER

5. Stellar activity

5.1. log RHK and Hα

The chromospheric emission from CaII H and K lines of HARPS-N spectra (log ) is directly provided by the HARPS-N pipeline (Lovis et al. 2011). We added the S-index measurements obtained from the HIRES dataset by Wright et al. (2004) (both datasets are calibrated with the Mt. Wilson activity survey, so no offset should be present) and converted into log by using the scaling relations by Noyes et al. (1984). We perform a tentative sinusoidal fit of the stellar activity cycle with the Levenberg–Marquardt fitting algorithm through the IDL package MPFIT2 (Fig. 7, a). According to the parameters of our fit, the length of this cycle is ~19 yr, but the huge gap in the data between 2003 and 2012 and an insufficient temporal coverage does not ensure the goodness of the fit. Since the long term correction does not lead to different results, we show the results only for the uncorrected time series. The upper left panel of Fig. 4 shows the periodogram of the . The periodicity at 40.2 d detected in the periodogram of the RV residuals is shown, located very close to one of the strongest peaks in the plot, having a statistical significance of 99.75%. Its third harmonics could be the responsible for the largest periodicity in the periodogram, around 10 d. In the lower left panel we compare the RV residuals with the values of : they show a moderate positive correlation with Spearman and Pearson coefficients of ~0.3 (see Table 3 for a summary of correlation coefficients), slightly lower than the values obtained by Lanza et al. (2016) in the solar case, 0.35 and 0.38, respectively. This is a remarkable indication of the physical origin of the 40.2 d periodicity, leading us to identify it with the rotation period, Prot, of the star. Furthermore, following Mamajek & Hillenbrand (2008) the expected rotation period for our target is 39.9 d, derived from the mean value of and the BV (0.738). This value, considering the stellar radius derived in Sect. 3, is also consistent with the adopted vsini for an edge-on inclination: larger vsini yields unphysical values for the inclination. We also analysed the periodograms of the measurements for separate seasons (Table 3). As in the case of the RVs residuals, the periodicity at 40.2 d is clearly present only in the first season while we only observe a periodicity of ~20 d (first harmonic of Prot) in the third one, also present in the RV residuals.

The HARPS-N spectra of HD 108874 were also analysed to extract the time series of the Hα index (following Gomes da Silva et al. 2011): the periodogram analysis reveals an excess of power distributed in a narrow envelope of peaks around 40 d, in agreement with the adopted Prot (upper right panel of Fig. 4). No significant correlation is found with the RV residuals (see Table 3). Despite the periodograms of the Hα and indices show similar periodicities, a very weak correlation is found between these two quantities (CP = 0.18, p-value = 0.16). This can be explained by the presence of plages on the stellar disc, revealed both by the Ca II and Hα lines, and filaments able to modify the emission of the Hα and to break the correlation between the time series, as demonstrated by Meunier & Delfosse (2009), Gomes da Silva et al. (2014), and by Scandariato et al. (2017), for early-M dwarfs. When the different timescales of these phenomena are modulated by the stellar rotation, their periodograms can actually show periodicities related to it.

thumbnail Fig. 4

Upper panels: periodograms of the (with S/N> 20) and of the Hα indices for HARPS-N data. The red solid line indicates the 40.2 d period, the red dotted lines represent its 1st, 2nd and the 3rd harmonics; the blue dotted line indicates the period of planet b. Lower panels: correlation with the RV residuals.

Open with DEXTER

Table 3

Summary of the values of Pearson (CP) and Spearman (ρ) correlation coefficients and corresponding p-values between RV residuals and all the activity indices investigated in this work.

5.2. Asymmetry indicators of the CCF

Stellar activity also results in a deformation of the line profile of the spectral lines, which can be quantified by several asymmetry indicators. We investigated these quantities that provide an independent evaluation of the stellar activity with respect to the chromospheric indices. A measurement of the CCF bisector velocity span (BVS) is directly provided by the HARPS-N pipeline, while a procedure presented in Lanza et al. (in prep.) estimates the values and the errors of ΔV (representing the RV shift produced by the asymmetry alone as defined in Nardetto et al. 2006 and reconsidered by Figueira et al. 2013) and of the quantity Vasy(mod)3, a modified version of the Vasy as defined by Figueira et al. (2013) but not affected by the RV shifts as in the case of the original definition. Periodograms of the asymmetry indices time series are presented in Fig. 5: all of them show a moderate amount of power around 40 d. Even in this case, there is no evidence of a clear correlation between the RV residuals and line profile indicators (Table 3).

thumbnail Fig. 5

Periodograms of the asymmetry indicators as derived by the HARPS-N pipeline (BVS) and Lanza et al. (in prep; ΔV, Vasy(mod)). Red solid and dotted lines represent the assumed Prot and its 1st, 2nd and the 3rd harmonics, respectively, while the blue dotted line indicates the period of planet b.

Open with DEXTER

Finally we analysed the FWHM of the CCF, initially affected by the HARPS-N defocusing which progressively enlarged it (left upper panel of Fig. 6). We tried to remove this effect by performing a polynomial fit, shown in the figure, before the focus correction. The resulting residuals (left lower panel) produce the periodogram in the right panel of Fig. 6. Besides a long-term periodicity close to the period of planet b, we found a peak around 20 d (Prot/ 2) and some power close to Prot, indicated on the figure. The correlation coefficients between RV residuals and the corrected FWHM in Table 3 (CP ~ ρ ~ 0.2), show a very weak linear correlation.

thumbnail Fig. 6

Left upper panel: original time series of the CCF FWHM, showing the trend due to the defocusing. Left lower panel: residuals of the CCF FWHM after the subtraction of the fit. Right panel: periodogram of the residual CCF FWHM. The red solid line represents the Prot, the dotted lines represent its 1st, 2nd and the 3rd harmonics.

Open with DEXTER

5.3. Activity modelling of raw time series

We modelled the contribution of the activity in the full RV dataset from HIRES, HRS and HARPS-N, following the approach of Boisse et al. (2011) as recently implemented in the open-source code PyORBIT4 (see Malavolta et al. 2016 for the details on the activity model, its implementation and the steps involved in parameter estimation). Due to the lack of simultaneous photometric data, we used the value of Prot from the activity indexes as prior and constrained it within ±0.5 d. We considered several combinations of harmonics for each dataset, and in all cases we found that the harmonics had RV amplitudes consistent with zero, so we decided to use only the sinusoids associated to Prot in the final fit. When including the activity model in the RV fit, the jitter parameter is reduced by 30% for the HIRES and HARPS-N datasets, while no improvement is visible in the (noisier) HRS dataset. The orbital parameters of the two planet companions are not affected by the activity correction since the RV modulation has a shorter time-scale with respect to their orbital periods, so the activity noise is averaged out during the fitting process. In Fig. 7 (panels b, c, d) we show the behaviour of the activity with time by comparing the fit of the activity cycle of the star described in Sect. 5.1 (a), the seasonal values of the RV jitter parameter when the activity model is not included in the fit (b), and the semi-amplitude and phase of the fitted RV sinusoid respectively (c, d). The jitter terms of the HIRES dataset show a gradual reduction, which may explain the negative slope of the corresponding activity cycle on the left side of panel a. The HRS jitter is probably dominated by instrumental errors, while HARPS-N shows a small amount of jitter, despite the correspondence at the maximum of the activity cycle. The amplitude of the sinusoid in the first season of HARPS-N data (c) is larger than that in the subsequent two seasons, explaining the excess of scatter in the time series (a).

thumbnail Fig. 7

a) Tentative fit of the stellar activity cycle for HARPS-N and HIRES data of log . Data points (grey dots) are plotted together with the mean values of each observing season (red dots). b) Fitted RV jitters for the separate observing seasons of HIRES, HRS and HARPS-N datasets, when activity model is not included in the data. c) Semi-amplitude of the fitted RV sinusoid. d) Phase of the fitted RV sinusoid.

Open with DEXTER

6. Refined orbital solution

To refine the orbital parameters of the two known planets around HD 108874, we modelled both the literature and our new HARPS-N data with two Keplerians, by including three RV zero points and three uncorrelated RV jitter terms for each dataset (HIRES, HRS, and HARPS-N). The jitter terms that account for possible RV scatter exceeding the nominal error bars were added in quadrature to the RV uncertainties, regardless of the origin of the jitter (stellar activity, instrumental effects, and/or a combination thereof). In addition to the offset and jitter terms, we fitted for the inferior conjunction times, the orbital periods, the RV semi-amplitudes, the orbital eccentricity, the argument of periastron and and (Ford 2005) of the two planets. Therefore, our model has 16 free parameters in total. We determined the posterior distributions of the model parameters with a Bayesian differential evolution Markov chain Monte Carlo (DE-MCMC) approach (Ter Braak 2006; Eastman et al. 2013). We ran 32 chains and used the same criteria as in Bonomo et al. (2014) and Desidera et al. (2014) for the removal of burn-in steps, convergence and good mixing of the chains. For all the parameters we considered uninformative priors. The values of fitted and derived system parameters and their 1σ uncertainties, which were computed as the medians and the 15.86% and 84.14% quantiles of their posterior distributions, are listed in Table 4. A general agreement is found between the literature and our results for the orbital parameters, except for the eccentricity of planet b. Our estimate, equal to 0.142 ± 0.011, is close to the value reported by Wr09, confirming that the planet is slightly more eccentric than the value derived by Wi09 (0.082 ± 0.021). A significant difference has been found for the orbital period of planet c: HARPS-N data revealed a larger value of the period that differs by 4σ from the value found by Wi09 (Pc = 1620 ± 24 d), and by 2σ from the one in Wr09 (Pc = 1680 ± 24 d). Our estimates rule out the supposed mean motion resonance (MMR) of 4:1, revealing that the system is close to a 9:2 MMR. Radial velocity residuals of the two-planets fit show no evidence of a long term trend, in agreement with the results presented in Bryan et al. (2016). The minimum masses of HD 108874 b and c are 1.255 and 1.094 MJ, with semi-major axes of 1.05 and 2.81 AU, respectively.

Table 4

Orbital parameters of HD 108874 b and c as derived from the present analysis.

7. System stability

We have investigated the long term stability of the system (by using the parameters in Table 4) first with a long term numerical integration of the planet orbits, obtaining a “nominal” solution, and then with a parametric exploration of the phase space around this solution. The direct N-body integration of the system over 1 Gyr, performed with SyMBA (the Symplectic Massive Body Algorithm; Duncan et al. 1998), confirms that it is stable with a quasi-periodic behaviour over that timescale. The exploration of the phase space around the nominal solution is performed with the Frequency Map Analysis (FMA; Laskar 1993; Šidlichovský & Nesvorný 1996; Marzari et al. 2003). This is a numerical tool for the detection of chaotic behaviour based on the analysis of the variation of the secular frequencies of the system. Its main advantage is in allowing users to identify unstable orbits with short term numerical integrations. The orbits of about 40 000 systems with orbital elements close to the nominal ones have been integrated with SyMBA over 5 Myr and their secular frequencies are analysed with the FMA to measure the diffusion speed in the phase space. Figure 8 shows the outcome of the FMA as a function of the initial values of the semi-major axes (a1 and a2) of the two planets. The initial eccentricities have been randomly selected around the nominal values, the semi-major axes are sampled in an interval given by the nominal value ±0.1 AU, while the initial inclinations are chosen between and . The unknown angles, that is, the initial mean anomalies and node longitudes, have been randomly chosen between and 360°, while the pericentre argument is taken in between the nominal value ±20°. The diffusion speed is measured as the dispersion of the main secular frequency of the system over running windows covering the integration timespan. The three large instability regions in Fig. 8 are related to the 5:1, 9:2 and 4:1 mean motion resonances from left to right, respectively.

thumbnail Fig. 8

Stability of the HD 108874 system. The values of the semi-major axes of the two planets are sampled around the nominal system (black square) according to the dispersion of the main secular frequency. The black box marks the 1σ uncertainties given in Table 4. The three yellow strips indicate the instability zones related to the 5:1, 9:2 and 4:1 MMR, respectively. A logarithm scale is used for the colour coding, where small values of the dispersion mean stability, while large values imply fast changes of the secular frequency and then chaotic evolution.

Open with DEXTER

Our nominal system (shown in the figure, along with the current uncertainty limits) lies in a stable zone located close to the 9:2 resonance. Close to a resonance, orbits may be chaotic while two planets trapped in MMR can be stable for long times depending on their location within the resonant region. A detailed exploration of the resonant behaviour would be needed if the planets were in MMR, however in our case the resonance is present only in the outer margin of the uncertainty box so that, statistically, this configuration is less probable than a non–resonant one. As in Wi09, we tested the possible existence of lower mass additional planets in the system. We performed a detailed investigation of the phase space between 0.5 and 10 AU within a full four-body model with the star, the two planets on their nominal orbits and a putative 1 M terrestrial planet. In Fig. 9 (left panel) we show the FMA of the region inside the orbits of the two nominal planets. Unstable systems over short timescales cannot be analysed and are represented by empty regions in the plots. By assuming 10-4 as a reference value between stable and unstable behaviour, a value tested with a few long term numerical integrations over 1 Gyr, two stable regions can be identified, one close to 0.1 AU and one extending from 0.25 to about 0.4 AU. Mean motion and secular resonances lead to chaotic evolution between 0.1 and 0.25 AU and beyond about 0.45 AU. In the region outside the two giant planets (right panel of Fig. 9), the stability is met around 8.5 AU and beyond, especially for low values of eccentricity.

thumbnail Fig. 9

Stability of a 1 M planet in the inner (left) and the outer regions (right) of the HD 108874 system as a function of the initial eccentricity and semi-major axis of the putative bodies. Colour scale ranges from black, representing stable regions, to yellow, which indicates higher degree of instability for planets with given separation and eccentricity. The empty (white) spaces indicate systems unstable over short timescales.

Open with DEXTER

8. The planetary system of HD 108874

8.1. Detection limits to additional planetary companions

On the basis of the HARPS-N dataset we evaluated the upper limits of the minimum mass for possible planetary companions as in Sozzetti et al. (2009), with a 99% confidence level, based on the F-test and χ2 statistics. Super Earths with Msini ≳ 5M in the HD 108874 system interior of 0.5 AU and objects with Msini ≳ 2M with orbital periods of a few days could have been detected. The grey areas in Fig. 10 represent the allowed regions for a stable Earth mass planet, as derived in Sect. 7, and the vertical dashed line indicates the orbital period of planet b.

thumbnail Fig. 10

Upper limits of the minimum mass of possible planetary companions to the HD 108874 system, evaluated with a 99.9% confidence level. Stable regions for an Earth mass planet are defined by grey areas; the vertical dashed line indicates the orbital period of planet b.

Open with DEXTER

The outcome of the dynamical analysis shows that no companions are allowed between the two giant planets, also excluded from HARPS-N data, at least for planets down to a few Neptune masses. The system is then dynamically full up to 7 AU, except for some stability strips in the inner region (Fig. 9, left panel). As reported in Sect. 7, a further planetary companion is potentially allowed in the outer region, beyond 11 AU: this hypothesis cannot be verified with the present instrumentation and the current observational baseline, at least for companions smaller than a brown dwarf. On the other hand, no evidence of a linear trend in the RV measurements is found (Sect. 6 and Bryan et al. 2016).

8.2. Constraining the system architecture

Our analysis indicates that, if present, planets with minimum mass down to 10 M should have been identified in the range between 1 and 200 d of period. The evidence from RV surveys (e.g. Mayor et al. 2011) and Kepler space-based photometry (e.g. Batalha et al. 2013) points towards a high frequency of close-in, packed systems of super Earths and Neptunes, but only in systems without any detected outer giant companion. Actually, from recent calculations (Izidoro et al. 2015) the populations of close-in super-Earth systems and Jupiter-like planets should be anti-correlated. Even if systems with hot Neptunes or hot super Earths with outer giant planetary companions do exist (Santos et al. 2016), systems hosting multiple super-Earths inside the orbits of giant planets are still uncovered. This configuration could be explained in terms of the role that giant planets may play in the terrestrial planet formation and evolution (e.g. Levison & Agnor 2003). Matsumura et al. (2013) performed simulations to investigate the fate of low-mass planets co-existing in the same planetary system with massive companions. Following their analysis, the architecture of the HD 108874 system could be the result of a quiet dynamical evolution of the two giant planets, which leaves no modifications to their original eccentricities (e< 0.3 for both of them and consistent with a planet migration mechanism, as confirmed in a study by Rein 2012) and orbital radii, but powerful enough to remove possible terrestrial planets through the ejection triggered by the secular perturbation or the merging with the star for crossing orbits with the inner giant planet. Jiang et al. (2015) studied the relation between the minimum mass and the period ratio for adjacent giant planet pairs observed in multiplanet systems, obtaining a clear correlation which is also confirmed for the two planets hosted by HD 108874. Their period–mass ratio falls well inside the level of scattering obtained with the calibrator pairs. The period ratio of the analysed system (~4.4) is expected to be typical for pairs of giant planets around not evolved low-mass stars (i.e. M< 1.4 M and log g> 4.0) when it is compared with the other systems considered, for example, in Fig. 1 in Sato et al. (2016). In their analysis they consider also the total minimum mass of the planet system: for HD 108874 it is 2.34 ± 0.26MJup, which appears to be larger with respect to the typical mass distribution.

8.3. The habitable zone of HD 108874

We evaluated the limits of the habitable zone (HZ) of our target by exploiting the analytical relation by Selsis et al. (2007). Our results, reported in Table 5, are estimated for the classical early Mars and recent Venus criteria described in Kopparapu et al. (2013) for the theoretical inner (runaway greenhouse) and outer limits with 50% cloudiness (with H2O and CO2 clouds respectively) and the extreme theoretical limits, with a 100% cloud cover.

Table 5

Boundaries of the habitable zone of HD 108874.

Since HD 108874 b is well inside the HZ we can consider the habitability of potential moons around it (Williams et al. 1997). Earth-mass moons revolving around Jupiter-mass planets have been shown to be dynamically stable for the lifetime of the solar system in systems where the stellar mass is larger than 0.15M (Barnes & O’Brien 2002). We considered this possibility in the case of HD 108874 b: if the total mass of the planet-satellite system, Mp + Ms, satisfies the relation Ms ≪ (Mp + Ms) ≪ M, being M the stellar mass, then PsatP p, where Psat is the satellite’s orbital period and P p is the circumstellar period of the planet-satellite system (Heller 2012). This condition can be translated in an upper limit for the orbital radius of the satellite around the giant planet that we evaluated in a range of masses between Titan and ten Earth-masses.

thumbnail Fig. 11

Upper limit of the orbital separation between HD 108874 b and an hypothetical satellite, evaluated for different values of the mass of the latter (smaller bodies of the solar system are depicted as reference).

Open with DEXTER

The result is shown in Fig. 11, where the position of a varying-mass satellite is identified with the name of the corresponding solar-system satellite or small planet and the orbital distance from the planet is indicated in Jupiter radii. For masses up to one Earth mass the upper limit of the orbital radius is ~53.95RJ, about 0.025 AU. In the case that HD 108874 b transits its star the presence of the moon can be detected by measuring the variation in the transit time (TTV) of the planet due to gravitational effects, as shown by Sartoretti & Schneider (1999). According to their Eq. (24), if we consider satellites of the size of the Earth, the TTVs is approximately eleven minutes. The uncertainty of our ephemeris (a few days), obtained with RV data only is not suitable to detect such a signal, but in the case of transit observations, for example, with the forthcoming CHEOPS satellite, better constraints are expected.

9. Conclusions

We have presented the analysis of the intensive RV monitoring of the star hosting planets HD 108874, with the HARPS-N spectrograph at TNG in the framework of the GAPS Programme. A significant periodicity of 40.2 d has been found in the RV residuals of the two known planets fit but after a full analysis of the stellar activity we conclude that it must be addressed to the rotation period of the star and not the presence of an additional low-mass planet to the system. This is an example of how the activity contribution to the RV is able to mimic a Keplerian modulation. We performed a refinement of the orbital parameters of the two giant planets, and the dynamical analysis shows that the system is in a stable configuration over 1 Gyr. Stable low-mass bodies are only allowed at small separation or very far from the star, even if HARPS-N data tend to exclude their presence.


Acknowledgments

This work was supported by INAF through the “Progetti Premiali” funding scheme of the Italian Ministry of Education, University, and Research. The authors acknowledge Dr. L. Bedin (INAF-OAPD) and Dr. G. Lodato (Università degli Studi di Milano) for their comments and suggestions. We thank the anonymous referee for her/his useful advices and comments.

References

Appendix A: Additional table

Table A.1

Time series of the HARPS-N spectra for HD 108874.

All Tables

Table 1

Stellar parameters of HD 108874.

Table 2

Summary of the datasets.

Table 3

Summary of the values of Pearson (CP) and Spearman (ρ) correlation coefficients and corresponding p-values between RV residuals and all the activity indices investigated in this work.

Table 4

Orbital parameters of HD 108874 b and c as derived from the present analysis.

Table 5

Boundaries of the habitable zone of HD 108874.

Table A.1

Time series of the HARPS-N spectra for HD 108874.

All Figures

thumbnail Fig. 1

Panel a): RV dataset for HD 108874. Green dots: HIRES at Keck; orange dots: HRS at HET; blue dots: HARPS-N at TNG; black solid line: global fit of the three sets. Panel b): RV time series for HARPS-N data only. Panel c): residuals after the subtraction of the two planets fit.

Open with DEXTER
In the text
thumbnail Fig. 2

Periodogram of RV residuals of HARPS-N dataset: both frequency and period domains are shown. Horizontal dashed lines represent the confidence level of the peaks. The inset shows the window function due to the temporal sampling of the data.

Open with DEXTER
In the text
thumbnail Fig. 3

Periodogram of HARPS-N RV residuals subdivided into each observing season. The red solid line marks the location of the 40.2 d periodicity, while the red dotted lines represent the 1st, 2nd and the 3rd harmonics of the main peak.

Open with DEXTER
In the text
thumbnail Fig. 4

Upper panels: periodograms of the (with S/N> 20) and of the Hα indices for HARPS-N data. The red solid line indicates the 40.2 d period, the red dotted lines represent its 1st, 2nd and the 3rd harmonics; the blue dotted line indicates the period of planet b. Lower panels: correlation with the RV residuals.

Open with DEXTER
In the text
thumbnail Fig. 5

Periodograms of the asymmetry indicators as derived by the HARPS-N pipeline (BVS) and Lanza et al. (in prep; ΔV, Vasy(mod)). Red solid and dotted lines represent the assumed Prot and its 1st, 2nd and the 3rd harmonics, respectively, while the blue dotted line indicates the period of planet b.

Open with DEXTER
In the text
thumbnail Fig. 6

Left upper panel: original time series of the CCF FWHM, showing the trend due to the defocusing. Left lower panel: residuals of the CCF FWHM after the subtraction of the fit. Right panel: periodogram of the residual CCF FWHM. The red solid line represents the Prot, the dotted lines represent its 1st, 2nd and the 3rd harmonics.

Open with DEXTER
In the text
thumbnail Fig. 7

a) Tentative fit of the stellar activity cycle for HARPS-N and HIRES data of log . Data points (grey dots) are plotted together with the mean values of each observing season (red dots). b) Fitted RV jitters for the separate observing seasons of HIRES, HRS and HARPS-N datasets, when activity model is not included in the data. c) Semi-amplitude of the fitted RV sinusoid. d) Phase of the fitted RV sinusoid.

Open with DEXTER
In the text
thumbnail Fig. 8

Stability of the HD 108874 system. The values of the semi-major axes of the two planets are sampled around the nominal system (black square) according to the dispersion of the main secular frequency. The black box marks the 1σ uncertainties given in Table 4. The three yellow strips indicate the instability zones related to the 5:1, 9:2 and 4:1 MMR, respectively. A logarithm scale is used for the colour coding, where small values of the dispersion mean stability, while large values imply fast changes of the secular frequency and then chaotic evolution.

Open with DEXTER
In the text
thumbnail Fig. 9

Stability of a 1 M planet in the inner (left) and the outer regions (right) of the HD 108874 system as a function of the initial eccentricity and semi-major axis of the putative bodies. Colour scale ranges from black, representing stable regions, to yellow, which indicates higher degree of instability for planets with given separation and eccentricity. The empty (white) spaces indicate systems unstable over short timescales.

Open with DEXTER
In the text
thumbnail Fig. 10

Upper limits of the minimum mass of possible planetary companions to the HD 108874 system, evaluated with a 99.9% confidence level. Stable regions for an Earth mass planet are defined by grey areas; the vertical dashed line indicates the orbital period of planet b.

Open with DEXTER
In the text
thumbnail Fig. 11

Upper limit of the orbital separation between HD 108874 b and an hypothetical satellite, evaluated for different values of the mass of the latter (smaller bodies of the solar system are depicted as reference).

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.