© ESO 2018

1 Introduction

At the end of their lives, low- and intermediate-mass stars, M ≲ 8 M_⊙, reach the asymptotic giant branch (AGB) and experience strong mass loss. The mass-loss process is thought to be the result of the action of stellar pulsations and strong convective motions that increase the density scale height of the atmosphere and allow the formation of dust, followed by radiation pressure acting on the newly-formed grains that drives the outflow (e.g. Höfner & Olofsson 2018). The details of the wind-driving mechanism in oxygen-rich AGB stars (with a carbon-to-oxygen number ratio lower than one) is particularly complex, because theoretical models indicate that the dust grains need to grow to tenths of micrometers in size to provide the required opacity through scattering (Höfner 2008). These large grains must have a low absorption opacity in order to remain relatively cold even close to the star and, hence, not be destroyed (Woitke 2006). Observations have shown that such large grains indeed occur (e.g. Norris et al. 2012; Ohnaka et al. 2016; Khouri et al. 2016a).

The sequence of dust species that condense in O-rich AGB stars and the composition of the grains that provide the scattering opacity are still debated (e.g. Gail et al. 2016; Gobrecht et al. 2016; Höfner et al. 2016). Nonetheless, based on the elemental abundance of dust-forming elements, the condensation of some type of silicates is thought to be necessary for the onset of the outflow. These silicates need to be iron-free in order to have low-absorption opacity (e.g. Bladh & Höfner 2012). Determining the composition of the grains seen in scattered light is difficult, because the scattering opacities are virtually independent of grain composition. Nonetheless, interferometric observations of dust emission and constraints from models of dust excess emission show that aluminium oxide grains exist very close to the star (at ~ 2 R_⋆, e.g. Zhao-Geisler et al. 2012; Karovicova et al. 2013; Khouri et al. 2015), while iron-bearing silicate grains are only seen at larger radii (at ≳ 5 R_⋆, e.g. Zhao-Geisler et al. 2012; Karovicova et al. 2013). The distribution of iron-free silicate grains is more difficult to determine from observations of dust emission, because these grains have low absorption cross sections, and therefore remain relatively cold and produceweak infrared emission.

To study the dust-formation and wind-driving processes and spatially resolve the wind-acceleration region, we observed Mira (o Ceti), the archetypal Mira-type variable, using the Atacama Large Millimetre Array (ALMA) and the Zurich imaging polarimeter (ZIMPOL) part of the spectro-polarimetric high-contrast exoplanet research (SPHERE) on the Very Large Telescope (VLT). The observations also resolve the stellar disc at visible and millimetre wavelengths. The ZIMPOL and ALMA observations were acquired 18 days apart, and Mira was at post-minimum light phase (φ ~ 0.7). Mira AB is a binary system consisting of an AGB star and a white-dwarf or K dwarf (e.g. Ireland et al. 2007) separated by roughly 90 au. The effect of the companion on the close environment of the AGB star is thought to be minimal (Mohamed & Podsiadlowski 2012). On larger scales, Mira B accretes material from the circumstellar envelope of Mira A and the slow, dense outflow from Mira A interacts with a faster, but thin, outflow from Mira B (e.g. Ramstedt et al. 2014). Throughout this paper we adopt a distance of 102 pc to the Mira system, which is the weighted average between the distance obtained by hipparcos (92 ± 10 pc, van Leeuwen 2007) and from period-luminosity relations (110 ± 9 pc, Haniff et al. 1995). The stellar radius of Mira varies strongly with wavelength, as demonstrated by our observations and by data at other wavelengths (see e.g. the discussion in Kamiński et al. 2016). For the discussion in this paper, we use the stellar radius measured from our ALMA observations, ~ 21 mas (see below), and a reference stellar radius of 15 mas obtained by Woodruff et al. (2009) from measurements at 1.25 μm at post-minimum stellar pulsation phase, φ = 0.7. We refer to this stellar radius measured in the infrared as $R_{\star }^{\text{IR}}$.

2 SPHERE observations

Mira was observed using filters NR, cnt748, and cnt820 by SPHERE/ZIMPOL on 27 November 2017 (ESO programme ID 0100.D-0737, PI: Khouri). The observations of the point-spread function (PSF) reference star (HD 12642) were carried out following the observations of Mira. HD 12642 has a diameter of ~ 2.3 mas at visible wavelengths (Bourges et al. 2017), and is, therefore, a suitable PSF reference. The data were reduced using the (SPHERE/ZIMPOL) SZ software package developed at the ETH (Eidgenössische Technische Hochschule, Zurich). The basic procedures are standard (bias frame subtraction, cosmic ray removal, and flat fielding; see e.g. Schmid et al. 2017) and essentially the same as those of the SPHERE reduction software provided by ESO. The images of Stokes Q, U, and I were produced using the expressions $$Q=\frac{Q_{+}-Q_{-}}{2}, U=\frac{U_{+}-U_{-}}{2} \text{and} I=\frac{Q_{+}+Q_{-}+U_{+}+U_{-}}{4},$$

where Q₊, Q₋, U₊, and U₋ are the images obtained with the half-wave plate oriented at 0°, 45°, 22.5°, and 67.5°, with respectto the north direction. The polarized intensity and the polarization degree are calculated using $$I_{\text{p}}^2=Q^2+U^2$$(1)

and $$p=I_{\text{p}}/I,$$(2)

respectively. In Table 1, we present the range of values of airmass and seeing during the observations as provided by ESO in the header of the raw fits files under the entries HIERARCH EZO TEL AIRM and HIERARCH ESO TEL IA FWHM, respectively. We also give the detector integration time (DIT), the number of exposures for a given orientation of the half-wave plate in a given cycle (NDIT), and the total exposure time, which is equal to 4 × DIT ×NDIT ×NQU, where NQU is the number of cycles for the acquisition of Q₊, Q₋, U₊, and  U₋ images. The detector saturated in the observations of HD 12642 using filter NR. Hence we are not able to compare models to the data obtained using this filter, as we do for filter cnt820 in Sect. 4.4. The adaptive optics system of SPHERE is expected to perform very well in the observational conditions and for bright sources such as Mira and HD 12642.

In the total intensity images, Mira A and Mira B can be seen in all filters, with the intensity of the compact companion decreasing relatively to that of the primary for longer wavelengths. In the images, the full-width at half maximum (FWHM) of the point-spread-function (PSF) reference star is 29.9 mas at 0.65 μm (with a small degree of saturation) and 29.5 mas at both 0.75 and 0.82 μm. The FWHM of Mira A is 42.1 mas at 0.65 μm, 43.9 mas at 0.75 μm, and 41.8 mas at 0.82 μm. The stellar disc of Mira A is, therefore, marginally resolved in the observations. We find Mira A to be asymmetric at 0.65 and 0.75 μm, extending in the north-east direction (see Fig. 1). At 0.82 μm, the stellar disc of Mira A does not appear significantly asymmetric.

The polarized intensity was computed from the observed images of the Stokes parameters Q and U (these images are shown in Fig. A.1). The data reveal structures whose relative strengths change considerably between filters. This can be seen in the images of the polarization degree shown in Fig. 1 for filters NR and cnt820. We do not discuss the observations using filter cnt748, because the images are consistent with those in filters NR and cnt820, but the signal-to-noise ratio in the images of polarized light is significantly lower. We have labelled the different dust structures observed in Fig. 1 to facilitate the discussion. Based on the spread of the values obtained for the individual exposures, the uncertainties on the reported polarization degree per pixel are ~ 0.02 at the brightest regions of structures α and β, ~ 0.01 at structure γ, and <0.01 closer to thecentral star. We note that the uncertainty on the average polarization degree of a given structure (which encompasses several pixels) would be significantly smaller than the values given above. The offsets in right ascension and declination given throughout the paper are with respect to the position of Mira A.

We detect polarized light close to Mira A, as expected, but also from around Mira B and from what seems to be a dust trail connecting the two stars. The degree of polarization peaks at ~ 3.2 % in NR, ~ 3.5 % in cnt748, and ~4.0% in cnt820. These are relatively low values when compared to other oxygen-rich AGB stars observed using ZIMPOL. This is particularly true when one considers that Mira A was at post-minimum light phase at the time of the ZIMPOL observationsand observations of other close-by AGB stars were obtained at a similar light phase. For instance, ZIMPOL observations of W Hya revealed a polarization degree peaking at 18% at minimum light phase and ~ 13% at pre-maximum light phase (Ohnaka et al. 2017), while observations of R Dor also show the peak of the polarization degree varying between ~7 and ~ 14% at two different epochs (Khouri et al. 2016a). Moreover, ZIMPOL monitoring of R Dor with observations at 11 epochs over eight months show the polarization degree to reach peak values >7% at all epochs (Khouri et al., in prep.). The polarization degree also peaks at larger distances from the star for Mira A (> 5 R_⋆) than for W Hya and R Dor (~2 R_⋆). One important difference between Mira A and the other two stars is that Mira A is a Mira-type variable, while W Hya and R Dor are classified as semi-regular pulsators.

