# The inner dust shell of Betelgeuse seen with high-angular-resolution polarimetry

X. Haubois, R. G. van Holstein, J. Milli, C. Pinte, A. López-Ariste, Ph.

Mathias, P. Kervella, G. Perrin, M. Montargès

► To cite this version:

X. Haubois, R. G. van Holstein, J. Milli, C. Pinte, A. López-Ariste, et al. The inner dust shell of Betelgeuse seen with high-angular-resolution polarimetry. Astronomy and Astrophysics - A&A, 2023, 679, 10.1051/0004-6361/202243458. insu-04473142

HAL Id: insu-04473142

https://insu.hal.science/insu-04473142

Submitted on 23 Feb 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# The inner dust shell of Betelgeuse seen with high-angular-resolution polarimetry*

X. Haubois<fnref n="1" />, R. G. van Holstein<fnref n="1" />, J. Milli<fnref n="2" />, C. Pinte<fnref n="2" />,<fnref n="3" />, A. López-Ariste<fnref n="4" />, Ph. Mathias<fnref n="5" />, P. Kervella<fnref n="6" />, G. Perrin<fnref n="6" />, and M. Montargès<fnref n="6" />

Received 3 March 2022 / Accepted 9 August 2023

## ABSTRACT

*Context.* The characteristics of the innermost layer of dust winds from red supergiants have not been identified. In 2019–2020, Betelgeuse exhibited an important dimming event that has been partially attributed to dust formation, highlighting the importance of understanding dust properties in the first stellar radii from the photosphere.

*Aims.* We aim to detect and characterize the inner dust environment of Betelgeuse at high spatial resolution.

*Methods.* We obtained SPHERE/ZIMPOL and SPHERE/IRDIS linear polarimetric observations from January 2019, before the dimming event, and compared them to a grid of synthetic radiative transfer models.

*Results.* We detect a structure that is relatively centro-symmetric with a 60 mas diameter (1.3–1.4 stellar diameter). We computed synthetic images using radiative transfer modeling assuming a spherical dust shell composed of MgSiO₃ grains. We find that most of the data are best reproduced with a dust shell whose outer radius is approximately 10 AU (i.e., ~2 stellar radii) and a maximum grain size in the 0.4–0.6 µm range. These results are close to the ones we obtained from 2013 NACO/SAMPOL data, indicating that the shell radius and grain size can show some stability for at least 6 yr despite morphological changes of the dust shell. The residuals after the subtraction of the best-fitting centro-symmetric model suggest complex asymmetric density structures and photospheric effects.

**Key words.** techniques: high angular resolution – supergiants – infrared: stars – stars: individual: Betelgeuse

## 1. Introduction

Dust is a vector and product of the mass-loss phenomenon at play in red supergiants (RSGs). This mechanism has an important role in the evolution of RSGs into type-II supernovae as well as in galactic chemical enrichment. The formation loci and properties of the first stage of a dust wind have not been established. Evidence for dust has been reported around the emblematic M2 type RSG Betelgeuse (α Orionis, HD39801) in its very close environment. Reconstructed images from speckle interferometry show an asymmetric structure at ~2–2.5 stellar radii (Roddier & Roddier 1985). More recent mid-infrared observations revealed a complex dust environment (Kervella et al. 2011) compatible with a silicate or alumina composition, confirming previous interferometric diagnostics from Verhoelst et al. (2006) and Perrin et al. (2007).

In December 2019, Betelgeuse underwent a dimming event that has been reported and analysed in various papers (e.g., Dupree et al. 2020; Montargès et al. 2021; Kravchenko et al. 2021; Taniguchi et al. 2022). However, it is not clear if this event can be entirely attributed to dust formation (Harper et al. 2020), and this has emphasized the importance of understanding the

conditions for dust formation when interpreting variability from RSGs.

Polarimetry at high angular resolution is a powerful technique for investigating circumstellar dust environments. In Haubois et al. (2019, hereafter H2019), we modeled a structure detected in linearly polarized light as a thin dust shell located only 0.5 stellar radii above Betelgeuse’s photosphere (corresponding to a dust shell radius of about 7.3 AU with a distance of 222 pc; Harper et al. 2017) from sparse aperture-masking polarimetric interferometry data (NACO/SAMPOL). This so-called inner dust shell (IDS) poses interesting constraints on the survival of dust grains. At that short distance from the photosphere, the local temperature can reach the 1500–2000 K range, implying dust material with high condensation and sublimation temperatures. From comparisons with SPHERE/ZIMPOL measurements, we show in H2019 that radiative transfer models favor a dust composition made of MgSiO₃ (enstatite) and Mg₂SiO₄ (forsterite) over Al₂O₃ (alumina). In an attempt to elaborate on the findings presented in H2019, we report here on SPHERE/ZIMPOL and SPHERE/IRDIS kinear polarimetric observations obtained January 18–19, 2019. After presenting the observations and their analysis in Sects. 2 and 3, respectively, we compare the images of polarized intensity with radiative transfer models in Sect. 4. Subsequently, we present the discussion and conclusions in Sects. 5 and 6, respectively.

<footnote n="1">¹ European Organisation for Astronomical Research in the Southern Hemisphere, Casilla 19001, Santiago 19, Chile
e-mail: xhaubois@eso.org</footnote>

<footnote n="2">² Univ. Grenoble Alpes, CNRS, IPAG, 38000 Grenoble, France</footnote>

<footnote n="3">³ School of Physics and Astronomy, Monash University, Clayton, VIC 3800, Australia</footnote>

<footnote n="4">⁴ IRAP, Université de Toulouse, CNRS, CNES, UPS, 14, Av. E. Belin. 31400 Toulouse, France</footnote>

<footnote n="5">⁵ IRAP, Université de Toulouse, CNRS, UPS, CNES, 57 avenue d'Azereix, 65000 Tarbes, France</footnote>

<footnote n="6">⁶ LESIA, Observatoire de Paris, Université PSL, CNRS, Sorbonne Université, Université de Paris, 5 place Jules Janssen, 92195 Meudon, France</footnote>

<footnote>* Based on SPHERE data obtained at the ESO VLT Melipal telescope (program 0102.D-0505).</footnote>
# The inner dust shell of Betelgeuse seen with high-angular-resolution polarimetry*

X. Haubois<fnref n="1" />, R. G. van Holstein<fnref n="1" />, J. Milli<fnref n="2" />, C. Pinte<fnref n="2" />,<fnref n="3" />, A. López-Ariste<fnref n="4" />, Ph. Mathias<fnref n="5" />, P. Kervella<fnref n="6" />, G. Perrin<fnref n="6" />, and M. Montargès<fnref n="6" />

Received 3 March 2022 / Accepted 9 August 2023

## ABSTRACT

*Context.* The characteristics of the innermost layer of dust winds from red supergiants have not been identified. In 2019–2020, Betelgeuse exhibited an important dimming event that has been partially attributed to dust formation, highlighting the importance of understanding dust properties in the first stellar radii from the photosphere.

*Aims.* We aim to detect and characterize the inner dust environment of Betelgeuse at high spatial resolution.

*Methods.* We obtained SPHERE/ZIMPOL and SPHERE/IRDIS linear polarimetric observations from January 2019, before the dimming event, and compared them to a grid of synthetic radiative transfer models.

*Results.* We detect a structure that is relatively centro-symmetric with a 60 mas diameter (1.3–1.4 stellar diameter). We computed synthetic images using radiative transfer modeling assuming a spherical dust shell composed of MgSiO₃ grains. We find that most of the data are best reproduced with a dust shell whose outer radius is approximately 10 AU (i.e., ~2 stellar radii) and a maximum grain size in the 0.4–0.6 µm range. These results are close to the ones we obtained from 2013 NACO/SAMPOL data, indicating that the shell radius and grain size can show some stability for at least 6 yr despite morphological changes of the dust shell. The residuals after the subtraction of the best-fitting centro-symmetric model suggest complex asymmetric density structures and photospheric effects.

**Key words.** techniques: high angular resolution – supergiants – infrared: stars – stars: individual: Betelgeuse

## 1. Introduction

Dust is a vector and product of the mass-loss phenomenon at play in red supergiants (RSGs). This mechanism has an important role in the evolution of RSGs into type-II supernovae as well as in galactic chemical enrichment. The formation loci and properties of the first stage of a dust wind have not been established. Evidence for dust has been reported around the emblematic M2 type RSG Betelgeuse (α Orionis, HD39801) in its very close environment. Reconstructed images from speckle interferometry show an asymmetric structure at ~2–2.5 stellar radii (Roddier & Roddier 1985). More recent mid-infrared observations revealed a complex dust environment (Kervella et al. 2011) compatible with a silicate or alumina composition, confirming previous interferometric diagnostics from Verhoelst et al. (2006) and Perrin et al. (2007).

In December 2019, Betelgeuse underwent a dimming event that has been reported and analysed in various papers (e.g., Dupree et al. 2020; Montargès et al. 2021; Kravchenko et al. 2021; Taniguchi et al. 2022). However, it is not clear if this event can be entirely attributed to dust formation (Harper et al. 2020), and this has emphasized the importance of understanding the

