Sudden extreme obscuration of a Sun-like main-sequence star: evolution of the circumstellar dust around ASASSN-21qj

In order to derive the fundamental stellar parameters of ASASSN-21qj, we applied the Bayesian inference code of del Burgo & Allende Prieto (2016, 2018) on the PARSEC v1.2S library of stellar evolution models (Bressan et al., 2012). This combination is an informed choice given the good statistical agreement of the dynamical masses of detached eclipsing binaries with their respective inferences, especially for main-sequence stars, where predicted masses are, on average, 4 % accurate (del Burgo & Allende Prieto, 2018).

We follow the method of del Burgo & Allende Prieto (2016), taking the absolute $G$ magnitude, $M_{G}$, the colour $G_{BP}-G_{RP}$ from Gaia DR3 (Gaia Collaboration et al., 2021), together with an abundance ratio of iron to hydrogen [Fe/H]= $-0.20~{}\pm~{}0.20$ dex. We obtained $M_{G}$ using the geometric distance for ASASSN-21qj of 556${}^{+3}_{-4}$ pc (Bailer-Jones et al., 2021).

The given metallicity is favoured by our overall prediction and matches the average value of Gaia DR3, although adopting a conservative figure for its uncertainty. Given the large distance to the star, we corrected the photometry from interstellar extinction by means of the Gaia DR3 estimates for the line-of-sight extinction in the $G$ band, $A_{G}$, and the reddening $E$($BP$-$RP$). Note this implementation is supported by the fact that our inferred value for the effective temperature is higher but comparable (within 200  K) to the Gaia DR3 estimate. The inputs for making our inference and the resulting fundamental parameters are presented in Table 1. They convey that ASASSN-21qj is a main sequence star whose fundamental parameters resembles those of the Sun.

The lightcurve at optical wavelengths including ASAS-SN $V$ and $g^{\prime}$ measurements (Shappee et al., 2014), LCOGT $g^{\prime}$ (this work), and ATLAS ‘cyan’ (653 nm) filter (Tonry et al., 2018; Heinze et al., 2018), and near-infrared wavelengths (NEOWISE 3.4 and 4.6 $\mu$m; NEOWISE Team, 2020) is presented in Figure 1, with a closeup of the multi-wavelength LCOGT observations given in Figure 2. The overall shape and duration of the dimming observed through the ASAS-SN/LCOGT and ATLAS observations are qualitatively similar, so we describe the event here with reference to the $g^{\prime}$ measurements from the former combined dataset. Around MJD 58200 a brief rise in the $g^{\prime}$ brightness is seen in ASAS-SN photometry, contemporaneously with the brightening at near-infrared wavelengths - this is marginal, but could be indicative of scattered light from the event that produced the dust responsible for the near-infrared excess and subsequent optical dimming. Prior to MJD 59000 no significant dimming of the star is observed in ASAS-SN monitoring, and the near-infrared excess remains constant. After that, a shallow decline in optical brightness is observed, as described in Rizzo Smith et al. (2021), and a near-simultaneous decrease in the near-infrared excess is observed. The optical brightness then rapidly declines around 59400 shortly before the onset of LCOGT monitoring. From MJD 59450 we record a variety of structure in the occultation event at timescales of hours to days, broadly consisting of a series of deep, quasiperiodic dimming events superimposed on a monotonically declining occultation.

The higher spatial resolution and sensitivity of LCOGT observations trace the dimming of below the 16 mag limit of ASAS-SN. At the onset of LCOGT observations we record a stellar magnitude of $g^{\prime}\simeq~{}$17 mag. Early in the evolution of the occultations two very deep dimming events were observed on MJDs 59476 and 59507 reducing the star to $g^{\prime}\simeq 20.3~{}$mag. The likelihood of two such deep events occurring within a single occultation being unconnected seems unlikely given the ongoing, decreasing overall dimming of the star in the intervening period. If we therefore assume that this is a recurrence of the same dense clump transiting the star, we infer an orbital period of 30.9 $\pm$ 0.5 days for the material responsible, equivalent to a semi-major axis of 0.19 $\pm$ 0.04 au.

