Evidence for localized onset of episodic mass loss in Mira

We report Multi-telescope interferometric measurements taken with the Interferometric Optical Telescope Array (IOTA) to provide imagery of the LPV Mira in the H-band. This wavelength region is well suited to studying mass loss given the low continuum opacity, which allows for emission to be observed over a very long path in the stellar atmosphere and envelope. The observed visibilities are consistent with a simple core-halo model to represent the central object and the extended molecular layers but, in addition, they demonstrate a substantial asymmetry. An analysis with image reconstruction software shows that the asymmetry is consistent with a localized absorbing patch. The observed opacity is tentatively associated with small dust grains, which will grow substantially during a multi-year ejection process. Spatial information along with a deduced dust content of the cloud, known mass loss rates, and ejection velocities provide evidence for the pulsational pumping of the extended molecular layers. The cloud may be understood as a spatially local zone of enhanced dust formation, very near to the pulsating halo. The observed mass loss could be provided by several such active regions around the star. This result provides an additional clue for better understanding the clumpiness of dust production in the atmosphere of AGB stars. It is compatible with scenarios where the combination of pulsation and convection play a key role in the process of mass loss.


Introduction
All stars of low to intermediate mass evolve through the AGB phase, with a dredging up of synthesized elements, and undergo extensive mass loss as long-period variables (LPVs). This return of heavy elements to the interstellar medium constitutes a critical step in the evolution of stars and galaxies. The LPV is a very complex phenomenon and our understanding of it remains incomplete. Both modeling and observations continue to guide our progress in this area.
LPVs vary in brightness cyclically, with periods on the order of a few hundred days. Photometric monitoring shows that these periodic variations are often slightly or even severely irregular. The improving capabilities at present for direct imagery of LPVs show that the LPV atmosphere is asymmetric, at least in some cases and potentially for all. Speckle and other imaging techniques have been used to obtain high angular resolution information about LPVs with ground-based telescopes and interferometers, as well as with HST. Scholz (2003) gives an extensive list of such observations prior to 2003. Creech-Eakman & Thompson (2009) summarize a seven-year study of 100 LPVs using the Palomar Testbed Interferometer and report extensive results for apparent sizes and their variation with position angle and time. Ragland et al. (2006), using the Interferometric Optical Telescope Array (IOTA), found that 1/3 of 56 nearby LPV stars surveyed showed asymmetry. The authors speculated that all of them would appear asymmetric if studied with sufficient resolution. This trend was confirmed by the observations of Wittkowski et al. (2011) at VLTI with AMBER which measured deviations from centro-symmetry on four Mira variables. Later, Wittkowski et al. (2016) used AMBER to observe largescale inhomogeneities that contribute a few percent of the total flux in six Mira variables. The detailed structure of the inhomogeneities show a variability on time scales of three months and above. This is compatible with recent 3D dynamical simulations of AGB stars which deal self-consistently with convection and pulsations (Freytag et al. 2017). Convection, together with pulsational shocks, contributes to the raising of material up to altitudes where dust can condense. This knowledge leads us to consider whether the spatial inhomogeneities revealed by observations and predicted by models could also lead to inhomogeneous zones of dust formation. Ohnaka et al. (2016) used visible polarimetric imaging observations of W Hya taken with VLT/SPHERE-ZIMPOL to detect three clumpy dust clouds at two stellar radii. The authors deduced large dust grains (0.4-0.5 µm) of Al 2 O 3 , Mg 2 SiO 4 , or MgSiO 3 . The grain size is consistent with the prediction of hydrodynamical models for mass loss driven by scattering due to Fe-free micron-sized grains (Höfner 2008). The detection of clumpy dust clouds close to the star lends support to dust formation induced by pulsation and large convective cells. Observations at similar spatial scales but at longer wavelengths in the millimeter range have also revealed asymmetries. Wong et al. (2016) detected bright SiO-emitting clumps with ALMA in Mira at distances around two stellar radii. Also with ALMA on Mira, Kamiński et al. (2017) revealed the anisotropy of the distribution of TiO 2 in the atmosphere of Mira up to 5.5 stellar radii. This may suggest inhomogeneous distributions for both dust and its precursors near the photosphere, which is compatible with clumpy nucleation sites.
The shapes suggest a structure that may eventually contribute to the form of the planetary nebula (PN) expected to follow the AGB phase. However, PNe almost certainly represent the shape of the outer envelope and beyond, arising perhaps at some distance from the asymmetrical physical phenomena occurring at yet deeper layers which give rise to them. In addition, part of the shape may be due to binarity. Measurements with higher angular resolution at wavelengths that can penetrate deeply into the envelope promise to stand as the best opportunity to tie together the most reliable part of the LPV theory -the unobservable interior processes -with the complexity of the observed exterior.
Here, we report the continuation of a series of LPV studies carried out at IOTA (Traub (1998), which offers excellent (u,v) coverage, which can be enhanced with multi-spectral measurements. In Perrin et al. (2004), we showed strong evidence in favor of the presence of a very deep, but surprisingly transparent, extended atmosphere that is analogous to, and perhaps similar to, the MOLsphere, which Tsuji ascribed to luminous, non-Mira stars. This interpretation was supported by other observations, for example Lacour et al. (2009) andLe Bouquin et al. (2009). Since dynamical models reproduce the characteristics of Tsuji's MOLsphere, without, however, an actual discrete layer (Bladh et al. (2015)), the historical term is no longer understood as descriptively correct. We follow our referee's recommendation to describe this region as "extended molecular layers." In our previous studies, we utilized pair-wise combination of telescopes and only reported the analysis based on fringe visibilities and spherical models. In this study, we use the IONIC beam combiner instrument of IOTA (Berger et al. 2003), which combines three telescopes simultaneously and measures closure phases. Squared visibilities and closure phases can be combined to reconstruct images. We present the observations and data reduction in Sect. 2. The data are first fit with a circularly symmetric model in Sect. 3. These data have allowed us to synthesize an image of the LPV Mira (omicron Ceti), presented in Sect. 4, which is sufficiently indicative of the asymmetric structure to motivate inclusion of an asymmetric component in our models and to derive some related parameters in Sect. 5. We follow with a discussion of how this model relates to other measurements and how it may fit in  the context of expected phenomena in LPV oscillation and mass loss in Sect. 6. Finally, we comment on the future of LPV imagery in Sect. 7.

Description of IOTA observations
The interferometric data presented herein were obtained using the IOTA interferometer (Traub et al. 2003), which was a long baseline interferometer operating at near-infrared wavelengths until 2006. It consisted of three 0.45 m telescopes movable among 17 stations along two orthogonal linear arms. IOTA could synthesize a total aperture size of 35×15 m, corresponding to an angular resolution of 10×23 milliarcseconds at 1.65 µm. Visibility and closure phase measurements were obtained using the integrated optics combiner IONIC (Berger et al. 2003); light from the three telescopes is focused into single-mode fibers and injected into the planar integrated optics (IO) device. Six IO couplers allow recombinations between each pair of telescopes. Fringe detection was done using a Rockwell PICNIC detector (Pedretti et al. 2004). The interference fringes were temporally-modulated on the detector by scanning piezo mirrors placed in two of the three beams of the interferometer.  The closure phases are plotted against the maximum spatial frequency of each triangle. The closure phases are compatible with 0(π) up to the first zero crossing of the squared visibilities near 50 cycles/arcsec. Beyond that spatial frequency, closure phases depart from 0(π) meaning the source is not centro-symmetric at scales smaller than the typical radius of the star.
Observations were carried out in the H band (1.5 µm ≤ λ ≤ 1.8 µm) divided into seven spectral channels. The science target observations were interleaved with identical observations of unresolved or partially resolved stars, used to calibrate the interferometer's instrumental response and effects of atmospheric seeing on the visibility amplitudes. The calibrator sources were chosen from the Borde et al. (2002) catalog. Mira was observed in October 2005 over nine nights and using five different configurations of the interferometer. Full observation information can be found in Table 1, including dates of observation, interferometer configurations, and visual phase computed from the Association Française des Observateurs d'Étoiles Variables (http://cdsarc.ustrasbg.fr/afoev/). The calibrators are listed in Table 2. Figure 1 shows the (u,v) coverage achieved during this observation run. The geometry of the IOTA interferometer and the position of the star on the sky constrained the frequency coverage, which is elongated and primarily located in the north-south direction.

Data reduction
Reduction of the IONIC visibility data was carried out using custom software by Monnier et al. (2006), which is similar in its main principles to those described by Coudé du Foresto et al. (1997). The premise is to measure the power spectrum of each interferogram (proportional to the target squared visibility, V 2 ) after correcting for intensity fluctuations and subtracting bias terms due to read noise, residual intensity fluctuations, and photon noise (Perrin 2003). Next a correction is applied by Table 3. Parameters of the star+shell model fit to the V 2 data. D T D shell T shell τ 21.1 mas 3400 K 40.00 mas 2000 K 0.70 the data pipeline for the variable flux ratios for each baseline by using a flux transfer matrix (Monnier 2001). Finally, raw squared visibilities are calibrated using the visibilities obtained by the same means on the calibrator stars. As in Lacour et al. (2008) less than a 1% uncertainty in visibility is associated with the calibrators. This corresponds to a 2% calibration error for V 2 . Therefore, we systematically added a 2% calibration error to all the calibrated squared visibilities in this paper. In order to measure the closure phase (CP), a real-time fringe tracking algorithm (Pedretti et al. 2005) ensured that interference occurred simultaneously for all baselines. We required that interferograms are detected for at least two of the three baselines in order to assure a good closure phase measurement. This technique, called baseline bootstrapping, allowed for precise visibility and closure phase measurements for a third baseline with very small coherence fringe amplitude. We followed the method of Baldwin et al. (1996) to calculate the complex triple amplitudes and derive the closure phases. Pair-wise combiners (such as IONIC) can have a large instrumental offset for the closure phase which can be calibrated by the closure phase of the calibrator stars.
The squared visibilities and the closure phases are presented in Fig. 2 and Fig. 3, respectively. The colors show different baseline orientations. Obviously the source is not circularly symmetric as the visibilities depend on azimuth. A zero crossing is visible in the visibility function near 50 cycles/arcsec, the signature of sharp edges for the main contributor to the flux. In addition, the closure phases depart from 0 modulo π beyond the first zero of the visibility function, showing a departure from centro-symmetry at scales that are smaller or comparable to that of the central object.

Model fitting with a circularly symmetric model
The V 2 data were first compared to a model star photosphere + partially transparent shell (Perrin et al. 2004) to derive the main parameters of the model: photosphere diameter, shell diameter, shell optical depth, shell temperature, and photosphere temperature. The degeneracy of the temperature fit is solved by fitting the flux. For this purpose, we used the phases from the AFOEV data, as shown in Table 1. For the average photometric phase of 0.43 the magnitudes interpolated from Whitelock et al. (2000) are H (-1.7 to -1.8) and K (-2.2 to -2.3). We note that the infrared fluxes are much less variable than the visible fluxes. The H and K bands gave weak constraints on the photospheric temperature, probably because the photospheric diameter is not uniquely determined. The flux gave a minimum value of 3000 K, and it was set to 3400 K, which is a rough average of the temperatures measured by Perrin et al. (2004) at earlier phases. This yields the global parameters for the model in Table 3. MiRA). As the (u,v) coverage was sparse, we used the same method as in Lacour et al. (2009) to regularize the image reconstruction process. The criterion minimized by MiRA in this section is a weighted sum of the squared distance to the data and of the squared distance to an initial image (SD2I in Appendix A). The initial image was built by fitting the data with a circularly symmetric star + shell model of Perrin et al. (2004), the parameters of which are given in Table 3. The total number of available data is 940 (705 squared visibilities and 235 closure phases). Images of 100x100 pixel were reconstructed. A hyper parameter µ controls the weight of the regularization in the reconstructed image with a total χ 2 reading: χ 2 tot = χ 2 data + µχ 2 rgl . χ 2 data is defined and discussed in Appendix A. The higher µ the higher the weight of the regularization and its impact on the reconstructed image. We tried several values of the hyper parameter and tested the impact on the image reconstructed with MiRA. We tried other reconstructions with MiRA or MACIM using various types of initial images and priors imposed by the regularization. They are discussed in Appendix A and the images are presented in Fig. A.1. For the level of information available at the angular resolution obtained with IOTA, all images share the same basic features: 1) a bright central object of size 20 mas elongated in the north-south direction; 2) surrounded by a fainter environment of 40 mas, roughly circular; 3) a dark patch in the environment that is west of the bright central object. From among the available images, we selected the one with the lowest hyper parameter, µ = 10 8 , as it is more constrained by the data relative to the others. The image is shown in Fig. 4.