conditions for dust formation when interpreting variability from RSGs.

Polarimetry at high angular resolution is a powerful technique for investigating circumstellar dust environments. In Haubois et al. (2019, hereafter H2019), we modeled a structure detected in linearly polarized light as a thin dust shell located only 0.5 stellar radii above Betelgeuse’s photosphere (corresponding to a dust shell radius of about 7.3 AU with a distance of 222 pc; Harper et al. 2017) from sparse aperture-masking polarimetric interferometry data (NACO/SAMPOL). This so-called inner dust shell (IDS) poses interesting constraints on the survival of dust grains. At that short distance from the photosphere, the local temperature can reach the 1500–2000 K range, implying dust material with high condensation and sublimation temperatures. From comparisons with SPHERE/ZIMPOL measurements, we show in H2019 that radiative transfer models favor a dust composition made of MgSiO₃ (enstatite) and Mg₂SiO₄ (forsterite) over Al₂O₃ (alumina). In an attempt to elaborate on the findings presented in H2019, we report here on SPHERE/ZIMPOL and SPHERE/IRDIS kinear polarimetric observations obtained January 18–19, 2019. After presenting the observations and their analysis in Sects. 2 and 3, respectively, we compare the images of polarized intensity with radiative transfer models in Sect. 4. Subsequently, we present the discussion and conclusions in Sects. 5 and 6, respectively.

<footnote n="1">¹ European Organisation for Astronomical Research in the Southern Hemisphere, Casilla 19001, Santiago 19, Chile
e-mail: xhaubois@eso.org</footnote>

<footnote n="2">² Univ. Grenoble Alpes, CNRS, IPAG, 38000 Grenoble, France</footnote>

<footnote n="3">³ School of Physics and Astronomy, Monash University, Clayton, VIC 3800, Australia</footnote>

<footnote n="4">⁴ IRAP, Université de Toulouse, CNRS, CNES, UPS, 14, Av. E. Belin. 31400 Toulouse, France</footnote>

<footnote n="5">⁵ IRAP, Université de Toulouse, CNRS, UPS, CNES, 57 avenue d'Azereix, 65000 Tarbes, France</footnote>

<footnote n="6">⁶ LESIA, Observatoire de Paris, Université PSL, CNRS, Sorbonne Université, Université de Paris, 5 place Jules Janssen, 92195 Meudon, France</footnote>

<footnote>* Based on SPHERE data obtained at the ESO VLT Melipal telescope (program 0102.D-0505).</footnote>
## 2. Observations

We observed Betelgeuse on January 17 and 18, 2019, with SPHERE (Beuzit et al. 2019) at the Very Large Telescope at visible wavelengths with ZIMPOL (Schmid et al. 2018) and at near-infrared wavelengths with IRDIS (Dohlen et al. 2008; de Boer et al. 2020; van Holstein et al. 2020). Data were taken in linear polarimetric modes in ten filters, including those of a point spread function (PSF) calibrator star, Phi² Orionis. The total intensity images are presented in Fig. 1, and Table 1 summarizes the observations. Atmospheric conditions were good and stable throughout the Betelgeuse and calibrator observations. However, we experienced a low-wind effect (Milli et al. 2018). This had the impact of decreasing the Strehl ratio (SR) as measured from the SPHERE images compared to the SPARTA SR measurements. The ZIMPOL and IRDIS data were reduced by the SPHERE Data Center (Delorme et al. 2017) and using the IRDAP<fnref n="1" /> pipeline (van Holstein et al. 2017, 2020), respectively. For ZIMPOL data, we estimated the instrumental polarization to be less than 0.1% for all filters by computing it over the halo of starlight in the images, a result consistent with the findings from Schmid et al. (2018). With IRDAP, we subtracted the instrumental polarization using the Mueller matrix model and then subtracted the stellar polarization from the halo of starlight. We followed the IRDAP convention in terms of orientation of polarization direction. From the *Q* and *U* images, we then computed images of the azimuthal Stokes parameters *Q*<sub>φ</sub> and *U*<sub>φ</sub> following the definitions of de Boer et al. (2020) as well as images of the linearly polarized intensity, *PI* = √(*Q*<sup>2</sup> + *U*<sup>2</sup>). The images of the Stokes parameters *Q* and *U* reveal “butterfly patterns” in all filters (see Fig. 2), as expected for centro-symmetric polarization created from scattering off a face-on structure such as a disk or a sphere.

In order to remove instrumental polarization (IP) effects, it is possible to use techniques such as minimizing the azimuthal Stokes parameter, *U*<sub>φ</sub>. (Avenhaus et al. 2014; de Boer et al. 2020). Because there should be no radial polarization, we applied a minimization technique to the negative *Q*<sub>φ</sub> signal in order to remove residual IP in our ZIMPOL data set. Due to the better spatial resolution delivered by ZIMPOL, there were enough pixels in the images to efficiently remove spurious polarization. This was not possible in the case of the IRDIS data set. We therefore used the *U*<sub>φ</sub> minimization technique in that case as the IP-generated butterfly pattern clearly dominated in these images.

In the V-band filter, we observed that the butterfly pattern had a different orientation from all the other filters. The reason for this could in principle be astrophysical as, for thick dusty stellar environments such as protoplanetary disks, multiple scattering can lead to non-azimuthal linear polarization effects (azimuthal Stokes *U*<sub>φ</sub> signal is non-null; see Canovas et al. 2015). However, in our case, the angle of linear polarization (AoLP) is azimuthal at all the other wavelengths from the visible to the near-infrared, so it is more likely that the cause for this V-band-specific polarization effect is of instrumental origin. Having corrected for instrumental polarization, we suggest that this issue originates from a so-far uncharacterized instrumental crosstalk.

The polarimetric design of ZIMPOL (Schmid et al. 2018) includes two half-wave plates (HWP2 and HWPZ) that together prevent the creation of instrumental crosstalk by SPHERE's

<figure><img src="image_1.png" /><figcaption>**Fig. 1.** Total intensity in the ten different filters. $\lambda$ is the central wavelength of each filter. North is up and east to the left. Units are W m$^{-2}$ sr$^{-1}$. The dashed orange and white circles mark angular sizes of 500 and 60 ms, respectively.</figcaption></figure>

<footnote n="1">¹ https://irdap.readthedocs.io (IRDIS Data reduction for Accurate Polarimetry).</footnote>
## 2. Observations

We observed Betelgeuse on January 17 and 18, 2019, with SPHERE (Beuzit et al. 2019) at the Very Large Telescope at visible wavelengths with ZIMPOL (Schmid et al. 2018) and at near-infrared wavelengths with IRDIS (Dohlen et al. 2008; de Boer et al. 2020; van Holstein et al. 2020). Data were taken in linear polarimetric modes in ten filters, including those of a point spread function (PSF) calibrator star, Phi² Orionis. The total intensity images are presented in Fig. 1, and Table 1 summarizes the observations. Atmospheric conditions were good and stable throughout the Betelgeuse and calibrator observations. However, we experienced a low-wind effect (Milli et al. 2018). This had the impact of decreasing the Strehl ratio (SR) as measured from the SPHERE images compared to the SPARTA SR measurements. The ZIMPOL and IRDIS data were reduced by the SPHERE Data Center (Delorme et al. 2017) and using the IRDAP<fnref n="1" /> pipeline (van Holstein et al. 2017, 2020), respectively. For ZIMPOL data, we estimated the instrumental polarization to be less than 0.1% for all filters by computing it over the halo of starlight in the images, a result consistent with the findings from Schmid et al. (2018). With IRDAP, we subtracted the instrumental polarization using the Mueller matrix model and then subtracted the stellar polarization from the halo of starlight. We followed the IRDAP convention in terms of orientation of polarization direction. From the *Q* and *U* images, we then computed images of the azimuthal Stokes parameters *Q*<sub>φ</sub> and *U*<sub>φ</sub> following the definitions of de Boer et al. (2020) as well as images of the linearly polarized intensity, *PI* = √(*Q*<sup>2</sup> + *U*<sup>2</sup>). The images of the Stokes parameters *Q* and *U* reveal “butterfly patterns” in all filters (see Fig. 2), as expected for centro-symmetric polarization created from scattering off a face-on structure such as a disk or a sphere.

In order to remove instrumental polarization (IP) effects, it is possible to use techniques such as minimizing the azimuthal Stokes parameter, *U*<sub>φ</sub>. (Avenhaus et al. 2014; de Boer et al. 2020). Because there should be no radial polarization, we applied a minimization technique to the negative *Q*<sub>φ</sub> signal in order to remove residual IP in our ZIMPOL data set. Due to the better spatial resolution delivered by ZIMPOL, there were enough pixels in the images to efficiently remove spurious polarization. This was not possible in the case of the IRDIS data set. We therefore used the *U*<sub>φ</sub> minimization technique in that case as the IP-generated butterfly pattern clearly dominated in these images.