Beam-shift effect. The ZIMPOL observations of Mira show signs of being affected by a beam-shift effect. This is a known problem that is caused by sub-pixel shifts between the images obtained at different orientations of the half-wave plate (Schmid et al. 2018). This leads to one negative and one positive lobe that appear in opposition in the images of the Stokes parameters Q and U, as can be seen in our data, particularly in the cnt820 images (see Fig. A.1). This beam-shift effect causes residuals to appear very close to the central star after polarimetric differential imaging. These residuals can be mistaken for a polarization signal, and need to be quantified when analysing polarized data in the innermost region.

We have examined the beam-shift effect in the cnt820 images. We find that a shift of ≲ 0.05 pixels (or 0.18 mas) in the north-east–south-west direction is needed to reproduce the lobes in the obtained Q images in filter cnt820, while a shift of similar magnitude roughly in the opposite direction reproduces the lobes seen in the U images in the same filter. These results are only valid for the observations we present since the beam-shift effect is highly dependent on many parameters (see discussion in Schmid et al. 2018). This shift produces artificial polarized light that reaches polarization degrees of ≾ 2% to the north-east and to the south-west of the central star up to ~70 mas from the centre. The polarization vectors produced by the beam-shift effect, in the specific case of the images in filter cnt820, are aligned (with respect to Mira A) tangentially in the north-east region and radially in the south-west region. Since the polarization vectors of single-scattered light from a central source are expected to be tangentially aligned, this would lead to constructive and destructive interactions with the polarized signal in the north-east and south-west, respectively. In this way, the polarization degree of Clump γ (see Fig. 1) is probably overestimated, especially in the cnt820 images. Moreover, the vectors of the polarized light we observe close to Mira A are aligned tangentially to the star in the north-east, but radially in the east (see Fig. 2). This is consistent with our estimation of the beam-shift effect in the cnt820 images.

Scattering of radiation by grains with sizes comparable to the wavelength of the incident radiation divided by 2π can cause polarization inversions for certain scattering angles (see e.g. Aronson et al. 2017). Such polarization inversions, or the presence of aligned dust grains (such as observed in e.g. VY CMa, Vlemmings et al. 2017b), could in principle be invoked to explain the radially aligned polarization vectors. The fact that we are able to reproduce the observed orientation and strength of the polarization vectors by introducing a beam-shift of the order of that seen in the Q and U images (Fig. A.1) shows that no astrophysical explanation needs to be invoked for the radially-aligned polarization vectors in the images we present.

The beam-shift effect appears to have a different outcome close to Mira B compared to Mira A at 0.82 μm (Fig. 2), although the same pattern and intensity of beam-shift effect is expected for both sources in the same image. We find two reasons for the observed difference. First, Mira B is relatively weak at 0.82 μm. Therefore, the polarization due to the beam-shift effect around Mira B is relatively weaker than around Mira A. This leads to polarization due to scattering of radiation off dust grains to dominate close to Mira B because the radiation field is dominated by Mira A. Second, at 0.82 μm the speckle ring at the adaptive-optics control radius (Schmid et al. 2018) appears on top of Mira B and significantlyaffects the total intensity and the polarization degree images. The brightest speckle at the adaptive-optics control radius can be seen in Fig. 2 as a dark region to the north of Mira B at 0.82 μm and to the west at 0.65 μm.

From our estimates of the residual polarization from the beam-shift effect and the polarization degree in the obtained images, we conclude that the intrinsic polarization degree is ≾ 2% in the innermost region around Mira A (r ≾ 84 mas). The beam-shift effect could in principle be removed by improving the alignment of the images during the reduction process (Schmid et al. 2018). However, since the effect is relatively small in the data we report and we can set a meaningful upper limit to the polarization degree close to the star, we do not attempt to improve the reduction.

Fig. 1 Polarization degree images obtained using ZIMPOL and filters NR (λC = 0.65 μm) and cnt820 (λC = 0.82 μm). The white contours show the total intensity at 80, 40, and 20% of the peak value. The red circle marks the position of the companion and the white ellipses on the lower right corners show the FWHM of the observed PSF reference star at each given filter. The offsets in declination and right ascension are given with respect to the position of Mira A.

Fig. 2 Polarization degree (shown by the colour map) and direction of the polarization vectors around Mira A (right panels) and Mira B (left panels). Upper and lower panels: image in the NR and the cnt820 filters, respectively. The red dots mark the positions of the two stars. The vectors are drawn when the polarization degree is higher than 1.5%. The offsets in declination and right ascension are given with respect to the position of Mira A.

3 ALMA observations

The ALMA observations were obtained on 9 November 2017 (project 2017.1.00191.S, PI: Khouri), about 18 days before the ZIMPOL observations. The observations were performed using four spectral windows (spws) of 1920 channels each. The velocity resolution was ~ 0.85 km s⁻¹ and the spws covered 329.25–331.1, 331.1–333.0, 341.35–343.25, and 343.2–345.1 GHz. The full calibration and imaging was performed in CASA 5.1.1. The initial gain calibration (bandpass, amplitude, absolute flux-density scale, and phase) was performed following the standard ALMA procedures. Observations of J0237+2848 were used for flux and bandpass calibrations and of J0217+0144 for gain calibration. We then imaged the line free channels and performed a self-calibration step using the stellar continuum. These phase solutions were also applied to the full dataset. The final image products were created using Briggs robust weighting (with robust parameter equal to 0.5), resulting in a beam with major and minor axes of 32.7 and 19.1 mas, respectively, and a position angle of 53.2 degrees. Root-mean squared (rms) noise levels of ~ 1.4 and ~ 2.0 mJy beam⁻¹ were achieved, respectively, in the final images of the higher and lower frequency spectral windows for a spectral resolution of ~ 1.7 km s⁻¹. In the continuum images, we find an rms noise of 0.04 mJy beam⁻¹.

By performing a fit of a uniform-disc source to the observed visibilities using the uvmultifit code (Martí-Vidal et al. 2014), we determine the size of Mira A to be (41.9 ± 0.1) × (42.5± 0.1) mas (at PA –83 ± 3 degrees) with a brightness temperature of 2100 K. Therefore, the stellar disc of Mira A is seen to be roughly circular at ~ 338 GHz at the time of the ALMA observations. We find a continuum flux density for Mira A of 250.6 ± 0.2 mJy. The residuals obtained after the subtraction of our uniform-disc model have maximum and minimum values of 1.9 and –2.8 mJy, respectively. We refer to the stellar radius at 338 GHz of 21.2 mas as $R_{\star }^{\text{338}\text{GHz}}$, that is $R_{\star }^{\text{338}\text{GHz}}=1.4 R_{\star }^{\text{IR}}$.

A spectrum extracted within an aperture of 50 mas in radius centred on the star reveals more than a 100 spectral lines of which approximately one third (consisting mainly of weak lines) has not yet been identified by us. The identified lines correspond to transitions from CO, ¹³CO, C¹⁸ O, SiO, ²⁹SiO, ³⁰SiO, Si¹⁷ O, ²⁹Si¹⁷O, SO, ³³SO, SO₂, TiO₂, AlO, PO, and H₂O. In this study, we focus on compact emission revealed by strong low-excitation lines of molecules such as SO and ¹³CO, on the vibrationally-excited (v = 1) rotational line CO J = 3−2 at 342.647 GHz (Gendriesch et al. 2009), and on lines of molecules relevant for dust formation (TiO₂ and AlO) incomparison to the images of polarized light. A comprehensive overview of the ALMA observations will be presented in a future paper. All the images of spectral line emission presented by us are stellar-continuum subtracted. The dark region at the stellar position is caused by molecular absorption against the stellar continuum.

In the ALMA data obtained by us, emission from transitions with upper-level excitation energy < 200 K arises from a region extending mostly in the northern hemisphere and reaching distances of up to ~ 200 mas in radius from the centre of Mira A (see Fig. 3). Two lobes can be identified in this emission region (see Fig. A.2), one extending to the east and blue-shifted and another extending to the north-west and slightly red-shifted with respect to the local standard of rest velocity of the centre of mass of the system, υ_LSR ~ 47 km s⁻¹ (e.g. Kamiński et al. 2017). Throughout this paper, we refer to this region (delimited by the red line in Fig. 3) as the compact molecular line emission region. As we argue in Sect. 4.1, this emission region is probably defined by a steep gas density decline at its outer edge.

