A dynamic view of V Hydrae Monitoring of a spectroscopic-binary AGB star with an alkaline jet ⋆

Context. The well studied carbon star V Hydrae is known to exhibit a complex asymmetric environment made of a dense equatorial wind and high-velocity outflows, hinting at its transition from the AGB phase to the asymmetric planetary nebula phase. In addition, V Hydrae also exhibits a long secondary period of 17 years in its light curve, suggesting the presence of a binary companion that could shape the circumstellar environment. Aims. In this paper, we aim to confirm the binary nature of V Hydrae by deriving its orbital parameters and investigating the e ff ect of the orbital motion on the circumbinary environment. Methods. In a first step, we used a radial-velocity monitoring performed with the HERMES spectrograph to disentangle the pulsation signal of the AGB from its orbital motion and to obtain the spectroscopic orbit. We combined the spectroscopic results with astrometric information to get the complete set of orbital parameters, including the system inclination. Next, we reported the time variations of the sodium and potassium resonance doublets. Finally, following the methods used for post-AGB stars, we carried out spatio-kinematic modelling of a conical jet to reproduce the observed spectral-line modulation. Results. We found the orbital solution of V Hydrae for a period of 17 years. We correlated the companion passage across the line of sight with the obscuration event and the blue-shifted absorption of alkaline resonant lines. Those variations were modelled by a conical jet emitted from the companion, whose opening angle is wide and whose sky-projected orientation is found to be consistent with the axis of the large-scale bipolar outflow previously detected in the radio-emission lines of CO. Conclusions. We show that the periodic variation seen for V Hydrae is likely to be due to orbital motion. The presence of a conical jet o ff ers a coherent model to explain the various features of V Hydrae environment.


Introduction
A star with an initial mass ranging from 0.8 to 8 M ⊙ will eventually evolve through the asymptotic giant branch (AGB) stage.During this phase, the star undergoes a process of significant mass loss (10 −7 to 10 −4 M ⊙ /yr) through a slow (5 to 20 km s −1 ) wind that enriches the interstellar medium (Habing & Olofsson 2004).When leaving the AGB phase, its effective temperature increases and goes through a short-lived transition of post-AGB and pre-planetary nebula (lasting about a few thousand years) ⋆ Based on observations made with the Mercator Telescope, operated on the island of La Palma by the Flemish Community, at the Spanish Observatorio del Roque de los Muchachos of the Instituto de Astrofísica de Canarias.Based on observations obtained with the HERMES spectrograph, which is supported by the Research Foundation -Flanders (FWO), Belgium, the Research Council of KU Leuven, Belgium, the Fonds National de la Recherche Scientifique (F.R.S.-FNRS), Belgium, the Royal Observatory of Belgium, the Observatoire de Genève, Switzerland and the Thüringer Landessternwarte Tautenburg, Germany.⋆⋆ Research fellow, FNRS, Belgium.
phase.The majority of planetary nebulae (PN) show a variety of aspherical (bi-or multi-polar) structures (Balick & Frank 2002).Those bipolar structures are expected to arise from binary interaction with a close-by companion (Jones & Boffin 2017).In the preceding phase, the post-AGB phase, binarity has been widely detected and studied (Oomen et al. 2018 and references therein).Common building blocks of post-AGB binaries have been identified (e.g.Kluska et al. 2022), consisting of a central binary surrounded by a stable circumbinary disk of gas and dust.They often exhibit a bipolar outflow which is likely launched by an accretion disc around the companion (e.g.Bollen et al. 2022 and references therein).Thanks to the near-infrared (NIR) interferometry, these building blocks can also be spatially resolved (e.g.Kluska et al. 2019;Corporaal et al. 2023).While these post-AGB binaries represent about a third of all known Galactic post-AGB stars (Kluska et al. 2022), they are poorly connected to the AGB progenitor sample, as evidence for binaries among these progenitors is lacking.Due to the intrinsic properties of AGB stars, detecting binaries is not a straightforward task.Given their photosphere is bright and dynamic, the detection of orbital motion is possible only if it exceeds the amplitude of the intrinsic pulsation of the stellar envelope (Jorissen 2004).
The late-type carbon star V Hydrae (V Hya) is among the best candidates to study the impact of a companion on the circumstellar environment of AGB stars.Secular photometric monitoring indicates that the star exhibits two photometric periods: a long period of 6160 ± 400 d with an amplitude of 3 mag and a short period of 530 ± 30 d with an amplitude of 1.5 mag (Knapp et al. 1999).While the 530 d period is commonly associated with the Mira pulsation of the AGB star, the origin of the 17 y period is still under debate.Several scenarios have been suggested, many of which invoke the presence of a binary companion: eclipse by a dusty cloud surrounding a companion (Knapp et al. 1999), extinction by the plasma ejection of the companion (Sahai et al. 2016), a triple system with an inner eccentric orbit of 8 y (Salas et al. 2019), along with possible episodic dust formation (Lloyd Evans 1997).In addition, the presence of a binary companion has been inferred on the basis of the near-UV (NUV) excesses (Sahai et al. 2008), detections of blue emission lines in the spectrum (Lloyd Evans 1991, Lloyd Evans 2009), and the exceptional structure of the large-scale environment of the star, as we describe below.
Several radio observations in the CO bands highlight the presence of an environment shaped by different velocity regimes: a slow (15 km s −1 ) dense equatorial wind, a medium-velocity (45 km s −1 ) spherical wind, and a fast (150-200 km s −1 ) bipolar outflow (Knapp et al. 1997).Kahane et al. (1996) mapped the fast wind component in the CO J=2-1 transition as a bipolar flow whose blue-shifted lobe is pointing toward the east with an inclination axis of about 30°.The presence of an eastward highvelocity knotty jet was also detected in the [S ii] emission by HST/STIS (Sahai et al. 2016) and its possible radio-equivalent by recent high-resolution ALMA observations together with the discovery of several concentric rings inside the dense equatorial structure (Sahai et al. 2022).All the reported observations confirm the large-scale complexity of the system, although no model (so far) has been able to assemble all the pieces of the puzzle and link the large-scale bipolar structure with the expected binary motion.
In this paper, we bring on new observations from the highresolution HERMES spectrograph (Raskin et al. 2011).The long-term spectral monitoring of V Hya provides unprecedented information on the system dynamics and new constraints on its binary motion.We report radial-velocity and spectral-line variations in phase with the long-term photometric period and interpret them in the framework of a binary system with a period of 17 y.The results are twofold: they provide the first spectroscopic evidence of a binary-AGB system and highlight the complex phenomena of the short-lived transition phase of V Hya, making it a case study for further investigation.This paper is structured as follows.The observations and data reduction are presented in Sect. 2. The determination of the orbital parameters is described in Sect.3. In Sect.4, we compare the spectroscopic signal with the visual photometry.In Sect.5, we analyze the temporal variation of selected spectral lines.The spatio-kinematic modeling of the Na i resonance doublet is introduced in Sect.6.In Sect.7, we discuss the possible origin of the long obscuration event and put our new spectroscopic results in the context of previous observations.In Sect.8, we outline our main findings.