The inner envelope
The image is generally consistent with the core-halo description of LPVs developed in Perrin et al. (2004). The two-component, circularly symmetric model correlates well with infrared spectroscopy, recorded in the same spectral region (Hinkle et al. 1984). The spectroscopy records molecular absorption against continuum emission from the effective photosphere. It is natural to associate the molecular absorption zone with the halo in Fig. 4 and to understand the core in the image as a representation of the deeper continuum-forming layer. This association is supported by the dynamical content of the spectra. The motions of the pulsating LPV atmosphere are followed in the Hinkle et al. (1984) study, showing in motions of gaseous CO the expansion, slowing, and falling back of the partially transparent halo. The radial velocity of this layer with respect to the Mira CM velocity (as inferred from the centroid of low excitation envelope emission lines) can be used to follow the halo expansion through a pulsational cycle. Using an average envelope speed of 7.5 km.s −1 Hinkle et al. (1984) for half of the 321-day period, the total extent of the pulsational motion is found to be 1.1 × 10 8 km. This corre- sponds well with the 1.4 × 10 8 km radius difference between the top of the 10 mas model halo and the photosphere (based on the Hipparcos distance of 2.8 × 10 15 km). Thus, to good approximation, the inner boundary of the imaged and modeled halo shows the H-band photosphere and the outer boundary of the halo maps, approximately, the extent of the periodically pulsating gaseous layers. The minimum observed CO excitation temperatures, ≈2700K, are in between the inferred model layer temperatures of 2000K and 3400K.
What remains for us to consider is the relation with the extended molecular layers described by Tsuji (2000). With a molecular excitation temperature of ≈ 1000K, it is, of course, outside the observed and modeled halo. The coupling of the extended layers to the star is not well observed. However, Hinkle & Lebzelter (2015) have shown in the LPV R Cas that gas at an intermediate temperature of 1300K shows a behavior that is intermediate with regard to the pulsating halo and the more stationary extended layers. The intermediate zone shows muchreduced mass motions, that are not coherent with the halo pulsation, and are consistent with "puffs" of mass loss, maintaining line-of-sight coherence over two pulsation cycles, without yet attaining escape velocity and propelling mass loss. While longer spectroscopic time series of additional targets are needed to confirm this finding, it is appears to be a direct observation of the process of populating the extended layers with material raised above the pulsating halo. Expecting a comparable phenomenon in the similar LPV Mira, we suggest that the evidence favors dust formation, or at least critical grain growth, which is likely to occur in either the halo or the extended layers, where a combination of low ambient temperature and relatively enhanced density favors dust formation.
A core-halo structure similar to Fig. 4 has been reported in interferometric imaging of the LPV T Lep (Le Bouquin et al. 2009). An important difference is that the T Lep star is approximately circularly symmetric. The Mira image is not symmetric. In the following, we consider what the asymmetry may have to tell us about the envelope structure and, possibly, about the mass-loss process.