Maximum recoverable scale. The lines we discuss in this study are ¹³CO J = 3−2, AlO N = 9−8, SO N_J = 8₈−7₇, CO v = 1, J = 3−2, and several lines of TiO₂. Among these, the emission regions of the ¹³CO and SO lines are expected to be more extended than the compact molecular line emission region we find. Such extended emission would be filtered out by the ALMA observations we report, which have a maximum recoverable scale of ~ 400 mas. The maximum extent of the structure we observe is <300 mas in the moment zero maps and ≾200 mas in the images of individual channels. This is significantly smaller than the maximum recoverable scale. The flux from lines with emission regions expected to be smaller than the maximum recoverable scale should be almost completely recovered. This is the case for AlO (Kamiński et al. 2016), CO v = 1 (Khouri et al. 2016b), and TiO₂ (Kamiński et al. 2017). To confirm this expectation, we performed simulations of observations of uniform discs by the array configuration used to observe Mira. We find that for a uniform disc with a diameter of 300 mas about 20% of the flux is lost when the visibilities are re-imaged, while for a 200 mas uniform disc model, only 5% of the flux is lost.

The J = 3−2 line of ¹³CO has an emission region much larger than the scales probed by the observations we report. Therefore, a large fraction of the smooth, large-scale emission is expected to be resolved out. This also affects the flux we estimate from the innermost region. To quantify how much flux density per beam is resolved out because of large-scale emission, we compared the flux we measure from the compact molecular line emission region to the observations obtained by Ramstedt et al. (2014). We imaged the ¹³CO J = 3−2 line (PI: Ramstedt) using Briggs robust weighting (with robust parameter equal to 0.5) and obtained images with a beam of 0.68 × 0.45 arcsec, and an area of 0.24 arcsec². The compact molecular line emission region (with an area of 0.046 arcsec²) is thus unresolved in these earlier observations. We find that the compact molecular line emission region accounts for 72% of the flux density in the central beam of the lower resolution image, while encompassing only 19% of the total area. Assuming that the resolved-outemission is smooth and that, therefore, the remaining 28% flux density is spread evenly over the area of 0.24 arcsec², we find that ~5% (28% × 0.046∕0.24) of the flux density in the central beam of the lower resolution observations is not recovered in the compact molecular line emission region. Therefore, the flux density measured in the compact molecular line emission region should have been about 7% larger. This emission would have been lost in our observations due to the filtering of large-scale emission.

Fig. 3 Comparison between a low-excitation SO NJ = 88−77 line observed using ALMA and the polarized light observed using ZIMPOL. The SO NJ = 88−77 integrated line intensity was normalized to the peak value in the image and is shown by the colour scale. The red contour shows the 10% level of the SO emission, and represents the edge of the emission region. The white contours show the 1.5, 2.3, and 3.1% levels of the polarization degree observed using ZIMPOL and filters NR, λC = 0.65 μm (top panel), and cnt820, λC = 0.82 μm (bottom panel). The white and red ellipses on the lower right corner show, respectively, the FWHM of the observed PSF reference star at the given filter in the visible and the beam of the ALMA observations. The offsets in declination and right ascension are given with respect to the position of Mira A.

4 Results

4.1 The distribution of gas and dust

First, we discuss the nature of the compact molecular line emission region seen in the ALMA images. We consider that the observed compact emission is most likely the consequence of a steep gas density decline at its outer edge. An alternative explanation could be that the low-excitation transitions are efficiently radiatively excited in this region or that there is a sudden change of the gas temperatureat the region’s outer edge. We consider these possibilities very unlikely, because lines of different molecules (SO, SO₂, ²⁹SiO, ¹³CO, AlO, TiO₂) with different upper-level energies show a very similar emission region. Moreover, the corresponding upper levels of the transitions have relatively low excitation energies and the relatively high gas densities we find (≳ 10⁸ cm⁻³, see Sect. 4.3) imply that these low-excitation transitions should be efficiently excited by collisions.

A clear correlation between the distributions of gas and dust can be seen when the ALMA and SPHERE data are compared. In Fig. 3, we show the moment zero map of the SO line and the polarization degrees at 0.65 and 0.82 μm. While Clump α and Structure β follow the edge of the compact molecular line emission region, Structure δ appears on top of the north-western lobe. The observed structures and the varying strength of polarized intensity is probably partially caused by different scattering angles in a somewhat hollow dust shell. In this case, the angle-dependence of the phase function and polarization efficiency can make the polarization degree vary significantly as a function of wavelength. Grain size variations between the different structures might also contribute to the observed differences in polarization degree. Since we do not know the geometry of this region along the line of sight, a detailed interpretation of the cause of the varying polarization degree is difficult.

Second, we discuss the possible causes for the existence of the compact molecular line emission region and the correlation between the distribution of gas and dust. Hydrodynamical models show that the propagation of a shock (e.g. Höfner et al. 2003; Ireland et al. 2011; Freytag et al. 2017) produces a density contrast at the shock front and causes dust formation to happen efficiently in the post-shocked gas (Freytag & Höfner 2008). These two features qualitatively match our observations but stellar-pulsation-induced shocks are expected to be important only very close to the star (up to a few stellar radii, Höfner & Olofsson 2018) before dissipating. Since the edge of the compact molecular line emission region extends up to $~12 R_{\star }^{\text{IR}}$, we consider it unlikely that the compact molecular line emission region is the consequence of pulsation-induced shocks.

Another possibility is that the mass-loss rate of Mira increased significantly recently. In this case, the edge of the molecular emission region would be the consequence of the density contrast between a previous lower mass-loss rate wind and a recent higher mass-loss rate wind. Interestingly, an X-ray outburst was seen towards Mira A by observations made in 2003 (Karovska et al. 2005). The authors speculated that the burst could have significant consequences for the mass-loss rate of Mira A. To cover a distance of ~ 18 au, to the edge of the compact molecular line emission region, over the period of 14 yr that separates the observationsreported by us from those by Karovska et al. (2005), gas would need to travel at a speed of ~ 6 km s⁻¹. This expansion speed is within the range expected for the inner regions of an AGB outflow and thus this is a plausible scenario. Finally, another possibility is that this high density region is caused by interactions with the binary companion, but the effect of the companion on gas very close to Mira A is expected to be small (Mohamed & Podsiadlowski 2012). Future observations that reveal the evolution of the system will be able to distinguish between these different possibilities.

Dust Clump γ is the one located closest to Mira A in the polarized-light images, located at a distance of $~3.5 R_{\star }^{\text{IR}}$ from the centre of Mira A. As discussed in Sect. 2, the polarized light detected in Clump γ is probably at least partially caused by the beam-shift effect. Nonetheless, the existence of this clump seems real, since it is seen in images in both filters. The amount of polarized light that originates from scattering of radiation on dust grains is, however, somewhat uncertain. The location of Clump γ is probably not only a projection effect because polarization through scattering would not be efficient if it were located too far out of the plane of the sky. Indeed, the image of the vibrationally-excited CO J = 3−2 line shows that Clump γ is most likely embedded in a region of very high gas density (see Fig. 4), with number density > 10¹⁰ cm⁻³ (see below).Interestingly, the stellar disc seems to extend in the direction of Clump γ at 0.65 μm (see Fig. 1) and 0.75 μm. A similar coincidence between an extension of the stellar photosphere and a dust clump has been reported for Betelgeuse (Kervella et al. 2018), which was interpreted to be caused by preferential mass loss from the pole of the star.

We note that the stellar continuum of Mira A measured by ALMA does not appear to be asymmetric. Therefore, while the images in visible light at 0.65 and 0.75 μm show an extension of Mira A toward Clump γ, those at 0.82 μm and in the sub-millimetre do not. Molecular opacity is expected to affect the observations at 0.65 and 0.75 μm more strongly than those at 0.82 μm, although the difference in source size between images at 0.82 μm and the other filters is small in the ZIMPOL observations (see Sect. 2). In the sub-millimetre, free electrons rather than molecules dominate the continuum opacity. Hence, a possible explanation for Mira A appearing extended only at 0.65 and 0.75 μm is that molecular opacity is large enough locally to affect the appearance of Mira A at these wavelengths but not enough to affect its appearance at 0.82 μm. In this scenario, the sub-millimetre opacity due to free electrons would also need to be small enough not to produce an effect on the 338 GHz continuum. This suggestion is speculative, given the data at hand, and more observations are needed to further investigate the cause of asymmetries of this type.