Spectroscopic data
V Hya has been monitored with the HERMES spectrograph mounted on the 1.2m Mercator telescope at La Palma (Raskin et al. 2011).A total of 130 spectra were taken on different dates between 2011 and 2023, with an average exposure time of 900 seconds.The high-resolution mode was used, providing a resolution of R=86,000 over a wavelength range from 3800 to 9000 Å.The dedicated python-based pipeline was used for the data reduction.The radial velocity (RV) was obtained by cross-correlating the object spectrum with a carbon-star template and fitting a Gaussian function to the mean line profile, as represented by the resulting cross-correlation function (CCF).The mask used to compute the CCF is built from a HERMES spectrum of the R-type carbon star BD +02 • 4338, Doppler-corrected to zero velocity from a correlation with the master Arcturus template.All lines with a depth getting below the relative level of 0.8 with respect to the local continuum were kept in the mask.This mask, covering the wavelength range 476.9 -654.7 nm, was selected because it provides the cleanest CCF for V Hya (deepest and least asymmetric).Thirty spectra were excluded from the orbital analysis because they bear the signature of shocks.Due to the presence of shocks linked to the stellar pulsation, the CCF profile at those epochs close to maximum light was distorted (asymmetric broadening or double-peaked; Alvarez et al. 2001).Table A.1 summarises the log of the observations, including the observing dates of each spectra, the computed RVs and their formal uncertainties (due to the quality of the CCF fit).An uncertainty of 0.07 km s −1 expressing the long-term stability of the spectrograph is added quadratically to the formal uncertainty of the CCF fit.This long-term accuracy is estimated from the average RV standard dispersion of a sample of RV standard stars monitored along with the science targets (see Fig. 10 in Raskin et al. 2011).

Photometric data
For more than a century, V Hya has been monitored by variablestar observers.We use the archival visual photometric data from the American Association of Variable Stars Observers, AAVSO (Kloppenborg, B. K. 2023).The overall dataset (9000 measurements) covers a time interval from 1913 to 2023, allowing us to monitor six 17-yr cycles.

Astrometric data
Astrometric measurements from both Hipparcos (Perryman et al. 1997) and Gaia (Gaia Collaboration et al. 2022) ESA missions have been used to constrain additional parameters that are not accessible from spectroscopic data alone, such as the orbital inclination.The difference between the instantaneous proper motion computed by the two missions and the average proper motion between the two reference epochs provides an estimate of the acceleration in the sky plane.That acceleration being non-null for V Hya (Brandt 2021) is a signature of a long-trend motion and gives additional constraints on the orbital solution.

Methods
In this section we describe the method used to retrieve the orbit from the RV curve.In a first step (hereafter, step 1), the con-tribution associated to the pulsation was fitted together with a Keplerian orbit.During a second step (hereafter, step 2), the cleaned (pulsation removed) RV curve was combined with the described astrometric data and fitted with the orvara 1 code based on a Markov chain Monte Carlo (MCMC) method (Brandt et al. 2021).

Disentangling the pulsation (step 1)
The RV curve exhibits two periodic variations: a short-period signal corresponding to the Mira pulsation (following the 530 d periodicity) as well as a long-term trend (see Fig. 1).The long trend can be attributed to the 17 y modulation observed in the light curve (Knapp et al. 1999) but, as this period is longer than the time spanned by the HERMES monitoring, an entire cycle is not visible.As mentioned in Sect.2.1, the errors in the RV from HERMES consist of two terms: the error intrinsic to the CCF fit (photon noise and template mismatch) and a 0.07 km s −1 uncertainty to account for a possible offset of the wavelength zero-point.To better constrain the fit, nine additional RV values obtained by Barnbaum et al. (1995) have been added to the RV curve obtained with the HERMES spectrograph.Fortunately, those data, falling on the previous orbital cycle, correspond to different orbital phases as compared to the HERMES data, allowing us to cover almost the whole orbit.To account for an eventual offset of the RV zero point reference between the two instruments, an additional 0.2 km s −1 uncertainty is added to the RV data from HERMES and represents the minimum value required to increase the relative statistical weight of the oldest RV data from Barnbaum et al. (1995) to a significant level, given their higher uncertainties (± 1.5 km s −1 ).
The RV curve was fitted by a combination of two functions: a Keplerian orbit and an asymmetric sine, sin a .The asymmetric sine is expected to mimic the S-curve observed in pulsating stars (Jorissen 2004) and is mathematically defined as: where ω 0 = 2π/P is the angular frequency, A is the signal amplitude, t 0 the time reference and Γ the tilt parameter satisfying |Γ| ≤ 1. Γ = 0 corresponds to the usual sine function of period, P, whereas Γ = +1 (resp.−1) corresponds to the saw-tooth function tilted to the left (right, resp.).
In the fitting process, the asymmetric sine period was fixed at 530 d and the Keplerian period was forced to remain in the range 6310±1220 days.Those values were obtained by computing the Lomb-Scargle periodogram of the light curve and taking the 3σ uncertainty interval for the long-period signal (see.Fig. B.1).In total, the fitting routine is composed of a nine-parameter function (3 for the asymmetric sine and 6 to describe the orbital motion).
After a first fitting, the contribution of the pulsation is subtracted from the RV and the cleaned signal is fitted with a Keplerian orbit alone using a χ 2 -minimisation.The method allows us to get precise orbital parameters at the price of imposing the period.The standard deviations of the orbital parameters are obtained by computing the square root of the diagonal elements of the covariance matrix.
1 https://github.com/t-brandt/orvara3.1.2.Orbital parameters with orvara (step 2) As a complementary analysis, meant to confirm our results and to obtain the full orbital solution (including the inclination), another independent method was used on the RV curve cleaned from the pulsation.The software package orvara, developed by Brandt et al. (2021), uses a parallel tempering MCMC method to simultaneously fit radial-velocity and astrometric measurements.Combining the measurements allows the companion mass and inclination for binary systems to be computed (Escorza & De Rosa 2023), which would otherwise be inaccessible with just the spectroscopic measurements.The code first fits the RV data then adds the absolute astrometric measurements: astrometric position (α, δ), proper motion (µ δ , µ α ), and Gaia DR3 parallax for each MCMC step.It allows us to simultaneously fit the full orbital parameters (period, P; eccentricity, e; semi-major axis, a; inclination, i; argument of the periastron of the visible star, ω; time of periastron passage ,T 0 ; angle of the ascending node Ω); the mass of the two stars; and an additional radial-velocity jitter, σ vr .Details on the numerical implementation can be found in the reference paper (Brandt et al. 2021).
For V Hya, the primary mass was set as a Gaussian prior of 1.9 ±1 M ⊙ and the Gaia DR3 parallax was used as prior, corresponding to a distance of 434±21 pc.As the parallax uncertainties for AGB stars is often underestimated, different prior values were tested, corresponding to the recomputed distances from Andriantsaralaza et al. (2022): 484 +107 −177 pc and 529 +248 −131 pc.However, no significant difference was found in the final orbital parameters and, therefore, the Gaia DR3 parallax has been adopted in the remainder of the paper.We used 15 MCMC temperatures and for each temperature, with 100 walkers run with 5 ×10 6 steps per walker.After a first try, the MCMC was oscillating between two orbital solutions: one with a high eccentricity for a period of 35 years and a mass ratio above 4 and another with a low eccentricity, a period of 17 years, and a smaller mass ratio.Such a high mass ratio would imply a companion that is more massive than V Hya initial mass, which is inconsistent with its evolutionary stage.We therefore re-started the MCMC with an additional prior on the eccentricity (e < 0.15) to favour the orbital solution with the shorter period, compatible with the observed modulation in the light curve.

