Chasing extreme planetary architectures I. HD

.


Introduction
Even though binaries were initially banished from exoplanet hunting campaigns, more than 200 planets 1 have been so far discovered in multiple systems orbiting one component or both components, referring to the S-type and P-type configuration respectively (see classification by Dvorak 1982).These systems have been initially discovered by accident like Gliese 86 Ab (Queloz et al. 2000) and γ Cep Ab (Hatzes et al. 2003), but more recently, dedicated surveys in radial velocity (RV), transit and imaging have been targeting them (e.g.Fontanive & Bardalez Gagliuffi 2021).Given that the formation of stars in multiple systems is a frequent byproduct of stellar formation, a natural question is to understand how the presence of a stellar companion can affect the planetary formation process.Planetary formation theories have been restricted for years to the case of a single star environment to understand the formation of our own Solar System (Mordasini et al. 2008), and more recently investigated for binary stars (e.g.Silsbee & Rafikov 2021).For the most frequent dynamical configuration observed for planet(s) in binaries, the S-type circumprimary one with a planet orbiting one component of the binary (generally the most massive one), models have predicted that the presence of a very close binary companion can hinder planetary formation around a primary star in several ways.The most problematic issue is probably not the disk truncation (see Jang-Condell 2015), but rather the dynamical excitation that could prevent mutual planetesimal accretion for objects in the 100 m-to-100 km size range (Thébault et al. 2006;Thebault 2011;Paardekooper et al. 2008;Fragner et al. 2011), hence obstructing the formation of a planet by core accretion, or ejecting the planet in unstable systems (e.g.Thebault & Haghighipour 2015).Dedicated observing campaigns suggest a transition at typical separation within 50-300 au indicating that short or intermediate-period binaries have statistically less chance to host planets or brown dwarf companions, compared to wide binaries that on the contrary have no influence on the architectures of planetary systems (Ziegler et al. 2021).This effect seems to be corroborated by the study of young stars for which short-separation (À 100 au) binaries have a lower probability of hosting circumstellar dust in the innermost few A&A proofs: manuscript no.chauvin-hd196885 Fig. 1.Overview of the binary star systems.The sizes of the circles representing the primary and secondary star indicate their masses according to the scale reported on the Bottom-Right.The planet masses are not to scale.The green zone marks the region of orbital stability for a coplanar configuration in the context of an S-type configuration.
au around each star, therefore with a depleted reservoir of solids for the formation of planets by core accretion (Duchêne 2010).
In the last decade, more intensive and dedicated surveys have been dedicated to search for S-type and P-type 2 exoplanetary systems using (and often combining) RV, transit, high-contrast and direct imaging observations.RV surveys, able to tackle multiple stellar signals at the precision required for detecting planets, were specifically carried out (Konacki 2005;Eggenberger & Udry 2007;Toyota et al. 2009;Desidera et al. 2010).Zucker & Mazeh (2002) initially showed that the planet properties in binaries were different to planets around single stars.Eggenberger et al. (2004) confirmed that the most massive short period planets were found in binary systems and that planets with short (P À 40 days) periods in binaries were likely to have low eccentricities.Tamuz et al. (2008) finally found that extremely high eccentric planets were all in binaries.
In parallel to Doppler surveys, high-contrast imaging campaigns have been able to access intermediate-separation (20 100 au) binaries where the companion influence is the most expected.Several deep imaging studies have revealed a number of additional close companions to known exoplanetary hosts (Patience et al. 2002;Luhman & Jayawardhana 2002;Chauvin et al. 2006;Eggenberger & Udry 2007;Law et al. 2014;Mugrauer & Ginski 2015;Baranec et al. 2016;Ginski et al. 2021).Eggenberger et al. (2007) carried out a dedicated survey to explore the multiplicity in two samples of stars with and without planets detected in RV in order to test the impact of duplicity on the giant planet occurrence.They found a lower corrected binary fraction for the planet-host sub-sample, particularly for physical separations shorter than 100 au (Eggenberger et al. 2008).
2 P-type orbits occur when the planet orbits both binary components, and generally referred as circumbinary planet For longer period (Á 100 au) binaries, catalogs compiled by Raghavan et al. (2006), Desidera & Barbieri (2007), and Bonavita & Desidera (2007) further updated by Bonavita & Desidera (2020) studied the effect of duplicity, and did not find any significant discrepancies in the planet frequency between binary and single stellar systems.More recently, the use of the Gaia offered the opportunity of a systematic search for (sub)stellar companions to exoplanet host (Kervella et al. 2019;Mugrauer 2019;Michel & Mugrauer 2021).Fontanive & Bardalez Gagliuffi (2021) recently derived an overall raw multiplicity rate of 23.2˘1.6% for hosts to exoplanets and brown dwarf companions, considering a volume-limited sample up to 200 pc crossmatched with the Gaia Data Release 2 catalog (Gaia Collaboration et al. 2016aCollaboration et al. , 2018)).
A complementary approach has been to target individual systems, and particularly the extreme ones that provide theoreticians with an ideal test bed for the theories of planetary formation, evolution and dynamical stability in a very hostile environment.Today, there are about 16 extreme exoplanetary systems (binaries with small semi-major axes of less than about 20 au hosting a planet; see Fig. 1).One emblematic case is γ Cep for which Benedict et al. (2018) managed to derive the mass and the orbital properties of the circumprimary planet using accurate astrometric measurements with HST.Orbiting at 1.94 au in a highly mutually inclined orbit (Φ " 70 ˘13 deg), the formation and survival of this 9.4 `0.7 ´1.1 M Jup exoplanet remains unresolved (e.g.Bazsó et al. 2017).Another interesting system is HD 196885 that Chauvin et al. (2006) discovered in the course of a VLT/CFHT deep imaging survey of exoplanet hosts.In this paper, we present the results of the long-term monitoring of this system combining RV, direct imaging, speckle and dual-field interferometry.The physical properties of the HD 196885 exoplanetary system are reported in Section 2. The observations and data reduction are presented in Section 3, together with the results in Section 4. In Section 5, we finally discuss their implication in terms of system architecture, formation and stability.

The HD 196885 exoplanetary system
HD 196885 is an F8V (V " 6.385 mag, H " 5.19 mag) star located at 34.00 ˘0.04 pc (Gaia EDR3, Gaia Collaboration et al. 2021). Based on CORALIE spectra, Sousa et al. (2006) derived a spectroscopic temperature, surface gravity, and metallicity of T eff " 6340 ˘39 K, log(g)" 4.46 ˘0.02 and [Fe/H]" 0.29 ˘0.05, respectively.The projected stellar rotational velocity v.sin(i) was estimated to 7.3 ˘1.5 km/s from ELODIE spectra (Correia et al. 2008).These results were supported by independent Lick observations reported by Fischer et al. (2009).The chromospheric activity level of logpR 1 HK q " ´5.01 (Wright et al. 2004) is relatively low when placed into context with the activity distribution of the solar neighborhood (Jenkins et al. 2011).Bolometric luminosity correction and evolutionary model predictions lead to an estimate of the luminosity and the mass of 2.4 L d and 1.3 M d , respectively.The corresponding stellar age derived from evolutionary tracks and from the activity level varies between 1.5 to 3.5 Gyr (see Correia et al. 2008;Fischer et al. 2009).Finally, we looked at wider physical companions within 3 arcmin from HD 196895 in Gaia EDR3 and found none.The companion BUP 1908 (192 arcseconds separation, position angle 8 ˝, V " 10.4 mag, TYC 1096-637-1) listed in the double star catalog, CCDM, is optical according to Gaia.As expected, HD 196885 shows a huge signature of binarity (with an S/N of 221) in the Hipparcos-Gaia proper motion anomaly cat-Table 1. Physical properties of the extreme planetary system HD 196885, the stellar host (Aa), the binary companion (B) and the circumprimary planet (Ab) based on the latest studies from Correia et al. (2008);Fischer et al. (2009); Chauvin et al. (2011);Gaia Collaboration et al. (2021).We consider below that the system is composed of the inner binary subsystem AaAb (composed by the pri-mary star Aa and the exoplanet Ab), and the outer binary subsystem AB (composed by the center of mass of the inner orbit A and the companion star B) alog (Kervella et al. 2022), but less constraining that the set of observations presented below.
In 2004, an RV variation indicating a planetary candidate was measured at Lick observatory and adjusted with a preliminary orbital period of P " 0.95 yr.The result was temporarily reported on the California Planet Search Exoplanet Web site, but finally withdrawn owing a significant residual drift in the orbital solution.Nevertheless, this star was included in a VLT/CFHT deep-imaging survey of stars hosting planets detected by RV observations (Chauvin et al. 2006).This led them to the discovery of a close (" 20 au) stellar binary companion to HD 196885, motivating dedicated follow-up observations with Adaptive-Optics (AO) imaging.The recent image of the system obtained by SPHERE at VLT is shown in Fig. 2. Following the CFHT discovery, astrometric monitoring and long-slit spectroscopic observations at VLT with NaCo enabled to resolve the orbital motion of the stellar companion and to confirm its low stellar-mass nature as an M1˘1V dwarf companion located at only 0.7 as (23 au in projected physical separation) and likely to be responsible for the trend seen in the Lick RV residuals (Chauvin et al. 2007).Using a double-Keplerian model for the binary star and the planet to adjust their ELODIE, CORALIE, and CORAVEL observations spread over 14 yr, Correia et al.

Radial velocities
In an attempt to revisit the physical and orbital properties of the HD 196885 system, we considered all data from past radial velocity surveys already published in Chauvin et al. (2011) for the primary host HD 196885 A: CORAVEL (  Correia et al. 2008;Fischer et al. 2009;Chauvin et al. 2011, and this work).SOPHIE data Table 2. Observation log of the new SPHERE, SOAR, and GRAVITY observations.The exposure time is given by Nexp, the number of exposures, multiplied by NDIT and DIT, which are the number and duration of integrations per exposure (in the case of HRCAM in ms).ω and τ 0 are the seeing and the coherence time respectively.For the modes, DBI refers to Dual-Band Imaging, SI to Speckle Interferometry and MED for medium spectral resolution interferometry.were acquired in High-resolution mode as part of the program dedicated to the follow-up of ELODIE candidates.Data from 2007 to June 2011 where reduced thanks to the SOPHIE pipeline (Bouchy et al. 2009).In June 2011, an update of the instrument was done.It consists mainly in the addition of octogonal fibers to light guiding from the telescope to the spectrograph (Bouchy et al. 2013).Data after this date were then reduced as in Heidari et al. (2022).

SPHERE adaptive-optics imaging
Following previous astrometric observations at VLT obtained with the NaCo instrument, the HD 196885 system has been monitored during the SpHere INfrared survey for Exoplanets (SHINE, Chauvin et al. 2017;Desidera et al. 2021;Langlois et al. 2021) at 2 different epochs on 5 June 2015 and 30 September 2017 (see Table 2) using the VLT/SPHERE high-contrast instrument (Beuzit et al. 2019).The observations were obtained with the modes IRDIFS-EXT that combine simultaneously the IRDIS (Dohlen et al. 2008) and IFS instruments (Claudi et al. 2008), although IRDIS observations were mainly used in the current study.The IRDIFS-EXT mode combines IRDIS in dualband imaging (DBI Vigan et al. 2010) mode with the K1K2 filter doublet λ K1 " 2.103 ˘0.102 µm, λ K2 " 2.255 ˘0.109 µm, and IFS in the YJH (0.97-1.66 µm) setting.Each observing sequence was performed with the field-tracking mode.The detail about the observations are reported in Table 2.Only the IRDIS observations at K1-band (with a lower background noise) were used to extract the relative astrometry of HD 196885 B relative to A.

SOAR speckle interferometry
Complementary astrometric measurements of HD 196885 have been obtained with the High-Resolution Camera (HRCAM) Tokovinin et al. (2021) at the SOAR 4.1-m telescope Tokovinin & Cantarutti (2008) 3 .While regular observations are made in the Stroemgren y, Cousins I, or Hα filters, all observations of HD 196885 were done with the I filter with a resolution of 36 mas.The data have been processed by a custom IDL pipeline developed by A. Tokovinin to extract the binary parameters (separation, position angle, and magnitude difference) by fitting the secondary peaks of the auto-correlation functions or the fringes from the power spectrum of the observations obtained by Fourier transform as explained in Tokovinin et al. (2010).HRCAM routinely delivers a precision of 1-3 mas in angular separation for objects brighter than V " 12. On every observing night "calibration binaries" are included, these are binaries with very well known orbits (grade 5 on the USNO orbit catalogue), from which the measurements are calibrated, leading to systematic errors less than 0.1 deg in position angle, and better than 0.2% in scale, i.e., smaller than the internal precision.The final precision of these measurements depends on a number of factors, but in this paper we adopt a uniform uncertainty of 3 mas, which is representative of all the measurements (Tokovinin 2018).The Observing log and the astrometric results are reported in Tables 2  and 3.   2.
These observations were obtained with the astrometric dualfield mode in off-axis: in the first exposure, the fringe tracker single-mode fiber is centered on the primary star HD 196885 A while the scientific single-mode fiber is centered on the stellar companion HD 196885 B. In a second exposure, the two fibers are swapped, to calibrate the metrology offset.This observation mode is different from the dual-field mode in on-axis, where the planet is bright enough to be directly observed (Gravity Collaboration et al. 2019).Here, the GRAVITY interferometer only observed the A and B components.The observations were preprocessed with the Public Release 1.5.0 (1 July 2021) of the ESO Gravity pipeline (Lapeyrere et al. 2014).The reduced data are then processed to extract astrometry using the GRAVITY python tools offered to the community 4 .The astrometric values are presented in Table 3.

Resolving the astrometric wobble due to the giant planet HD 196885 Ab
The main objective of our GRAVITY program was to measure the astrometric wobble in the separation of a binary system related to presence of the giant planet, and in such a circumstance be able to derive the orbital properties of the binary companion and the giant planets, but also to confirm the planetary nature of HD 196885 Ab.To date, such observation have only been performed in the context of the γ Cep system by Benedict et al. (2018).Combining all the available measurements obtained from radial velocity, adaptive-optics imaging, speckle interferometry, and dual-field interferometry, we have revisited the system properties using a Bayes-based orbital elements estimation inspired in a previous study by our team to triple hierarchical stellar systems (Villegas et al. 2021).

Code description and set up
The orbital parameters estimation code has been developed in a Bayesian Markov-Chain Monte Carlo-based framework, in the context of the so-called hierarchical approximation to combine astrometric with spectroscopic measurements for multiple systems.The current implementation of this code 5 is an adaptation of the Bayesian inference algorithm for binary stars developed by Videla et al. (2022), which is based on the No-U-Turn Sampler (NUTS) algorithm (Hoffman et al. 2014).The NUTS algorithm is an MCMC method that avoids the random-walk behavior and the sensitivity to correlated parameters of other MCMC methods commonly used in astronomy (such as the Metropolis-Hastings within Gibbs sampler (Ford 2005), the Parallel Tempering sampler (Gregory 2005(Gregory , 2011)), the Affine Invariant MCMC Ensemble sampler (Hou et al. 2012), and the Differential Evolution Markov Chain sampler (Nelson et al. 2013)) by incorporating first-order gradient information of the parameters space to guide the sampling steps with an adaptive criterion for determining their lengths.The NUTS algorithm qualities allow for performing a quick and efficient estimation of the target distribution, specially in the case of the high-dimensionality problem in hand.We estimate the posterior distribution of the joint orbital parameter space of a triple hierarchical system that we adapted to explore the planetary architecture in S-type binaries.The parameter space consists of the orbital parameters of the inner binary subsystem AaAb (composed by the primary star Aa and the exoplanet Ab), and the orbital parameters of the outer binary subsystem AB (composed by the center of mass of the inner orbit A and the companion star B), assuming of course the same value of the parallax π for both subsystems.Therefore, the parameter space is defined by the vector ϑ " tP AaAb , t PAaAb , e AaAb , a AaAb , ω AaAb , Ω AaAb , i AaAb , q AaAb u Y 4 https://version-lesia.obspm.fr:/repos/DRS_gravity/gravi_tools3 5 The original software is available on GitHub: BinaryStars codebase: https://github.com/mvidela31/BinaryStarsunder a 3-Clause BSD License.Its evolution used in this study will soon be released.
tP AB , t PAB , e AB , ω AB , Ω AB , i AB , q AB , π, V 0 u6 .We use uniform priors on all the orbital parameters according to its physical valid range (e.g., q AaAb " U r0,1s ), except on the system parallax, where we use a Gaussian prior π " Np29.408, 0.027 2 q mas (according to the value and standard deviation reported by Gaia DR3, Gaia Collaboration et al. 2016b).Based on the spectral types reported for the primary and secondary stellar components reported in Table 1, we also used truncated Gaussian priors on the stellar masses M Aa " Np1.245, 0.067 2 q M d and M B " Np0.485, 0.059 2 q M d within the interval rµ ´2σ, µ `2σs (since it accounts " 95% of the Gaussian probability density function), using the recent calibration provided by Abushattal et al. (2020) 7 .We run 10000 iterations on each of four independent Markov chains, discarding the first half for the warm-up phase.Each chain was initialized from an starting point determined through the quasi-Newton method L-BFGS (Liu & Nocedal 1989).The code was implemented using the probabilistic programming language Stan (Carpenter et al. 2017).

4.2.
A massive Jupiter on a mutually inclined orbit relative to the stellar companion?
The results of the orbital fitting exploration led to the identification of nine for families of solutions (considering a truncated Gaussian prior for the mass of HD 196885 A) that are reported in Table A. Their log-posterior density (up to a constant) distributions is also shown in Fig A .1.In Table 4 and in Fig. 3 (in blue), we report the most likely solutions of the family Sol. 2 which are compared with the observed relative astrometric and radial velocity measurements of the system.The posterior solutions of each orbital and physical parameters of the exoplanetary companion HD 196885 Ab, and the ones of the stellar companion B are shown in Fig. 4. The long-term monitoring over almost 40 yr combined with 15 yr of direct imaging set a relatively tight constraints on the orbital properties of the stellar companion.The solutions are very similar to the ones already derived by Chauvin et al. (2011) indicating a moderately eccentric (e AB " 0.417) orbit, a 70 yr period, an inclination of 120.43 deg, and a time to periastron in 1982.87 yr.However, the addition of SOAR, SPHERE, and GRAVITY astrometric measurements from the last decade enables us to nail down previous uncertainties.
For the planet HD 196885 Ab, the gain of almost two orders of magnitude in astrometric accuracy obtained by GRAVITY at VLTI (30 ´50 µas precision) compared to what was regularly achieved by NaCo and SPHERE at VLT (1 ´2 mas precision) completely opens a new window to detect and characterize the periodic wobble caused by the S-type planet on the angular separation between both stellar components A and B. The four GRAVITY measurements obtained in August 2019, May 2021, August 2021, and October 2021 covering less than 60% of the orbital period of HD 196885 Ab, allow to unambiguously reject all families of solutions associated with stellar and brown dwarf masses for Ab, and to pinpoint one specific family of retrograde solutions of Jupiter-like planets with low masses (M Ab " 3.39 M Jup ), with eccentricity (e AaAb " 0.44), and in- We note that other islands of solutions remain possible in the planetary mass regime but for different orbital configurations: retrograde and prograde solutions with higher mutual inclination (Φ " 60 ´80 deg or Φ " 110 ´120 deg, respectively) for low masses for Ab (M Ab " 2 ´4 M Jup ) and moderate eccentricities, as well as one family of prograde solution for higher masses (M Ab " 12´14 M Jup ), moderate eccentricity and mutual inclination (Φ " 130 deg).

Origin and fate
Our orbital determination clearly suggests that the 3-body system HD 196885 is presumably non-coplanar.Figure 4 shows the current favored family of retrograde solutions for which mutual inclination Φ lies between typical values of 27.6 deg (first quartile) and 37.8 deg (third quartile) (note that the MAP value is equal to 24.4 deg).This non negligible inclination makes it a good candidate for a von Zipele Kozai Lidov-Lidov mechanism (or resonance) Kozai (1962) mainly affecting the inner orbit.This dynamical behavior is common in non-planar multiple systems (Beust & Dutrey 2006;Krymolowski & Mazeh 1999;Ford et al. 2000).It is characterized by coupled large amplitude eccentricity and inclination oscillations mainly affecting the inner orbit.Regarding the relatively large individual masses involved here, this could lead the system to instability, which justifies a dedicated stability study.Table 4. Summary of the orbital parameters derived for HD 196885 Ab and HD 196885 B: the semi-major axis (a), the eccentricity (e), the inclination (i), the longitude of ascending node (Ω), the argument of periastron (ω), the time of periastron passage (t P ).We report the maximum a posteriori estimates (the most probable solution), with the first and third quartiles (as an uncertainty measure) as subscripts and superscripts, respectively.Note that the MAP values are not necessarily equal to the second quartile (the median).To do it, we used the HJS integrator (Beust 2003), a symplectic package specifically designed to study hierarchical stellar systems.The same package was successfully indeed used  4 over 10 8 yr.

Parameters
Left-Right: Semi-major axes, eccentricities, and mutual inclination.In the Left and Central plots, the orange curves stand for the inner orbit and the blue one for the outer orbit, respectively.
by Beust & Dutrey (2006) to investigate the multiple system GG Tauri.Here we performed an integration over 10 8 yr starting from the orbital solution from Table 4.The time step was fixed to 0.05 yr, which is a very conservative assumption, as this represents 1/70 of the smaller orbital period (AaAb orbit), while tests by Beust (2003) show that with only 1/20, the total energy is preserved to 10 ´6 relative fraction.But as we fear that the system might be close to instability given the high individual masses of the various bodies involved, this appeared to be a reasonable assumption.Figure 5 shows the result of this computation, taking the orbital fitting solution from Table 4 as starting point.The evolution is shown over 10 8 yr.This configuration appears stable over the whole integration, as can be seen from the stability of the semi-major axes.Nonetheless the eccentricity and mutual inclination plots suggest von Zipele Kozai Lidov cycles that affect the inner orbit, while the outer orbit is only little affected.This von Zipele Kozai Lidov resonance remains moderate here, due to the rather small mutual inclination.Hence the stability of the triple system is not affected.The formation mechanism of such a system is unclear.We currently do not know if planets in such tight binaries follow similar or different formation processes than the ones orbiting single stars.For instance, Cadman et al. (2022) have investigated the role of binary companions triggering fragmentation in selfgravitating discs.They find that binarity triggers disk instability for intermediate-close (100-400 au) binaries while the presence of a binary companion within 100 au has a strong negative role against the occurrence of instability.For core accretion, in the specific case of HD 196885 b, Thebault (2011) conclude that the perturbed dynamical environment is strongly hostile to planetesimal accretion in the region where the planet is located today (note that this study considered a coplanar case).Finally, Marzari & Barbieri (2007) propose that the present dynamical configuration may not reflect the formation process since the binary orbit may have changed in the past once the planet formed.Disrupted triple system or single stellar encounter could be alternative processes.Consequently, this formation and dynamical picture must be further re-investigated in the light of upcoming astrometric observations and more detailed simulations.First, the solution depicted in Table 4 may not be unique, and a dedicated follow-up with GRAVITY should enable to rapidly secure the most probable family of solutions.Due to the complexity of the parameter space, other configurations might fit the data as well.Second, other effects could need to be taken into account such as a combination of von Zipele Kozai Lidov cycles and tides (e.g.Fabrycky & Tremaine 2007;Correia et al. 2011;Beust et al. 2012).

Conclusion
In this study, we have revisited the global architecture of the extreme planetary system HD 196885.This system is composed of a tight F8V-M1V binary hosting a S-type circumprimary planet around the A component.Based on the combination of new measurements in radial velocity, speckle interferometry, highcontrast imaging and high-precision astrometry with interferometry, we could narrow down the orbital properties of both HD 196885 B and HD 196885 Ab.The binary companion orbits on an inclined and retrograde plan (i AB " 120.43 deg) with a semi-major axis of 19.78 au, and an eccentricity of 0.415 tightly constrained by almost 40 yrs and 20 yrs of radial and astrometric monitoring, respectively.For the radial velocity planet HD 196885 Ab, we unambiguously confirm its mass (M Ab " 3.39 M Jup ) and planetary nature, rejecting all families of pole-on solutions corresponding to stellar and brown dwarf masses.The most favored island of orbital solutions found is associated with moderate eccentricity (e AaAb " 0.46), and moderate-inclination close to 143.04 deg.This results points toward a moderate mutual inclination (Φ " 24.36 deg) between the orbital plans of B and Ab.This configuration seems to be dynamically stable over time, although the eccentricity and mutual inclination variations clearly reveal the presence of moderate von Zipele Kozai Lidov cycles that may affect the inner planet.Further observations in the coming years will allow to reject remaining islands of orbital solutions for Ab.Further dynamical simulation should also enable to incorporate additional effects as tides to assess the inner orbit dynamical evolution and stability.

Fig. 2 .
Fig. 2. SPHERE/IRDIS observations of HD 196885 AB obtained in September, 2017.The early-M1V dwarf companion HD 196885 B is resolved at K1-band at a separation of 604 mas and a position angle of 45.7 deg.

Fig. 3 .
Fig. 3. Top-Left: Astrometric position of HD 196885 B relative to A combining SOAR, NaCo, SPHERE and GRAVITY observations from August 2005 to August, 2021.Specific zoomed plots on the NaCo astrometric measurements obtained between August, 2005 and August, 2009 (Middle-Center), the SPHERE and SOAR ones from June, 2015 to May, 2018 (Top-Right), and the GRAVITY ones from May, 2019 to August, 2021 (Middle-Right) are shown.Bottom-Left: Radial velocity measurements of HD 196885 A from CORAVEL, ELODIE, CORALIE, Lick and SOPHIE covering 40 yr of monitoring.Top-Right: Zoom on the radial velocity signal induced by HD 196885 Ab.The results of the orbital fitting exploration are reported with the thick blue curves.

3. 4 .
VLTI/GRAVITY dual-field interferometry HD 196885 was also monitored in the context of monitoring program with GRAVITY (Gravity Collaboration et al. 2017) and the UTs in MED resolution at K-band.Over two years, we have obtained a total of four observations on 15 August 2019, 29 May 2021, 25 August 2021, and 25 September 2021.The weather conditions ranged from good to very good, and the observing log of the observations is also reported in Table

Fig. 4 .
Fig. 4. Marginal posterior distributions of the physical and orbital parameters of HD 196885 Ab including the period, periastron passage, eccentricity, the semi-major axis... (two first rows).On the two last rows are reported the marginal posterior distributions for HD 196885 B. The maximum a posteriori (the most probable sample of the posterior distribution) are reported in red.

Fig. 5 .
Fig. 5. Dynamical evolution of HD 196885 system as computed with the HJS integrator starting from the orbital solution of Table4 over 10 8 yr.Left-Right: Semi-major axes, eccentricities, and mutual inclination.In the Left and Central plots, the orange curves stand for the inner orbit and the blue one for the outer orbit, respectively.

Table 3 .
Chauvin et al. 2011)of HD 196885 B relative to A extracted from NaCo (seeChauvin et al. 2011), SOAR, IRDIS, and GRAVITY (refered as GRAV. in the Table)observations.For GRAVITY, the Pearson's coefficient ρ quantify the correlation between the ∆α and ∆δ uncertainties.