In the V-band filter, we observed that the butterfly pattern had a different orientation from all the other filters. The reason for this could in principle be astrophysical as, for thick dusty stellar environments such as protoplanetary disks, multiple scattering can lead to non-azimuthal linear polarization effects (azimuthal Stokes *U*<sub>φ</sub> signal is non-null; see Canovas et al. 2015). However, in our case, the angle of linear polarization (AoLP) is azimuthal at all the other wavelengths from the visible to the near-infrared, so it is more likely that the cause for this V-band-specific polarization effect is of instrumental origin. Having corrected for instrumental polarization, we suggest that this issue originates from a so-far uncharacterized instrumental crosstalk.

The polarimetric design of ZIMPOL (Schmid et al. 2018) includes two half-wave plates (HWP2 and HWPZ) that together prevent the creation of instrumental crosstalk by SPHERE's

<figure><img src="image_2.png" /><figcaption>**Fig. 1.** Total intensity in the ten different filters. $\lambda$ is the central wavelength of each filter. North is up and east to the left. Units are W m$^{-2}$ sr$^{-1}$. The dashed orange and white circles mark angular sizes of 500 and 60 ms, respectively.</figcaption></figure>

<footnote n="1">¹ https://irdap.readthedocs.io (IRDIS Data reduction for Accurate Polarimetry).</footnote>
<table><caption>**Table 1.** Central wavelength, DIT × NDIT (number of DITs), seeing, and airmass for the Betelgeuse exposures as well as the full width at half maximum (FWHM) and SRs measured from the Phi² Ori images.</caption><thead><tr><th>Filter</th><th>Central<br>wavelength (nm)</th><th>DIT (s) × NDIT</th><th>Seeing (′′)</th><th>Airmass</th><th>FWHM<br>(mas)</th><th>Strehl<br>ratio (%)</th></tr></thead><tbody><tr><td colspan="7"><strong>ZIMPOL (January 18, 2019)</strong></td></tr><tr><td>V</td><td>554</td><td>1.1 × 20</td><td>0.42–0.48</td><td>1.21–1.23</td><td>25.4</td><td>28.5</td></tr><tr><td>Hα_CONT</td><td>644.9</td><td>1.1 × 20</td><td>0.44–0.59</td><td>1.27–1.32</td><td>21.1</td><td>44.3</td></tr><tr><td>Hα_NB</td><td>656.5</td><td>1.1 × 20</td><td>0.44–0.59</td><td>1.27–1.32</td><td>21.2</td><td>46.3</td></tr><tr><td>TiO_717</td><td>716.8</td><td>1.1 × 20</td><td>0.44–0.56</td><td>1.24–1.27</td><td>22.5</td><td>50.2</td></tr><tr><td>KI</td><td>770.2</td><td>1.1 × 20</td><td>0.45–0.66</td><td>1.32–1.37</td><td>23.7</td><td>48.4</td></tr><tr><td>748cont</td><td>747.4</td><td>1.1 × 20</td><td>0.45–0.66</td><td>1.32–1.37</td><td>23.5</td><td>51.3</td></tr><tr><td>820cont</td><td>817.3</td><td>1.1 × 20</td><td>0.44–0.56</td><td>1.24–1.27</td><td>24.2</td><td>58.6</td></tr><tr><td colspan="7"><strong>IRDIS (January 17, 2019)</strong></td></tr><tr><td>N_HeI</td><td>1085</td><td>0.84 × 30</td><td>0.53–0.78</td><td>1.19–1.10</td><td>41.7</td><td>73.8</td></tr><tr><td>N_CntJ</td><td>1213</td><td>0.84 × 30</td><td>0.55–0.70</td><td>1.18–1.19</td><td>42.1</td><td>78.4</td></tr><tr><td>N_PaB</td><td>1283</td><td>0.84 × 30</td><td>0.52–0.70</td><td>1.20–1.22</td><td>39.2</td><td>80.5</td></tr></tbody></table>

**Notes.** Filters parameters can be retrieved at www.eso.org/sci/facilities/paranal/instruments/sphere/inst/filters.html

image derotator (K-mirror) when the derotator moves in field-tracking (P2) mode. HWP2 rotates the polarization directions to be measured by ZIMPOL such that they are parallel and perpendicular to the plane of reflection of the derotator. As a result, the derotator does not create instrumental crosstalk. HWPZ, located within ZIMPOL, then rotates the to-be-measured polarization directions back to their original orientation. However, before the light reaches HWPZ, it is reflected off several silver-coated mirrors. Because the polarization directions are now not aligned parallel or perpendicular to the planes of reflection of these mirrors, these mirrors could create crosstalk, predominantly at short wavelengths. To validate this scenario, it is necessary to measure Mueller matrices in the V band and in other filters for comparison. This instrumental task is still pending.

To correct the orientation of the butterfly pattern of the *Q* and *U* images in the *V* band, we referred to procedures outlined in Appendix C of Avenhaus et al. (2014) and in de Boer et al. (2020) for the case of the correction of instrumental polarization. We first set the AoLP to be azimuthally oriented. We did this by finding the angular offset for the linear polarization that minimizes the sum of the squared signal in the image of the azimuthal Stokes parameter, *U*<sub>φ</sub>. After fitting the offset, we converted the images back to *Q* and *U*, and the butterfly patterns then presented a similar orientation as in the other filters. In practice, we performed this step in the same iteration loop as the *Q*<sub>φ</sub> minimization designed to remove residual instrumental polarization as described above. The comparison between the data before and after the *U*<sub>φ</sub> and *Q*<sub>φ</sub> minimizations is presented in Fig. 2. However, even though the orientation of the butterfly pattern was corrected by this technique, the polarized intensity cannot be fully retrieved as we do not know the amount of polarization signal that is lost due to the crosstalk. A modeling of this instrumental loss remains beyond the scope of the present work. We therefore disregard V-band filter data in the rest of the paper.

For the ZIMPOL observations, the PSF reference star Phi² Ori (stellar parameters from Bordé et al. 2002, and references therein) could not be used to calibrate the photometry of our images as atmospheric transparency conditions were too unstable (the so-called THICK category, transparency variations above 20%) between Betelgeuse and this calibrator (see

Fig. 3). However, the IRDIS data points (wavelengths above 1 µm) were taken the night before, when the atmosphere was more transparent (so-called CLEAR conditions). Those near-infrared Betelgeuse measurements were calibrated with Phi² Ori measurements and fitted by an ATLAS9 stellar atmosphere model (Castelli & Kurucz 2003) via the pysynphot package (STScI Development Team 2013). We used an effective temperature of 3600 K, a solar metallicity (M/H) and a log *g* of 0.0 (Levesque & Massey 2020; Taniguchi et al. 2022), and a Cardelli-Clayton-Mathis (CCM) extinction law (Cardelli et al. 1989) with *R*<sub>v</sub> = 3.1 and the extinction parameter, *A*<sub>v</sub>, as our only adjusted parameter. We find *A*<sub>v</sub> = 0.58 ± 0.10, which is in agreement with the range of values reported in Fig. 2 of Taniguchi et al. (2022) for the same period of time. We then used this curve to calibrate all our SPHERE photometric measurements in the ZIMPOL filters. Uncertainties were computed as the standard deviation of the total intensity from cubes of reduced frames, which corresponds to a duration of ~15 min for each filter.

## 3. Analysis