Results
The orbital parameters found with the two methods (step 1 and step 2) are listed in Table 1.Despite their different implementation and assumptions, both methods produce a low-eccentricity orbit (compatible with a circular nature) with a period consistent with the light-curve modulation.The orvara solution has larger uncertainties, but it is expected to be more robust and does not require to impose the period.Below, we discuss the values obtained and their sensitivity regarding the astrometric measurement.

Spectroscopic fit (step 1)
The pulsation was fitted with an asymmetric sine amplitude of 5.18 km s −1 .Integrating the curve over half a period to obtain the corresponding radial distance covered leads to a value of ∆R = 0.25 au.For a stellar radius of 2.5 au (see Sect. the mass of the primary.The inclination angle cannot be obtained with the spectroscopic signal alone.It requires additional information from astrometric measurements and an estimation of the primary mass.

Astrometric-spectrometric orbit (step 2)
The right panel of Table 1 gives the orbital parameters obtained when simultaneously fitting the RV curve and astrometric data.The χ 2 red is the overall reduced χ 2 , taken as the sum of the χ 2 for the Hipparcos proper motion (χ 2 H ), the Gaia proper motion (χ 2 G ), and the long-term Hipparcos-Gaia proper motion (χ 2 HG ).The different χ 2 can be found in the central column of Table 2. Figure 2 displays the best fit for the RV curve and the proper motions from Hipparcos and Gaia.The corner plot of the MCMC chain is displayed in Appendix B (Fig. C.1).The best-fitting value for the inclination is 37.7 ± 2.2°.This lies within the range of previous estimations: about 30°for Knapp et al. (1997), and45°for Sahai et al. (2022).

Sensitivity with respect to the Hipparcos data
In Fig. 2, the data point corresponding to the proper motion in right ascension at Hipparcos epoch, µ H α * , appears as an outlier in   2021), as the 60/40 weighted average between the instantaneous proper motion for the original and the second reduction (Perryman et al. 1997, van Leeuwen 2008).For V Hya, the values in right ascension between the two reductions (-14.21±1.41mas/yr and -11.02±1.14mas/yr, resp.) are not compatible within their uncertainties (as it is for the declination component, 2.72±1.14mas/yr and 2.29±1.11mas/yr resp.).This discrepancy can be explained by the variability-induced mover (VIM) flag associated with V Hya for the first reduction.Using an improved chromaticity correction adapted for red stars, Pourbaix et al. (2003) did not confirm the VIM solution; whereas in the second reduction, its recomputed solution is a standard five-parameter solution that does not bear any VIM processing.It indicates the uncertain nature of the µ H α * datapoint.We thus decided to keep the 60/40 mean value computed in the catalog of acceleration (Brandt 2021), but to assess the robustness of the obtained orbit regardless the µ H α * , a sensitivity analysis was performed.The µ H α * value was artificially offset by +4, +2, -2, and -4 mas/yr and, for each new value, orvara was run again to fit the RV curve and the new astrometric indicators, using the same MCMC strategy described in Sect.3.1.2,but without any prior on the eccentricity.The results are the following: for the positive offset all orbital parameters stay unchanged, except for a change in the inclination sign.The solution with an inclination of i=180-38=142°(and Ω = 180+20°) is favoured over the solution with i=38°(and Ω = 180-20°).It corresponds to a µ δ signal in anti-phase, compared to the best-fitting curve from Fig. 2 (right panel), resulting in a change in the orbit motion direction with respect to the previous solution.For the negative offset, the obtained orbital solution is bimodal: either the previous estimated solution (period of 17 yr with e < 0.1, i = 38°and a = 11.2 au) or the solution with a period that is twice larger (and, thus, a mass ratio greater than 5), implying a companion mass higher than 4 M ⊙ .This solution is more probable for higher shift of the µ α * parameter but corresponds to an important increase of the χ 2 (see Table 2).As discussed in Sect.3.1.2,such a massive companion is considered as unlikely.The fit of the PM and RV curves for the four artificial µ H α * are displayed in the appendix, together with their respective corner plot.We conclude that the obtained orbital parameters are robust  against the unreliable µ H α * value.Therefore, we took the orbital parameters (i, Ω, M 1 /M 2 , a) derived from the astrometry as indicative, rather than precise, and considered the orbital motion direction as undetermined: either counterclockwise rotation of the companion around the primary on the sky (i = 40°solution) or clockwise rotation (i = 180-40°solution) Such a degeneracy could be lifted with additional constraints from astrometric data (forthcoming Gaia DR4 catalog) or through interferometric observations (Planquart et al. 2024).