Fig. 4 Comparison between the CO v = 1, J = 3−2 and stacked TiO2 lines observed using ALMA and the polarized light observed using ZIMPOL. The colour map shows the stacked emission from the nine unblended TiO2 lines (see Table 2) normalized to the peak value. The 2 and 3% levels of the polarization degree at 0.65 and 0.82 μm are shown by the dashed white and orange lines, respectively. The yellow solid contour shows the 10% level of the CO v = 1, J = 3−2 emission. The white and yellow ellipses on the lower right corner show, respectively, the FWHM of the observed PSF reference star in the NR filter (which is very similar to that using filter cnt820) and the beam of the ALMA observations. The offsets in declination and right ascension are given with respect to the position of Mira A.

4.2 Dust around Mira B and the dust trail

The dusty Trail ϵ marked in Fig. 1 seems to connect the two stars. It is not obvious from the ZIMPOL images whether this traces an accretion flow to Mira B, the interaction between the outflows of the two stars, or some other process. We note that the ALMA observations do not reveal any structure related to this in the many spectral lines observed. A different bridge-like component was seen in X-rays using Chandra and in ultraviolet light using the Hubble Space Telescope (Karovska et al. 2005). While the dust trail we report bends through the southern hemisphere, the bridge in the X-rays connects the two sources through a roughly straight line. It is not clear whether this is because of time variability of the gas and dust density distribution or if these different observations trace distinct and different structures.

The polarized light detected towards the companion at 0.65 and 0.82 μm is most likely produced by scattering of radiation off dust grains. Ireland et al. (2007) found evidence of a disc surrounding Mira B, with an outer radius of ~ 10 au, using observations of thermal dust emission at wavelengths between 8 and 18 μm. The scattered light we see might be produced by the structure identified by them. Although Mira A outshines Mira B at visible wavelengths, radiation from the compact object still dominates the radiation field in its immediate vicinity. From the total intensity images, we find that Mira A is brighter than Mira B by factor of ~ 45 and ~ 1200 at 0.65 and 0.82 μm, respectively. By considering the separation between the two stars of ~90 au (Ireland et al. 2007, for an assumed distance of 107 pc), we estimate that Mira B dominates the radiation field at 0.65 and 0.82 μm up to distances of ~13 au and ~2.5 au in the direction of Mira A, respectively. This is reflected in the orientation of the polarization vectors around Mira B, which are mostly tangential to a circle centred on Mira B at 0.65 μm and to a circle centred on Mira A at 0.82 μm (see Fig. 2). The polarization vectors appear tangentially aligned to Mira B at 0.65 μm over a region of ~ 60 mas in radius around the companion star. This corresponds to ~6 au and is smaller by a factor of two than the estimate above. We note that determining the direction of the polarization vectors in that region is not straightforward. For instance, the polarization vectors to the west of Mira B in the image in filter NR are radially aligned with respect to Mira B, which we attribute to the beam-shift effect (see Sect. 2).

4.3 Gas model

Emission from the N = 9−8 transition of AlO and from several transitions of TiO₂ (see Table 2) is detected in the ALMA observations. These molecules contain elements thought to be important for dust formation (aluminium and titanium) and an analysis of the observed emission in comparison to the polarized light is particularly interesting. Doing a similar analysis for Si from the data at hand is not straightforward, because emission from the ²⁹SiO J = 8−7 line is very optically thick, and the emission becomes insensitive to column density in this regime. Observations of lower abundance isotopologues of SiO might be the best way to constrain the SiO abundance distribution. To estimate the total column density of gas, we use the lines v = 1, J = 3−2 of ¹²CO and v = 0, J = 3−2 of ¹³CO. We originally intended to also use the v = 0, J = 3−2 line of C¹⁸O, but our data show that it seems to be weakly masing, as inferred from the spatially compact and spectrally narrow emission observed.

4.3.1 CO v = 1, J = 3−2

We modelled the vibrationally excited CO J = 3−2 line with the same code and approach as used by Khouri et al. (2016b) and Vlemmings et al. (2017a). In the code, local thermodynamical equilibrium is assumed to calculate the level populations in a shell around the central star. We assume a CO abundance of 4 × 10⁻⁴ relative to H₂. The velocity, density, and temperature of the gas are defined by power laws and spherical symmetry is assumed. After the level populations are obtained, images at different velocities are calculated by solving the radiative transfer equation at several lines of sight spread radially from the centre of the star to the edge of the molecular shell. We convolved the model images with the beam of the ALMA observations obtained by us and extracted model spectra from different apertures to compare with the observations.

The model star was assumed to be a uniform disc when projected on the plane of the sky. The stellar radius and temperature were set to $R_{\star }^{\text{338}\text{GHz}}=2.15$ au and T_⋆ = 2100 K based on our fit to the continuum visibilities. We adopt a temperature profile of the type $$T_{\text{kin}}(r)=T_{\circ }{(R_{\star }^{\text{338} \text{GHz}}/r)}^{ϵ}.$$

An initial input for the values of T° and ϵ was obtained using the brightness temperature of the ¹³CO line at the υ_LSR of the system. We note that the observed brightness temperature might be affected by resolved-out flux. The flux loss is larger for the intensity at υ_LSR than for the integrated flux (discussed in Sect. 3), because at the υ_LSR the emission region is larger. We estimate the resolved-out flux in the ¹³CO line using the observations of Ramstedt et al. (2014) and find that 10% of the flux density per beam, at υ_LSR, from the compact molecular line emission region is not recovered by our observations. Therefore, we calculated brightness temperature maps correcting for this flux loss, and we find the brightness temperature profile to be described by T° = 900 K and ϵ = 0.4. The outer radius of the modelled region (of 4 $R_{\star }^{\text{338} \text{GHz}}$, or 84 mas) was defined based on the observed emission region of the CO v = 1, J = 3−2 line.

The velocity and density profiles are defined by $$\upsilon (r)={\upsilon }_{\circ }+({\upsilon }_{\text{f}}-{\upsilon }_{\circ }){[(r-R_{\star }^{\text{338} \text{GHz}})/(R_{\text{out}}-R_{\star }^{\text{338} \text{GHz}})]}^{\beta }$$

and $$n(r)=n_{\circ }\times {(R_{\star }^{\text{338} \text{GHz}}/r)}^{\eta },$$

respectively. The free parameters in the model that were varied to fit the data are the stochastic velocity, υ_sto, the velocity at $R_{\star }^{\text{338} \text{GHz}}$, υ°, the velocity at R_out, υ_f, the exponent of the velocity law, β, the density at $R_{\star }^{\text{338} \text{GHz}}$, n°, the exponent of the density law, η, the temperature at $R_{\star }^{\text{338} \text{GHz}}$, T°, and the exponent of the temperature law, ϵ. The velocity of Mira with respect to the local standard of rest, υ_LSR, was also adjusted to fit the observations. We find υ_LSR = 47.7 ± 0.1 km s⁻¹.

Our fitting procedure was carried out as follows. We first searched for the region of the parameter space that provided the best fits to the data by eye. Then, we calculated grids of models to refine the best models and to estimate uncertainties on the obtained parameters. The best models were chosen by minimizing the ${\chi }^2\equiv {\sum }_{i=1}^n{(S_i^{\text{obs}}-S_i^{\text{mod}})}^2/{\sigma }_i^2$, where $S_i^{\text{obs}}$ and $S_i^{\text{mod}}$ are the observed and modelled flux density values and σ_i is the standard deviation within a given observed spectrum. In the grid calculations, we varied all the free parameters, T°, ϵ, υ_sto, υ°, υ_f, β, n°, η, and υ_LSR. The grid was progressively refined until uncertainties on the parameters were well defined. We calculated close to 2 × 10⁴ models spanning the ranges given in Table 3. The fits had 54 degrees of freedom (DF), which is the number of fitted spectral channels minus the number of free parameters. Our best fit reached a reduced-χ² = χ²∕DF ~ 4.0. This shows that a more complex model (non-spherically symmetric) is required to fully reproduce the exquisite ALMA images of the molecular gas emission close to Mira. Nonetheless, we are able to constrain the free parameters in the context of our simplified model.

The best fitting model found by us is shown in Fig. 5 and the corresponding parameters and uncertainties are given in Table 3. In our best model, the gas density decreases from ${2.5}_{-0.5}^{+0.8}\times {10}^{11} {\text{cm}}^{-3}$ at $r=R_{\star }^{\text{338} \text{GHz}}$ (21 mas) to (1.0 ± 0.5) × 10¹⁰ cm⁻³ at $r=4 R_{\star }^{\text{338} \text{GHz}}$ (84 mas). The model requires infalling material close to the star to reproduce the observed red-shifted absorption (see Fig. 5). As expected, the free parameters are not all independent. For instance, lower values of T° require lower values of ϵ and vice-versa. The same is true for n° and η. We find a mass in the vibrationally excited CO region of (3.8 ± 1.3) × 10⁻⁴M_⊙.