The final reduced images of Stokes *Q* and *U* and the linearly polarized intensities, PI, are shown in Fig. 2. The PI images exhibit an overall ring-like structure, compatible with a spherical dust shell seen in polarized light. It is noteworthy that PI drops in the center of images, which is expected for dusty environments observed with scattering angles close to 0°. This is also expected because the positive and negative *Q* and *U* signals cancel each other out in this part of the image due to the finite spatial resolution. The AoLP is oriented azimuthally in all filters. Most of the signal is included in a 60 mas radius and peaks at a distance of 30–40 mas from the center (Fig. 4). Along this ring, there is a high azimuthal variation that is consistent over the visible filters. In the near-infrared, the azimuthal profile is much smoother owing to the lower spatial resolution of the images. To illustrate this, we identified four azimuthal sectors of PI values included between a distance of 22 and 60 mas. We plot their variation as a function of azimuth and wavelength in Fig. 5. The highest relative azimuthal variation (with respect to the average value over the same filter) is about 80% and is observed in
<table><caption>**Table 1.** Central wavelength, DIT × NDIT (number of DITs), seeing, and airmass for the Betelgeuse exposures as well as the full width at half maximum (FWHM) and SRs measured from the Phi² Ori images.</caption><thead><tr><th>Filter</th><th>Central wavelength (nm)</th><th>DIT (s) × NDIT</th><th>Seeing (′′)</th><th>Airmass</th><th>FWHM (mas)</th><th>Strehl ratio (%)</th></tr></thead><tbody><tr><td colspan="7">ZIMPOL (January 18, 2019)</td></tr><tr><td>V</td><td>554</td><td>1.1 × 20</td><td>0.42–0.48</td><td>1.21–1.23</td><td>25.4</td><td>28.5</td></tr><tr><td>Hα_CONT</td><td>644.9</td><td>1.1 × 20</td><td>0.44–0.59</td><td>1.27–1.32</td><td>21.1</td><td>44.3</td></tr><tr><td>Hα_NB</td><td>656.5</td><td>1.1 × 20</td><td>0.44–0.59</td><td>1.27–1.32</td><td>21.2</td><td>46.3</td></tr><tr><td>TiO_717</td><td>716.8</td><td>1.1 × 20</td><td>0.44–0.56</td><td>1.24–1.27</td><td>22.5</td><td>50.2</td></tr><tr><td>KI</td><td>770.2</td><td>1.1 × 20</td><td>0.45–0.66</td><td>1.32–1.37</td><td>23.7</td><td>48.4</td></tr><tr><td>748cont</td><td>747.4</td><td>1.1 × 20</td><td>0.45–0.66</td><td>1.32–1.37</td><td>23.5</td><td>51.3</td></tr><tr><td>820cont</td><td>817.3</td><td>1.1 × 20</td><td>0.44–0.56</td><td>1.24–1.27</td><td>24.2</td><td>58.6</td></tr><tr><td colspan="7">IRDIS (January 17, 2019)</td></tr><tr><td>N_HeI</td><td>1085</td><td>0.84 × 30</td><td>0.53–0.78</td><td>1.19–1.10</td><td>41.7</td><td>73.8</td></tr><tr><td>N_CntJ</td><td>1213</td><td>0.84 × 30</td><td>0.55–0.70</td><td>1.18–1.19</td><td>42.1</td><td>78.4</td></tr><tr><td>N_PaB</td><td>1283</td><td>0.84 × 30</td><td>0.52–0.70</td><td>1.20–1.22</td><td>39.2</td><td>80.5</td></tr></tbody></table>

**Notes.** Filters parameters can be retrieved at www.eso.org/sci/facilities/paranal/instruments/sphere/inst/filters.html

image derotator (K-mirror) when the derotator moves in field-tracking (P2) mode. HWP2 rotates the polarization directions to be measured by ZIMPOL such that they are parallel and perpendicular to the plane of reflection of the derotator. As a result, the derotator does not create instrumental crosstalk. HWPZ, located within ZIMPOL, then rotates the to-be-measured polarization directions back to their original orientation. However, before the light reaches HWPZ, it is reflected off several silver-coated mirrors. Because the polarization directions are now not aligned parallel or perpendicular to the planes of reflection of these mirrors, these mirrors could create crosstalk, predominantly at short wavelengths. To validate this scenario, it is necessary to measure Mueller matrices in the V band and in other filters for comparison. This instrumental task is still pending.

To correct the orientation of the butterfly pattern of the *Q* and *U* images in the *V* band, we referred to procedures outlined in Appendix C of Avenhaus et al. (2014) and in de Boer et al. (2020) for the case of the correction of instrumental polarization. We first set the AoLP to be azimuthally oriented. We did this by finding the angular offset for the linear polarization that minimizes the sum of the squared signal in the image of the azimuthal Stokes parameter, *U*<sub>φ</sub>. After fitting the offset, we converted the images back to *Q* and *U*, and the butterfly patterns then presented a similar orientation as in the other filters. In practice, we performed this step in the same iteration loop as the *Q*<sub>φ</sub> minimization designed to remove residual instrumental polarization as described above. The comparison between the data before and after the *U*<sub>φ</sub> and *Q*<sub>φ</sub> minimizations is presented in Fig. 2. However, even though the orientation of the butterfly pattern was corrected by this technique, the polarized intensity cannot be fully retrieved as we do not know the amount of polarization signal that is lost due to the crosstalk. A modeling of this instrumental loss remains beyond the scope of the present work. We therefore disregard V-band filter data in the rest of the paper.

For the ZIMPOL observations, the PSF reference star Phi² Ori (stellar parameters from Bordé et al. 2002, and references therein) could not be used to calibrate the photometry of our images as atmospheric transparency conditions were too unstable (the so-called THICK category, transparency variations above 20%) between Betelgeuse and this calibrator (see

Fig. 3). However, the IRDIS data points (wavelengths above 1 µm) were taken the night before, when the atmosphere was more transparent (so-called CLEAR conditions). Those near-infrared Betelgeuse measurements were calibrated with Phi² Ori measurements and fitted by an ATLAS9 stellar atmosphere model (Castelli & Kurucz 2003) via the pysynphot package (STScI Development Team 2013). We used an effective temperature of 3600 K, a solar metallicity (M/H) and a log *g* of 0.0 (Levesque & Massey 2020; Taniguchi et al. 2022), and a Cardelli-Clayton-Mathis (CCM) extinction law (Cardelli et al. 1989) with *R*<sub>v</sub> = 3.1 and the extinction parameter, *A*<sub>v</sub>, as our only adjusted parameter. We find *A*<sub>v</sub> = 0.58 ± 0.10, which is in agreement with the range of values reported in Fig. 2 of Taniguchi et al. (2022) for the same period of time. We then used this curve to calibrate all our SPHERE photometric measurements in the ZIMPOL filters. Uncertainties were computed as the standard deviation of the total intensity from cubes of reduced frames, which corresponds to a duration of ~15 min for each filter.

## 3. Analysis

The final reduced images of Stokes *Q* and *U* and the linearly polarized intensities, PI, are shown in Fig. 2. The PI images exhibit an overall ring-like structure, compatible with a spherical dust shell seen in polarized light. It is noteworthy that PI drops in the center of images, which is expected for dusty environments observed with scattering angles close to 0°. This is also expected because the positive and negative *Q* and *U* signals cancel each other out in this part of the image due to the finite spatial resolution. The AoLP is oriented azimuthally in all filters. Most of the signal is included in a 60 mas radius and peaks at a distance of 30–40 mas from the center (Fig. 4). Along this ring, there is a high azimuthal variation that is consistent over the visible filters. In the near-infrared, the azimuthal profile is much smoother owing to the lower spatial resolution of the images. To illustrate this, we identified four azimuthal sectors of PI values included between a distance of 22 and 60 mas. We plot their variation as a function of azimuth and wavelength in Fig. 5. The highest relative azimuthal variation (with respect to the average value over the same filter) is about 80% and is observed in
<figure><img src="image_2.png" /><figcaption>**Fig. 2.** Stokes *Q* and *U* and polarized intensity images for the ten filters. The left and right parts correspond to the data before and after *Q*<sub>φ</sub> and *U*<sub>φ</sub> minimizations, respectively. Red and blue circles represent the photospheric radius (22 mas; Montargès et al. 2016) and the position of the IDS in H2019 (32 mas), respectively. On the polarized intensity map, angles of polarization are marked as orange segments. North is up, and east is to the left. Units are W m<sup>-2</sup> sr<sup>-1</sup>.</figcaption></figure>
<figure><img src="image_3.png" /><figcaption>**Fig. 2.** Stokes *Q* and *U* and polarized intensity images for the ten filters. The left and right parts correspond to the data before and after *Q*<sub>φ</sub> and *U*<sub>φ</sub> minimizations, respectively. Red and blue circles represent the photospheric radius (22 mas; Montargès et al. 2016) and the position of the IDS in H2019 (32 mas), respectively. On the polarized intensity map, angles of polarization are marked as orange segments. North is up, and east is to the left. Units are W m<sup>-2</sup> sr<sup>-1</sup>.</figcaption></figure>
<figure><img src="image_3.png" /><figcaption>**Fig. 3.** Synthetic flux density (ATLAS9 model with CCM reddening law) and derived flux density from SPHERE observations. The black dots represent the ATLAS9 average flux densities in SPHERE filters. As the conditions did not allow us to trust the ZIMPOL (visible) photometric calibrator measurements, we used the synthetic model fitted in the near-infrared to perform the photometric calibration.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Fig. 4.** Normalized vertical cut of the polarized intensities in each filter.</figcaption></figure>

the Ha-NB filter. The PI mean values vary from sector to sector, with a similar behavior at all wavelength, which indicates azimuthal density variations. This also suggests the components of the polarizing structure are similar, which in terms of dust population could point to a uniform composition and a similar grain size distribution for the four sectors.

## 4. Centro-symmetrical modeling

In order to derive physical parameters from our observations, we compared the SPHERE images to MCFOST (Pinte et al. 2006, 2009) radiative transfer simulations (similarly to H2019). We convolved radiative transfer models with a PSF that was extracted from SPHERE observations of Phi² Ori in the same instrumental and atmospheric stability conditions. We chose not to model the degree of linear polarization as the intensity at the IDS location is highly contaminated by the stellar flux. For the same reason, we did not model the total intensity images

<figure><img src="image_5.png" /></figure>