Nature of the companion
The derived mass function (Table 1) allows us to constrain the nature of the companion, hereafter referred to as 'V Hya B'.With the inclination derived at step 2, V Hya B is expected to be more massive than V Hya A. With the masses of the two components being highly correlated (see corner plot in Fig. C.1), no precise mass can be retrieved for the companion without imposing V Hya A mass.Taking a primary mass of 1 M ⊙ as lower limit, the companion mass would be of 1.9 M ⊙ .Such a mass is too high to correspond to a white dwarf and therefore the companion is likely to belong to the main sequence phase.A main sequence companion with a mass above 2.0 M ⊙ would corresponds to an A-type star.The main-sequence nature of V Hya B was previously inferred from the measured UV excess by GALEX (Sahai et al. 2008).Such UV excess would require a hot main sequence star.
To test the nature of the companion obtained from our orbital analysis with the UV excess, the synthetic spectra of the two stars were superimposed on the system photometry from UV to near-IR (NIR).The V Hya B spectrum was modelled as an A0 star with a temperature of 9950 K from Kurucz models (Castelli & Kurucz 2003) and a radius 1.5 R ⊙ corresponding to a mass of 2 M ⊙ .The AGB star was modelled by a black-body of 2700 K (inferred temperature from Knapp et al. 1999), as no synthetic spectra for carbon stars exist for λ < 300 nm.The SED is displayed in Fig. 3. To account for the photometric variability induced by the stellar pulsation, the uncertainties of the photometric data points have been set to 2 mag in the optical and 1 mag in the near-infrared.The spectra are reddened by combining the effect of two sources: the interstellar extinction and the extinction due to the circumstellar envelope (CSE) of the AGB star.
The visual interstellar extinction was estimated to be A V = 0.09 mag (Lallement et al. 2019), and the extinction curve of Gordon et al. ( 2009) was adopted.The interstellar contribution is of negligible effect compared to the effect of the circumstellar dust.The effect of the dust extinction is obtained assuming a CSE composition made of amorphous carbon, AmC (Rouleau & Martin 1991).The dust grains are treated as spherical, with their radius a distributed following the standard MNR distribution (dn/da = a −3.5 , with a min = 0.005 µm and a max = 0.25 µm; Mathis et al. 1977).The total extinction curve is displayed in the inset of Fig. 3, and the total A V is 3.90 mag.The additional 3.54 mag corresponds to a column density for AmC dust grains of the order of 10 12 cm −2 .
The cold nature of the primary does not make it possible to reproduce the observed flux in the UV, whose expected flux in the far-UV (FUV) filter lies several orders of magnitude below its measured value.On the contrary the A0 spectrum, dominating the overall flux for wavelengths below 300 nm, reproduces perfectly the GALEX photometric data in the NUV and FUV filters.The SED confirms the assumption of Sahai et al. (2008) who attributes the UV excess of V Hya to the presence of a hot companion.It has to be stressed that the black body radiation curve overestimates the flux at short wavelengths for red stars such as carbon stars, and explains the observed discrepancy between the model and the bluest optical point (Johnson B) in Fig. 3.Under these circumstances, we could expect to detect signatures of the hot companion even up to 450 nm, as suggested by Lloyd Evans (1991).

Comparison with the light-curve modulations
In this section, we compare the RV periodic motion with the light-curve variation.3: Spectral energy distribution of the V Hya system from UV to NIR (colored square symbols) and model spectra (black: V Hya as a black body, blue: companion as an A0-type star).The diamonds represent the expected FUV and NUV fluxes from the AGB star (black) and the hot companion (blue).The inset displays the adopted extinction curve from the dust (interstellar + circumstellar).

The pulsation cycle
To highlight the RV signal associated to the pulsation, the contribution from the Keplerian orbit was removed.For the light-curve variation, the 17 yr signal was fitted by taking a 530-d rolling mean, thus, the signal associated to the Mira variation can be obtained by removing the obtained rolling mean.
In Fig. 4, the RV and light curve signals cleaned from the long variation are plotted as a function of the pulsation phase.The two signals are found to be in quadrature, with the RV signal being in advance with respect to the AAVSO light curve (∆ϕ =π 2 ).Such a quadrature delay is expected and reflects the change of the AGB photospheric radius: as the star undergoes expansion of its envelope, its brightness decreases (descending part of the light curve).The expansion leads to blue-shifted RVs for the observer (lower lobe of the RV signal).The same phase delay was already found for the CS Mira star R CMi (see Fig. 5 from Jorissen 2004).

The 17-yr cycle
In Fig. 5, the RV solution (extrapolated to early cycles) is overplotted over more than a century of photometric observations, covering six dimming events.Their duration and depth vary from cycle to cycle, with a mean amplitude of 3 mag and a mean duration of 30% of the cycle.The shape of the dimming events is symmetric: the ingress and egress have the same slope, with a value that varies from 0.001 to 0.003 mag/d, depending on the cycle.
The two long-period signals are again found to be in quadrature.But, in contrast with the short-period signal associated to the pulsation, the RV signal is phase-delayed (∆ϕ = π 2 ) with respect to the light curve, hinting at a different mechanism than intrinsic variations to explained this phase lag.This time, it can easily be explained in the framework of binary motion: the minimum light found by the AAVSO monitoring always corresponds to the superior conjunction, when the primary, namely, V Hya A, is behind its companion and its RVs switch from red-shifted to blue-shifted.A direct consequence is that the obscuration causing the deep minima in the light curve occurs during the companion passage in front of V Hya.
Due to the particular shape of the eclipse, with a long duration and a deep obscuration, and the system inclination, a simple eclipsing scenario by a stellar body has to be rejected.Instead, a more extended and opaque object attached to the companion is needed.The inclination angle of about 35−40°(see Table 1) implies that the eclipsing object must have an important extension in the direction perpendicular to the orbital plane.A possible explanation could be that the primary star is eclipsed by a vertical dense outflow.

