© ESO 2020

1. Introduction

Near the end of their lives, all stars in the Universe initially more massive than 0.8 M_⊙ lose 40−80% of their mass via a stellar wind. For stars less than 8 M_⊙, this phase of intense mass-loss is known as the asymptotic giant Branch (AGB) phase. The mass-loss process close to the star is complex (e.g. Bladh et al. 2015) but, beyond a few stellar radii, the radial force on an AGB wind is thought to be due to radiation pressure on solid-state dust grains formed in the cool and extended atmospheres of the AGB star (Höfner & Olofsson 2018). It has typical mass-loss rate values between 10⁻⁸ and 10⁻⁴M_⊙ yr⁻¹, and velocities between 5 and 20 km s⁻¹ (Ramstedt et al. 2008). For mass-loss rates above ∼10⁻⁷M_⊙ yr⁻¹, the associated time scale for stars to shed their envelope is much shorter than the nuclear burning time scale, implying that mass loss dominates and/or determines all further evolution (McDonald et al. 2018).

Historically, these winds were believed to possess a uniformly spherical geometry (Neri et al. 1998). Hence, they were also modelled as such, and most stellar wind properties presented in the literature (such as wind mass-loss rates) were derived under this assumption. However, recent observations have revealed that in fact most, if not all, AGB nebulae harbour a high degree of morphological complexity. Among these structural curiousities are, for example, bipolar structures, arcs, shells, clumps, spirals, tori, bubbles, and rotating discs (Decin et al. 2012, 2015, 2020; Bowers & Johnston 1990; Cox et al. 2012; Maercker et al. 2012; Balick et al. 2013; Ramstedt et al. 2014; Kervella et al. 2015, 2016). A range of dynamical processes that operate on different length and time scales have been proposed to explain some of these structures (Woitke 2006; Cox et al. 2012; van Marle et al. 2014; Chen et al. 2016). However, the frequent occurrence of cylindrical morphologies among planetary nebulae (PNe), their presence among the post-AGB PN predecessors (Bujarrabal et al. 2015), and the recent first detections of comparable morphologies around AGB stars (Kervella et al. 2016; Doan et al. 2017; Homan et al. 2018a) are steadily shifting the scientific consensus towards binary activity as the dominant shaping mechanism.

Most high-resolution observations of the circumstellar environments of AGB stars exhibit strong indications that binarity is the source of the spectrum of perturbations of the circumstellar gas and dust distributions. Furthermore, it has been shown that the misinterpretation of companion-perturbed AGB winds as smooth nebulae affects derived mass-loss rates (Homan et al. 2015; Decin et al. 2019). Hence, it is important to perform detailed observations, analyses, and modelling of the targets that exhibit deviations from spherical symmetry in order to obtain a better grasp on the physical and chemical conditions that dominate the wind-launching and wind-companion interaction zones.

Mass is one of the decisive evolutionary properties of stars in general. The rate of mass transport away from the stellar mantle by these stellar winds, which varies as a function of age, governs the change of the stellar mass with time, which in turn has decisive repercussions on the evolutionary path of the star and nucleosynthesis in its core. Hence, a proper grasp on stellar, and by extension galactic chemical, evolution cannot be achieved without a detailed understanding of the wind physics over a star’s life cycle.

It is in this context that we present the latest results on the circumstellar environment of the M-type AGB star EP Aquarii. For an overview of the general CSE properties we refer to Homan et al. (2018b), Hoai et al. (2019) and Tuan-Anh et al. (2019), who discuss the CSE gas properties based on ALMA cycle 4 observations of the CSE at a resolution of ∼0.2″. We provide a short summary of the main results: EP Aquarii is an M-type, semi-regular variable star (spectral type M8III), and is currently in the AGB evolutionary phase. It has a pulsation period of 55 days (Lebzelter & Hinkle 2002), and HIPPARCOS measured its distance to the Earth to be approximately 135 parsec (Lebzelter & Hinkle 2002). It possesses an effective temperature of T_* ≃ 3200 K, a mass of M_* ≃ 1.7 M_⊙, and a luminosity is ∼4800 L_⊙ (Dumm & Schild 1998). It is known to have a radial velocity relative to the local standard of rest (LSR) of v_* = −34 km s⁻¹, a mass-loss rate of ∼1.2 × 10⁻⁷ M_⊙ yr⁻¹, and a terminal wind velocity of ∼11 km s⁻¹ (Nhung et al. 2015). High spatial resolution observations in the last few years have revealed the circumstellar environment of EP Aqr to be extremely complex. ALMA observations of CO emission reveal that the wind of EP Aqr possesses a pronounced face-on equatorial density enhancement (EDE) that harbours a spiral. This EDE is encased in two, nearly pole-on outflow hemispheres that could be considered bi-conical, and that exhibit a multitude of nearly circular filamentary structures. In the inner wind, the CO emission patterns resemble the wind-Roche-lobe-overflow hydrodynamical models (Mohamed & Podsiadlowski 2012). The SiO emission within the same dataset reveals a peculiar void of emission about 0.45″ to the west of the continuum brightness peak (CBP) position, located at the base of the spiral, and coincident with the position of the companion in the obovementioned hydrodynamical models. As a consequence, the void has been suggested to be a local environment generated by a companion. Finally, the SiO emisison also exhibited enhanced gas dynamics that has been interpreted as high-velocity jets.

In the current work we adopt a distance of 135 pc (Lebzelter & Hinkle 2002) in order to retain consistency with the work by Homan et al. (2018b). However, the literature also presents an alternatitve distance measurement made by Gaia of 124 pc (Brown 2018). All analyses presented in this work scale with the distance. Hence, if the Gaia distance is more accurate, this difference implies a reduction of all quantitative analyses by about 8%.

In this first paper, we present the latest cycle 6 band 6 ALMA data of the wind of EP Aqr at a resolution that is finer by a factor ten, that has been combined with the older dataset. In addition, we also present VLT/SPHERE data on the target, showing the distribution of dust in the inner wind. Subsequent papers will present the three-dimensional hydrodynamical modelling of these data, and the three-dimensional radiative post-processing of the model. The paper is organised as follows: In Sect. 2 we give an overview of the ALMA observations and how they have been prepared for analysis, and continue with a similar report on the SPHERE/ZIMPOL data in Sect. 3. Subsequently, we discuss the nature of the continuum emission in Sect. 4. We proceed with an elaborate description of the ALMA molecular emission in Sect. 5. These observational results are discussed in Sect. 6. Finally, we summarise our findings in Sect. 7.

2. ALMA data acquisition and reduction

2.1. Observations

EP Aquarii was observed by ALMA in Band 6 for project code 2018.1.00750.S (PI W. Homan) using the extended configuration C43-10, hereafter referred to as TM0. These observations complement data taken in three configurations at lower resolutions, described by Homan et al. (2018b), and aspects of data reduction common to all observations are described there. EP Aqr was observed at RA 21:46:31.88046, Dec −02:12:45.5471. The antenna configuration consisted of 46 antennae, and covered baseline lengths between 83 and 16196 m, which are sensitive to angular scales up to 0.48″. As the phase reference source we used J2156−0037 (position RA 21:56:14.758, Dec 00:37:04.594), and the QSO targets J2232+1143 and J2253+1608 were used as bandpass calibrators. The target:phase reference switching cycle was about 80 s, with about 50 s on target, for a total of 99 min on target, and the angular separation was 3°. J2206−0031 was observed as a check source at 22:06:43.282598 −00.31.02.49604, 246.5° from the reference source. During observations, the precipitable water vapour (PWV) was between 1 and 2.6 mm.

Listed in Table 1 are the four spectral windows (spw) that were observed, each with central frequencies at 217.128, 220.422, 230.562 and 235.176 GHz. Each spw has a width of 0.468692 GHz, divided into 1920 channels, resulting in a velocity channel width of ∼0.32 km s⁻¹, depending on frequency.

2.2. Data reduction and imaging

The data were calibrated using the ALMA pipeline. We split out the pipeline-calibrated target data, adjusted to the Local Standard of Rest (v_LSR) in the direction of the target, with one set of data at full spectral resolution and one averaged to 32 channels per spw, for continuum. Scanning the visibility spectra, we selected the line-free channels and used these to image the continuum, which were combined over all epochs. We then used this image as a starting model for self-calibration of phase and then amplitude. The corrections from the continuum self-calibration were applied also to the line data. The best TM0 continuum image has a noise rms level of 0.032 mJy, a 21 × 23 mas² beam size at a PA (PA) of −212°.

After continuum-subtraction, we imaged 100−150 channels of each spw, corresponding to the maximum velocity extent of the previously detected molecular emission, and used a pixel size of 0003 and a total image size of 6″. The rms noise, given in Table 1, was measured in the channels devoid of signal.

2.3. Combination of multi-array data

The combination of the lower-resolution data is described in Homan et al. (2018b). We followed the same procedure to add the extended-configuration data to the combined dataset. At each epoch, the pointing position was adjusted for the proper motion of EP Aqr, based on the HIPPARCOS data available at proposal submission, (24.98 ± 0.75, 19.54 ± 0.33) mas yr⁻¹ (van Leeuwen 2007).