The large step between the continuum brightness temperature and the gas kinetic temperature at the stellar radius is a somewhat inconsistent feature of our model. We explored making the gas kinetic temperature at $R_{\star }^{\text{338} \text{GHz}}$ (2.15 au) equal to the temperature of the model star (2100 K) while increasing the value of the exponent of the temperature power law, so that the temperature would drop to ~500 K at $r~4 R_{\star }^{\text{338} \text{GHz}}$. We find that models with such a temperature profile are unable to simultaneously reproduce the absorption towards the star and the emission from around it, because the absorption becomes too weak, in a relative sense. This shows that our model requires fairly cold gas, T ~ 1000 K, already very close to the star. One way to make our model more consistent would be to introduce a break in the temperature power law, so that the temperature profile would be steep very close to the star and more shallow after reaching T ~ 1000 K. We have not explored this possibility, however, to avoid increasing the number of free parameters in the model.

Our model assumes a sharp transition between the region where continuum opacity dominates and the region where CO v = 1, J = 3−2 line opacity dominates, while in reality the transition could be gradual. We note, however, that we obtain a good fit to the stellar disc of Mira at 338 GHz using a uniform disc, and, therefore, there is no evidence for a slow radial decrease of continuum opacity at the time of the observations. Hence, the large difference between the stellar brightness temperature and the gas kinetic temperature at $R_{\star }^{\text{338} \text{GHz}}$ indicates a sharp decline of the kinetic temperature from ~2100 to ~ 1000 K at the radius of the 338 GHz continuum, $R_{\star }^{\text{338} \text{GHz}}$. This result seems in qualitative agreement with the hydrodynamical models of Ireland et al. (2011), since at some phases a sharp temperature decline of approximately the required magnitude is observed at ~2 au.

The total mass we find in the CO v = 1, J = 3−2 line emission region, (3.8 ± 1.3) × 10⁻⁴M_⊙, is almost two orders of magnitude larger than that found by Khouri et al. (2016b), ~ 6 × 10⁻⁶ M_⊙, by modelling spatially unresolved observations of the CO v = 1, J = 3−2 line towardsMira at roughly the same pulsation phase and using a similar approach to the one in this work. Although variability might explain, at least partially, the different masses found, the lack of spatial information in their study also plays an important role. Since the star was unresolved in the observations used by the above authors, they assumed a stellar brightness temperature of 3000 K to derive the stellar radius at sub-millimetre wavelengths. This stellar brightness temperature is clearly too high based on the observations we report. This leads to an underestimation of the stellar radius by a factor of approximately two. The model by these authors was also simpler than ours, since the authors assumed constant excitation temperature and gas density over the line-formation region. The values of the excitation temperature they considered range from 1000 to 2000 K, which are higher than those we find. This shows how spatial information provides crucial constraints for this type of study.

The gas density we find at r = 21 mas, ${2.5}_{-0.5}^{+0.8}\times {10}^{11} {\text{cm}}^{-3}$, is almost two orders of magnitude lower than that obtained by Wong et al. (2016), 10¹³ cm⁻³, at the same region by calculating models to fit ALMA observations of H₂ O and SiO line emission acquired during science verification time in the 2014 ALMA Long Baseline Campaign. At r = 84 mas our result, (1.0 ± 0.5) × 10¹⁰ cm⁻³, is more comparable to theirs, ~5 × 10¹⁰ cm⁻³. This could be because of variability, since gas with radial speeds of ~ 8 km s⁻¹ would be able to travel through this region in the period of three years that separates the observations we report from those of Wong et al. (2016).

The hydrodynamical models calculated by Ireland et al. (2011) for the parameters of Mira produce gas densities that are very variable through the different pulsation cycles. The values presented vary between approximately 10¹² and 2.5 × 10¹⁰ cm⁻³ at 2.15 au and 10⁹ and 3 × 10⁷ cm⁻³ at 5 au (which is the maximum radius for which gas densities are presented). While the density we find at 2.15 au ($r=R_{\star }^{\text{338} \text{GHz}}$), $({2.5}_{-0.5}^{+0.8})\times {10}^{11} {\text{cm}}^{-3}$ is in agreement with the values reported by Ireland et al. (2011), we find higher densities at 5 au, (3.0 ± 1.0) × 10¹⁰ cm⁻³, than those their models predict.

Fig. 5 Best model fit to the CO v = 1, J = 3−2 line observed using ALMA (dashed black line) and preferred model after constraints from the 13CO J = 3−2 line are considered (dot-dashed red line). The parameters of these two models are given in Table 3. Upper panels: observed spectra of the CO v = 1, J = 3−2 line (solidblue line) and residuals (thin blue line with error bars) obtained by subtracting the best-fit model from the observed spectra. The spectra were extracted from circular regions with radii of 20 mas (left panel), 50 mas (middle panel), and 100 mas (right panel). Lower panels: gas temperature (left panel), density (middle panel), and velocity profiles (right panel) of the two models. The velocity is positive for material flowing away from the stellar surface.

4.3.2 Further constraints with the ¹³CO, J = 3−2 line

Now we compare models to the observed ¹³CO v = 0, J = 3−2 line based on the results from the CO v = 1, J = 3−2 line calculations. The ¹³CO v = 0, J = 3−2 models were computed assuming local thermodynamical equilibrium (LTE), as for the CO v = 1, J = 3−2 line. We extracted the ¹³CO v = 0, J = 3−2 line spectrum using circular apertures with radii of 20 and 84 mas ($4 R_{\star }^{\text{338} \text{GHz}}$). The area of the vibrationally excited CO region corresponds to ~40% of the total area of the compact molecular line emission region, while this region accounts for approximately 50% of the flux we recover in the ¹³CO line.

We calculated models considering the best-fit models, within the uncertainties given in Table 3, and assuming a ¹²C-to-¹³C ratio of 10 ± 3 (Hinkle et al. 2016). We find that we obtain reasonably good fits when using models with masses at the lower end of the one-σ uncertainty interval we obtain from the CO v = 1, J = 3−2 line, of (3.8 ± 1.3) × 10⁻⁴ M_⊙. These models have steeperdensity profiles and shallower temperature profiles than the best-fitting model from our grid but, nonetheless, still provide a reasonable fit to the CO v = 1, J = 3−2 observations (see Fig. 5). Our fit improves by considering a ¹²C-to-¹³C ratio of 13, at the edge of the one-σ uncertainty interval given by Hinkle et al. (2016). The preferred model fit, taking into account the ¹³CO line, is compared to the data in Fig. 6 and its parameters are given in Table 3. Based on the model predictions for the ¹³CO J = 3−2 line, we conclude that the mass in the vibrationally excited CO region is ~ 2.5 × 10⁻⁴ M_⊙. This corresponds to a total H₂ column density of 3.4 × 10²⁴ cm⁻² in this innerregion.

4.3.3 Depletion of dust-forming elements

Having constrained the total amount of gas, we are now able to study the abundances of molecules that bear dust-forming elements, TiO₂ and AlO. Of the 26 lines of TiO₂ that fall within the frequency range of the ALMA observations (based on measurements by Kania et al. 2011), about 20 can be identified in the spectrum. Of these, most are blended with stronger or equally strong lines. We identified nine lines for which the integrated line flux could be determined reliably. We extracted the fluxes using a circular aperture with $4 R_{\star }^{338 \text{GHz}}$ radius (see Table 2). For TiO₂ lines with E_u > 160 K, the vibrationally excited CO region accounts for ≳70% of the total flux from the compact molecular line emission region, while for the three lines with the lowest excitation energy this fraction is between 50 and 60%. From the rotational diagram shown in Fig. 7, we derive an excitation temperature of 580 ± 65 K and a column density of (2.7 ± 0.7) × 10¹⁶ cm⁻². We note that the excitation temperature obtained is comparable to the lowest values obtained in the CO model at $r~4 R_{\star }^{\text{338} \text{GHz}}$. The column density we find is a factor of 2.5 larger than the value of N_{TiO_2} = (1.15 ± 0.12) × 10¹⁶ cm⁻² obtained by Kamiński et al. (2017) using ALMA observations acquired in 2015. These authors also found a significantly lower excitation temperature of 174 ± 7 K. Their results are based on a region similar to the vibrationally excited CO region we consider. This difference could be explained by variability of the source, or by uncertainties and differences between the two analyses. Interestingly, for TiO, Kamiński et al. find an excitation temperature of 491 ± 96 K, which is comparable to what we find for TiO₂. Using the H₂ column density obtained from the CO model, our results imply a TiO₂-to-H₂ ratio of ~ 8 × 10⁻⁹.