<figure><img src="image_6.png" /><figcaption>**Fig. 5.** Variation in the PI in the azimuthal direction. Top: azimuthal variation in the PI between 22 and 60 mas. Bottom: image at 717 nm, showing the four sectors.</figcaption></figure>
<figure><img src="image_4.png" /><figcaption>**Fig. 3.** Synthetic flux density (ATLAS9 model with CCM reddening law) and derived flux density from SPHERE observations. The black dots represent the ATLAS9 average flux densities in SPHERE filters. As the conditions did not allow us to trust the ZIMPOL (visible) photometric calibrator measurements, we used the synthetic model fitted in the near-infrared to perform the photometric calibration.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>**Fig. 4.** Normalized vertical cut of the polarized intensities in each filter.</figcaption></figure>

the Ha-NB filter. The PI mean values vary from sector to sector, with a similar behavior at all wavelength, which indicates azimuthal density variations. This also suggests the components of the polarizing structure are similar, which in terms of dust population could point to a uniform composition and a similar grain size distribution for the four sectors.

## 4. Centro-symmetrical modeling

In order to derive physical parameters from our observations, we compared the SPHERE images to MCFOST (Pinte et al. 2006, 2009) radiative transfer simulations (similarly to H2019). We convolved radiative transfer models with a PSF that was extracted from SPHERE observations of Phi² Ori in the same instrumental and atmospheric stability conditions. We chose not to model the degree of linear polarization as the intensity at the IDS location is highly contaminated by the stellar flux. For the same reason, we did not model the total intensity images

<figure><img src="image_6.png" /></figure>

<figure><img src="image_7.png" /><figcaption>**Fig. 5.** Variation in the PI in the azimuthal direction. Top: azimuthal variation in the PI between 22 and 60 mas. Bottom: image at 717 nm, showing the four sectors.</figcaption></figure>
either, as an accurate modeling of the PSF would be very challenging in our case. Our modeling approach is to reproduce the absolute values of the PI signal from the circular structure whose average radius is located at 30 mas. As a first approximation, we modeled the shape seen in the PI images as a spherical shell composed of spherical dust grains. The main parameters of the model are:

- a distance of 222 pc (Harper et al. 2017);

- an effective stellar temperature of 3600 K (Levesque & Massey 2020);

- a stellar radius of 1000 $R_\odot$ (roughly 22 mas photospheric diameter at 222 pc; Montargès et al. 2016);

- a stellar mass of 15 $M_\odot$, compatible with mass ranges published in the literature, for example Luo et al. (2022);

- an internal and external shell radii ($R_{in}$, $R_{out}$);

- a dust mass ($M_d$);

- a size distribution for spherical dust grains between a minimum and a maximum value ($a_{min}$ and $a_{max}$). We assumed a power law with two values for its index ($a_{exp} = 0$ and the commonly adopted value from interstellar medium studies: $a_{exp} = -3.5$; Mathis et al. 1977) and a growth scenario where small grains are always present in the distribution: we set $a_{min}$ to 0.05 µm;

- dust composition: we chose enstatite (MgSiO$_3$) as it is the dust species favored by modeling done in H2019. The optical constants of the enstatite are derived from Jäger et al. (2003) and taken from the Jena University Database<fnref n="2" />. We used the Mie theory to compute optical properties.

As the shell is optically thin, $M_d$ acts as a scaling factor to the absolute value of the polarized intensity. It does not influence the geometrical extent of the modeled dust shell. We fixed its value to $M_d = 10^{-12} M_\odot$. We ran a grid of MCFOST models to explore the remaining parameters, $R_{in}$, $R_{out}$, and $a_{max}$:

- $a_{max}$ ranged from 0.1 to 1.1 µm with steps of 0.1 µm;

- $R_{in}$ ranged from 5 AU (roughly the stellar radius) to 15 AU, with a step of 1 AU;

- $R_{out}$ was set to always be higher than $R_{in}$ by at least 1 AU and up to a maximum of 11 AU (2.25 stellar radii), with a step of 1 AU.

We then computed a $\chi^2$ value for each average value per sector compared to the ones extracted from the SPHERE images. We took the median of the $\chi^2$ values across the nine filters and plotted it as a function of the grid parameters. Finally, similarly to an analysis presented in Sect. 4.4 of Pinte et al. (2008), we computed a probability for each point of the parameter space (i.e., one model) as $P = \exp^{-\chi^2_{median}/2}$. All probabilities were normalised so that the sum of the probabilities of all models over the entire grid is equal to 1. These results are presented in Fig. 6 and were obtained by marginalising (i.e., summing) the probabilities of all models where one parameter is fixed successively to its different values, that is to say, the probabilities of our 3D parameter space are projected successively on each of the dimensions. Potential correlations between parameters are entirely accounted for with this approach, and the error bars extracted from the probability curves take the interplay between parameters into account.

The models that yield the highest probabilities are presented in Table 2. Whatever the sector, there is an overall decreasing probability for increasing $R_{in}$ values. This statement is valid to a lesser extent for the $R_{out}$ distribution of probability because

<table><caption>**Table 2.** Best model parameters by sectors for all models and models whose temperature is below 1500 K.</caption><thead><tr><th>Sector</th><th>R<sub>in</sub><br>(AU)</th><th>R<sub>out</sub><br>(AU)</th><th>a<sub>exp</sub></th><th>a<sub>max</sub><br>(µm)</th></tr></thead><tbody><tr><td>S1 – all</td><td>6</td><td>8</td><td>-3.5</td><td>0.5</td></tr><tr><td>S1 – T < 1500 K</td><td>8</td><td>12</td><td>-3.5</td><td>0.6</td></tr><tr><td>S2 – all</td><td>6</td><td>8</td><td>0</td><td>0.6</td></tr><tr><td>S2 – T < 1500 K</td><td>8</td><td>12</td><td>-3.5</td><td>0.6</td></tr><tr><td>S3 – all</td><td>6</td><td>8</td><td>0.</td><td>0.6</td></tr><tr><td>S3 – T < 1500 K</td><td>9</td><td>10</td><td>-3.5</td><td>0.6</td></tr><tr><td>S4 – all</td><td>6</td><td>8</td><td>-3.5</td><td>0.5</td></tr><tr><td>S4 – T < 1500 K</td><td>8</td><td>9</td><td>-3.5</td><td>0.6</td></tr></tbody></table>

the 1–6 AU range presents important variations, which points to some degree of indetermination for this parameter. The best models in terms of probability have the same extension for all sectors: between 6 and 8 AU (1.2 and 1.6 stellar radius, approximatively). The distribution of $a_{max}$ peaks between 0.4 and 0.6 µm, and the best models all have values of 0.5 or 0.6 µm. The probability distribution of $R_{in}$ indicates that the dust shell would start right above the photosphere, where the presence of dust is hard to justify given the high temperatures (see Sect. 5.1). A $R_{in}$ value of 5 AU corresponds approximately to the full width at half maximum of the PSF that we achieved in these observations (Table 1). Therefore, the estimation of $R_{in}$ is probably sensitive to convolution effects.

# 5. Discussion

## 5.1. Dust temperature

Despite the relative agreement with observations, the range of dust grain temperature is uncomfortably high for some MgSiO$_3$ models. Figure 7 shows the temperature as a function of distance to the star. Models with distances to the star center of less than 8 AU (~1.6 stellar radius) imply grain dust temperatures higher than 1500 K, which is approximatively the MgSiO$_3$ sublimation temperature. These models are consequently not entirely plausible even if they produce a reasonable fit to the SPHERE observations.

Discarding these, the remaining models lead to the distribution of probabilities plotted in blue in Fig. 6. The corresponding best models are summarized in Table 2. The inner radius, $R_{in}$, is set to start from 8 AU due to the temperature constraint, and this value corresponds to the $R_{in}$ for the best models except for sector 3 (9 AU). As noted above, the outer dust shell radius, $R_{out}$, is not clearly constrained. For models with a dust temperature lower than 1500 K, this is particularly true for sectors 1 and 2, where the best models have $R_{out} = 12$ AU. It is noteworthy that the second-best models for sectors 1 and 2 are the best models for 3 and 4, with $R_{out} = 10$ and 9 AU, respectively. All these models systematically have the same maximum size for the grain distribution ($a_{max}$): 0.6 µm. Absolute normalized values of residuals (or relative difference) for these models are shown as 1D profiles in Fig. 8. The highest absolute values of residuals are reached for the 0.820 µm filter and longer central wavelength filters. This suggests a more complex dust population than what was simulated in our grid of enstatite dust models. Below this wavelength, residuals are below 70% for all sectors. In order to illustrate the departure from spherical symmetry, we show

<footnote n="2">² http://www.astro.uni-jena.de/Laboratory/OCDB/amsilicates.html</footnote>
either, as an accurate modeling of the PSF would be very challenging in our case. Our modeling approach is to reproduce the absolute values of the PI signal from the circular structure whose average radius is located at 30 mas. As a first approximation, we modeled the shape seen in the PI images as a spherical shell composed of spherical dust grains. The main parameters of the model are:

- a distance of 222 pc (Harper et al. 2017);