Interpretation of the image
The bright core identified with the photosphere shows an elongated structure. A number of studies have reported strong azimuth dependence in the apparent size of LPV stars. Various models have been proposed, ranging from underlying non-radial oscillation to less dramatic but equally puzzling phenomena such as an azimuthally dependent envelope shape. While non-radial oscillation has not been ruled out for Mira by the new image, the strong asymmetry that is a primary feature of the observation indicates a more complex appearance than a simple elliptical shape of the photosphere or envelope. This is confirmed in the next section with the technical argument that the G. Perrin et al.: Evidence for localized onset of episodic mass loss in Mira data are hard to fit solely with an elongated structure and, in fact, an asymmetrical brightness is required.
A natural question is to ask whether this differential effect represents a dimming or a brightening. This can be answered with reasonable confidence by referring to the AAVSO light curve for the last 15 years reproduced in Fig 5. This shows that the 2005 maximum of Mira was one of the faintest in recent years and more than a magnitude fainter than the second following maximum. At the same time, the image appears to show that most of the stellar disk is rather uniformly bright. We tentatively identify the asymmetry with a differential darkening phenomenon.
A related question considers whether the dimming represents a zone of reduced brightness on the photosphere or a region of increased opacity in the envelope. If the dimming is due to a dark region on the star, a possible explanation would be a long-lived convective structure. Convection is expected to be important in LPVs and present at the same time as pulsation. We note that pulsation does not transport energy effectively and convection is still required for energy transport (Wood 2007). The relative apparent brightness differential across the photosphere (approximately 6×), which implies a similar differential in total surface emittance, also seems too large in the H band for a temperature differential. The dark region covers a significant fraction of the hemisphere and even though large convective cells are predicted, this scale seems unprecedented for a convective effect. The 3D hydrodynamical models of AGB stars combine the effects of pulsation and convection, providing evidence that convection plays a role to produce surface brightness asymmetries (Freytag et al. 2017). These simulations show that dark patches cannot be present alone without bright patches. The spatial convective signature is better characterized as localized, high-contrast bright spots associated with hot, rising material, rather than localized dark spots of falling, cool material. Since the present data set do not demonstrate or require a bright patch, we take this as an additional confirmation for a dark region.
The darkening could be due to an increased opacity in the halo or above. This appears more consistent with the apparent large spatial scale -the dark spot appears to significantly obscure a region of the halo. It is also consistent with previous evidence for irregularities in the near circumstellar region of LPVs, including Mira (Lopez et al. 1997). Therefore, we suggest that the Mira image is showing the effects of a localized enhancement in opacity above the halo. This could be consistent with localization within the intermediate layer feeding the extended molecular layers or with that material itself, with coherent structure lasting over several pulsational periods, but not indefinitely consistent with the observed shell expansion velocity ≈few km/sec - Hinkle & Lebzelter (2015). It is reasonable to find that the cause is absorption by dust located in or near the extended molecular layers as the dust condensation radius is expected in this region, as discussed in, for example, Reid & Menten 1997, Perrin et al. 2004, and Wittkowski et al. 2007. In the next section, an estimate of the mass of the clump with respect to the net mass loss rate leads to another plausible explanation for a location near or somewhat interior to this region, which Tsuji (2000) describe as a MOLsphere.
Previous images of LPVs with large telescope aperture masking (Woodruff et al. 2008) show a rich variety of information about the upper layers of LPVs, but they do not probe asymmetries in the apparent shape of the photosphere well owing to the lack of (u,v) sampling beyond the first visibility null. Multi-telescope interferometry with COAST has reported some significant non-zero closure phases, but did not distinguish the spatial location of the responsible asymmetries (Burns et al. 1998). Hence the inference that asymmetry and inhomogeneity in the circumstellar dust opacity is consistent with previous imaging work, but not well constrained by it. The objective of the next section is to study what constraints on the opacity distribution are provided by the Mira interferometric image.