Emission from the AlO N = 9−8 line peaks very close to the star. As shown in Fig. 8, roughly 60% of the emission of the compact molecular line emission region arises from the vibrationally excited CO region. Assuming LTE excitation and an excitation temperature of 580 ± 65 K, as found for TiO₂, we find that an AlO column density of (7 ± 2) × 10¹⁵ cm⁻² is needed to reproduce the observed emission from the vibrationally excited CO region ($r\le 4 R_{\star }^{\text{338GHz}}$). This translates to an AlO-to-H₂ ratio of ~2 × 10⁻⁹. Assuming a higher excitation temperature of 1000 K would cause the column density and AlO-to-H₂ ratio to increase only by a factor of two. The column density we find is comparable to that reported by Kamiński et al. (2016), ~ 5 × 10¹⁵ cm⁻², for an excitation temperature of 329 ± 51 K.

The Ti and Al abundances with respect to H₂ for a gas with solar composition (Asplund et al. 2009) and in which hydrogen is fully molecular are 1.8 × 10⁻⁷ and 5.6 × 10⁻⁶, respectively. Therefore, the TiO₂ molecules account for ~ 4.5% of the expected titanium and the AlO molecules account for <0.1% of the expected aluminium. Kamiński et al. (2017) found a 3.4 times higher column density of TiO than of TiO₂ towards Mira. If we assume the same ratio is valid at the epoch of our observations, our results indicate that these two molecules could account for close to a fifth of all titanium. Hence, our results indicate that a significant fraction of titanium is locked in other carriers. We note that the lack of TiO lines in our data and the relatively large uncertainties make a firm conclusion on the combined abundances of TiO and TiO₂ not possible, however. Moreover, observations by Kamiński et al. (2017) revealed that atomic Ti also exists in the extended atmosphere of Mira A, but the authors were unable to derive the Ti column density. Hence, we conclude that the titanium not found in TiO and TiO₂ molecules may be in atomic form or in dust grains. Determining the amount of free Ti atoms is therefore necessary for measuring the depletion of this element. The low abundance of AlO makes a scenario possible in which aluminium is very efficiently depleted into dust, but the uncertainties on the abundance of other aluminium-bearing species also hampers a firm conclusion. Nonetheless, thanks to the relatively larger solar abundance of aluminium, we are able to constrain the amount of Al atoms that can be locked in grains of different sizes by studying the polarized light observations with ZIMPOL (see below). An alternative explanation for the low estimated column density of AlO could be that the assumption of LTE is not valid. We consider this to be unlikely, however, because the high gas densities we derive in the region (> 10⁹ cm⁻³) imply that the level populations of AlO are thermalized (Kamiński et al. 2016).

Fig. 6 Comparison between the observed 13CO J = 3−2 line and predictions from the models that provide good fits to the CO v = 1, J = 3−2 line. The dashed black line shows the best-fitting model to only the CO v = 1, J = 3−2 line, and the dot-dashed red line shows the preferred model after the analysis of the 13CO J = 3−2 line. The models are the same as shown in Fig. 5 and with the parameters given in Table 3. Upper and lower panels: spectra extracted using apertures with radii equal to 20 and 84 mas, respectively.

Fig. 7 Rotational diagram obtained using the unblended TiO2 lines given in Table 2.

Fig. 8 Comparison between the AlO N = 9−8 and the CO v = 1, J = 3−2 lines observed using ALMA. The colour map shows the AlO N = 9−8 line in Jy km s−1 beam−1. The yellow contours show emission in the CO v = 1, J = 3−2 line at 80, 60, 40, 20, and 10% levels of the peak emission. The red contour shows the 10% level of the SO emission (same as in Fig. 3). The offsets in declination and right ascension are given with respect to the position of Mira A.

4.3.4 The aluminium budget

Gobrecht et al. (2016) calculated the abundances of aluminium-bearing species and of aluminium oxide dust particles considering a shock-induced chemistry in the extended atmosphere of a model AGB star. The model was tailored to fit the star IK Tau, with mass-loss rate and gas densities in the inner wind at least an order of magnitude larger than those of Mira. They found that a much larger fraction of the aluminium can be in AlOH than in AlO, and that AlOH can account for ~ 30% of all aluminium at a few stellar radii. Observations of Mira have suggested that AlOH is located only very close to the stellar photosphere, however, at temperatures of ~2000 K (Kamiński et al. 2016). Although the authors were unable to constrain the column density of AlOH, the abundance does not seem to be as high as expected based on the models of Gobrecht et al. (2016). Decin et al. (2017) observed two AGB stars, one with a higher mass-loss rate than Mira (and the same one modelled by Gobrecht et al. 2016, IK Tau) and another with a mass-loss rate more comparable to that of Mira (R Dor). They found that the abundances of AlO and AlOH are much lower than the expected aluminium abundance and that these molecules together account for only 2% of the total aluminium. This does not support the high abundance of AlOH predicted by Gobrecht et al. (2016). Moreover, the abundance of AlO in R Dor was found to be almost two orders of magnitude larger than that of AlOH.

The emission from AlOH gas in the vibrationally excited CO region can be estimated based on our gas model. Assuming solar composition, we find that even if only ~ 10% of aluminium is in AlOH molecules (or an AlOH abundance ~2.8 × 10⁻⁷), AlOH emission would produce a strong, optically thick N = 11−10 line at 346.155 GHz with a flux density ~1 Jy. The signal estimated for only this small region is about an order of magnitude larger than reported by Kamiński et al. (2016) within the beam of the Atacama Pathfinder Experiment (APEX) telescope. Therefore, we conclude that AlOH most likely does not account for a significant fraction of aluminium atoms in the close environment of Mira A. The aluminium budget is, therefore, still very incomplete. As discussed by Kamiński et al. (2016), aluminium is expected to be in atomic form at temperatures ≳ 2200 K, which is considerably higher than the gas temperatures we derive. In broad terms, the search for other aluminium-bearing molecules has not revealed clear candidates for the bulk of aluminium (Kamiński et al. 2016; Decin et al. 2017).

4.4 Aluminium oxide dust grains and molecular clusters

A very likely explanation of the low column density of AlO is that aluminium atoms are very efficiently depleted into dust. Mid-infrared interferometric observations of dust thermal emission towards other Mira variables show aluminium-bearing dust to form at radii ≲ 2 R_⋆ (e.g. Zhao-Geisler et al. 2012; Karovicova et al. 2013), models for dust emission suggest aluminium oxide grains must exist in the extended atmosphere (Khouri et al. 2015), and dust condensation models predict aluminium oxide dust to form efficiently from the gas phase (Gobrecht et al. 2016; Höfner et al. 2016). Therefore, in this section we study the polarized light observed using ZIMPOL in order to constrain the presence and properties of aluminium oxide dust grains that are expected to exist close to the star.

Theoretical models for dust nucleation and growth find that aluminium oxide grains can form efficiently around O-rich AGB stars. Gobrecht et al. (2016) predict a significant fraction of Al atoms (~ 30%) to be found in aluminium oxide dust grains (Al₂ O₃), provided that the particles can withstand the passage of shocks. Moreover, they found that such grains grow to sizes between 0.1 and 0.3 μm within 2 R_⋆, with R_⋆ = 1.7 au in their models. Höfner et al. (2016) calculated the growth of amorphous aluminium oxide dust on top of seed dust particles considering a small range of input parameters of model O-rich stars that produce mass-loss rates of the same order as that of Mira. They found that more than 50% of the aluminium condenses into dust grains, if the grains have low absorption opacity. The particles reach sizes ≳ 0.1 μm at ~ 1.7 R_⋆, with the stellar radius in their models being equal to either 1.4 or 1.8 au for different models. To constrain the size of possible aluminium oxide dust particles located in the vibrationally excited CO region, we calculated continuum radiative-transfer models using MCMax (Min et al. 2009) in its spherically symmetric mode and with a similar approach to that of Khouri et al. (2016a).

4.4.1 The lattice structure and optical constants of Al₂ O₃

The optical constants of the considered solid have a major effect on the results of the calculations, because they directly determine the scattering and absorption properties of the model grains. We discuss the optical properties of Al₂ O₃ solids based on the complex index of refraction, m(λ) = n(λ) + iκ(λ). The real and imaginary parts of the index of refraction mostly affect the scattering and absorption properties of the grains, respectively. In this work, we are mainly interested in the value of the real part of the index of refraction, n, at 0.82 μm.