Time series analysis of spectral lines
As discussed in the introduction, bipolar outflows have been observed and studied in post-AGB binaries and their signature can be clearly detected in optical spectra.Bollen et al. (2022) showed that the Hα line in the spectra of post-AGB binaries with jets shows a variable profile whose time-series can be modeled to obtain different jet parameters.In this section, we report on the temporal variation of particular spectral lines: the sodium and potassium doublets, Hα, and the C 2 (0,0) vibrational band, and we connect them with the stellar and orbital star periodicity.

Sodium doublet
The most interesting features are found in the sodium doublet Na i D-lines.In Fig. 6, the Na i D-profiles at two different orbital phases are plotted, together with a reference spectrum.The phase 0.5 corresponds to the light minimum (superior conjunction, V Hya behind its companion) and the phase 0 corresponds to the inferior conjunction.A strong and variable blue-shifted absorption is present.This absorption is visible in both lines of the doublet and is most blue-shifted around phases 0 and less  blue-shifted around phases 0.5.The profiles of the two sodium lines do not seem to differ.
To highlight the variation of the Na i D-line profile, the spectra are plotted as a function of the orbital phase and interpolated between the phases to obtain a smooth representation of the variation.This is shown in Fig. 7.A sample of spectra at different orbital phases can be found in Fig. D.1.Since the HERMES monitoring duration is shorter than the orbital phase, the phases between 0.05 and 0.35 are not yet covered and are left blank in the figures.The 17 year periodicity of the Na i modulation was confirmed by prior spectroscopic monitoring of Lloyd Evans & Carter (1993).Indeed, a weakening of Na i absorption followed by an increasing of its emission core is reported at the beginning of the fading event of 1993 (see red rectangle in Fig. 5).The large absorption band seen in Fig. 6 is also visible in Fig. 7 and reaches its maximum blue-shifted value, -150 km s −1 at phase 0.5, and its minimum at phase 0. At phase 0.5, an emission (represented in red) is also present on the right side of the blueshifted absorption.The combination of blue-shifted absorption and emission is reminiscent of P-Cygni profiles.The periodicity of the absorption pattern seems to be the same as the orbital one but it is phase-shifted.More precisely, they are in quadrature as the absorption band is at its maximum blue-shifted value at maximum light when the star velocity equals the centre-of-mass velocity.This is easily seen in Fig. 7, thanks to the sine black line representing the orbital motion.This phase shift is quite puzzling since one would expect it to be in phase if the absorption was due to the primary star or in anti-phase if it was due to the companion motion.Thus this feature cannot be related to the motion of the stars but must have a different origin.A possible explanation is the presence in the binary system of a conical outflow with a velocity gradient.As the system goes along its orbital motion, the cone obstructs V Hya light with a varying angle and therefore a varying absorption velocity.
We can see that even though the shapes are different (Na i D2 exhibits a double-peak emission whereas Na i D1 only a leftcentred peak), their equivalent widths are quite similar.The ratio between D2 and D1 equivalent widths range from 0.8 to 1, for a spectral window of 2 Å.Those values being well below 2 (the expected values for the case of an optically-thin media, and ratio of the oscillator strength), indicating that the media is therefore optically thick.

Potassium doublet
Peculiarities are also observed in another alkaline metal: the potassium doublet K i D-lines at λ = 7664.89Å and 7698.96Å.A large blue-shifted absorption band is present and shares the same phase dependency as the sodium doublet, being most blueshifted at phase 0.5.The absorption is narrower, especially at phase 0 where the bandwidth reaches about 20 km s −1 .Pseudoemission lines following the orbital variation are visible in the red-shifted part of the spectrum (in red in Fig. 7).Their evolution follows the orbital trend, indicating that those lines are of photospheric origin, contrary to the alkaline resonance lines whose variations are phase-shifted, hinting at a different origin above the photosphere.

Swan band
The Swan bands are vibrational bands of the C 2 molecule.Their emission is temperature-dependent and begins to develop above 2450 K (King 1948).The most intense one is the (0,0) transition at 5165 Å.As V Hya is a cool carbon star, only emission in the (0,0) band is observed in its spectrum.In Fig. 7, the C 2 profile seems to change with the orbital phase.The line profile around λ = 5165 Å is quite different between superior and inferior conjunctions.Around phase 0.5 a pseudo-emission seems to be present but this feature is absent around phase 0. Such an abrupt emission of C 2 was also reported during the obscuration event of 1993 (Lloyd Evans & Carter 1993, Lloyd Evans 1997) confirming the relation between the C 2 emission and the star drop of brightness.

Balmer lines
No amplitude modulation correlated with the orbital motion is found in the Balmer lines.Only the Hα emission is strong enough to be observed with a high signal-to-noise ratio (bottom right panel of Fig. 7) and mainly shows radial-velocity variation that follows the orbital trend with a faint amplitude modulation due to the Mira pulsation.No blue-shifted absorption as for the Na i or K i lines is seen.Compared to other carbon stars with a Mira pulsation, the intensity of Hα for the star V Hya is quite reduced overall, even outside the eclipse phases, reaching at most five times the continuum.As previously reported by Lloyd Evans (1991), the Hγ line is found in absorption and could be associated to the hot companion signature.Indeed, no absorption is found in the Balmer lines of the reference carbon star X Tra.

Spatio-kinematic modelling of the outflow
As discussed in Sect 5.1, the presence of varying blue-shifted absorption in the sodium doublet can be interpreted as the signature of a high-velocity outflow attached to the companion star.Bollen et al. (2021) studied the jet imprints on Balmer lines for a sample of post-AGB binaries for a variety of jet configurations.Their modelling approach consists in fitting the model spectral lines created by a parametric jet model to the phase-resolved Hα profile using an MCMC routine.
To deduce the jet structure (geometry, velocity, and density profile) we used the same spatio-kinematic modelling engine2 and adapted it for the phase-resolved sodium line at 5895.92 Å. Hereafter, we summarise the overall fitting routine adopted to fit V Hya sodium spectral line over the range 0 to -250 km s −1 .We refer to the reference paper (Bollen et al. 2019) for a more general description of spatio-kinematic jet structures and their numerical implementation.Among the jet models available, we restrict ourselves to the stellar jet model, which is similar to the model developed by Thomas et al. (2011) for the Red Rectangle Nebula.It has been shown that the fitting routine provide similar results for different configurations, hence we choose the most generic model containing the smallest number of free parameters and assumptions about the launching site of the outflow.The stellar jet model and its different parameters are described in Appendix E.1.As the absorption feature (believed to be the jet imprint) is present at all observed phases which cover nearly 70% of the orbit, the jet opening angle is expected to be larger than the inclination angle of the system.Therefore, no observed spectrum can be used as template for the unattenuated (background) spectrum in the fitting process.Instead, we used the high-resolution spectrum of the carbon star X Tra as template for the unattenuated spectrum.The spectrum was modulated in RV to reproduce V Hya orbital motion and in flux at to reproduce the 530 d intrinsic variability.The dynamic spectrum is shown in Appendix E.2.
The procedure is meant to quantitatively model the absorption of the stellar contribution through the jet.However, in the case of V Hya, an extra emission pattern is seen in the sodium profile and needs to be modelled and implemented on top of the Bollen et al. (2019) method.The modulated emission pattern is not found in other carbon-star spectra, suggesting it to be the result of the star-jet interaction rather than of the star photosphere.We interpreted the blue-shifted emission as re-scattering of V Hya light in our line of sight.Additional details on the implementation and normalisation of those two effects are de-Table 3: Jet parameters, range and best-fit values.For the definition of the parameters, see Appendix E.1.