Dust patch model
The image in Fig. 4 clearly shows asymmetries that require a more complex model than the one used in Sect. 3. We have tried several approaches to explain the data with a central source elongated in the north-south direction and various offsets with respect to the more circularly symmetric shell. These did not lead to satisfactory fits as, in particular, they were inconsistent with the closure phases.
Looking at the reconstructed image, however, it appears generally compatible with the simple circularly symmetric model of Sect. 3 in the sense that it is globally symmetric at low spatial frequency scale and with characteristic sizes of 20 and 40 mas for the star and shell around it. With the model temperatures fixed as described above, the inner and outer diameters and the optical depth of the shell are essentially given by the first visibility lobe and the zero crossing around 50 cycles/arcsec. The parameters of the star+shell model are therefore those of Table 3. The star+shell model has five free parameters. It defines a featureless surface brightness distribution B 0 (x, y). After determining the minimal impact of adding high spatial resolution features to the star+shell model, we decided to separately fit the parameters of Table 3 and the parameters of the more complex model described below. This method with two separate steps led to satisfactory results.
To fit the closure phases and account for the azimuth dependency of the squared visibilities, we considered a dark patch to the west combined with the star+shell model. We tried a single patch as increasing the number of the patches made the fit more complex but did not bring improvement.
Experimenting with several patch options, a smooth edge better matches the data. We used two types of smoothing functions: exponential and Gaussian. Technically, the star+shell model is multiplied by an attenuation function: ( The dark patch model then has four parameters: 1) the patch offset α and δ with respect to the star center; 2) the attenuation strength a (between 0 and 1); 3) the patch FWHM. The total number of parameters of the star+shell+patch model is 9 and the corresponding surface brightness distribution is B(x, y) with B(x, y) = B 0 (x, y)× f (x, y).
The closure phases outside the ±150 • range were found to have a diverging impact on the fit as they are very sensitive to noise. We therefore decided to crop them to search for the most robust solution. The best parameters we found are listed in Table 4. The error bars were derived by varying the χ 2 by 1. It is clear, however, that these are lower limits since the model does not perfectly account for the data. We did not consider it useful to refine error bars as the goal is to obtain orders of magnitude for physical parameters of the dust clump later in Sect. 5.3 and as we used two different models for the patch which give an idea of the scatter for these parameters as reported in Table 5.
The surface brightness distributions, residues, and closure phase fits are presented in Fig. 6. The fit to the closure phases is quite good, even at high spatial frequencies, with residuals smaller than 25 • . The shape of the closure phases is consistent with a smooth spot as a uniform spot did not provide such nice fit. The V 2 residuals at high visibility are quite large. We investigated this issue by fitting the low spatial frequency data (below 40 cycles/arcsec) with a uniform disk model with an azimuth-dependent diameter (10 • bins). The amplitude of the scatter remained the same albeit more symmetric around 0. We concluded that the scatter is real and does not depend on the model we chose to interpret the data. Surface brightness distribution, V 2 residues, closure phases (red), model closure phases (black), and closure phases residues. The closure phase and closure phase residues are displayed as a function of the maximum spatial frequency of the triangle of baselines.
The two models yield qualitatively similar results and although the χ 2 for the Gaussian model is a bit lower, both were retained for further consideration. The total χ 2 may look large but the two branches of models capture most of the visibility and closure phase features that are in very good agreement with the reconstructed image so that we consider that it is a satisfactory description of the object.