Based on this potential period, additional observations were scheduled around the next couple of times this clump was expected to occult the star at MJDs 59537 and 59569. For the first of these, the expected occultation was a day late, demonstrating the clump was not periodic, and at most quasiperiodic. Whilst later than predicted, a relatively deep occultation did occur reducing the star to $g^{\prime}>16~{}$mag, from a baseline of 14 mag either side of this dip. From this we concluded that the structure of the clump was evolving, becoming larger and more diffuse over time. At the second predicted epoch we found a delay of 4 days from the predicted time. The clump had significantly broadened only dimming the star to $g^{\prime}>15.5~{}$mag with two distinct dips within the overall event. The overall structure had changed to span several days with an extended tail characteristic of a comet-like transit before the star resumed its pre-dip brightness. From this we concluded that the original spatially compact clump could be fragmenting into several distinct sub-clumps, or involve multiple bodies that could now be discerned due to their increasing separation. Following this, another two deep occultations are recorded in the sparser $g^{\prime}$ monitoring following MJD 59600; these follow the established pattern of being increasingly extended in time and shallower in dimming, although the measurements are too sparse (a cadence of every other day) to identify any substructure.

The occultation of ASASSN-21qj can be described as the conflation of two effects, namely a broad, low level extinction covering the star between MJDs 59000 and 60000 that is interspersed with multiple short duration events lasting a few days to a few tens of days. The light curve is shown in Figure 3 with the broad and narrow dimming events denoted by triangles. We seek to disentangle these two effects on the stellar brightness in order to determine the distribution of depths for narrow events that are superimposed on the broader overall dimming.

We model the broad extinction as an asymmetric Gaussian scaled to the shallow ingress and egress at the wings of the extinction event. The shape of this broad profile is described by a date of peak dimming $t_{\rm peak}$, peak depth $A_{\rm peak}$, and standard deviations before and after the peak $\sigma_{\rm prior}$ and $\sigma_{\rm post}$. The broad, shallow component of the obscuration has a maximum depth of 1 magnitude at MJD 59450 and standard deviations of 50 days before the peak, and 100 days after the peak. This asymmetric shape is characteristic of transiting dust clouds associated with exocomets.

Narrow, deep events in the time series are identified by searching for data points for which the preceding and subsequent data points are brighter by 3-$\sigma$, these events are modelled as symmetric Gaussians in similar fashion to the broad component. The peak depths of the narrow events are determined with reference to the level of broader extinction occurring at the same time. We identify 35 narrow dimming events in the $g^{\prime}$ light curve. During high cadence observations we see strong variability in the extinction within a single night’s photometry. We are thus insensitive to rapid changes in the structure of the cloud mainly comprising smaller dimming events. Likewise, the adopted method to identify narrow dimming events does not track broader, shallow events as can be seen by the lack of dimming events identified in the lightcurve during egress from overall extinction beyond MJD 59800. The estimated number of occultation events is therefore a lower limit to the total number.

We also examine the complete light curve using the Lomb Scargle periodogram as implemented in astropy to search for significant periods. Combining the $g^{\prime}$ measurements from LCOGT and ASAS-SN we obtain a peak period at 19.6 days, with a false alarm probability $<1\%$. Examining the data sets separately, we find peak periods of 29.9 days for the ASAS-SN observations, and 41.6 days for the LCOGT observations, both with false alarm probabilities $<1\%$. The LCOGT periodogram is noisy with several substantial aliases due to the irregular sampling and shorter time span of the observations compared to the ASAS-SN periodogram. Taking the peak period from the combined periodogram is therefore suspect. At the very least we can state there is good evidence for periodic behaviour between 19 and 42 days present in the light curve of ASASSN-21qj, consistent with the detected infrared excess.

The presence of a 31 day quasi-periodicity was previously inferred based on two deep dimming events early in the occultation (see Figures 2 and 3). Combined with the speculative periodic behaviour seen here and the presence of near-infrared excess indicative of close-in dust, we can estimate the size and structure of dust clumps responsible for the narrow, deep events within the occultation. Assuming the dust lies on circular orbits with periods between 19 and 42 days (equivalent to a semi-major axes from 0.14 to 0.23 au), the transiting velocity of the dust would be 74 to 60 km/s. The two deepest, narrowest features have total durations of 5.4 days, equivalent to a cloud extent along the orbit of 0.23 to 0.19 au. Given the near complete extinction of the stellar flux during these events, the cloud must also have a substantial vertical extent, with a scale height around 0.03 to completely cover the stellar disc.