- an effective stellar temperature of 3600 K (Levesque & Massey 2020);

- a stellar radius of 1000 $R_\odot$ (roughly 22 mas photospheric diameter at 222 pc; Montargès et al. 2016);

- a stellar mass of 15 $M_\odot$, compatible with mass ranges published in the literature, for example Luo et al. (2022);

- an internal and external shell radii ($R_{in}$, $R_{out}$);

- a dust mass ($M_d$);

- a size distribution for spherical dust grains between a minimum and a maximum value ($a_{min}$ and $a_{max}$). We assumed a power law with two values for its index ($a_{exp} = 0$ and the commonly adopted value from interstellar medium studies: $a_{exp} = -3.5$; Mathis et al. 1977) and a growth scenario where small grains are always present in the distribution: we set $a_{min}$ to 0.05 µm;

- dust composition: we chose enstatite (MgSiO$_3$) as it is the dust species favored by modeling done in H2019. The optical constants of the enstatite are derived from Jäger et al. (2003) and taken from the Jena University Database$^2$. We used the Mie theory to compute optical properties.

As the shell is optically thin, $M_d$ acts as a scaling factor to the absolute value of the polarized intensity. It does not influence the geometrical extent of the modeled dust shell. We fixed its value to $M_d = 10^{-12} M_\odot$. We ran a grid of MCFOST models to explore the remaining parameters, $R_{in}$, $R_{out}$, and $a_{max}$:

- $a_{max}$ ranged from 0.1 to 1.1 µm with steps of 0.1 µm;

- $R_{in}$ ranged from 5 AU (roughly the stellar radius) to 15 AU, with a step of 1 AU;

- $R_{out}$ was set to always be higher than $R_{in}$ by at least 1 AU and up to a maximum of 11 AU (2.25 stellar radii), with a step of 1 AU.

We then computed a $\chi^2$ value for each average value per sector compared to the ones extracted from the SPHERE images. We took the median of the $\chi^2$ values across the nine filters and plotted it as a function of the grid parameters. Finally, similarly to an analysis presented in Sect. 4.4 of Pinte et al. (2008), we computed a probability for each point of the parameter space (i.e., one model) as $P = \exp^{-\chi^2_{median}/2}$. All probabilities were normalised so that the sum of the probabilities of all models over the entire grid is equal to 1. These results are presented in Fig. 6 and were obtained by marginalising (i.e., summing) the probabilities of all models where one parameter is fixed successively to its different values, that is to say, the probabilities of our 3D parameter space are projected successively on each of the dimensions. Potential correlations between parameters are entirely accounted for with this approach, and the error bars extracted from the probability curves take the interplay between parameters into account.

The models that yield the highest probabilities are presented in Table 2. Whatever the sector, there is an overall decreasing probability for increasing $R_{in}$ values. This statement is valid to a lesser extent for the $R_{out}$ distribution of probability because

<table><caption>**Table 2.** Best model parameters by sectors for all models and models whose temperature is below 1500 K.</caption><thead><tr><th>Sector</th><th>R<sub>in</sub><br>(AU)</th><th>R<sub>out</sub><br>(AU)</th><th>a<sub>exp</sub></th><th>a<sub>max</sub><br>(µm)</th></tr></thead><tbody><tr><td>S1 – all</td><td>6</td><td>8</td><td>-3.5</td><td>0.5</td></tr><tr><td>S1 – T < 1500 K</td><td>8</td><td>12</td><td>-3.5</td><td>0.6</td></tr><tr><td>S2 – all</td><td>6</td><td>8</td><td>0</td><td>0.6</td></tr><tr><td>S2 – T < 1500 K</td><td>8</td><td>12</td><td>-3.5</td><td>0.6</td></tr><tr><td>S3 – all</td><td>6</td><td>8</td><td>0.</td><td>0.6</td></tr><tr><td>S3 – T < 1500 K</td><td>9</td><td>10</td><td>-3.5</td><td>0.6</td></tr><tr><td>S4 – all</td><td>6</td><td>8</td><td>-3.5</td><td>0.5</td></tr><tr><td>S4 – T < 1500 K</td><td>8</td><td>9</td><td>-3.5</td><td>0.6</td></tr></tbody></table>

the 1–6 AU range presents important variations, which points to some degree of indetermination for this parameter. The best models in terms of probability have the same extension for all sectors: between 6 and 8 AU (1.2 and 1.6 stellar radius, approximatively). The distribution of $a_{max}$ peaks between 0.4 and 0.6 µm, and the best models all have values of 0.5 or 0.6 µm. The probability distribution of $R_{in}$ indicates that the dust shell would start right above the photosphere, where the presence of dust is hard to justify given the high temperatures (see Sect. 5.1). A $R_{in}$ value of 5 AU corresponds approximately to the full width at half maximum of the PSF that we achieved in these observations (Table 1). Therefore, the estimation of $R_{in}$ is probably sensitive to convolution effects.

# 5. Discussion

## 5.1. Dust temperature

Despite the relative agreement with observations, the range of dust grain temperature is uncomfortably high for some MgSiO$_3$ models. Figure 7 shows the temperature as a function of distance to the star. Models with distances to the star center of less than 8 AU (~1.6 stellar radius) imply grain dust temperatures higher than 1500 K, which is approximatively the MgSiO$_3$ sublimation temperature. These models are consequently not entirely plausible even if they produce a reasonable fit to the SPHERE observations.

Discarding these, the remaining models lead to the distribution of probabilities plotted in blue in Fig. 6. The corresponding best models are summarized in Table 2. The inner radius, $R_{in}$, is set to start from 8 AU due to the temperature constraint, and this value corresponds to the $R_{in}$ for the best models except for sector 3 (9 AU). As noted above, the outer dust shell radius, $R_{out}$, is not clearly constrained. For models with a dust temperature lower than 1500 K, this is particularly true for sectors 1 and 2, where the best models have $R_{out} = 12$ AU. It is noteworthy that the second-best models for sectors 1 and 2 are the best models for 3 and 4, with $R_{out} = 10$ and 9 AU, respectively. All these models systematically have the same maximum size for the grain distribution ($a_{max}$): 0.6 µm. Absolute normalized values of residuals (or relative difference) for these models are shown as 1D profiles in Fig. 8. The highest absolute values of residuals are reached for the 0.820 µm filter and longer central wavelength filters. This suggests a more complex dust population than what was simulated in our grid of enstatite dust models. Below this wavelength, residuals are below 70% for all sectors. In order to illustrate the departure from spherical symmetry, we show

<footnote n="2">2 http://www.astro.uni-jena.de/Laboratory/OCDB/amsilicates.html</footnote>
<figure><img src="image_7.png" /><figcaption>**Fig. 6.** Probabilities of the various parameters for two values of *a*<sub>exp</sub> per sector for the whole range of models (black) and for models whose dust temperatures do not exceed 1500 K (blue). The triangles represent the parameters of the best model for each case.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**Fig. 7.** Dust grain temperature as a function of distance for the best models of the four sectors. The red segment of the curves marks the zone where the grain dust density is located.</figcaption></figure>

observation, model, and residual images for the best model with *T* < 1500K, corresponding to sector 2 in Fig. 9. We believe it is possible to reproduce the azimuthal variations with an appropriate 3D density structure. However, this task is beyond the scope of the present paper.

Different structures are present in the first stellar radii: a luke-
warm chromosphere (Lim et al. 1998; O’Gorman et al. 2017)
whose temperature at 1.7 stellar radius is about 3000 K (Fig. 4 in
O’Gorman et al. 2020), a 2000–2500 K MOL sphere at about
1.3 stellar radius (Perrin et al. 2004; Tsuji 2006; Ohnaka et al.
2011; Montargès et al. 2014), and the inner dust environment
that we focus on in the present paper. O’Gorman et al. (2020,
Sect. 4.1) propose that the coexistence of different temperatures
can be reconciled if the different components are distributed
in small inhomogeneities, smaller than the angular resolution.
For an IDS at 1.5 stellar radius, the local temperature needs to
allow dust survival; these inhomogeneous dusty zones would
be preferably found above the coolest areas of the convective
surface.

<figure><img src="image_9.png" /><figcaption>Fig. 8. Residuals per sector as defined in Fig. 5.</figcaption></figure>

The radial location and morphology of the azimuthal IDS
structures reported here are similar to the features of the ²⁸SiO
emission line intensity map presented in Fig. 5 of Kervella et al.
(2018) from 2016 ALMA observations. They are both included
<figure><img src="image_8.png" /><figcaption>**Fig. 6.** Probabilities of the various parameters for two values of $a_{exp}$ per sector for the whole range of models (black) and for models whose dust temperatures do not exceed 1500 K (blue). The triangles represent the parameters of the best model for each case.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>**Fig. 7.** Dust grain temperature as a function of distance for the best models of the four sectors. The red segment of the curves marks the zone where the grain dust density is located.</figcaption></figure>