The index of refraction differs between Al₂ O₃ solids with different lattice structure. For amorphous Al₂ O₃ at short wavelengths, Suh (2016) calculated the index of refraction based on experimental data obtained by Begemann et al. (1997) for wavelengths between 7.8 and 500 μm. For crystalline Al₂ O₃, Harman et al. (1994) and Koike et al. (1995) measured the index of refraction of α-Al₂O₃ and γ-Al₂O₃, respectively. Koike et al. (1995) studied two samples of γ-Al₂O₃, one commercially available and the other obtained as a combustion product. Edlou et al. (1993) measured the optical constants of an Al₂ O₃ coating deposited by reactive magnetron sputtering. The lattice structure of the coating was not specified by the authors. A summary of the relevant values obtained from the literature is given in Table 4.

The large spread in values of κ reported for Al₂ O₃ solids leads to very uncertain absorption properties of Al₂ O₃ grains, while the small spread in the value of n suggests some confidence in the value of the scattering properties. If the value of κ is relatively too large, the equilibrium temperature of the grains might be too high for them to remain amorphous, or even form at all, close to the star. Höfner et al. (2016) found that aluminium oxide dust can form very efficiently if the grains have low absorption cross sections (corresponding to low values of κ). The authors found that if the value of κ for aluminium oxide grains is similar to that of Mg₂ SiO₄ (~ 10⁻⁴), aluminium would be fully depleted at r ~ 2.1 au (1.5 $R_{\star }^{\text{IR}}$). If the same applies to aluminium oxide grains around Mira A, full condensation of aluminium could happen already below the millimetre and visible stellar photosphere.

4.4.2 Model calculations and results

The grain opacities were calculated using the hollow-spheres approximation (Min et al. 2003). We used the complex refractive index obtained by Suh (2016) in our model calculations, since grains with these properties scatter less efficiently. The scattering opacity we find at 0.82 μm using the data from Suh (2016) is ~25% smaller than that obtained using the optical constants provided by Koike et al. (1995). MCMax calculates the scattering of radiation considering the full angle-dependent Mueller matrix and produces the images of the Stokes parameters I, Q, and U. We convolved these images with the observed PSF reference to obtain images at the same resolution as the observed ones (I^conv, Q^conv, and U^conv). Then, we constructed the images of the polarized intensity and the polarization degree using the convolved images and Eqs. (1) and (2). We only model the observations at 0.82 μm.

Our modelling focuses on the region traced by the vibrationally excited CO line, with an outer radius of $5.6 R_{\star }^{\text{IR}}=4 R_{\star }^{\text{338} \text{GHz}}=8.6$ au. We consider a solar abundance of aluminium of 2.8 × 10⁻⁶ relative to H (Asplund et al. 2009), which implies the same value for the maximum abundance of Al₂ O₃ relative to H₂. Based on the gas model presented in Sect. 4.3.1, the solar abundance of Al, and the mass of an Al₂ O₃ molecule (1.7 × 10⁻²² g), we write the radial density profile of Al₂ O₃ dust grains as $${\text{ρ}}_{\circ }^{{\text{Al}}_{\text{2}}{\text{O}}_{\text{3}}}(r)\approx 4.8\times {10}^{-28} n_{\circ } f {(1.4\times R_{\star }^{\text{IR}}/r)}^{\eta } \text{g} {\text{cm}}^{-3},$$

for $R_{\text{in}}^{\text{Dust}}<r<R_{\text{out}}^{\text{Dust}}$. In this expression, f is the fraction of aluminium that condenses into Al₂ O₃ dust, $R_{\text{in}}^{\text{Dust}}$ is the inner radius of the dust envelope (see below), and $R_{\text{out}}^{\text{Dust}}=5.6 R_{\star }^{\text{IR}}$ is the outer radius of all model dust envelopes. The values of n° and η are adopted from the gas model.

We calculated models to fit the cnt820 images assuming that condensation of aluminium happens at $r~22.5 \text{mas}=1.5 R_{\star }^{\text{IR}}$ or at $r~\text{}37.5 \text{mas}=2.5 R_{\star }^{\text{IR}}$, based, respectively, on theoretical results and on interferometric mid-infrared observations of other Mira variables, as discussed above. The star was approximated by a black body with a temperature of 2500 K. This assumption does not strongly affect the calculations of the scattered light. The free parameters in our dust models are f and the size of the dust grains. We consider an upper limit to the polarization degree of 2% in the vibrationally excited CO region (see discussion in Sect. 2).

For both condensation radii considered, our models show that only < 10, < 1, < 1, and ~ 100% of the available aluminium atoms can be locked in 1.0, 0.3, 0.1, and 0.02 μm aluminium oxide grains, respectively. Otherwise, a stronger polarized-light signal at visible wavelengths would be expected. This result is valid for all the values of n given in Table 4, and is therefore independent of the assumed lattice structure. Al₂ O₃ grains with sizes of 0.1 μm could account for >10% of the aluminium atoms in the vibrationally excited CO region only if the value of n of such grains was <1.05. This is significantly smaller than what the experimental results suggest for Al₂ O₃ solids with different lattice structures. We find that a relatively larger population (~ 10%) of grains with sizes ~1 μm could exist, although the opacity per unit mass they provide is larger than that provided by ~ 0.1 μm grains. The lower polarization degree is caused by the more complex scattering phase function of the 1.0 μm grains. We consider it unlikely that a population of grains of ~1 μm exists around Mira without a significant amount of grains with size between 0.02 and 0.5 μm. Moreover, even if ~ 1 μm-sized grains do exist, these would still only account for a small fraction of the expected amount of aluminium atoms close to Mira A. We conclude that if aluminium atoms in the vibrationally excited CO region are efficiently locked in aluminium oxide grains, these must have relatively small sizes ≲ 0.02 μm. This is not consistent with the larger grains (~0.1 μm) found in models for Al₂ O₃ dust formation and growth in the extended atmospheres of AGB stars. We consider it unlikely that a significant fraction of aluminium atoms are locked in ~0.1 μm grains of a dust species other than aluminium oxide (such as in a grain with an Al₂ O₃ core and a silicate mantle). This is because aluminium would account for a smaller fraction of the mass of a single grain in that case and, therefore, significant depletion of aluminium would lead to a larger amount of dust. Such dust species could therefore only exist if they were very transparent at visible wavelengths. Optical constants of Mg₂ SiO₄ and MgSiO₃ (Jäger et al. 2003) indicate that these silicates have scattering properties that are similar to those of Al₂ O₃.

Decin et al. (2017) suggested the bulk of aluminium atoms could be locked in (Al₂O₃)_n clusters, which would have sizes consistent with the limits we set. These clusters might still produce an infrared excess similar to that expected from amorphous aluminium oxide grains, and therefore be detected in interferometric observations in the mid-infrared. Although experiments have shown that large aluminium oxide clusters (with more than 34 Al₂O₃ molecules) have a similar spectral signature to that of amorphous aluminium oxide (van Heijnsbergen et al. 2003), the opacity per mass of such clusters is not known.

4.4.3 The effect of Al₂O₃ opacity at mid-infrared wavelengths

Considering either the optical constants given by Koike et al. (1995) for crystalline aluminium oxide or Begemann et al. (1997) for amorphous aluminium oxide, we find that the inner part of the envelope would become optically thick at the aluminium oxide resonance at ~11 μm for full aluminium condensation. The opacity due to the efficient formation of large Al₂O₃ clusters cannot be determined based on the data available in the literature at the moment. This would cause the size of Mira A to be roughly two times larger at 11 μm (45 mas in radius in our models) than at millimetre wavelengths. Weiner et al. (2003) observed Mira A between 1999 and 2001 with the Infrared Spatial Interferometer (ISI) at ~ 11 μm and found a stellar radius varying between 23 and 31 mas. This is larger than the radius measured in the near-infrared by roughly a factor of two, and therefore comparable to the millimetre and visible sizes we report (~ 21 mas). Hence, this measurement does not support the existence of a shell of aluminium-oxide dust that is optically thick at ~ 11 μm. Models by Weiner (2004) and Ohnaka et al. (2006) explained the observed 11 μm size of Mira A with an optically thick layer of hot (2200 K) water, although some discrepancies between model and observations remained. Ohnaka et al. (2006) argued that dust formation on the upper atmosphere could also contribute significantly to the opacity at 11 μm.