The positions of the continuum peaks at each epoch are given in Table 2. The extended configuration peak is offset by (40, 18) mas from the pointing position. Position errors due to noise are < 1 mas and the check source’s apparent position was displaced by ∼4 mas, mainly due to errors in phase transfer from the reference source. This suggests that the offset of (63, 9) mas noted by Homan et al. (2018b) was not entirely due to noise (although the uncertainty was greater at the lower resolution), but there was a small error in the input positions and/or proper motions. The recently published Gaia DR2 (Brown 2018), gives an epoch J2000 position of 21:46:31.84675 −02:12:45.9017, with an uncertainty of 0.46 and 0.39 mas, respectively, and with proper motions of 28.204 ± 0.763, 18.794 ± 0.682 mas yr⁻¹ in right ascension and declination, respectively. This predicts a position at epoch 2019.5 of 21:46:31.8834 (±0015) and −02:12:45.535 (±00133). Our 2019 position is 21:46:31.8832 −02:12:45.529, an offset of (0.003″, 0.006″), which is within the Gaia uncertainties.

The relative error in the pointing positions between 2016 and 2019 due to the difference between HIPPARCOS and Gaia proper motions is < 10 mas, which is less than the position uncertainties of all but the most recent data. Therefore, we did not adjust the positions, but combined the data aligned with the phase centre of the earliest observation. This proved sufficiently accurate to avoid imaging artefacts, but the 2019 TM0 position given above should be used for astrometry of the stellar position.

The continuum peak flux densities for all epochs are very similar at resolutions < few 100 mas, with a value of ∼16.6 ± 0.5 mJy beam⁻¹, and a maximum change of 6%. This is consistent with the ALMA flux scale accuracy of 7%. The larger-scale continuum emission varies up to 20% (see Table 2). These variations are probably intrinsic to the immediate stellar environment, and not calibration errors (see Sect. 4.1). Stellar variability will affect the immediate vicinity and radiatively excited lines, especially near the star, while for lines formed outside a few stellar radii the effects will be diluted and probably differ in different directions in an inhomogenous envelope. We therefore did not attempt to modify the line flux density scales.

We combined the line-free continuum from all epochs, with equal weight for each dataset. The parameters of the combined image are given in Table 2. We combined the full spectral resolution data for the two spectral windows with very similar coverage, around 217 and 230 GHz. The TM0 data were given a weight of 0.2 relative to the combined 2016 data, in order to give sufficient weight to the shorter spacings to detect extended emission (with imageable on scales up to 15″) and avoid deconvolution artefacts. With a pixel size of 5 mas, we imaged the channels around the SiO and CO lines using multi-scale clean, and the primary beam correction was applied. The rms noise, given in Table 3, was measured in quiet channels.

3. SPHERE data acquisition and reduction

We observed EP Aqr using SPHERE/ZIMPOL (Beuzit et al. 2019; Schmid et al. 2018) in a Service Mode observing run (0102.D-0006(A)) in four filters, V, CntHα, BHα and Cnt748, on August 17th 2018, along with HD 1921 which was used as a point-spread function (PSF) calibrator star. The raw data was reduced using the ESO reflex data pipeline v0.24.0 (Freudling et al. 2013). The total observed intensity (I), polarised flux () and degree of linear polarisation () were computed from the Stokes parameter maps, with Q and U the Stokes parameters. The reduction products are shown in Fig. 1, which focuses on the innermost 0.3″ × 0.3″ of the total field of view. At larger length-scales, the images are dominated by noise, as shown in Fig. A.3. The observations using the V and Cnt748 filters are excluded from the analysis as the CCD is saturated for the science target. The PSF calibrator star is shown in Fig. A.2. A summary of the technical details of the observations can be found in Table 4.

Fig. 1. Overview of the SPHERE/ZIMPOL observations of the wind of EP Aqr. Top row: CntHα filter, bottom row: BHα filter. First column: total observed intensity, second column: polarised flux, third column: degree of linear polarisation (DoLP), fourth column: orientation of the polarisation vector. The beam size is shown in the bottom left of the leftmost panels. Both the intensity and polarised flux were normalized to the maximum intensity of each filter, thus preserving the degree of linear polarisation. The green circle in the centre represents the radius of the stellar photosphere as determined by the Stefan-Boltzmann relation. The white dashed line in the DoLP panels represents the orientation of the symmetry axis of the system as determined by Homan et al. (2018b).

4. Continuum data description

In this section we present the continuum data obtained by the ALMA continuum subtraction procedure, and the SPHERE filters that show no saturation of the CCD.

4.1. ALMA continuum emission

In Fig. 2 the ALMA continuum emission of the latest data combination is shown. The continuum peak is 16.32 mJy beam⁻¹, with a total flux within the 3σ_rms (σ_rms = 2.83 × 10⁻⁵ Jy beam⁻¹) contour of 23.5 mJy, with a radius of approximately 40 milliarcseconds (mas). The emission is predominantly centrally condensed, with a slight north-east to south-west elongation, and with 2 symmetric filaments around the 135° position angle (PA, measured north-to-east) axis appearing in the faintest emission. These are sidelobe remnants of the combination, probably present due to the variability of the target between epochs of observation, which would limit the continuum dynamic range in the combined image. The same features can be recognised in the dirty beam shown in Fig. A.1. The resolution and signal-to-noise are good enough to deconvolve with tha ALMA beam, permitting the analysis of the geometry of the central source. Fitting it with a 2D Gaussian results in a measured size of 14.0 ± 0.6 mas × 11.5 ± 0.7 mas, with a PA of 54 ± 8°, confirming its slight elongation.

Fig. 2. Continuum emission of the combined ALMA data on the EP Aqr system. Contours are drawn at 3, 6, 12, 24, 48 and 96 times the continuum rms noise value (σrms = 2.83 × 10−5 Jy beam−1). The offset position of the target with respect to the telescope pointing (ICRS 21:46:31.8475, −02:12:45.930) is due to transverse motion over the field of view. The ALMA beam is shown in the bottom right corner.

Table 2 shows the results of fitting a two-dimensional Gaussian component to the emission. The consistent increase in the Gaussian flux peak with larger beam areas suggests that ∼15% of the flux is on scales between 0.1 and 1. However, there are also indications of variability since the integrated flux density over the Gaussian at 150 mas resolution of the 2016 data is 17.89 mJy while in 2019, at 22 mas resolution, the integrated flux density is 21.35 mJy, amounting to a ∼20% increase. In order to allow for irregularities we also measured the total flux within the 3 σ_rms contour for the combined emission observed in 2016 (Homan et al. 2018b, ACA+TM2+TM1) and the new TM0 data. The 2016 (lower resolution) total flux was 19.3 mJy within 0.7″ from the peak, whilst the 2019 total flux is 24.6 mJy within 0.1″ from the peak. This is again a ∼20% increase, greater than the expected combined flux scale and stochiastic uncertainties of ∼10%.

Assuming that the star has a brightness of 4800 L_⊙ and a temperature of 3236 K (Winters et al. 2007), Homan et al. (2018b) estimated the continuum stellar flux contribution (by assuming black-body emission in a bandwidth 1.6 GHz around the rest frequency of 230.5 GHz) to be 22.7 mJy. It is known that the observable extent of an AGB star increases with wavelength (Reid & Menten 1997; Matthews et al. 2018). So, for EP Aqr, we can assume photospheric radius R_⋆ of 166 R_⊙, or about 6 mas at 135 pc, as estimated by Dumm & Schild (1998). The 230 GHz optically thick surface is unlikely to be more than ∼2 R_⋆ (e.g. Vlemmings et al. 2019) so even with uncertainties it has an upper limit of about 20 mas radius. Therefore, it seems likely that the radio atmosphere of the star in fact contributes only 15−18 mJy, hence additional, sporadic contributions to emission on larger scales are implied. These may be associated with companion interaction and/or a specific dust production events, and/or to variability of the radio atmosphere.

4.2. SPHERE observations

Figure 1 shows the results of the unsaturated SPHERE/ZIMPOL observation of the inner wind of EP Aqr at length-scales in the order of 0.1″. The intensity images, which have a resolution of ∼27 milliarcseconds, exhibit a strong centrally condensed signal at the stellar position. The linearly polarised flux map appears similar, though it possesses some instrumental/reduction artefacts that split the otherwise predominantly round central feature in half vertically. This central artefact is seen on all ZIMPOL polarimetric images when the star is not resolved, including the PSF (see Fig. A.2). It does not question the real signal. The maps reveal two patches with an elevated degree of linearly polarisation localised to the north-west and the south-east of the stellar position. The presence of instrumental/reduction artefacts in the polarised flux maps and in the maps of the degree of linear polarisation (DoLP) suggests that the central portion of the data cannot be trusted. These bright DoLP regions, which are quite similar in shape and amplitude, exhibit maxima around ∼15% of the total intensity, but are on average closer to ∼9% of the total captured light. Furthermore, they appear to trace sections of a ring centred on the stellar position, but decrease in intensity and width as they approach the north-west to south-east axis. At their widest point, they measure about 0.09″, and both exhibit some substructure in the form of slight local DoLP enhancements of the order of a few percent, and with typical angular sizes of 0.05″. Despite these small variations, the overall distribution appears predominently smooth, and not clumpy. The locations where the DoLP drops down to the detection limits appear to be well aligned with the orientation of the symmetry axis of the wind of EP Aqr as determined by Homan et al. (2018b), and which is overplotted on the DoLP plots in Fig. 1 as a white dashed line. We elaborate on this coincidence in Sect. 6.3 in the discussion.