Jet parameter
1.5 − 2.5 2.19± 0.16 scribed in Appendix E.3.In summary, the general fitting procedure is a three-step process: first, the synthetic model spectrum is used to fit the stellar-jet model to the phase-dependent observations and to obtain the best-fitting absorbing jet.Then, from the obtained spatio-kinematic parameters, the additional source term is modelled.Lastly, the jet parameters are fine-tuned by fitting the absorbing jet on an updated spectrum consisting in the synthetic model and the modulated source term.

Best-fit model
The jet parameters are listed in Table 3, together with their bestfitting value and their initial range.In the table, i represents the inclination of the orbital plane, α is the jet opening angle, θ cav is the angle of the jet cavity, θ tilt is the angle between the jet axis and the normal to the orbital plane, 0 is the gas velocity at the jet axis, α is the velocity at the cone edge, p v and p d are the power indices governing the density and velocity laws (see Appendix E.1), c τ is the scaling factor (see Appendix E.1), and R 1 is the radius of the primary star.Due to a strong correlation between the maximal jet opening angle, α, and the system's minimal inclination, we constrained the inclination angle, i, to its value obtained from the orbital analysis (see Table 1).For the definitions of all the other parameters, we refer to Appendix E.1.The velocity and density structure across the jet is displayed in Fig. 9.The velocity follows a decreasing and nearly-linear trend starting from 0 at the jet axis to α at the jet edge.The maximum deprojected velocity, 0 , is below 200 km s −1 , although such a value is well above the velocity of the wind feeding the AGB circumstellar envelope.
The normalised density profile is maximum at the edge of the cone and minimum along its axis, with a cavity angle θ cav .Such a hollow cone is also encountered for post-AGB binaries (Bollen et al. 2019) and reveals the presence of a cavity along the jet axis.The presence of a cavity in the V Hya system was already put forward to explain the launching mechanism of the [S ii] jets observed by Sahai et al. (2016) and Scibelli et al. (2019).The observed and modelled dynamic spectra of the Na i D1 line is displayed in Fig. 8.It can be seen that the model reproduces the temporal modulation of the blue-shifted absorption quite well.The obtained jet geometry and orbital configuration are sketched in Fig. 10.

Limitation of the model
The above-described analysis allows us to conclude that a highvelocity outflow launched by a companion orbiting V Hya is able to satisfactorily reproduce the observed spectral modulation.The choice of such a modelling tool was motivated by previous large-scale observations of high-velocity outflows (Knapp et al. 1997, Hirano et al. 2004) and by the close proximity of V Hya evolutionary stage to the post-AGB phase where those structures are widely found.However, no direct evidence of the companion motion is seen in the spectrum as its flux contribution is expected to negligible above 450 nm (see Sect. 3.4), and the exact launch site of the jet remains unknown.Further analyses are needed to constrain the temperature and the absolute density structure through a full 3D radiative transfer code.Those ingredients and the putative properties of the companion (see Sect. 3.4), are needed to fully characterise the physical mechanism feeding and launching the jet.

Discussion
In this section, we discuss the possible origin of the long-period obscuration event in the framework of the binary model obtained.Alternative explanations not invoking a companion are then discussed.Finally, we confront the obtained model of the V Hya system with previous large-scale observations to build a multi-scale consistent picture of the system.

Obscuration by dust
As highlighted by the orbital analysis in Sect.4, the dimming episodes of about 3 mag occurring every 17 years are associated with the superior conjunction and the high-velocity absorption of alkaline resonance lines.These observations can be explained in the framework of an absorbing gaseous jet passing in front of the line of sight.However, a gaseous jet alone cannot be responsible for the visual dimming event, as such a broad-band attenuation requires particles with a gray absorbing power, such as dust grains.A light curve modulation induced by the variation of the dust grains column density through the line of sight with the orbital phase would require a non-uniform distribution of particles above the orbital plane, due to the inclination well below 90 degrees.For a mean magnitude drop of 3 mag, the observed modulation can only be recovered by assuming the highest dust density at the jet centre, with a dust column density that is three times denser at the superior conjunction when the line of sight crosses the axis of the cone.The dust distribution through the jet polar angle would be opposite to what was found for the gas density in the jet, better described by a hollow cone at its centre.It raises the question of the physical mechanism behind the dust and gas distributions.Only a better characterisation of the absolute gas density and temperature using a full 3D radiative transfer code would allow us to determine more precisely the dust formation limit, the mass-to-dust ratio, and the exact dust spatial distribution in the system and understand their possible coupling with the gas.However, such a complex investigation is out of the scope of the current paper, as we are focused on the system dynamics.Interestingly, some post-AGB systems also exhibit a long-term modulation of their light curve, but of a smaller amplitude, (<0.5 mag) associated with the radial-velocity and atomicline variations.In RV-Tauri stars, those obscurations occur at the inferior conjunction, the opposite situation as encountered here, explained by the particular orbital configuration of the system where the line of sight grazes the circumbinary disc due to the high orbital inclination (Manick et al. 2019).A few other post-AGB stars with a bipolar morphology show small-amplitude obscuration events together with radial-velocity and Hα variation: For instance, V501 Lup exhibits obscurations of about 0.5 mag in the visible and near-infrared which share the orbital periodicity (2727±26 days) and could be the result of dust absorption and scattering by the nebula (Manick et al. 2021).For 89 Her, a long-term obscuration with a periodicity of 5000 d is reported together with modulation of Hα, Na i D, and Ba ii lines.With the the long period of the light curve being larger than the orbital one (289 days, Bollen et al. 2021), two scenarios have been suggested: an orbital wobbling of the circumbinary disc inducing a misalignement with the hour-glass structure or recurrent massejection episodes (Gangi et al. 2021).