observation, model, and residual images for the best model with $T < 1500$K, corresponding to sector 2 in Fig. 9. We believe it is possible to reproduce the azimuthal variations with an appropriate 3D density structure. However, this task is beyond the scope of the present paper.

Different structures are present in the first stellar radii: a luke-warm chromosphere (Lim et al. 1998; O’Gorman et al. 2017) whose temperature at 1.7 stellar radius is about 3000 K (Fig. 4 in O’Gorman et al. 2020), a 2000–2500 K MOL sphere at about 1.3 stellar radius (Perrin et al. 2004; Tsuji 2006; Ohnaka et al. 2011; Montargès et al. 2014), and the inner dust environment that we focus on in the present paper. O’Gorman et al. (2020, Sect. 4.1) propose that the coexistence of different temperatures can be reconciled if the different components are distributed in small inhomogeneities, smaller than the angular resolution. For an IDS at 1.5 stellar radius, the local temperature needs to allow dust survival; these inhomogeneous dusty zones would be preferably found above the coolest areas of the convective surface.

<figure><img src="image_10.png" /><figcaption>**Fig. 8.** Residuals per sector as defined in Fig. 5.</figcaption></figure>

The radial location and morphology of the azimuthal IDS structures reported here are similar to the features of the $^{28}$SiO emission line intensity map presented in Fig. 5 of Kervella et al. (2018) from 2016 ALMA observations. They are both included
<figure><img src="image_10.png" /><figcaption>**Fig. 9.** Polarized intensity images. Left: observed images in polarized intensity. Middle: best model with dust temperatures below 1500 K for sector 2. Right: residuals (observation minus model).</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>**Fig. 10.** Evolution of the IDS in polarized intensity. Images are not deconvolved, and the convention is the same as in Fig. 2.</figcaption></figure>

within the first 3 stellar radii (~15 AU, ~66 mas) and present arc-shaped structures. The exact connection between the structures on the two maps is not clear as condensation times and the dynamics at play are not known. Nevertheless, having such a distribution of gaseous SiO in the immediate circumstellar environment supports an at least partial silicate composition for the IDS, which is also compatible with predictions of the thermodynamical calculations of dust condensation in Betelgeuse (Gupta & Sahijpal 2020).

## 5.2. Mass loss

The amount of dust mass reported here, 1e-12 $M_\odot$, is significantly lower than the IDS mass of 1.5e-10 $M_\odot$ reported in H2019 and the dust mass-loss rate of 2e-9 $M_\odot$ yr$^{-1}$ reported in Josselin & Plez (2007). As noted above, it is possible that the actual IDS mass needed to reproduce the asymmetric features is much higher. It should also be noted that we have not probed models involving nonspherical dust particles and that results presented here could also vary with a larger exploration of dust properties.

It is also possible that the IDS exhibits variability with a minimum timescale of years. This is supported by the morphological variation in the IDS between observations reported in Kervella et al. (2016, however, the non-deconvolved PI images were not published) and the observations reported here, which are shown in Fig. 10. In both images, most of the polarized flux is enclosed within a ~70 mas radius and the maximum is found ~30 mas from the center. However, the spatial distributions of the intensity at the two epochs differ significantly in both radial and azimuthal directions. In the 2015 image, there is a significant polarized flux emission within the photospheric (red) radius, whereas in 2019 this area is characterized by a drop, creating a globally more ring-like structure than in 2015.

## 6. Conclusion

We have presented polarimetric high-spatial-resolution images of the Betelgeuse IDS obtained with SPHERE/ZIMPOL and SPHERE/IRDIS. Within our modeling hypotheses, the polarized signal arises from stellar light scattered by a centro-symmetric dust shell that extends up to several AU above the photosphere. Following H2019, we assumed dust grains are composed of $MgSiO_3$ and that the derived maximum dust grain size is in the 0.4–0.6 µm range, 0.6 µm for the best-fitting model. Though our modeling is limited, the description it gives of the polarized structure around Betelgeuse in terms of grain size distribution did not show any great variation between 2013 (reported in
<figure><img src="image_11.png" /><figcaption>**Fig. 9.** Polarized intensity images. Left: observed images in polarized intensity. Middle: best model with dust temperatures below 1500 K for sector 2. Right: residuals (observation minus model).</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>**Fig. 10.** Evolution of the IDS in polarized intensity. Images are not deconvolved, and the convention is the same as in Fig. 2.</figcaption></figure>

within the first 3 stellar radii (~15 AU, ~66 mas) and present arc-shaped structures. The exact connection between the structures on the two maps is not clear as condensation times and the dynamics at play are not known. Nevertheless, having such a distribution of gaseous SiO in the immediate circumstellar environment supports an at least partial silicate composition for the IDS, which is also compatible with predictions of the thermodynamical calculations of dust condensation in Betelgeuse (Gupta & Sahijpal 2020).

## 5.2. Mass loss

The amount of dust mass reported here, 1e-12 $M_\odot$, is significantly lower than the IDS mass of 1.5e-10 $M_\odot$ reported in H2019 and the dust mass-loss rate of 2e-9 $M_\odot$ yr$^{-1}$ reported in Josselin & Plez (2007). As noted above, it is possible that the actual IDS mass needed to reproduce the asymmetric features is much higher. It should also be noted that we have not probed models involving nonspherical dust particles and that results presented here could also vary with a larger exploration of dust properties.

It is also possible that the IDS exhibits variability with a minimum timescale of years. This is supported by the morphological variation in the IDS between observations reported in Kervella et al. (2016, however, the non-deconvolved PI images were not published) and the observations reported here, which are shown in Fig. 10. In both images, most of the polarized flux is enclosed within a ~70 mas radius and the maximum is found ~30 mas from the center. However, the spatial distributions of the intensity at the two epochs differ significantly in both radial and azimuthal directions. In the 2015 image, there is a significant polarized flux emission within the photospheric (red) radius, whereas in 2019 this area is characterized by a drop, creating a globally more ring-like structure than in 2015.

## 6. Conclusion

We have presented polarimetric high-spatial-resolution images of the Betelgeuse IDS obtained with SPHERE/ZIMPOL and SPHERE/IRDIS. Within our modeling hypotheses, the polarized signal arises from stellar light scattered by a centro-symmetric dust shell that extends up to several AUs above the photosphere. Following H2019, we assumed dust grains are composed of $MgSiO_3$ and that the derived maximum dust grain size is in the 0.4–0.6 µm range, 0.6 µm for the best-fitting model. Though our modeling is limited, the description it gives of the polarized structure around Betelgeuse in terms of grain size distribution did not show any great variation between 2013 (reported in
H2019) and 2019. It is possible that this points to a relatively sta-
ble dust condensation radius and dust characteristics over time.
Our centro-symmetric model leads to large azimuthal residuals
that may be the results of density variations, for instance due to
asymmetric outflows, projection effects, and/or illumination
asymmetries. In order to further constrain the mechanisms at
play in this IDS, further work should include a modeling of these
effects as well as a monitoring of the IDS, whose morphology is
believed to vary over several months to 1 yr. This should allow
the links between the IDS and surface phenomena, which vary
on similar timescales, to be revealed.

*Acknowledgements.* This work has made use of the SPHERE Data Centre, jointly operated by OSUG/IPAG (Grenoble), PYTHEAS/LAM/CESAM (Marseille), OCA/Lagrange (Nice), Observatoire de Paris/LESIA (Paris), and Observatoire de Lyon.

References

Avenhaus, H., Quanz, S. P., Schmid, H. M., et al. 2014, ApJ, 781, 87
Beuzit, J. L., Vigan, A., Mouillet, D., et al., 2019, A&A, 631, A155
Bordé, P., Coudé du Foresto, V., Chagnon, G., et al. 2002, A&A, 393, 183
Canovas, H., Ménard, F., de Boer, J., et al. 2015, A&A, 582, A7
Cardelli, J. A., Clayton, G. C., & Mathis, J. S. 1989, ApJ, 345, 245
Castelli, F., & Kurucz, R. L. 2003, Modell. Stellar Atmos., 210, A20
Cohen, M., Walker, R. G., Carter, B., et al. 1999, AJ, 117, 1864
de Boer, J., Langlois, M., van Holstein, R. G., et al. 2020, A&A, 633, A63
Delorme, P., Meunier, N., Albert, D., et al. 2017, SF2A-2017: Proceedings of the Annual meeting of the French Society of Astronomy and Astrophysics
Dohlen, K., Langlois, M., Saisse, M. et al. 2008, SPIE Conf. Ser., 7014, 1266
Dupree, A. K., Strassmeier, K. G., Matthews, L. D., et al. 2020, ApJ, 899, 68
Gail, H.-P., Tamanai, A., Pucci, A., et al. 2020, A&A, 644, A139
Gupta, A., & Sahijpal, S. 2020, MNRAS, 496, L122
Harper, G. M., Brown, A., Guinan, E. F., et al. 2017, AJ, 154, 11
Harper, G. M., Guinan, E. F., Wasatonic, R., et al. 2020, ApJ, 905, 34

Haubois, X., Norris, B., Tuthill, P. G., et al. 2019, A&A, 628, A101
Jäger, C., Dorschner, J., Mutschke, H., Posch, T., & Henning, T. 2003, A&A, 408, 193
Josselin, E., & Plez, B. 2007, A&A, 469, 671
Kervella, P., Perrin, G., Chiavassa, A., et al. 2011, A&A, 531, A17
Kervella, P., Lagadec, E., Montargès, M., et al. 2016, A&A, 585, A28
Kervella, P., Decin, L., Richards, A. M. S., et al. 2018, A&A, 609, A67
Kravchenko, K., Jorissen, A., Van Eck, S., et al. 2021, A&A, 650, A17
Levesque, E. M., & Massey, P. 2020, ApJ, 891, L37
Lim, J., Carilli, C. L., White, S. M., et al. 1998, Nature, 392, 575
Luo, T., Umeda, H., Yoshida, T., et al. 2022, ApJ, 927, 115
Mathis, J. S., Rummel, W., & Nordsieck, K. H. 1977, ApJ, 217, 425
Milli, J., Kasper, M., Bourget, P., et al. 2018, Proc. SPIE, 10703, 107032A
Montargès, M., Kervella, P., Perrin, G., et al. 2014, A&A, 572, A17
Montargès, M., Kervella, P., Perrin, G., et al. 2016, A&A, 588, A130
Montargès, M., Cannon, E., Lagadec, E., et al. 2021, Nature, 594, 365
Nozawa, T., Maeda, K., Kozasa, T., et al. 2011, ApJ, 736, 45
O’Gorman, E., Kervella, P., Harper, G. M., et al. 2017, A&A, 602, A10
O’Gorman, E., Harper, G. M., Ohnaka K., et al. 2020, A&A, 638, A65
Ohnaka, K., Weigelt, G., Millour, F., et al. 2011, A&A, 529, A163
Perrin, G., Ridgway, S. T., Coudé du Foresto, V., et al. 2004, A&A, 418, 675
Perrin, G., Verhoelst, T., Ridgway, S. T., et al. 2007, A&A, 474, 599
Pinte, C., Ménard, F., Duchêne, G., & Bastien, P. 2006, A&A, 459, 797
Pinte, C., Padgett, D. L., Ménard, F., et al. 2008, A&A, 489, 633
Pinte, C., Harries, T. J., Min, M., et al. 2009, A&A, 498, 967
Roddier, F., & Roddier, C. 1985, ApJ, 295, L21
Safonov, B., Dodin, A., Burlak, M., et al. 2020, ArXiv e-prints [arXiv:2005.05215]
Schmid, H. M., Bazzon, A., Roelfsema, R., et al. 2018, A&A, 619, A9
STScI Development Team 2013, Astrophysics Source Code Library, [record ascl:1303.023]
Taniguchi, D., Yamazaki, K., & Uno, S. 2022, Nat. Astron., 6, 930
Tsuji, T. 2006, ApJ, 645, 1448
van Holstein, R. G., Snik, F., Girard, J. H., et al. 2017, SPIE Conf. Ser., 10400, 1040015
van Holstein, R. G., Girard, J. H., de Boer, J., et al. 2020, A&A, 633, A64
Verhoelst, T., Decin, L., van Malderen, R., et al. 2006, A&A, 447, 311
H2019) and 2019. It is possible that this points to a relatively sta-
ble dust condensation radius and dust characteristics over time.
Our centro-symmetric model leads to large azimuthal residuals
that may be the results of density variations, for instance due
to asymmetric outflows, projection effects, and/or illumination
asymmetries. In order to further constrain the mechanisms at
play in this IDS, further work should include a modeling of these
effects as well as a monitoring of the IDS, whose morphology is
believed to vary over several months to 1 yr. This should allow
the links between the IDS and surface phenomena, which vary
on similar timescales, to be revealed.