5. ALMA molecular emission

In this section we describe the observational results for the molecular lines for which we successfully combined the 2016 data with the latest observations. The 2016 data did not contain the ¹³CO and SO₂ spectral lines observed in the 2019 campaign, so we do not include these in this article for the sake of consistency. All mentions of blue-shift and red-shift are with respect to the stellar velocity, which is assumed to be v_LSR = −34 km s⁻¹. For optimal appraisal of the fine detail in the data figures we recommend that all images presented in this paper are displayed on screen.

5.1. ¹²CO J = 3−2 emission

In this section we describe the nature of the ¹²CO J = 3−2 emission. In addition to the initial analyses performed by Homan et al. (2018b), the increased angular resolution allows for a better visualisation of the finer details in the observed emission patterns. Figure 3 shows the line channel maps of the combined data, and Figs. A.4 and A.5 show the spectral line and the moment 0 map of the combined data, respectively.

Previous work by Homan et al. (2018b), Hoai et al. (2019) and Tuan-Anh et al. (2019) has shown that the emission features in the high-resolution interferometric data on EP Aqr are caused by a face-on, dense equatorial density enhancement (EDE), symmetrically encased by two separate wind hemispheres, which are the sources of the central peak and the broad plateau in the spectral line, respectively. The dense EDE was shown to harbour a counter-clockwise oriented spiral shape, and its narrow spectral extent of ∼3 km s⁻¹ was interpreted as indication that its thickness was quite limited. In the centre of the EDE, they detected highly complex emission pattern at the base of the spiral. The channel maps, zoomed-in on the inner 4″ × 4″ of the ¹²CO emission, and spectrally confined to the central peak in the emission line (Fig. A.4), are shown in Fig. 4. At the stellar velocity a large, nearly complete ring of emission can be seen surrounding the continuum brightness peak (CPB) position. The CBP position is not in the centre of the circle, but rather to its west. To the east of the CBP position, a central arc can be seen extending from north to south, which grazes the CBP position, and appears to divide the circular feature into a western and an eastern half. Due to the curved nature of this central arc and the presence of the bright central emission centred on the CBP position, the low-emission parts appear as crescent-shaped, one larger one to the east of the CBP position, and a smaller one to the west. The correlation between the central, northwards curving arc, and the larger circular feature could be interpreted as the first winding of a clockwise spiral. However, this is opposite to the orientation of the larger spiral that Homan et al. (2018b) discovered in the EDE, and can therefore not straightforwardly be recognised as the base of the spiral in the EDE. The overall shape of the emission in EP Aqr’s inner wind was tentatively interpreted by Homan et al. (2018b) as the spectral manifestation of morphology associated with the wind-Roche-lobe-overflow scenario proposed by Mohamed & Podsiadlowski (2012), and is visualised in their manuscript’s Fig. 3, though their image must be flipped vertically to be comparable with the situation in EP Aqr’s wind.

With the combined datacubes, the emission distribution of the inner part of the EDE of EP Aqr is much better resolved, and can thus be investigated in much finer detail. The large, nearly circular arc covering the south, east, and part of the north of the central channel maps can now be seen to be composed of a smooth background covered in higher-emission clumps with sizes up to a few beams at most. It has a radius of approximately 1.25″, and its width, measured from the 18 × σ_rms contour, varies between 0.3″ and 0.45″. Furthermore, the eastern crescent appears to be quite perturbed, and actually not very crescent-shaped at all. Only the western crescent-shaped void, which is only visible in the blue-shifted channels, remains distinctly recognisable. Its axis of symmetry lies east-to-west, and on this axis the void is centred on a position located approximately 0.4″ to the west of the CBP. The void is present in all blue-shifted channels between 0 and −1.5 km s⁻¹, beyond which it blends with it surroundings. This velocity width corresponds exactly with the spectral half-width of the EDE.

The red-shifted channels (1 to 3 km s⁻¹ range) reveal an interesting small emission feature to the north and west of the CBP. At 3−2 km s⁻¹, this feature is clustered predominantly around 0.4″ to the west of the CBP position. As the channel velocity gradually progresses towards the stellar velocity, the emission cluster becomes more elongated, and traces an arc-like course over the north of the CBP position, and to the east. And though barely resolved, this arc remains clearly distinct from the brightness peak of the CO maps. Then, around 1 km s⁻¹, the arc blends with the additional emission features that start to dominate the inner wind of EP Aqr. Both the blue-shifted and red-shifted emission show some very interesting behaviour around 0.4″ to the west of the CBP position, which coincides with the location of the significant emission void in the SiO emission discovered by Homan et al. (2018b). We elaborate on the peculiarity of this position in Sect. 6.1 in the discussion.

Fig. 3. Continuum-subtracted channel maps of the 12CO emission. The labelled velocities have been corrected for v*. The central observed frequency is the average over the observation epochs. The contours are drawn at 3, 12, and 24 times the rms noise value outside the line (σrms = 6.86 × 10−4 Jy beam−1). Length scales are indicated in the bottom left panel. The maps are centred on the continuum peak position, which is indicated by the yellow star symbol. The beam size is illustrated in the bottom-right corner of the bottom-left panel (0.087″ × 0.079″).