Origin of the long-period variability
For the sake of completeness, it ought to be mentioned that alternative scenarios not involving a binary system can be invoked to explain the obscuration of variable stars: a long-period secondary pulsation period that is twelve times longer than the fundamental tone of 530 days or an episodic dust condensation from the star.Below, we describe both of them and argue why we consider them less likely.For the first one, as discussed in Sect.4, the phase shift between the radial velocity and light curves is found to be +90°for the long cycle.Therefore, the usual scenario of radial pulsation induced by the so-called kappa-mechanism (associated instead with a -90°phase shift) cannot be responsible for the 17 yr periodic fading.Besides, this long period would make V Hya an outlier in the usual period-luminosity diagram, presented by Wood (2000), being too long to belong to the Dsequence ("long-secondary periods").
The second interpretation, namely, episodic dust formation, is the usual mechanism put forward to account for the abrupt fading by several magnitudes of R Coronae Borealis (R CrB) variables.Those stars exhibit light curves with important dimming events separated by several years (Clayton 1996).As these obscurations are correlated with intrinsic dust formation, their occurrence is far from periodic and their fading episodes show an abrupt decrease followed by a gradual restoration of the star brightness.The only exception is the sub-class of DY-Per objects, labeled as cold hydrogen-deficient R CrB stars, whose light curves exhibit a rather regular periodicity in their fading.Considering the regularity of V Hya light curve over a century of observations and previous hints in favour of the binary hypothesis (Hirano et al. 2004, Sahai et al. 2008, Sahai et al. 2016) we consider intrinsic variability scenarios to be less likely to explain the long-term periodicity of the star.

Connection to large-scale bipolar outflows
The characterisation of the orbital parameters (including the system inclination) associated with the detection of the highvelocity jet allows us to link the different building blocks of the system; in particular, the jet that is supposed to be attached to the companion, with its sky-projected orientation.A sketch of the system is displayed in Fig. 10.As discussed, the inclination found is consistent with the previous inclination estimates from the bipolar outflow measurement.Additionally, the ascending node being close to 180°(see Table 1), the jet axis (perpendicular to the orbital plane) points towards the east.Its direction is consistent with previous measurements, at a larger scale, of the high-velocity blue-shifted [S ii] emission and the high-velocity parabolic outflow found in the radio-emission of CO molecules (Sahai et al. 2016(Sahai et al. , 2022)).The line of nodes and jet axis angle have an uncertainty of ±10°due to the inclination degeneracy (see Sect. 3.1.2).
Figure 11 superimposes the sky-projected jet axis and the system's line of node on the high-velocity emission obtained from ALMA observations (PI: Sahai R., Project ID: 2018.1.01113.S).Both geometry and orientation between the two structures match.This result is to be underlined as the two sets of observations, having entirely independent approaches, provide a complementary set of evidence for the system morphology.On the one hand, in our study, the conical geometry and orientation of the jet is deduced by fitting phase-dependent blue-shifted absorption of optical resonance lines and compare them with the system orbit.On the other hand, the bipolar jet Fig. 11: Sky-projected system orientation compared to radio observation.The contour levels is a stereogram view of the 12 CO J = 2−1 high-velocity emission.The contour levels are drawn at 6, 8, and 10 mJy/beam for two maps equally shifted by ±100 km s −1 from the centre-of-mass velocity.The blue-shifted contours are displayed in blue, and the red-shifted in red.The dashed green line represents the orbital line of node, the black arrow the jet axis and the dashed black line sketches the cone edges.
is mapped using a one-shot observation in the radio-emission of 12 CO J = 2−1 at larger scale.
This compatibility strongly suggests that they may share a common origin: the modeled conical jet could be seen as the innermost part of a large-scale parabolic structure.This scenario does demonstrate the link between the dynamical binary interaction (as inferred from our spectral monitoring) and the morphology of the large-scale bipolar nebula.

Conclusion
Thanks to an unprecedented twelve-year-long spectroscopic monitoring of V Hya, we performed a temporal study of the dynamical processes at the origin of its two periodic variations.The main results obtained in this paper are two-fold.
In the first step, the radial-velocity curve revealed a short variation of 530 days superimposed on a long trend.Combining the RV data with secular visual photometry and astrometric acceleration, we disentangled the signal associated with the Miralike pulsation of 530 days from the 17 year cycle.The latter is interpreted as the signature of a Keplerian orbit for a period equal to the long photometric cycle of 6351 days, a low eccentricity (e = 0.02 ± 0.02), and an inclination of 37°±2°, compatible with a main sequence companion.
In a second step, the time-resolved spectroscopy unveiled the presence of a high-velocity gas stream in the system.The resonant lines of sodium and potassium display phase-dependent absorption, blue-shifted up to -150 km s −1 at the epoch corresponding to the visual minimum, and superior conjunction.We interpret the gaseous jet as the imprint of an interacting binary system.A spatio-kinematic modeling of the jet launched from the companion orbiting around V Hya is fitted to reproduce the observed spectral-line modulation.The obtained jet model is described as a hollow cone with a wide opening angle and a velocity gradient along its radial direction.The highest velocities are oriented along the jet axis, whose projection on the sky plane is pointing toward the east.In the meantime, the jet geometry obtained (velocity and direction) is consistent with previous large-scale observations of the high-velocity bullets found in [S ii] lines by Sahai et al. (2016), and high-velocity bipolar flow from the molecular lines of CO (Knapp et al. 1997, Hirano et al. 2004, Sahai et al. 2022).
It is, to our knowledge, the first case where the temporal variation of the light curve, high-velocity gas, and the AGB star velocity are reported to be correlated.Our proposed scenariosimilar to what is found in many post-AGB binaries (e.g.Bollen et al. 2022 and references therein), where the evolved star orbits around the jet launched by its main sequence companionexplains the reported observations in a common framework.To conclude, this multi-epoch study reveals the complex behavior of V Hya, which can be seen as a benchmark case to better understand the shaping mechanism induced by a binary system at the origin of bipolar planetary nebulae.
Further studies reproducing the radiative transfer inside the circumstellar and circumbinary environments of V Hya system could help to constrain the connection between the binary interaction and the large-scale bipolar morphology.Additionally, enlarging the presented analysis to a larger sample of AGB stars with known companion or bipolar morphology, for instance, TU Tau (Gray et al. 2023), R Crt (Khouri et al. 2020), UV Aur (Herbig 2009), RR UMi (Dettmar & Gieseking 1983) would be needed to obtain a statistical distribution of the orbital parameters of such systems.performing a first-order ray-tracing of a grid of scattering points evenly spaced in the cone.The final source function is obtained by summing the contribution of each scattered photon depending on their final projected velocity and weighted by their probability to be re-scattered in the direction of the observer.This probability only depends on the deflection angle between the photon initial direction and the line of sight.
The orbital modulation, β, accounts for the normalisation correction due to the different origins of two emitting sources unresolved by the instrument.The scattering source being emitted from an extended region offset from the star, it does not follow V Hya orbital motion and is therefore not affected by the periodic obscuration event.The spectra displayed in the analysis are always normalised by averaging over a given spectral window dominated by photospheric pseudo-emission.The obtained window-normalised spectra can therefore be compared, independently of their flux modulation.The side effect of such a normalisation is that any extrinsic emission, unaffected by V Hya obscuration, is artificially modulated by a factor inversely proportional to the stellar flux modulation: appearing brighter when the star fades and nearly vanishing at the opposite phase.In practice, such extrinsic emission will vary in strength with the orbital phase, ϕ, according to a correction factor, expected to be related to the light-curve modulation, equal to: where ∆m ϕ is the visual light-curve modulation expressed in magnitude.