Interpretation of the attenuation function
We concluded that the dark patch is probably located above the pulsation/convective zone and is accounted for by absorption from a patch of material. It therefore equally absorbs the radiation from the star and from the inner envelope. We neglect the H-band emission by the patch as it is at a higher, cooler elevation. The attenuation function can be locally equated to an optical depth effect: f (x, y) = e −τ(x,y) .
G. Perrin et al.: Evidence for localized onset of episodic mass loss in Mira We can therefore derive the local optical depth from the ratio of the obscured to unobscured surface brightnesses. The local optical depth can then be converted into a mass of absorbing material. Let ρ [g.cm −3 ] be the density of the absorbing material, κ [g −1 .cm 2 ] the opacity and z(x, y) the local depth of absorbing material. In the following relation: where ρz(x, y) is the mass column density in g.cm −2 . The integration of this quantity over the stellar disk, therefore, yields the total mass m of absorbing material, m: The integral can be evaluating using the surface brightnesses of the star+shell and star+shell+patch models: We can also evaluate the average attenuation by the patch with the ratio of the integrated surface brightnesses: The results are given in Table 5. Uncertainties were derived with the Monte-Carlo technique using the values of Table 4. The two patch models yield similar results. The patch causes an average attenuation of the H-band flux of 0.925±0.007 (minimum relative surface brightnesses of 2 × 10 −3 and 2 × 10 −1 for the exponential and Gaussian respectively, compatible with the dynamic range of the reconstructed image). We use the average value of (6.31 ± 0.72) × 10 25 cm 2 for the mass-opacity product (optical depth integrated over the star surface), denoted as m.κ H from this point forward.

Interpretation of the attenuation as a dust opacity effect
As discussed in Sect. 4, our favored hypothesis to explain the lower surface brightness at the spot is dimming by a dust patch. Al 2 O 3 is a good candidate for dust formation at high temperature in the molecular shell of Mira stars (see e.g. Perrin et al. (2015)). But the condensation of Al 2 O 3 at distances and concentrations implied by observations requires high transparency of the grains in the visual and near-IR region to avoid destruction by radiative heating (Höfner et al. 2016). We therefore consider silicate grains to explain the obscuration. Amorphous Mg 2 SiO 4 grains are the main candidate for driving the winds of O-rich AGBs according to Bladh & Höfner (2012). We interpret the dimming in H measured with IOTA/IONIC and in V from the AAVSO as being due to this single dust component. Within the scope of the single-grain type hypothesis, the relative V/H extinction constrains the dust grain characteristics. We have decided to consider two types of magnesium-silicate dust grains: MgSiO 3 and Mg 2 SiO 4 . As a remark, we note that Khouri et al. (2018) did their modeling with aluminum oxides, adding that magnesium-silicate dust grains have scattering properties that are similar to those of Al 2 O 3 .
Opacities can be obtained from the literature for specific grain sizes and types. We performed calculations of the opacity of dust in the Mie approximation using the Fortran code, bhmie, from Bohren & Huffman (1983), which was modified by Bruce Draine (http://scatterlib.wikidot.com/mie). Optical constants were taken from Jäger et al. (2003). The numerical values are listed at the Jena group web page: http://www.astro.uni-jena.de/Laboratory/OCDB/amsilicates.html (OCDB database).
The extinction efficiency Q ext was computed at 0.55 and 1.65 µm wavelengths. In the Mie approximation, the extinction efficiency can be translated to an opacity via the equation: where a is the radius of the dust sphere and ρ is the density, all in CGS units. The typical density of 3 g.cm −3 from the OCDB database was used. The opacities were computed for radii ranging from 0.5 nm to 100 µm to explore the different regimes of the Mie scattering. The plots of Fig 7 show the opacities for the two species of dust as a function of radius at 0.55 and 1.65 µm as well as the κ V /κ H ratio. The opacities are quite similar for the two species of dust except for radii below 0.01 µm.
The AAVSO data of Fig 5 show that the two V-band maxima before and after the median date of the observation (JD2453654) are brighter, respectively, by a 0.4 and 0.0 magnitude. Assuming the difference in maximum magnitude is entirely due to the obscuration of a persisting dust cloud, the V-band data give an upper limit on the attenuation due to dust of 0.4 magnitude. We made the following assumptions to compute the flux attenuation by the dust patch in V: 1) the stellar brightness in the visible is dominated by the brightness of the molecular shell, the photosphere is not visible: this is supported by the large optical depth of the shell at visible wavelengths derived both from observations (diameters in the visible are much larger, see for example the first observational evidence by Labeyrie et al. (1977)) and by the 3D hydrodynamical simulations; 2) the diameter of the molecular shell is the same as in H; 3) the dark patch characteristics are the same as at H except for the optical depth.
With these assumptions, the maximum global flux attenuation can be computed as a function of the strength parameter a in Eq. 1. The maximum attenuation is 5% with the Gaussian model, whereas it reaches 10% with the exponential model. In consequence, we used the exponential model to explore the variations of the attenuation and of the m.κ V product versus a. The maximum for the geometrical characteristics of the patch derived in the H band is an attenuation of 10% of the flux corresponding to a difference in magnitude of ∆V = 0.11. This is compatible with the maximum brightness decrease with respect to the two neighboring V-band maxima. This cannot, however, be the only source of the dimming in V as the AAVSO data show cycle-to-cycle variations up to a magnitude or more. Unless several such patches are episodically produced. Here, we acknowledge the simplicity of our approach as the reasons for the cycle-to-cycle V magnitude variations of Mira are likely to be multiple and not strictly due to dust patches occurring at the surface of the star. Our simple hypothesis, therefore, provides an upper value for the opacity of the dust patch in V. Following this approach, it yields a maximum mass-opacity product at visual wavelengths of m.κ V = 1.49 10 26 cm 2 and, thus, a maximum ratio for the dust opacity between V and H of: 49 × 10 26 (6.31 ± 0.72) × 10 25 2.36 ± 0.27.
The relative fraction of V/H extinction assumed does not strongly impact the deductions concerning dust mass. In any event, the inferred V-band obscuration cannot explain the large variations of the maximum visual magnitude shown by the AAVSO measurements.
According to Fig. 7, this is compatible with both MgSiO 3 and Mg 2 SiO 4 dust. The particle radii ranges are shown in Table 6, along with the implied total mass of the dust in the patch, and also the total Article number, page 7 of 12page.12 A&A proofs: manuscript no. Mira_GP_21_08_2020 mass of the associated gas -assuming a gas/dust ratio of 89.5 (Knapp 1985).