Fig. 4. Identical channel maps to Fig. 3, but with contours drawn at 12, 18, 24, 30, 36, 42, and 60 times the rms noise value outside the line (σrms = 6.86 × 10−4 Jy beam−1, and zoomed in on the central 4″ × 4″ to highlight the complexity of the inner wind region.

The brightness peak of the CO emission can be seen to shift position with respect to the CBP position as a function of velocity, as is illustrated in Fig. 5. For the most red-shifted velocities (∼3 km s⁻¹), the CO brightness peak is located approximately 2 mas to the west of the CBP position. and it gradually shifts eastwards as velocities become more blue-shifted. Around stellar velocity it centres on the CBP position, and at higher blue-shifted speeds (∼ − 2 km s⁻¹) the CO brightness peak is distinctly positioned 5 mas to the east of the CBP position. This spectral progression is precisely in agreement with the findings of Homan et al. (2018b) who discovered what appeared to be gas in (rigid) rotation from their analysis of the SO₂ emission within a radius of 2″. In Fig. 6 the position shift of the CO brightness peak as a function of projected speed is captured in a moment 1 map. The map is indeed comparable to the moment 1 map of SO₂ revealed by Homan et al. (2018b), but due to the increased resolution it can be seen to be a factor 2 more compact. The white zero-velocity transition zone between the blue- and red-shifted emission, which has been interpreted as the orientation of the rotation axis of the gas in the plane of the sky, has a PA of ∼145°.

Fig. 5. Identical to Fig. 4, but zoomed in on the central 1″ × 1″ to highlight the position change of the CO emission brightness peak with respect to the continuum brightness peak position (star symbol).

Fig. 6. Moment 1 map of the inner CO emission for the −3 km s−1 to 3 km s−1 velocity range, and centred on the continuum brightness peak position. Only the signal above 24 times the noise rms is shown. The beam is shown in green in the bottom right corner. The black dashed line shows the approximate orientation of the rotation axis, with a PA of 145°.

These data also highlight the large amount of sub-structure that is present in the blue-shifted and red-shifted outflow hemispheres of this AGB stellar wind, which is manifested by the wings of the spectral line. The structure observed here consist of distinct filaments, which appear to trace segments of nearly complete circular shells in the maps. An example of such a nearly complete circle can be seen around velocities of 5.2−5.8 km s⁻¹ (see Fig. 3). And while not perfectly identical, such structures can clearly be recognised both in the red-shifted and the blue-shifted line wing, and have smaller radii at larger projected velocities.

As was first noted by Homan et al. (2018b), these arcs exhibit a slight deviation from concentricity, with the red-shifted emission pulling slightly to the south-east, and the blue-shifted emission to the north-west. The arcs in the blue- and red-shifted line wings also differ quite substantially in number and size. Though generally difficult to differentiate from each other due to the weak, clumpy emission, we have counted them by eye by summing up the generally circular features that appear in the emission. At the 3 times σ_rms level only one shell-like feature can be discerned. The other shells appear in higher signal-to-noise regions, and must therefore be appraised at the 6 times the σ_rms level. These contours reveal 6 distinct shells for the red-shifted line wing (see Fig. A.6). In totality, the outflow hemispheres thus exhibit 7 separate shells in the part of the outflow that is moving away. These arcs have typical widths that increase with radius. Assuming a constant velocity of 11 km s⁻¹ in this part of the wind, we can calculate the deprojected radii of the most obvious shells in the data. We found arc widths of approximately 1.5″ for the largest ring (at 3σ_rms level) located at a distance of 7″, of 0.4″ for rings with a radius of 3″, down to 0.3″ for rings with a radius of 1.2″. Though the arcs have an irregular outline at the enclosing contour level, the distance between the arcs appears to vary between 0.5″ and 0.8″, except for the two outermost ones, which appear somewhat more detached (> 2″) form the bulk of the emission. In the blue-shifted wings of the outflow the arcs are much more difficult to differentiate, with only 2 recognisable arc-like features at the 3σ_rms level, and no other clearly distinguishable arc-like features at higher rms noise levels. We comment on the peculiarity and potential origin of these arc-like features in Sect. 6.4.

5.2. ²⁸SiO J = 5−4 emission

The channel maps of the combined array data of the SiO emission can be found in Fig. 7. As discussed by Homan et al. (2018b), the emission is comparatively smooth and spherical compared to the CO emission. This is supported by the shapes of the emission line and moment 0 map, shown in Figs. A.7 and A.9, respectively. This is to be expected for the SiO molecule, since it is highly abundant in the inner winds of AGB stars and has a strong dipole moment. This implies that the stellar infrared radiation dominates the molecular energy level excitation close to the star. Hence, contrary to the CO molecule that has a very low electric dipole moment, any collisional effects, which act as a proxy for local density, only have a limited influence on the population of the SiO energy levels with the result that density variations are not expected to be prevalent in SiO emission. The SiO spectral line does not exhibit any unusual features besides a slight but important asymmetry, which is caused by decreased emission in the negative portion of velocity-space. This can be attributed to the emission void first detected in Homan et al. (2018b). The higher-resolution data now permit a better analysis of the properties of this void, for which we show zoomed-in channel maps in Fig. A.10.

Fig. 7. Continuum-subtracted channel maps of the SiO emission. The labelled velocities have been corrected for v*. The central observed frequency is the average over the observation epochs. The contours are drawn at 3, 6, 12, 24, 48, and 96 times the rms noise value outside the line (σrms = 1.5 × 10−3 Jy beam−1). Length scales are indicated in the bottom left panel. The maps are centred on the continuum peak position, which is indicated by the yellow star symbol. The beam size is illustrated in the bottom-right corner of the bottom-left panel (0.085″ × 0.079″).

The void, which is clearly resolved in the combined data, has an elongated shape, with a vertical (north-to-south) size of approximately 0.4″, and a horizontal size of approximately 0.2″. Its lowest emission point drops down well below the 3σ_rms detection limit, and is located approximately 0.4″ directly west of the continuum peak position, which at a distance of 135 parsec corresponds to a physical distance of ∼55 AU. The void appears precisely at the stellar velocity in the EDE, and remains a dominant spectral feature in the data down to velocities of −6 km s⁻¹, which are well inside the wind hemispheres. Beyond this velocity it blends with the surrounding low emission. The deep penetration of this void into the blue-shifted wind hemisphere implies that the void is in fact more like a hollow tube, which originates in the EDE (see Sect. 5.1, and Homan et al. 2018b; Hoai et al. 2019; Tuan-Anh et al. 2019). This tube extends nearly perpendicular to the plane of the EDE, penetrating the blue-shifted portion of the outflow while not penetrating the red-shifted part at all. In combination with the suggestion by Homan et al. (2018b) that this void may be the manifestation of a local environment created by a photodissociating radiation field caused by the presence of a companion, this asymmetric tube-like behaviour yet again provides evidence that the outflow of EP Aqr consists of two wind hemispheres, separated by a dense and flat EDE. With these improved data we can use the size of the void to estimate the local conditions needed to explain the observations. We defer further discussion to Sect. 6.1.

These new data also enable the identification of some structures that were not detected in the 2016 data. In particular, the red-shifted channels around 6 km s⁻¹ reveal a distinct arc or outline of a shell-segment to the east of the emission zone, at the level of the 3σ_rms contour. The shell-segment has a radius of approximately 2.6″, a thickness of approximately 0.6″, and extends most notably to the east and south-east. As mentioned above, detecting density variations so clearly is quite peculiar considering the dominant excitation mechanism for SiO is IR pumping. It is unclear why this particular arc emerges from the SiO emission.

For the sake of completion we show the position-velocity diagram of the SiO emission along the PA = 150° axis in Fig. A.8, constructed using slit width equal to the beam size. It possesses the typical shape associated with spherical emission in a predominantly radial outflow. Thanks to the high-resolution of the data, even the blue wing absorption over the stellar surface van be discerned. However, of main interest is that the data clearly shows an extension beyond what is typically associated with the wind speed of EP Aqr, as first discovered by Homan et al. (2018b). These features are asymmetric in velocity-space, we measure their maxima to extend to approximately −18 and +22 km s⁻¹. In later work, Tuan-Anh et al. (2019) interpreted these feautres as narrow high-velocity polar streams of gas that originate from within a radius of 25 AU from the stellar source. If so, then these polar streams could be generated by a tertiary companion that resides within 10 AU from the AGB star. We elaborate on this thought in Sect. 6.4.

In addition, even though the SiO emission distribution appears relatively smooth in the channel maps, the radial emission profile shown in Fig. 8 reveals that in fact the SiO emission does possess some low-amplitude substructure. In this figure, we show how the emission varies on an east-to-west oriented axis with a width of one beam, and the emission is averaged over the beam width. The axis runs through the CBP position, for the channel at stellar velocity. This is the channel through the bulk of the EDE surrounding EP Aqr, and is also the first channel to exhibit the void that appears in the right (western) side of the profile, highlighted by a green tag. The emission on this side is clearly weaker than its eastern counterpart, as can be seen in the difference between the blue and grey dashed profiles in Fig. 8. The eastern profile, which is presumably the least affected by the presence of the void, shows monotonically decreasing emission following what appears to be two distinct regimes. The emission in the inner portion of the wind can be approximated by a power-law with a slope of −0.77. Around 1.8″, this regime abruptly shifts to a much steeper radial slope, which can be approximated by a power-law with a slope of −2.6. This behaviour can be reasonably attributed to a radiative transfer effect, with the outer regions probing optically thin SiO emission, and the slope shift representing the transition towards optically thick SiO emission at smaller radii. Under this assumption, we can to first order state that if the SiO abundance and level excitation is approximately constant within this region, the density drops on average as r^−2.6 in the EDE of EP Aqr. Detailed three-dimensional radiative transfer modelling of discs around AGB stars has shown that such steep density profiles could be typical in the context of AGB EDE’s (Homan et al. 2017, 2018a). However, these objects were found to be very compact and to show signs of differential rotation. They are therefore quite different from the EDE in the wind of EP Aqr. Thus, we are not certain whether the same conclusions apply in this context. At this point, however, we cannot conceive any other mechanisms that would explain such a steep slope.

Fig. 8. Top: SiO emission profile of the central channel, extracted from east to west through the center. The intensity axis has units of Jy beam−1. The eastern (left) profile is overplotted as a grey dashed line over the western (right) profiles to highlight the asymmetry. The SiO void is highlighted by the green tag. The purple tags indicate the locations of the equidistant bumps in the inner portion of the SiO profile. Bottom: eastern profile in log-scale. The overplotted red dotted line exhibits a shift in regime between a slope of −0.77 for the inner 1.8″, and −2.6 for the outer part.

Finally, we note the presence of 4 small consecutive bumps that occur roughly every 0.45″ present a first order deviation from a smoothly decaying profile. These bumps do not have an equivalent counterpart in the CO emission (see Fig. A.11). We elaborate on what may be the origin these bumps in Sect. 6.4.

6. Discussion

6.1. Evidence for companion in the wind

The strongest evidence for the presence of a companion is the overall morphology of the wind. As has been extensively analysed by Homan et al. (2018b) and Hoai et al. (2019), the resolved emission patterns can be best explained by a nearly face-on, spiral harbouring equatorial density enhancement (EDE), encased on both sides by less dense wind hemispheres that could be considered a bi-conical outflow. For slowly rotating giants such as AGB stars, such global features are impossible to explain without invoking the presence of a companion. This companion will bend the otherwise radial wind streaming lines towards the orbital plane of the system, generating an equatorial overdense region with an additional narrow spiral generated from its gravity wake (Kim 2011; Liu et al. 2017). As such, the global morphology of the wind is shaped. There is, however, additional evidence in the observed emission patterns that a companion is present in the wind.

Further evidence can be found at the base of the spiral, which shows a distinct hole in the emission that is especially pronounced in SiO. This so-called emission void, which is located about 0.4″ to the west of the continuum brightness peak position, was first discovered by Homan et al. (2018b), and we describe its properties in Sect. 5.2. For a CO inter-spiral-arm distance of about 3″ − 4″, the dynamical crossing time is ∼200−250 years for an assumed wind speed of ∼10 km s⁻¹. This crossing time corresponds with the orbital period of an object located at approximately 50 AU, or 0.37″ from the central star. Considering the uncertainty on the wind speed in EP Aqr’s EDE, which so far has not been properly measured due to its nearly face-on orientation, this back-of-the-envelope calculation shows the striking similarity between the derived location of the tentative companion and the position of the SiO emission void.

Additionally, close inspection of the distribution of CO emission (whose distribution predominantly follows the wind density) in the inner wind of EP Aqr reveals a clear filament of gas originating a the location of the SiO void at ∼3 km s⁻¹, and curving northwards over the CBP position at lower speeds (see yellow arc in Fig. A.12). This feature was first detected by Homan et al. (2018b), but as rightly pointed out by Tuan-Anh et al. (2019), an error was made in the sign of the considered velocity-channel. This CO emission filament appears only in the red-shifted channels, as described in Sect. 5.1. Comparison with hydrodynamical models such as the ones by Mohamed & Podsiadlowski (2012) and Liu et al. (2017) reveals that, for companion masses comparable to the mass-losing star’s, this could be the formation of the spiral shock as a result of interaction with the companion. These hydrodynamical models predict that a large portion of the matter lost by the AGB star leaves the system in two distinct streams that flow along the Lagrange points L1 and L3. The stream along L3 is typically associated with reflex motion of the mass-losing AGB star around the centre of mass of the binary system (henceforth reflex-spiral), while the flow along L1 is caused by the direct attraction of the outflow material by the companion (Liu et al. 2017). The latter is launched into the CSE at augmented speeds after its closest interaction with the companion, and creates what is typically referred to as the wake-spiral. Considering the orientation of the larger-scale CO spiral discussed by Homan et al. (2018b), the companion must be rotating around the AGB star in a clockwise fashion. This implies that the position of the filament and its arc-like spectral evolution around the CBP position, while remaining in the red-shifted portion of the cube, indeed fit the dynamics of the gas that probes the formation of the wake-spiral. By extension of the argument, the filament that defines the southern edge of the western crescent could then be the reflex-spiral (see red arc in Fig. A.12). We present a schematic diagram of the proposed gas flow pattern in the inner wind of EP Aqr in Fig. 9. These features, which emerge from the higher-resolution data, strengthen Homan et al. (2018b)’s initial hypothesis that the detected shapes are associated with the hydrodynamical wind shaping caused by the presence of a companion at 0.4″ to the west of the CBP position.

Fig. 9. Schematic diagram of the wind-companion interaction zone. The two red dots represent the AGB star and the companion, separated by 0.4″. The pale outline in the central grey (less dense) zone represents the local density enhancement as a consequence of the interaction hydrodynamics. The blue (orange) dotted arrow line represents the flow along L1 (L3), which combine into the larger scale spiral observed in the EDE. The three arrows on top and bottom represent the (bi-conical) AGB outflow. We have indicated the dust ring observed with SPHERE/ZIMPOL in fuchsia.

There exists a spatial and spectral coincidence between the SiO emission void (see Fig. 7) and a similar, but spatially somewhat larger crescent-shaped void in CO (see Fig. A.12). As mentioned in Sect. 5.1, Homan et al. (2018b) suggested that the crescent shaped decreased emission zones could be the morphological manifestation of wind-Roche-lobe-overflow (WRLOF; Mohamed & Podsiadlowski 2012), which pointed them towards invoking the presence of a companion. The current data reveals that the location of the centre of the crescent shape in CO and the centre of the void in SiO are approximately equal (see Fig. A.12). This could suggest that (1) both the void in SiO and the crescent in CO may have the same origin, and/or (2) that the WRLOF scenario alone may not the best model to explain EP Aqr’s inner wind morphology. This could imply that the observed shapes may be the result of a different type of wind-companion hydrodynamical interaction, combined with emission-dampening physics such as cooling or photodissociation. It must be noted, however, that there are similar, albeit smaller, emission gaps between arcs at many PAs in the inner wind region in the −2 to 2 km s⁻¹ velocity range. Hence, the discussed crescent-shaped void could also simply be yet another such gap between the complex sequence of arcs and clumps, and its coincidence with the SiO gap merely coincidental.

6.2. The nature of the companion

In this section we will estimate the properties of the companion under the assumption that the observed emission hole in SiO is a consequence of photodissociation. Homan et al. (2018b) estimated that if the companion was a main-sequence star, then its non-detection in the continuum yielded an upper limit on its mass of 0.1 M_⊙. However, it is questionable whether such an M6 dwarf can generate the required pull on the gas to shape the wind to the extent that is observed (El Mellah et al. 2020). Furthermore, with a typical temperature of 2800 K and luminosity of ∼10⁻³ L_⊙ (Siess et al. 2000), its capability to sustain a flux of sufficient UV photons to explain the emission void observed in SiO and CO, either photospherically or via the creation of a hot accretion disc, is unlikely.

With the current high-resolution data the SiO hole is resolved, and hence we can estimate the size and temperature of the photodissociating blackbody by relating the predicted UV photodissociation time scale with the dynamical time scale associated with the wind of EP Aqr crossing the SiO emission void. Assuming the companion is in the centre of the void, half of the void’s width (∼27 AU/2) divided by the wind speed (∼11 km s⁻¹) corresponds to approximately 5.8 years.

In order to calculate the photodissociation rate of SiO, we have made use of the cross-sections made available via the Leiden database for photodissociation¹ which shows that it is photodissociated through five absorption lines at wavelengths between 100 and 140 nm (Heays et al. 2017). Using this data, we have calculated the UV photodissociation rate β(r_sc) of SiO for 21 evenly spaced black body temperatures between 5000 K and 15 000 K, at an arbitrary scaling radius of r_sc which is assumed to be 50 times the radius of the black body. The results of this calculation can be found in Table A.1. The actual photodissociation rate at a distance r from the black body can then be calculated by accounting for geometric dilution and extinction through dust, captured in the following expression

(1)

where κ_d is the dust opacity and ρ_d the dust density at r. These will be assumed to both be constant, and have values of κ_d = 7 cm² g⁻¹ (Demyk et al. 2017) and ρ_d equal to the gas density ρ_g times an assumed dust-to-gas mass ratio of 200. For this exercise, ρ_g is calculated using mass conservation. However, because we currently have no proper constraints on the density contrast between the EDE and the outflow, we use the ratio of the integrated fluxes between the CO line peak and the plateau as proxy for the density contrast between the EDE and the outflow hemispheres. This is an absolute lower limit on the actual value, but considering the unknowns this is currently the best that can be done. The photodissociation time scale is then given by β(r)⁻¹. Setting this photodissociation time scale equal to the dynamical void crossing time scale produces a space of solutions for black body temperature and radius combinations. The solution curve is visualised in Fig. 10.

Fig. 10. Blackbody temperature and radius combinations that would open-up a hole of the observed size by photodissociation of SiO. The blue curve is calculated for a wind speed of 11 km s−1. The translucent blue zone represents the space of solutions for a maximum deviation of 90% on the assumed wind speed. The red hatched zone represents the white dwarf radii found in the literature.

The back-of-the-envelope nature of the above calculation implies that the results in Fig. 10 should be interpreted within an order of magnitude margin. To provide an awareness on the uncertainties associated with this solution curve, we illustrate the effect of an incorrect assessment of the wind velocity, and hence the dynamic crossing time scale at the location of the companion, which was assumed to be 11 km s⁻¹. The shaded zone in Fig. 10 represents the shift of the solution curve under a 90% increase and decrease in the assumed wind speed. This is an indicative uncertainty, and by no means does it fully represent the actual error in the curve of solutions. Nevertheless, the generated space of solutions provides an interesting spectrum of companion candidates.

In the low companion temperature regime, and assuming the uncertainties on the black-body radii are substantial, the space of solutions could be consistent with a solar-type star. In the high companion temperature regime only an extremely compact object would be able to explain the observed size of the SiO hole. However, such objects are known to be highly magnetic, and much hotter than the associated temperature on the curve. Hence, such objects can be ruled out. The current analysis indicates that the space of solutions is distinctly consistent with white dwarf (WD) radii (Tremblay et al. 2017; Parsons et al. 2017; Joyce et al. 2018; Romero et al. 2019) in the 8000−6000 K temperature range. The evolution time of any more massive WD progenitor (few Gyr) to reach the WD state, and the typical WD cooling time (few Gyr) to reach such low temperatures is indeed consistent with typical evolutionary lifetime of 1 M_⊙ star. We used the WD cooling curves given by Renedo et al. (2010) to obtain the cooling time for different WD masses. The results are given in Fig. 11. The hatched zone represents the cooling time versus WD mass combinations that are consistent with the observations. This result can be used to estimate the mass-range of the tentative WD in order to be able to explain to observations. To this end, we use the combined main-sequence (MS) mass-luminosity relation and MS nuclear time scale relation to estimate the typical MS lifetime of stars with different masses. For a 1.7 M_⊙ star such as EP Aqr, this yields a lifetime of 2.65 Gyr, which is the hard upper limit on the WD cooling time (represented by the black vertical line). The actual maximal WD cooling time would then be the lifetime of EP Aqr, minus the lifetime of the WD progenitor. This maximal WD cooling time is represented by the vertical dashed lines: each line representing the WD cooling upper limit for WDs with the respective colour-coded progenitor masses. Finally, combining this with the WD initial-final-mass relation (Cummings et al. 2018), which are plotted as the horizontal dotted lines, we find that the space of solutions is only consistent for WD masses between ∼0.65 and 0.8 M_⊙, which would have had progenitor masses between 2.1 and 3.3 M_⊙.

Fig. 11. WD mass as a function of WD cooling time. The hatched region between the two broad grey curves represents the range of WD temperatures as derived from Fig. 10. The black vertical line is the main sequence lifetime of a 1.7 M⊙ star, and therefore the WD cooling time upper limit. The zone to the left of the vertical dashed lines are the maximal cooling times for WD’s with respective colour-coded progenitor masses. The horizontal dotted lines represent the WD inital-final-mass relation for the respective colour-coded progenitor masses (5.0 M⊙ is outside of the y-range of the plot). The blue and green dot represent WD masses that would have had enough time to cool to reach a temperature between 6000 K and 8000 K.

Hot WDs are typically detected through their strong UV signatures. Montez et al. (2017) show a GALEX spectrum of EP Aqr, which shows no definitive WD signatures inbetween the dominant photospheric emission. This is not surprising considering the strong observational bias associated with observing even much closer, nearly naked WDs, which coincidentally also creates the misconception that WDs are always very hot. Indeed, studies of volume-limited observational samples do confirm that WD temperatures of 5000−6000 K are in fact the norm (Hollands et al. 2018), even if they are never detected. Hence, it is uncertain whether any observational campaign in the near future, targetting the detection of EP Aqr’s companion will be successful.

For the sake of completion, we note that the space of solutions does not appear to be consistent with the presence of a hot accretion disc around the companion. Such objects have predicted sizes of up to hundreds of solar radii and temperatures up to ten thousand degrees Kelvin (Saladino et al. 2018). However, following the trend in Fig. 10 an object of such size would have to have temperatures lower than 1000 K in order to explain the data. Additionally, such discs are not expected to arise for companions with such large orbital radii (55 AU), since their accretion mechanism would lie well in the Bondi-Hoyle-Lyttleton regime (Hoyle & Lyttleton 1939; Bondi & Hoyle 1944; Bondi 1952), isotropic accretion from wind material within the Hill radius, and would hence be rather inefficient. As discussed by Homan et al. (2018b) a companion mass upper limit of 0.1 M_⊙ was derived to explain the non-detection of a main-sequence companion in the continuum. This points the conclusion of our analysis on the nature of the void-inducing companion towards one of two scenarios: Either the extinction along the line-of-sight in the analysis by Homan et al. (2018b) was severely underestimated, in which case the companion could be a brighter main-sequence solar-like star, or the extinction was well estimated, in which case the companion could likely be a white dwarf.

As a final note, we emphasise that the above discussion is based on the distinct assumption that the companion is hot. The possibility of course exists that this assumption is false, and that the companion is of a very different nature. For example, the emission void could also be an extremely cool, dusty and opaque local environment. Such an object would suppresses the radiative excitation of molecules such as SiO locally and, via shadowing of the stellar radiation field, in a tube-like region deep into the radially trailing wind hemisphere. Such an enormous cool, dusty and opaque dust clump has been detected in the inner wind of the red supergiant VY CMa (O’Gorman et al. 2015), and its origins remain mysterious. For EP Aqr too it remains unclear how such a dust clump would accumulate around the companion. AGB stars are known to occasionally undergo eruptions of a convective cell on their stellar surface, which could generate local conditions that are ideal for the formation of copious amounts of dust. However, it seems unlikely that such a clump of dust would have survived the journey to the location of the companion considering the intricate hydrodynamical streaming lines caused by the wind-companion interaction. The local growth of the clump could have occurred via simple Bondi-Holyle-Lyttleton accretion, but it is not clear what mechanism would have cooled the compressing gas sufficiently to explain the observation.

6.3. Shape of the polarised light from SPHERE

The patterns that emerge in the DoLP maps in the SPHERE data, shown in in Fig. 1, can be explained by three scenarios:

Scenario 1. It is possible that the observed features are generated by directionally confined mass-ejections, caused by the occasional eruptions of the convective cells that dominate the stellar surface. These eruptions are known to generate shocks throughout the loosely bound outer layers of the star that create conditions that are ideal for dust formation (Freytag et al. 2017). The symmetry of the signal in the DoLP maps would hence be a mere coincidence.

Scenario 2. This same symmetry, however, suggests that the mechanism that determines the dust distribution throughout the inner wind may be less random than the above mentioned scenario suggests. If the dynamical time scale of dust propagation through the inner wind is smaller than the pulsation time scale that fuels the production of dust, and if the dust-producing shocks are randomly distributed in time and space, then the dust may accumulate in a shell surrounding the AGB star. Such features have been invoked in order to explain SED dust excesses (Khouri et al. 2015; Van de Sande et al. 2018). Accounting for the degree of linear polarisation as a function of scattering angle given by the Rayleigh scattering model, such a shell would manifest as a ring in the DoLP maps. This scenario qualitatively explains the data rather well. We note that the DoLP signal becomes dominated by noise beyond ∼150 mas, implying that the width of the DoLP signal is by no means a reliable representation for the width of the dust shell. Furthermore, the gaps in the DoLP maps would imply that the shell is incomplete, and that some mechanism has prevented efficient dust formation along the north-west to south-east axis.

Scenario 3. However, though not perfect, the striking coincidence between the orientation of the symmetry axis of the system (as previously determined by Homan et al. 2018b) and the patches of decreased DoLP suggests yet a third scenario. Scattering theory predicts that forwards and backwards scattering produces no linearly polarised light, implying that most of the detected scattered light from the pole-on, nearest outflow hemisphere of EP Aqr only contributes a negligible amount to the DoLP signal. The same can be said for the farthest outflow hemisphere. Hence, most of the detected DoLP signal probably originates from the nearly face-on equatorial density enhancement (EDE), in which the dust would form a ring around the AGB star. For a perfect face-on orientation this ring of dust would manifest as a ring in the DoLP maps. However, EP Aqr having a known non-negligible system inclination between 4° and 18° (Homan et al. 2018b) the DoLP is expected to decrease along the system’s symmetry axis, exactly as observed in the SPHERE data.

Because most physical and geometrical details of EP Aqr’s wind are unknown, and because there exist strong degeneracies between some of these properties (such as e.g. dust density and thickness of the dust ring), the current data can by no means be used to provide a quantitative differentiation between the above mentioned scenarios. Nevertheless, considering the morphology of EP Aqr’s wind the third scenario provides a satisfying, consistent, and qualitative explanation for the detected shape in the DoLP signal. Furthermore, under the assumption of its validity, it provides additional constraints on the system properties. For one, attributing the decreased DoLP to inclination implies that the EDE morphology reaches all the way down to the inner wind regions that are located approximately 0.1″ away from the stellar surface. The reason for this is not entirely clear if we assume that the tentative companion that is responsible for the SiO void (located at a physical distance of 55 AU) is also responsible for the presence of the EDE. If massive enough, then this companion could in principle bend the otherwise radial streaming lines of the equatorial portion of the wind even more towards the orbital plane, but it is unclear to what degree this can happen for a companion at such a substantial distance, and how important this effect can be near the mass-losing star. Whatever the equatorial wind flattening mechanism may be, the SPHERE data suggest that it may already be operational at distances as close as 0.1″ from the stellar surface. Still under the assumption that the dust that contributes to the observed DoLP maps is located in an EDE surrounding the AGB star, we can use the DoLP signal to provide additional constraints on the symmetry axis of the system and the inclination angle at these length-scales.

In Fig. 12 we show how the DoLP varies as a function of PA in the plane of the sky. The (blue) curve was constructed by averaging the DoLP over a sequence of radial axes originating from the intensity brightness peak position, which we assume to be the centre of the AGB star. To confine the analysis to the region in the field of view which is not dominated by noise (see Fig. A.3) we have set an outer radius determined by the shortest segment (over all radial axes) for which the signal remained above 2 times the noise σ_rms, plotted in grey in Fig. 12. The blue curve has a periodic aspect: It exhibits two distinct troughs with a ⟨DoLP⟩ of around 0.03 and two distinct peaks with a ⟨DoLP⟩ around 0.06. Each of these pairs are approximately separated by an angle of π radians, they are in geometrically polar opposition. We can estimate the PA of the symmetry axis of the system using the angles where the ⟨DoLP⟩ troughs are located, which is at a PA of 138.5°. Considering the crudeness of the measurement, this is in remarkable agreement with the result of Homan et al. (2018b) who used the offset between the centres of the blue-shifted and red-shifted wind hemispheres to measure the PA for the symmetry axis of the wind of EP Aqr to be 150°, and with the PA of the rotation axis of EP Aqr’s the SO₂ emission (Homan et al. 2018b), and of the CO emission shown in Fig. 6.

Fig. 12. Radially averaged DoLP as a function of angle, shown in blue. The 2σrms noise level is shown in grey. The dashed line in the DoLP image in the top right corner is a visual representation of the outer radius beyond which the DoLP does not contribute to the profile.

The peak-to-trough ⟨DoLP⟩ ratio provides a crude estimate for the inclination angle of the system. The peaks are quite stable around a ⟨DoLP⟩ of 0.06, while the troughs vary quite substantially between 0.035 and 0.024. This translates to a decrease in peak ⟨DoLP⟩ between 60% and 40%, respectively. Assuming the dust is confined to a cylindrical and flat EDE, then back-of-the-envelope Rayleigh scattering theory predicts scattering angles between ∼35° and ∼45°. This inclination angle is much higher than the previously determined inclination angle of the system of 4°–18° derived from the CO analysis by Homan et al. (2018b). This difference could point to a fundamental morphological difference between the inner and the outer EDE, but in all probability it can simply be attributed to the crudeness of the inclination estimate. A more reliable estimate for the inclination of the dust disc will be derived from the modelling of the data that will appear in a future work.

6.4. Tentative tertiary companion

The structure observed in the line wings of the CO emission are evidence of an active shaping mechanism which not only affects the global shape of the outflow by forming a spiral-harbouring EDE, but also appears to induce smaller-scale perturbations that manifest as arcs in the outflow hemispheres. In AGB binaries, such arcs are typically attributed to the reflex motion of the AGB star around the centre of mass of the binary system, which produces a spiral shock that due to the geometry of the system can extend to far above and below the orbital plane. For a constant outflow velocity of 11 km s⁻¹, the typical period of the wobbling that would be necessary to generate periodic features at 0.5″ to 0.8″ intervals would be between 30 and 45 years. These could never have been generated by the orbital motions of the AGB binary due to the tentative companion that is thought to be responsible for the SiO emission void.

Alternatively, the presence of an actively accreting companion in close orbit may be able to explain the concentric arcs. Tuan-Anh et al. (2019) proposed that the augmented speeds observed in the inner wind of EP Aqr (see Fig. A.8) could be the manifestation of narrow polar jets that are launched less than ∼25 AU away from the central star. They correctly conclude that the tentative companion responsible for the SiO void cannot be the engine that launches these collimated outflows, since it is located too far away from the AGB star. However, a closer companion at a distance of ∼10 AU could be a viable candidate to explain these features. Such an object would have an orbital period of ∼25 years, which lies quite close to the dynamical crossing time between the observed arcs. This estimate was made under the assumption that the wind speed within 2″ has already reached its terminal value. Any decrease in wind speed would increase the period estimate, bringing it even closer to the dynamical crossing time between two consecutive shells. Collimated streams of gas, ejected nearly perpendicular to the orbital plane, would also form nearly-concentric shell-like features in the wind hemispheres as time progresses due to the mostly radial radiation pressure.

In addition, the radial SiO emission profile shown in Fig. 8 exhibits a sequence of 4 consecutive bumps to the east, with a seemingly constant separation of 0.45″. It is by no means clear what the origin of these bumps is. However, dynamical arguments can be used to quickly rule out pulsational shocks or the companion in the void. Hence, we propose that these bumps could also be the signature of spiral shock waves created by the interaction of the wind with a much closer companion. This is supported by the fact that the bumps appear to shift outwards when the PA of the profile is progressively increased, as shown in Fig. 13. Assuming the wind has a speed of about 11 km s⁻¹ in the inner wind, the distance between the consecutive bumps would correspond with a period of ∼25 years.

Fig. 13. Stacked plot of three different SiO emission profiles taken at different PAs (north to east), vertically shifted by an arbitrary factor to highlight the progression of the bumps. The top profile represents the eastern wing of the profile shown in Fig. 8. The middle profile is taken at a PA of 45°, and the bottom profile at a PA of 90°. The dots highlight the bumps to show their outwards shift.

Approximately 10% of solar-type main-sequence stars are found in triple/quadruple systems (Moe & Di Stefano 2017). Furthermore, for the currently determined orbital separation of ∼10 AU (0.07″) the likelihood for stars with an inital mass greater than 1.5 M_⊙ (Fulton et al. 2019) (such as EP Aqr Dumm & Schild 1998) of having a > 1 M_Jup planetary companion is higher than 40%. The presence of such a closer companion could also provide the necessary angular momentum to explain the gas in rotation in the innermost regions of the wind, as was revealed by the moment 1 map of CO shown in Fig. 6, and a similar map of the SO₂ emission shown by Homan et al. (2018b). And finally, it may also provide a plausible explanation for the reason why the SPHERE data appear to indicate that the EDE of EP Aqr reaches all the way down to 0.1″ distances from the central star, while the most obvious companion candidate resides at 0.4″.

Triple systems are known to be notoriously unstable. However, a multitude of binary systems have been discovered that have planets in stable S-type² orbits (Thebault & Haghighipour 2015). It has been shown that a companion star in a wide binary (> 50 AU) has only a little impact on the formation or distribution of S-type planets, which has been shown to follow the same trend as for single stars (Roell et al. 2012). It is not clear whether such a system remains stable throughout the evolution of the binary, nor whether there exists a set of progenitor orbital configurations that would result in the configuration we observe today, especially considering that EP Aqr has most likely experienced the complete evolution of the secondary from a main-sequence star to the WD phase (see Sect. 6.2).

7. Summary

In this manuscript we present the latest high-resolution (0.025″) band 6 ALMA data investigating the wind of the M-type AGB star EP Aquarii, supported by additional SPHERE/ZIMPOL observations of the inner wind though the V, CntHα, BHα, and Cnt748 filters. These are follow-up observations of previous investigations, which have shown the wind of EP Aquarii to harbour a multitude of morphological complexities: a face-on equatorial density enhancement (EDE) that contains a spiral, a pole-on bi-conical outflow, a highly complex inner wind morphology that resembles the wind-Roche-lobe-overflow models in the literature, a distinct emission void in SiO approximately 0.45″ to the west of the continuum brightness peak (CBP) position, a signal that is indicative of gas in rotation within the inner 0.2″ of the wind, and thin polar high-velocity gas streams. All these features hinted towards the presence of a companion in EP Aqr’s wind. The new and better resolved data strongly support this hypothesis: The now resolved SiO emission void has been measured to reside precisely 0.4″ to the west of the CBP position. This position also appears to be at the centre of the complex emission patters that are observed in CO that strongly resemble the typical hydrodynamical gas flow patterns associated with the interaction of a mass-losing AGB star with a companion. In particular, the flow along the L1 and L3 Lagrange points can be identified, before they coalesce into the spiral that manifests at larger length-scales. These flow patterns are consistent with a companion located at 0.4″ to the west of the CBP position. In addition, the central brightness peak of the CO emission also exhibits the signal of rotation discovered in earlier work, but now shown to be confined to the inner 0.1″. Assuming the emission void that emerges in the SiO maps is caused by photodissociation, we tentatively calculated the properties of a black body that could provide sufficient UV photons to explain the size of the void. The space of solutions is particularly consistent with a white dwarf, with a mass between 0.65 and 0.8 M_⊙, though a solar-like companion could not be definitively excluded. The SiO radial profile exhibits four periodic consecutive bumps that cannot be explained by the presence of the above mentioned companion, nor stellar surface pulsations. We tentatively propose that these could be caused by a tertiary companion that that resides as close as 10 AU from the AGB star. The SPHERE data reveal a ring of polarised light around the AGB star, with a radius of ∼0.1″. Two symmetrically situated reductions in the degree of linear polarisation signal suggest that the dust responsible for the signal is confined to an inclined plane, with an undetermined thickness, and an inclination of around 45°. These are consistent with the previously determined wind morphology, and show that the EDE in the wind persists to radii below the position of the primary companion candidate, supporting the suggestion that a tertiary companion could already be moulding the wind at smaller radii away from the AGB star.

Acknowledgments

We would like to express our sincere thanks to Dr. Mark Hollands (Warwick Uni., UK) for his valuable contributions on white dwarf properties. We would also like to thank Dr. E. Lagadec and Dr. J. Milli for their contributions to the processing and troubleshooting of the SPHERE data. W.H., L.D. acknowledge support from the ERC consolidator grant 646758 AEROSOL. E.C. acknowledges funding from the Ku Leuven C1 grant MAESTRO C16/171007. TJM is grateful to the STFC for support via grant number ST/P000321/1. This project has received funding from the European Union’s Horizon 2020 research and innovation program under the Marie Skłodowska-Curie Grant agreement No. 665501 with the research Foundation Flanders (FWO) ([PEGASUS]² Marie Curie fellowship 12U2717N awarded to M.M.). This paper makes use of the following ALMA data: ADS/JAO.ALMA#2011.0.00277.S and ADS/JAO.ALMA#2018.1.00750.S. ALMA is a partnership of ESO (representing its member states), NSF (USA), and NINS (Japan), together with NRC (Canada), NSC and ASIAA (Taiwan), and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NR.A.O, and NAOJ. Based on observations collected at the European Southern Observatory under ESO programme(s) 0102.D-0006(A). This work has made use of the SPHERE Data Centre, jointly operated by OSUG/IPAG (Grenoble), PYTHEAS/LAM/CeSAM (Marseille), OCA/Lagrange (Nice), Observatoire de Paris/LESIA (Paris), and Observatoire de Lyon/CRAL, and supported by a grant from Labex OSUG@2020 (Investissements d’avenir – ANR10 LABX56).

References

Appendix A: Additional material

Fig. A.1. ALMA dirty beam of the combined data.

Fig. A.2. Same as Fig. 1 but for the calibrator star HD 1921.

Fig. A.3. Degree of polarised flux observed through the BHa filter, and associated rms noise map on the right. It is clearly visible that the field of view becomes dominated by noise beyond the features described in Fig. 1.

Fig. A.4. Spectral line of the CO emission for the combined data. Different colours correspond to different circular extraction aperture diameters, as indicated by the plot legend.

Fig. A.5. Moment 0 map of the CO emission for the combined data. Contours are plotted et 3, 6, 12, 24, 48, and 96 times the noise rms (σrms = 6.86 × 10−4 Jy beam−1) of the signal-free edge of the map. The ALMA beam size is shown in the bottom right corner.

Fig. A.6. CO channel in the red wing exhibiting the circular arc-like features at the 6 times the rms noise level.

Fig. A.7. Spectral line of the SiO emission for the combined data. Different colours correspond to different circular extraction aperture diameters, as indicated by the plot legend.

Fig. A.8. Narrow-slit SiO position-velocity diagram through the CBP position. The red dashed lines are drawn at 12 km s−1. The green dashed lines show the maximum speed of the compact SiO emission, at −18 and +22 km s−1.

Fig. A.9. Moment 0 map of the SiO emission for the combined data. Contours are plotted et 3, 6, 12, 24, 48, and 96 times the noise rms (σrms = 1.5 × 10−3 Jy beam−1) of the signal-free edge of the map. The ALMA beam size is shown in the bottom right corner.

Fig. A.10. Identical SiO channel maps to Fig. 7, but zoomed in on the central 2″ × 2″ to highlight the resolved void west of the peak brightness position.

Fig. A.11. CO emission profile of the central channel, extracted from east to west through the center. The coloured tags indicate the location of the bumps (purple) and the void (green) highlighted in Fig. 8.

Fig. A.12. Identical channel maps to Fig. 4, slightly more zoomed in. The white cross on the position of the centre of the SiO emission void to highlight the coincidence between its position and the features in the CO emission. The yellow arc to the north in the red-shifted channels highlights the emission that is associated with the launching of the wake-spiral. The red arc to the south in the blue-shifted challes highlights the emission associated with the reflex-spiral.

All Tables

All Figures

Fig. 1. Overview of the SPHERE/ZIMPOL observations of the wind of EP Aqr. Top row: CntHα filter, bottom row: BHα filter. First column: total observed intensity, second column: polarised flux, third column: degree of linear polarisation (DoLP), fourth column: orientation of the polarisation vector. The beam size is shown in the bottom left of the leftmost panels. Both the intensity and polarised flux were normalized to the maximum intensity of each filter, thus preserving the degree of linear polarisation. The green circle in the centre represents the radius of the stellar photosphere as determined by the Stefan-Boltzmann relation. The white dashed line in the DoLP panels represents the orientation of the symmetry axis of the system as determined by Homan et al. (2018b).

In the text

Fig. 2. Continuum emission of the combined ALMA data on the EP Aqr system. Contours are drawn at 3, 6, 12, 24, 48 and 96 times the continuum rms noise value (σrms = 2.83 × 10−5 Jy beam−1). The offset position of the target with respect to the telescope pointing (ICRS 21:46:31.8475, −02:12:45.930) is due to transverse motion over the field of view. The ALMA beam is shown in the bottom right corner.

In the text

Fig. 3. Continuum-subtracted channel maps of the 12CO emission. The labelled velocities have been corrected for v*. The central observed frequency is the average over the observation epochs. The contours are drawn at 3, 12, and 24 times the rms noise value outside the line (σrms = 6.86 × 10−4 Jy beam−1). Length scales are indicated in the bottom left panel. The maps are centred on the continuum peak position, which is indicated by the yellow star symbol. The beam size is illustrated in the bottom-right corner of the bottom-left panel (0.087″ × 0.079″).

In the text

	Fig. 4. Identical channel maps to Fig. 3, but with contours drawn at 12, 18, 24, 30, 36, 42, and 60 times the rms noise value outside the line (σrms = 6.86 × 10−4 Jy beam−1, and zoomed in on the central 4″ × 4″ to highlight the complexity of the inner wind region.
In the text

	Fig. 5. Identical to Fig. 4, but zoomed in on the central 1″ × 1″ to highlight the position change of the CO emission brightness peak with respect to the continuum brightness peak position (star symbol).
In the text

	Fig. 6. Moment 1 map of the inner CO emission for the −3 km s−1 to 3 km s−1 velocity range, and centred on the continuum brightness peak position. Only the signal above 24 times the noise rms is shown. The beam is shown in green in the bottom right corner. The black dashed line shows the approximate orientation of the rotation axis, with a PA of 145°.
In the text

Fig. 7. Continuum-subtracted channel maps of the SiO emission. The labelled velocities have been corrected for v*. The central observed frequency is the average over the observation epochs. The contours are drawn at 3, 6, 12, 24, 48, and 96 times the rms noise value outside the line (σrms = 1.5 × 10−3 Jy beam−1). Length scales are indicated in the bottom left panel. The maps are centred on the continuum peak position, which is indicated by the yellow star symbol. The beam size is illustrated in the bottom-right corner of the bottom-left panel (0.085″ × 0.079″).

In the text

Fig. 8. Top: SiO emission profile of the central channel, extracted from east to west through the center. The intensity axis has units of Jy beam−1. The eastern (left) profile is overplotted as a grey dashed line over the western (right) profiles to highlight the asymmetry. The SiO void is highlighted by the green tag. The purple tags indicate the locations of the equidistant bumps in the inner portion of the SiO profile. Bottom: eastern profile in log-scale. The overplotted red dotted line exhibits a shift in regime between a slope of −0.77 for the inner 1.8″, and −2.6 for the outer part.

In the text

Fig. 9. Schematic diagram of the wind-companion interaction zone. The two red dots represent the AGB star and the companion, separated by 0.4″. The pale outline in the central grey (less dense) zone represents the local density enhancement as a consequence of the interaction hydrodynamics. The blue (orange) dotted arrow line represents the flow along L1 (L3), which combine into the larger scale spiral observed in the EDE. The three arrows on top and bottom represent the (bi-conical) AGB outflow. We have indicated the dust ring observed with SPHERE/ZIMPOL in fuchsia.

In the text

	Fig. 10. Blackbody temperature and radius combinations that would open-up a hole of the observed size by photodissociation of SiO. The blue curve is calculated for a wind speed of 11 km s−1. The translucent blue zone represents the space of solutions for a maximum deviation of 90% on the assumed wind speed. The red hatched zone represents the white dwarf radii found in the literature.
In the text

Fig. 11. WD mass as a function of WD cooling time. The hatched region between the two broad grey curves represents the range of WD temperatures as derived from Fig. 10. The black vertical line is the main sequence lifetime of a 1.7 M⊙ star, and therefore the WD cooling time upper limit. The zone to the left of the vertical dashed lines are the maximal cooling times for WD’s with respective colour-coded progenitor masses. The horizontal dotted lines represent the WD inital-final-mass relation for the respective colour-coded progenitor masses (5.0 M⊙ is outside of the y-range of the plot). The blue and green dot represent WD masses that would have had enough time to cool to reach a temperature between 6000 K and 8000 K.

In the text

	Fig. 12. Radially averaged DoLP as a function of angle, shown in blue. The 2σrms noise level is shown in grey. The dashed line in the DoLP image in the top right corner is a visual representation of the outer radius beyond which the DoLP does not contribute to the profile.
In the text

Fig. 13. Stacked plot of three different SiO emission profiles taken at different PAs (north to east), vertically shifted by an arbitrary factor to highlight the progression of the bumps. The top profile represents the eastern wing of the profile shown in Fig. 8. The middle profile is taken at a PA of 45°, and the bottom profile at a PA of 90°. The dots highlight the bumps to show their outwards shift.

In the text

	Fig. A.1. ALMA dirty beam of the combined data.
In the text

	Fig. A.2. Same as Fig. 1 but for the calibrator star HD 1921.
In the text

	Fig. A.3. Degree of polarised flux observed through the BHa filter, and associated rms noise map on the right. It is clearly visible that the field of view becomes dominated by noise beyond the features described in Fig. 1.
In the text

	Fig. A.4. Spectral line of the CO emission for the combined data. Different colours correspond to different circular extraction aperture diameters, as indicated by the plot legend.
In the text

	Fig. A.5. Moment 0 map of the CO emission for the combined data. Contours are plotted et 3, 6, 12, 24, 48, and 96 times the noise rms (σrms = 6.86 × 10−4 Jy beam−1) of the signal-free edge of the map. The ALMA beam size is shown in the bottom right corner.
In the text

	Fig. A.6. CO channel in the red wing exhibiting the circular arc-like features at the 6 times the rms noise level.
In the text

	Fig. A.7. Spectral line of the SiO emission for the combined data. Different colours correspond to different circular extraction aperture diameters, as indicated by the plot legend.
In the text

	Fig. A.8. Narrow-slit SiO position-velocity diagram through the CBP position. The red dashed lines are drawn at 12 km s−1. The green dashed lines show the maximum speed of the compact SiO emission, at −18 and +22 km s−1.
In the text

	Fig. A.9. Moment 0 map of the SiO emission for the combined data. Contours are plotted et 3, 6, 12, 24, 48, and 96 times the noise rms (σrms = 1.5 × 10−3 Jy beam−1) of the signal-free edge of the map. The ALMA beam size is shown in the bottom right corner.
In the text

	Fig. A.10. Identical SiO channel maps to Fig. 7, but zoomed in on the central 2″ × 2″ to highlight the resolved void west of the peak brightness position.
In the text

	Fig. A.11. CO emission profile of the central channel, extracted from east to west through the center. The coloured tags indicate the location of the bumps (purple) and the void (green) highlighted in Fig. 8.
In the text

Fig. A.12. Identical channel maps to Fig. 4, slightly more zoomed in. The white cross on the position of the centre of the SiO emission void to highlight the coincidence between its position and the features in the CO emission. The yellow arc to the north in the red-shifted channels highlights the emission that is associated with the launching of the wake-spiral. The red arc to the south in the blue-shifted challes highlights the emission associated with the reflex-spiral.

In the text