*Acknowledgements.* This work has made use of the SPHERE Data Centre, jointly operated by OSUG/IPAG (Grenoble), PYTHEAS/LAM/CESAM (Marseille), OCA/Lagrange (Nice), Observatoire de Paris/LESIA (Paris), and Observatoire de Lyon.

References

Avenhaus, H., Quanz, S. P., Schmid, H. M., et al. 2014, ApJ, 781, 87

Beuzit, J. L., Vigan, A., Mouillet, D., et al., 2019, A&A, 631, A155

Bordé, P., Coudé du Foresto, V., Chagnon, G., et al. 2002, A&A, 393, 183

Canovas, H., Ménard, F., de Boer, J., et al. 2015, A&A, 582, A7

Cardelli, J. A., Clayton, G. C., & Mathis, J. S. 1989, ApJ, 345, 245

Castelli, F., & Kurucz, R. L. 2003, Modell. Stellar Atmos., 210, A20

Cohen, M., Walker, R. G., Carter, B., et al. 1999, AJ, 117, 1864

de Boer, J., Langlois, M., van Holstein, R. G., et al. 2020, A&A, 633, A63

Delorme, P., Meunier, N., Albert, D., et al. 2017, SF2A-2017: Proceedings of the Annual meeting of the French Society of Astronomy and Astrophysics

Dohlen, K., Langlois, M., Saisse, M. et al. 2008, SPIE Conf. Ser., 7014, 1266

Dupree, A. K., Strassmeier, K. G., Matthews, L. D., et al. 2020, ApJ, 899, 68

Gail, H.-P., Tamanai, A., Pucci, A., et al. 2020, A&A, 644, A139

Gupta, A., & Sahijpal, S. 2020, MNRAS, 496, L122

Harper, G. M., Brown, A., Guinan, E. F., et al. 2017, AJ, 154, 11

Harper, G. M., Guinan, E. F., Wasatonic, R., et al. 2020, ApJ, 905, 34

Haubois, X., Norris, B., Tuthill, P. G., et al. 2019, A&A, 628, A101

Jäger, C., Dorschner, J., Mutschke, H., Posch, T., & Henning, T. 2003, A&A, 408, 193

Josselin, E., & Plez, B. 2007, A&A, 469, 671

Kervella, P., Perrin, G., Chiavassa, A., et al. 2011, A&A, 531, A17

Kervella, P., Lagadec, E., Montargès, M., et al. 2016, A&A, 585, A28

Kervella, P., Decin, L., Richards, A. M. S., et al. 2018, A&A, 609, A67

Kravchenko, K., Jorissen, A., Van Eck, S., et al. 2021, A&A, 650, A17

Levesque, E. M., & Massey, P. 2020, ApJ, 891, L37

Lim, J., Carilli, C. L., White, S. M., et al. 1998, Nature, 392, 575

Luo, T., Umeda, H., Yoshida, T., et al. 2022, ApJ, 927, 115

Mathis, J. S., Rummel, W., & Nordsieck, K. H. 1977, ApJ, 217, 425

Milli, J., Kasper, M., Bourget, P., et al. 2018, Proc. SPIE, 10703, 107032A

Montargès, M., Kervella, P., Perrin, G., et al. 2014, A&A, 572, A17

Montargès, M., Kervella, P., Perrin, G., et al. 2016, A&A, 588, A130

Montargès, M., Cannon, E., Lagadec, E., et al. 2021, Nature, 594, 365

Nozawa, T., Maeda, K., Kozasa, T., et al. 2011, ApJ, 736, 45

O'Gorman, E., Kervella, P., Harper, G. M., et al. 2017, A&A, 602, A10

O'Gorman, E., Harper, G. M., Ohnaka K., et al. 2020, A&A, 638, A65

Ohnaka, K., Weigelt, G., Millour, F., et al. 2011, A&A, 529, A163

Perrin, G., Ridgway, S. T., Coudé du Foresto, V., et al. 2004, A&A, 418, 675

Perrin, G., Verhoelst, T., Ridgway, S. T., et al. 2007, A&A, 474, 599

Pinte, C., Ménard, F., Duchêne, G., & Bastien, P. 2006, A&A, 459, 797

Pinte, C., Padgett, D. L., Ménard, F., et al. 2008, A&A, 489, 633

Pinte, C., Harries, T. J., Min, M., et al. 2009, A&A, 498, 967

Roddier, F., & Roddier, C. 1985, ApJ, 295, L21

Safonov, B., Dodin, A., Burlak, M., et al. 2020, ArXiv e-prints
[arXiv:2005.05215]

Schmid, H. M., Bazzon, A., Roelfsema, R., et al. 2018, A&A, 619, A9

STScI Development Team 2013, Astrophysics Source Code Library, [record ascl:1303.023]

Taniguchi, D., Yamazaki, K., & Uno, S. 2022, Nat. Astron., 6, 930

Tsuji, T. 2006, ApJ, 645, 1448

van Holstein, R. G., Snik, F., Girard, J. H., et al. 2017, SPIE Conf. Ser., 10400, 1040015

van Holstein, R. G., Girard, J. H., de Boer, J., et al. 2020, A&A, 633, A64

Verhoelst, T., Decin, L., van Malderen, R., et al. 2006, A&A, 447, 311