Discussion
In the previous sections, we modeled the IOTA/IONIC data on Mira in the H band using a patch of absorbing material obscuring the corehalo, that is, the modeled H-band photosphere and pulsating envelope. We show, based on a combination of these data with visual photometric data, that this material is compatible with silicate dust. As summarized in Table 6, the inferred dust mass associated with the absorbing material is strongly dependent on the dust physical characteristics. Our data set is limited and the characteristics of the dust remain rather speculative in this analysis and further interpretation must rely on indirect evidence. Since micron size grains are commonly inferred from models of the dust at several AU from an LPV (Ohnaka et al. 2016), such a grain size might be favored for the observed patch. However, from Table 6, the associated dust/gas mass would be only ≈ 10 −12 M . This is much less than the known yearly mass loss rate for Mira of 3 × 10 −7 M /yr (Ryde et al. 2000). In order for such a patch to contribute significantly to the observed mass loss rate, it would be necessary to suggest an implausibly large number of patches or high frequency and speed of ejection. So we tentatively reject the large grain option for the observed patch.
There is a precedent for considering small grains. Recent studies of the mira W Hya by Ohnaka et al. (2017) show a shift toward smaller grain sizes near photometric minimum (0.54), compared to their earlier study at 0.77 (Ohnaka et al. (2016)). This is not surprising since the minimum brightness corresponds approximately to maximum expansion and greatest atmospheric cooling. Dynamical models reported by Höfner et al. (2016), Aronson et al. (2017) and Liljegren et al. (2017) indeed show dust formation near the brightness minimum. The observations reported here are also near minimum.
From Table 6, a small grain size would be consistent with a reservoir of gas/dust mass in the patch ≈ 10 times larger than the annual stellar mass loss if the gas/dust ratio can be applied at short distance from the photosphere and comparable to or a bit below the yearly mass loss rate otherwise. Therefore, it is reasonable to consider the possibil-ity that the patch can be linked with the mass ejection process and what this implies for the ejection efficiency.
Only a single patch was observed. However, from the model suggesting that the patch is near but outside the halo and from the knowledge that approximately 1/3 of the progeny of LPVs show an asymmetry (Bobrowsky 2007, for Mira, here, attributed to dust patches), it can be estimated from geometry and assumed random distribution that at any given time there may be on the order of five such patches around Mira, of which zero or one would typically be in the line of sight to the halo at any given time. Thus, the elevation of a few percent of each of the mass of each of these "puffs" of material over each pulsation cycle, on average, could support the observed Mira mass loss rate. Chandler et al. (2007) reported their detection of episodic mass loss around Mira at 11.15 µm with the ISI interferometer. Our detection in H band is compatible with their finding.
With low envelope expansion velocities, compared to the escape velocity, and being on the same order as those observed in the extended molecular layers, puffs of rising dust-enhanced material would have ample time for grain growth during the elevation to the vicinity of 2+ stellar radii, where much larger grain sizes (with much reduced opacity) are found by modeling spatially resolved polarimetric observations: Norris et al. (2012) have detected 0.3 µm grains of Mg 2 SiO 4 or MgSiO 3 for W Hya, R Dor and R Leo at 2 stellar radii or less; Ohnaka et al. (2016) detected 0.4-0.5 µm grains of Al 2 O 3 , Mg 2 SiO 4 or MgSiO 3 between ≈2-3 stellar radii in clumpy dust clouds for W Hya; Adam & Ohnaka (2019) detected 0.1 µm grains of Al 2 O 3 , Mg 2 SiO 4 , or MgSiO 3 in a clumpy shell around IK Tau between 2 and 5 stellar radii; Khouri et al. (2016) measured the asymmetric photosphere of R Dor attributed to scattering by dust. Finally, Khouri et al. (2018) observed clumpy structures in visible polarized light around Mira at larger distance than our measurements and within ≈8 photospheric radii. They used aluminum oxide grains to model their data and concluded that aluminum could be in dust phase if the grain size is less than 0.02 µm. All other existing observations, apart from our own and those of Khouri et al. (2018), detect a few clouds of larger dust particles that account for a good fraction of the mass loss. This is a discrepancy that is to be resolved with more measurements. The greater grain size at higher elevation satisfies the pre-Article number, page 8 of 12page.12 G. Perrin et al.: Evidence for localized onset of episodic mass loss in Mira ferred conditions to drive mass loss by radiation pressure on the grains (Höfner 2008), that is, 0.36-0.66 µm at 2-3 R .
These findings give rise to a few more questions. The first two regard the nature of the mechanism that could explain the existence of such a clump and how long can it exist in the atmosphere or above. A binary companion sufficiently close to Mira A could reduce the local effective gravity enough to enhance mass loss and it could also account for the symmetric structures commonly observed in planetary nebulae. The known companion Mira B is too far away to provide a gravity reduction by more than 1%. A much closer companion would be needed, even a massive planet, but there is no evidence for another very close companion at this point. We find the dust clump to be on the side opposite to Mira B as Chandler et al. (2007) did for the excess of dust emission measured with ISI. Khouri et al. (2018) did not find a systematic asymmetry at larger distance. They found a dust trail connecting Mira A and Mira B; however, there is no clear evidence for an influence of Mira B on the location of sites for dust production.
The size of the clump is larger than convective cells, such as those simulated in 3D hydrodynamical codes. However, the collective impact of complex convective structure could produce larger scale irregularities. The coupling between pulsations and convection could be invoked to explain such a clump of dust. A pioneer study of this coupling was performed by Anand & Michalitsanos (1976) with a mechanical analog of a giant or a supergiant star. Their conclusion was that strong interactions can be expected with asymmetrical fluctuations of the surface brightness. Kiss et al. (2000) invoked this mechanism as a possible explanation for period and amplitude changes revealed in photometric monitoring. A single pulsation with effective convective coupling in the vicinity of a single or few large convective cells could well enhance the efficiency of mass ejection into the envelope, generating an over-density. This over-density would presumably follow a ballistic trajectory, with much of the material falling back toward the photosphere, where it would encounter the next levitating shock. It is easy to suppose that this pattern could preserve the over-density for some time, until it is reduced by leakage up, especially into the extended molecular layers, with a timeframe that is consistently similar to the pulsational period. Freytag & Höfner (2008) performed 3D numerical simulations of an AGB star combining convection, pulsation, and dust formation which confirm the Anand & Michalitsanos (1976) predictions. In their simulations, the bright convective cells provide a strongly asymmetric brightness distribution at and above levels they describe to be the top of the convective envelope which we believe corresponds approximately to what we have referred to as the core, or the depth to which near-IR imagery can penetrate. Asymmetries in the density distribution appear most strongly in the highest levels, which are detected in the near-IR through their partial transparency and molecular absorption, which vary through the pulsational period. Freytag & Höfner (2008) note a typical lifetime for convective structure of years, compared with a pulsational time of one year. Something that could prove particularly significant for the upper envelope opacity, Freytag & Höfner (2008) predicted localized dust formation in regions of enhanced density and depressed temperature which occur naturally in the complex flows and shocks. This has been confirmed by Höfner & Freytag (2019) with more recent 3D simulations of atmospheric shock waves generated by large-scale convective flows and pulsations. They find that the clumpy environment of dust can form around M-type AGB stars. This is clearly compatible with what we report here and with the result of Ohnaka et al. (2016); Adam & Ohnaka (2019); Norris et al. (2012). This is also compatible with the dimming of the visual photometry at the time of our observations and during the few cycles that came before.