Colour-colour plots ($g^{\prime}-r^{\prime}$ vs $r^{\prime}-i^{\prime}$ and $g^{\prime}-r^{\prime}$ vs $r^{\prime}-z^{\prime}$) are used to determine the reddening vectors associated with the source extinction from the LCOGT time series photometry. We use the computer code optool (Dominik et al., 2021) to calculate the $Q_{\rm ext}(=Q_{\rm abs}+Q_{\rm sca})$ values for a range of dust grain compositions and two size distributions. We used the in-built optical constants from optool to generate the $Q_{\rm ext}$ values for materials including astronomical silicate (Draine, 2003), various crystalline and amorphous iron-magnesium silicates (Dorschner et al., 1995; Jaeger et al., 1998; Fabian et al., 2001; Suto et al., 2006), carbon-bearing species (Draine & Malhotra, 1993; Henning & Stognienko, 1996), quartz (Kitamura et al., 2007), and iron particles (Henning & Stognienko, 1996).

We consider two grain size distributions: 1) a power law ($dN=a^{-q}da$) between a minimum (free) and maximum grain size (fixed as 1 mm, typical for debris discs) with a fixed size distribution exponent of 3.5 appropriate for a collisional cascade (Dohnanyi, 1969) and 2) a log-normal with a width of 0.5 and a free peak size, considering contributions from sizes between 0.01 and 10 $\mu$m. In each case, the free parameter was allowed to vary between 0.1 and 1.1 $\mu$m in steps of 0.01 $\mu$m when calculating $Q_{\rm ext}$, determined from the complex optical constants assuming Mie theory. For a Sun-like star the typical blowout size for a dust grain is $\simeq$0.5 $\mu$m (Krivov, 2010), so we include transient dust with this range of grain sizes. The optical depth $\tau$ and extinction of a dust cloud which fully covers the stellar disc are then calculated from the $Q_{\rm ext}$ values for each combination of composition and size distribution, which are then compared to the measured reddening vectors by least-squares fit to identify those which best match the observations.

The colour-colour plots for some of the materials tested in this analysis are presented in Figures 4 and 5, along with the best-fit material and size from the modelling in Figure 6. We find adequate fits to the reddening curves for several materials with both power law and log-normal size distributions. For both size distributions, a wide variety of grain sizes and materials can be matched to $g^{\prime}-r^{\prime}$ vs $r^{\prime}-i^{\prime}$ colour, but only a handful provide a good match to both $g^{\prime}-r^{\prime}$ vs $r^{\prime}-i^{\prime}$ and $g^{\prime}-r^{\prime}$ vs $r^{\prime}-z^{\prime}$ simultaneously. For both size distributions the best fit composition was amorphous pyroxene (Mg_1.6Fe_0.4Si₂O₆). For the power law size distribution the best fit minimum grain size was 0.29${}^{+0.01}_{-0.18}~{}\mu$m, whereas for the log-normal size distribution the best fit peak size was 0.88${}^{+0.05}_{-0.03}~{}\mu$m. If we favour the power law distribution, the presence of small, sub-blowout dust grains indicates ongoing collisional grinding of larger bodies. The dissipation of the occulting screen of material during the event between dips would also be consistent with the rapid removal of these grains.

The data points associated with the two deepest occultation measurements $g^{\prime}\simeq 20~{}$mag, on MJDs 59476 and 59507, do not follow the same reddening vector as the majority of the data points. If we assume the grains in those two deep occultations have the same size distribution and composition as the remainder of the cloud, we can infer a minimum grain size of $\geq 1.2~{}\mu$m for the power law distribution, limited by the wavelength range from which the reddening vector is derived. This suggests some degree of grain size segregation within the occulting cloud. As the overall level of dimming decreases with time, this is consistent with the idea that smaller grains will spread more quickly from an impact or breakup, either through kinetic displacement or subsequent spreading by radiation forces, whereas larger grains may remain more tightly associated with their point of formation.