Fig. 1 :
Fig. 1: Radial velocities of V Hya.Top: RV curve as a function of time.The full line is the Keplerian orbit and the dashed line is the composite signal (orbit and pulsation).Middle: Phase-folded cleaned orbit of V Hya.Bottom Residual values with the ±3σ from the mean is represented as the grey-shaded region.
fitted model.The instantaneous values of the proper motion at Hipparcos epoch are computed by Brandt (

Fig. 2 :
Fig. 2: orvara results: RV curve (left), proper motions in right ascension (centre) and declination (right).In all plots, the black thick line represents the best-fitting orbit, while 40 other well-fitting orbits are included and colour-coded as a function of the companion mass.
Fig.3: Spectral energy distribution of the V Hya system from UV to NIR (colored square symbols) and model spectra (black: V Hya as a black body, blue: companion as an A0-type star).The diamonds represent the expected FUV and NUV fluxes from the AGB star (black) and the hot companion (blue).The inset displays the adopted extinction curve from the dust (interstellar + circumstellar).

Fig. 4 :
Fig.4: RV curve cleaned from the Keplerian motion, the AAVSO light curve with the 530 d rolling-mean signal removed, and their best-fitting functions phase-folded with respect to the pulsation period.The red dots correspond to the RV data affected by shocks and not used in the analysis (see Sect. 2).
Fig.5: AAVSO lightcurve of V Hya recorded since 1913, displaying six evenly spaced dimming events.The blue curve is the extrapolated spectroscopic orbit.The red regions are the observing epochs of sodium absorption (see Sect. 5.1) from the HERMES monitoring (filled rectangle), and from the observing dates of LloydEvans & Carter (1993) (dashed rectangle).-200-100 0 -200 -100 0 RV[km/s]

Fig. 6 :
Fig. 6: Na i doublet near superior conjunction (phase 0.46, top) and near inferior conjunction (phase 0.98, bottom).The dashed line is the spectrum of the reference carbon star X Tra used as template for the modelling.The vertical dashed lines indicate the centres of the two components of the doublet.

FluxFig. 7 :
Fig. 7: Dynamic spectra of the sodium doublet (top row), potassium doublet (middle row), C 2 (0,0) (bottom left), and Hα (bottom right) of V Hya as a function of the orbital phase.In each panel, the observed phases are shown twice to guide the eyes, time runs down.The colour represents the pseudo-continuum-normalised fluxes taken as the median value of the spectral window.The black curve represents the primary motion while the vertical line marks the systemic velocity, and the horizontal lines on the left indicate the observed phases.

Fig. 8 :Fig. 9 :
Fig. 8: Observed dynamic spectrum of the Na i D1 line (left) and the model spectrum obtained by the spatio-kinematic fitting of a conical jet (right).

Fig. 10 :
Fig. 10: Geometric representation to scale of the binary system where V Hya A is represented as an orange filled circle and V Hya B as a black dot (to be visible, its radius was increased).V Hya A mass was set to 2 M ⊙ , corresponding to a mass of 2.5 M ⊙ for V Hya B. Upper left panel: System depicted edge-on at the superior conjunction.The black line is oriented along the line of sight.The green cross indicates the centre of mass.Lower left panel: System viewed face-on, with the Roche equipotential in black line and the Roche lobe of the two stars in dashed gray.Right panel: System projected on the sky-plane at the superior conjunction for the counterclockwise solution, assuming a parallax of 2.31 mas.The black line is the best projected orbit, the dashed black line represents the line of nodes.The effect of the mass-ratio uncertainty on the orbital separation is illustrated by the colored ellipse.

Fig. E. 3 :
Fig. E.3: Illustration of the stellar light propagating through the cone: the light passing through the gray shaded area is attenuated by the jet cone (see the blue line tracing the ray).The light emitted in another direction but passing through the rest of the cone can be re-scattered or re-emitted into the line of sight (red arrow).
Fig. E.4: Normalisation effect on the Na i D1 line profile for spectra taken at two different phases.(a) The normalisation is obtained by dividing each spectrum by its mean over a 100 Å window.This procedure is incorrect (see text).(b) The normalisation between the two spectra takes into account the relative flux ratio of 0.33 between the two phases corresponding to a drop of 1.31 mag.The RV zero point is set to λ = 5895.92Å.

Table 1 :
Orbital parameters from V Hya for the two methods described in Sect.3.1.Left panel: parameters from step 1.Right panel: parameters from step 2, parameters with asterisk (*) are discussed in Sect.3.3.DoF stands for 'degrees of freedom'.