Conclusion and outlook
The Mira image we obtained in H band at high spatial resolution offers supporting evidence for a possible outcome of pulsation-convection coupling and offers a possible direction for pursuing the evidence further.
Confirmation should be possible with a series of interferometric images obtained with spectral information and polarimetry, possibly combined with visible resolved polarimetry, as already shown in some works cited in this paper. The three requirements include: excellent image quality, a combination between spectral and spatial resolution, and a suitable time series with a few month time interval over a few periods, that is, over several years. The quality of the images can be addressed in two ways.
First, for the higher spatial frequency baselines, good (u,v) coverage is essential. This is possible at present only at the VLTI, where the large selection of baselines could be optimized for a target list, although the need for rather contemporaneous UV coverage may be a challenge. Second, as is clearly shown by the aperture masking and other results mentioned above, these stars have a lot of low spatial frequency power. It is essential to measure this correctly in order to properly interpret the high frequency measurements. This can be provided in the visible or nearinfrared with extreme adaptive optics as shown by Ohnaka et al. (2016) with SPHERE or by aperture masking like the SAM mode of SPHERE on the VLT (Cheetham et al. 2016) and with MICADO on the ELT in the future (Lacour et al. 2014), or with short-exposure imaging with adaptive optics in a moderate speckle mode (Kervella et al. 2009). Midinfrared wavelengths can be reached without adaptative optics (Poncelet et al. 2007) in this same mode. Near-IR polarimatric measurements are needed to simultaneously constrain scattering and absorption effects of dust. Spectral resolution is required to constrain temperatures and to separate the continuum (or pseudo-continuum) from the regions of the spectrum dominated by molecular opacity. Spectral resolutions as high as 4000 can be obtained with GRAVITY on the VLTI (Gravity Collaboration et al. 2017) in the K band. With such data, we could envisage the detection of a correlation between the brightness distribution on the photosphere and the evolving opacity in the envelope over subsequent years. The Freytag & Höfner (2008) models specifically predict episodes of dust formation in the envelope, as close as 1-2 stellar radii from the surface. The combination of continuum and molecular opacities will help distinguish between dust and other continuous opacities. Some preliminary results on correlations between the visual light curve and the diameter measured in different band passes in the K band were recently published by Wittkowski et al. (2018) using GRAVITY. In addition, these data could be combined with mid-infrared four-telescope observations with the MATISSE instrument on VLTI (Lopez et al. 2014) to further characterize the properties of the dust. MACIM uses a simulated annealing algorithm with the Metropolis sampler. For the kind of considered data, the objective function is not convex and the minimum of J(x) over Ω may not be unique. The result produced by MiRA and MACIM thus depend on the initial image, x init , and on the regularization settings, µ and J rgl (x).
For the reconstructions with MiRA, we considered three different regularizations. The first one, which measures the squared distance of the image to the initial image is given by: The second one, favors compact images: with θ j the 2-dimensional angular position of pixel j relative to the center of the field of view, γ > 0 the full width at half maximum (FWHM) of the assumed prior image (θ j and γ being in the same units). Attempting to minimize J Compact (x) for x ∈ Ω imposes that the brightest parts of the sought image be mostly concentrated around the center of the field of view. Hence this regularization corresponds to a "compactness" prior (Thiébaut & Young 2017). The third one is given by:  where ∇ j x denotes the 2D spatial gradient of x at pixel j and τ > 0 is a selected threshold. This kind of regularization penalizes the image gradients and thus imposes a smooth image. The penalization is however less severe for the image gradients larger than τ in magnitude than for the image gradients smaller than τ in magnitude which helps preserving the edges in the image. Hence this regularization imposes an edge-preserving smoothness prior.
For MACIM, the regularizer J Dark Energy (x) is a dark interaction energy regularizer. This regularizer is the sum of all pixel boundaries with zero flux on either side of the pixel boundary. It gives higher weight to large regions of dark space in between regions of flux.
In addition to the image presented in Fig. 4 and reconstructed with the SD2I regularization and an intial image obtained by fitting a circularly symmetric star+ shell model to the V 2 data, we reconstructed four other images with MiRA and one with MACIM. These images are presented in Fig. A.1 together with the image of Fig. 4. The details for each image are listed in Table A.1. Images of 100x100 pixels were reconstructed with MiRA while 40x40 pixels were used for MACIM. The pixel size is 1 mas for MiRA and 1.44 mas for MACIM. As explained in Sect. 4, different µ values were tried for the reconstruction of image 1 with MiRA and the best result was obtained for µ = 10 8 . For images 3 to 6, also reconstructed with MiRA, the value of µ was adjusted to minimize the χ 2 . The parameter τ is in units of intensity per pixel and therefore depends on the pixel size and on the level of the flux normalization (Renard et al. 2011). For the image reconstructed with MACIM, µ was set to 1.
The reduced χ 2 are given Table A.1. The number of degrees of freedom is the number of independent parameters subtracted to the number of data (940, i.e. 705 V 2 and 235 closure phases). The number of independent parameters is not as easy to define as in classical χ 2 fitting. The maximum number of independent parameters is the number of pixels, e.g. 10,000 for the images reconstructed with MiRA. This number is larger than the number of data but not all pixels are independent because 1) they oversample the angular resolution elements in the image and 2) the regularization induces correlations between pixels. What matters here are the relative values to compare the images, the absolute values of the reduced χ 2 are therefore only indicative. We have chosen to define the reduced χ 2 as: χ 2 r = χ 2 /#data = χ 2 /940. Our definition of χ 2 r is thus an underestimation of the reduced χ 2 . The errors provided by the pipeline have been modified for image 2 reconstructed with MACIM: a 5% relative error was applied to the V 2 data with an additive default floor of 0.0001 and a minimum error of 0.5 • was used for the closure phase data. Raw error bars were used for the images reconstructed with MiRA (1 and 4 to 6). We have computed the reduced χ 2 in two different ways: for χ 2 r,1 we compared the reconstructed image with the V 2 and closure phase data using the raw error bars while for χ 2 r,2 the modified error bars were used. Values are given in Table A.1. The lower values of the χ 2 r,2 estimator are a good indicator that the error bars produced by the pipeline are underestimated as noted in Sect. 5.1. The reduced χ 2 r,2 are all close to 1 meaning all images are reconstructed with the same quality with respect to the data.
In addition we have reconstructed the dirty beam with MiRA ( Fig. A.2). All the visibilities of the sampled (u, v) points were set to 1 for this image. The white ellipse on the image gives the FWHM of the central peak and gives the scale of the resolution in the reconstructed images.
In spite of the aforementioned issue of the solution not being unique, the reconstructions in Fig. A.1 clearly show that the result does not strictly depend on the given algorithm that was used. We see no major influence of the initial image on the final image, in particular, there is no influence coming from a circular or absence of symmetry. The Compact regularization with MiRA makes the brightest peak appear a bit larger compared to the other reconstructions but does not change the image substantially. The features in the image cannot be explained by the structures of the dirty beam. Overall, the reconstructed images are very similar and exhibit the same remarkable structures with a 20 mas north-south elongated central object surrounded by a fainter 40 mas circular environment and a dark patch in the environment that is west of the bright central peak. We conclude that the image reconstruction process is quite robust despite the limited (u, v) coverage. Image 1 featured in Sect. 4 is fairly representative of all the reconstructed images.