More recent measurements with mid-infrared interferometric instrument (MIDI) show an asymmetric and relatively larger star between 10 and 14 μm, with radii as large as 70 mas (Paladini, priv. comm.). This would be comparable to what we expect from an optically thick aluminium oxide shell. The complexity of the source at the time of the MIDI observations makes an objective conclusion not possible. Moreover, the distribution of gas close to the star at the time of the ISI and MIDI observations was most likely very different from that shown by the observations we report. New spatially resolved observations at ~ 11 μm, possibly together with coordinated ALMA and SPHERE observations, will help further constrain the presence and the properties of aluminium oxide grains in the close vicinity of Mira A.

5 Summary

We imaged molecular line emission and polarized light, produced by scattering of radiation on dust grains, around the Mira AB system using ALMA and SPHERE/ZIMPOL. The images in polarized light show dust around Mira A, but also around Mira B and a dust trail connecting the two stars. We find a region of high gas densities around Mira A, which seems to be delimited by a steep density decline at its outer edge. The dust grains that produce the observed scattered light cluster at the edge of this region. If the proposed density decline is caused by the propagation of a shock front, this would be in qualitative agreement with models that predict dust formation to happen very efficiently in the post-shocked gas. The radius at which the density drop is seen ($~12 R_{\star }^{\text{IR}}$) is larger than the distances at which pulsation-induced shocks dissipate in hydrodynamical models, however. Alternative possibilities are that the mass-loss rate of Mira has recently increased creating an inner region with high gas densities, or that interactions with the binary companion enhance the density in this region. We model the CO v = 1, J = 3−2 line to obtain the radial profiles of the gas density, temperature, and velocity very close to the star. The gas model was refined by comparing different models to the ¹³CO J = 3−2 line data. We calculate the abundances of TiO₂ and AlO and find that these account for only 4.5% of Ti and < 0.1% of Al expected for a gas of solar elemental abundance. The lack of TiO lines and observations of other potential Ti carriers (such as atomic Ti) in our dataset makes a firm conclusion on the Ti depletion difficult, but we estimate that TiO and TiO₂ account together for ~20% of the expected Ti. It is even harder to determine the aluminium gas-phase budget, but we are able to constrain the properties and sizes of aluminium oxide dust grains that might be an important reservoir of aluminium atoms. We find that grains located within 84 mas ($5.6 R_{\star }^{\text{IR}}$) from the star must have an average size ≲ 0.02 μm to account for a significant fraction of Al. If Al₂ O₃ grains with sizes ~0.1 μm are present, these must account for less than 1% of the aluminium atoms. We also find that opacity due to efficient aluminium oxide dust formation could cause the star to appear very large at wavelengths close to the 11 μm feature. Recent interferometric observations in the mid-infrared suggest that Mira A might indeed be very large at these wavelengths, but a detailed analysis is hampered by the complexity of the source and the time lag between the mid-infrared and the ZIMPOL and ALMA observations that we report. Quasi-simultaneous observations in visible scattered light, infrared thermal dust emission, and molecular line emission will help better constrain the aluminium budget in Mira and other O-rich AGB stars.

Acknowledgements

T.K., W.V., H.O., and M.M. acknowledge support from the Swedish Research Council. E.D.B. acknowledges support from the Swedish National Space Agency. This paper makes use of the following ALMA data: ADS/JAO.ALMA#2017.1.00191.S and ADS/JAO.ALMA#2012.1.00524.S. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada), NSC and ASIAA (Tai- wan), and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ.

Appendix A Additional figures

Fig. A.1 Images of Stokes parameters Q and U in filters NRand cnt820. The colour maps show the observed counts toward the Mira AB system and their range differs to highlight the beam-shift effect (upper panels) and the polarized-light signal (lower panels). The offsets in declination and right ascension are given with respect to the position of Mira A.

Fig. A.2 Velocity-channel maps of the SO NJ = 88−77 transition at 344.253 GHz given in Jy beam−1. Each map is the average of emission within a channel 1.9 km s−1 wide centred on the given local-standard-of-rest velocity. The offsets in declination and right ascension are given with respect to the position of Mira A.

References

All Tables

All Figures

Fig. 1 Polarization degree images obtained using ZIMPOL and filters NR (λC = 0.65 μm) and cnt820 (λC = 0.82 μm). The white contours show the total intensity at 80, 40, and 20% of the peak value. The red circle marks the position of the companion and the white ellipses on the lower right corners show the FWHM of the observed PSF reference star at each given filter. The offsets in declination and right ascension are given with respect to the position of Mira A.

In the text

Fig. 2 Polarization degree (shown by the colour map) and direction of the polarization vectors around Mira A (right panels) and Mira B (left panels). Upper and lower panels: image in the NR and the cnt820 filters, respectively. The red dots mark the positions of the two stars. The vectors are drawn when the polarization degree is higher than 1.5%. The offsets in declination and right ascension are given with respect to the position of Mira A.

In the text

Fig. 3 Comparison between a low-excitation SO NJ = 88−77 line observed using ALMA and the polarized light observed using ZIMPOL. The SO NJ = 88−77 integrated line intensity was normalized to the peak value in the image and is shown by the colour scale. The red contour shows the 10% level of the SO emission, and represents the edge of the emission region. The white contours show the 1.5, 2.3, and 3.1% levels of the polarization degree observed using ZIMPOL and filters NR, λC = 0.65 μm (top panel), and cnt820, λC = 0.82 μm (bottom panel). The white and red ellipses on the lower right corner show, respectively, the FWHM of the observed PSF reference star at the given filter in the visible and the beam of the ALMA observations. The offsets in declination and right ascension are given with respect to the position of Mira A.

In the text

Fig. 4 Comparison between the CO v = 1, J = 3−2 and stacked TiO2 lines observed using ALMA and the polarized light observed using ZIMPOL. The colour map shows the stacked emission from the nine unblended TiO2 lines (see Table 2) normalized to the peak value. The 2 and 3% levels of the polarization degree at 0.65 and 0.82 μm are shown by the dashed white and orange lines, respectively. The yellow solid contour shows the 10% level of the CO v = 1, J = 3−2 emission. The white and yellow ellipses on the lower right corner show, respectively, the FWHM of the observed PSF reference star in the NR filter (which is very similar to that using filter cnt820) and the beam of the ALMA observations. The offsets in declination and right ascension are given with respect to the position of Mira A.

In the text

Fig. 5 Best model fit to the CO v = 1, J = 3−2 line observed using ALMA (dashed black line) and preferred model after constraints from the 13CO J = 3−2 line are considered (dot-dashed red line). The parameters of these two models are given in Table 3. Upper panels: observed spectra of the CO v = 1, J = 3−2 line (solidblue line) and residuals (thin blue line with error bars) obtained by subtracting the best-fit model from the observed spectra. The spectra were extracted from circular regions with radii of 20 mas (left panel), 50 mas (middle panel), and 100 mas (right panel). Lower panels: gas temperature (left panel), density (middle panel), and velocity profiles (right panel) of the two models. The velocity is positive for material flowing away from the stellar surface.

In the text

Fig. 6 Comparison between the observed 13CO J = 3−2 line and predictions from the models that provide good fits to the CO v = 1, J = 3−2 line. The dashed black line shows the best-fitting model to only the CO v = 1, J = 3−2 line, and the dot-dashed red line shows the preferred model after the analysis of the 13CO J = 3−2 line. The models are the same as shown in Fig. 5 and with the parameters given in Table 3. Upper and lower panels: spectra extracted using apertures with radii equal to 20 and 84 mas, respectively.

In the text

	Fig. 7 Rotational diagram obtained using the unblended TiO2 lines given in Table 2.
In the text

Fig. 8 Comparison between the AlO N = 9−8 and the CO v = 1, J = 3−2 lines observed using ALMA. The colour map shows the AlO N = 9−8 line in Jy km s−1 beam−1. The yellow contours show emission in the CO v = 1, J = 3−2 line at 80, 60, 40, 20, and 10% levels of the peak emission. The red contour shows the 10% level of the SO emission (same as in Fig. 3). The offsets in declination and right ascension are given with respect to the position of Mira A.

In the text

	Fig. A.1 Images of Stokes parameters Q and U in filters NRand cnt820. The colour maps show the observed counts toward the Mira AB system and their range differs to highlight the beam-shift effect (upper panels) and the polarized-light signal (lower panels). The offsets in declination and right ascension are given with respect to the position of Mira A.
In the text

	Fig. A.2 Velocity-channel maps of the SO NJ = 88−77 transition at 344.253 GHz given in Jy beam−1. Each map is the average of emission within a channel 1.9 km s−1 wide centred on the given local-standard-of-rest velocity. The offsets in declination and right ascension are given with respect to the position of Mira A.
In the text