Based on the inferred power law size distribution, the dimming associated with the deepest occultations requires an intervening dust mass of $3.5\times 10^{-12}~{}M_{\oplus}$. Assuming that the dust cloud lies on a circular orbit with a period of 30.9 days (semi-major axis of 0.19 au), based on the very deep dips early in the evolution of the event, the integrated occluding dust mass was at least 1.50 $~{}\pm~{}0.04\times 10^{-9}~{}M_{\oplus}$ (comparable to the total mass of Hale-Bopp) during the initial deep occultations. From the degree of dimming and the shape of the occultation, the dust cloud responsible must extend several stellar radii both vertically and horizontally, making this a lower limit to the total mass in the system.

Since the dust cloud is predicted to lie close to the star, we searched archival infrared measurements for evidence of excess emission from the system, indicative of a substantial circumstellar disc. We find no evidence of excess at near- or mid-infrared wavelengths in archival observations spanning $>$10 years, including Spitzer/IRAC measurements from GLIMPSE360 (Whitney et al., 2008) and the AllWISE catalogue (Wright et al., 2010). The absence of significant infrared excess, along with its inferred large age, rules against ASASSN-21qj being surrounded by a misaligned dense, optically thick circumstellar disc (e.g. Bredall et al., 2020; Ansdell et al., 2020), strengthening the case for an exocometary interpretation for the origin of its variability.

However, the NEOWISE observations reveal an intriguing brightening of the source several years prior to the onset of the occultation events we set out to monitor. This near-infrared brightening is first seen in measurements around MJD 58200, shortly after the optical brightening seen in the ASAS-SN photometry around MJD 58100. As of the most recent NEOWISE epoch around MJD 59900, the observed infrared excess has dissipated. The same material could therefore be responsible for both the reflected and emitted light seen in the time series observations such that whatever event caused the occultations we now observe was captured serendipitously. This near-infrared brightening is not seen in the NEOWISE measurements of any other exocomet host star, making the presence of a transient significant infrared excess unique to ASASSN-21qj. None of the other stars in close proximity to ASASSN-21qj (within 1$\arcmin$) exhibit a similar pattern in NEOWISE photometry, from which we rule out an instrumental origin. Assuming the excess emission to be intrinsic to ASASSN-21qj, we derive dust temperatures by fitting a blackbody curve to the averaged NEOWISE photometry at each epoch after subtraction of a scaled stellar photosphere model. The excess is consistent with temperatures in the range 1300 to 700 K, equivalent to blackbody distances of $\simeq$ 0.05 to 0.15 au from the star. Accounting for the absorption and emission properties of dust with the minimum size and composition derived previously (i.e. 0.29 $\mu$m amorphous pyroxene), we would expect distances of 0.12 to 0.24 au at the temperatures inferred from the infrared excess (1300${}^{+2700}_{-1000}$ to 700${}^{+100}_{-85}$ K).

An upper limit to the dust mass responsible for the unusual events reported here can then be inferred from the near-infrared excess. Radiative transfer modelling using Hyperion (Robitaille, 2011), adopting the previously calculated power law size distribution, dust composition, and spatial distribution, requires $\simeq 2\times 10^{-6}~{}M_{\oplus}$ of dust to replicate the excess.

Hot dust exozodiacal dust around main sequence stars has been identified in around 10 to 20% of systems, correlated with the presence of an outer cool debris disc (Absil et al., 2006, 2013; Ertel et al., 2014, 2016; Absil et al., 2021). The dust masses associated with these hot exozodis amounts to a few $10^{-9}~{}M_{\oplus}$ (Kirchschlager et al., 2017). Around very young stars, near-infrared excess and variability thought to be associated with terrestrial planet formation has been observed, with inferred dust masses around $10^{-4}~{}M_{\oplus}$ (Meng et al., 2015; Su et al., 2019, 2020, 2022). Although ASASSN-21qj is an old star, the mass of dust seen here is more comparable to the young planet-forming systems.