11institutetext: Institute for Solar Physics, Dept. of Astronomy, Stockholm
University, AlbaNova University Centre, SE-10691 Stockholm, Sweden
22institutetext: Laboratory for Atmospheric and Space Physics, University of
Colorado, 3665 Discovery Drive, Boulder, CO 80303, USA
22email<EMAIL_ADDRESS>33institutetext: Department of Physics,
Catholic University of America, 620 Michigan Avenue, Northeast, Washington, DC
20064, USA 44institutetext: NASA Goddard Space Flight Center, Heliophysics
Sciences Division, 8800 Greenbelt Road, Greenbelt, MD 20771, USA

# On the radiative losses in the chromosphere during a C-class flare

Rahul Yadav 1122 J. de la Cruz Rodríguez 11 Graham S. Kerr 3344 C. J. Díaz
Baso 11 Jorrit Leenaarts 11

(Draft: compiled on at UT)

###### Abstract

Context. Solar flares release an enormous amount of energy ($\sim$1032 erg) in
the corona. A substantial fraction of this energy is transported to the lower
atmosphere, which results in chromospheric heating. The mechanisms that
transport energy to the lower solar atmosphere during a flare are still not
fully understood.

Aims. We aim to estimate the temporal evolution of the radiative losses in the
chromosphere at the footpoints of a C-class flare, in order to set
observational constraints on the electron beam parameters of a RADYN flare
simulation.

Methods. We estimated the radiative losses from hydrogen, and singly ionized
Ca and Mg using semi-empirical model atmospheres, which were inferred from
multi-line inversion of observed Stokes profiles obtained with the CRISP and
CHROMIS instruments at the Swedish 1-m Solar Telescope. The radiative losses
were computed taking into account the effect of partial redistribution and
non-local thermodynamic equilibrium. To estimate the integrated radiative
losses in the chromosphere the net cooling rates were integrated between the
temperature minimum and the height where the temperature reaches 10 kK. We
have also compared our time series of radiative losses with those from RADYN
flare simulations.

Results. We have obtained a high spatial-resolution map of integrated
radiative losses around the flare peak time. The stratification of the net
cooling rate suggests that the Ca IR triplet lines are responsible for most of
the radiative losses in the flaring atmosphere. During the flare peak time,
the contribution from Ca ii H & K and Mg ii h & k lines are strong and
comparable to the Ca IR triplet ($\sim$32 kW m-2). Since our flare is a
relatively weak event the chromosphere is not heated above 11 kK, which in
turn yields a subdued Ly$\alpha$ contribution ($\sim$7 kW m-2) in the selected
limits of the chromosphere. The temporal evolution of total integrated
radiative losses exhibits sharply-rising losses (0.4 kW m-2 s-1) and a
relatively slow decay (0.23 kW m-2 s-1). The maximum value of total radiative
losses is reached around the flare peak time, and can go up to 175 kW m-2 for
a single pixel located at footpoint. After a small parameter study, we find
the best model-data consistency in terms of the amplitude of radiative losses
and the overall atmospheric structure with a RADYN flare simulation in the
injected energy flux of $5\times 10^{10}$ erg s-1 cm-2.

###### Key Words.:

Sun: chromosphere – Sun: heating – Sun: flares

##  1 Introduction

Solar flares are the most energetic events on the Sun. A significant amount of
energy, up to 1032 erg, is released during an intense solar flare. A typical
flare is characterized by a rapid increase in emission over a wide range of
the electromagnetic spectrum, which affects the entire solar atmosphere, i.e.,
the photosphere, the chromosphere and, the corona (Fletcher et al., 2011). In
the standard flare model (Carmichael, 1964; Sturrock, 1966; Hirayama, 1974;
Kopp & Pneuman, 1976), the fundamental mechanism of solar flares is believed
to be magnetic reconnection in the corona, converting the magnetic energy into
kinetic and thermal energy. Following energy release in the corona, a
substantial amount of energy is also transported to the lower solar
atmosphere, resulting in intense plasma heating and ionisation, ultimately
leading to a substantial increase in the local radiative output in flare
ribbons/footpoints. Observations have demonstrated that during flares a major
amount of the energy budget is radiated in the chromosphere (Benz, 2017),
making understanding the flaring chromosphere key to understanding flare
energetics.

There is unambiguous evidence of the presence of energetic particles in solar
flare footpoints, and a large body of work points to the prominent role of
directed beams of non-thermal electrons as a significant mechanism by which
energy is transported from corona to chromosphere (Brown, 1971; Hudson, 1972;
Kontar et al., 2011; Holman et al., 2011). However, it is not known to what
extent other energy transport mechanisms play a role. Examples include direct
in-situ heating at loop tops, which is conducted through to the lower
atmosphere, the acceleration of protons and heavier ions alongside electrons,
or the presence of flare-induced magnetohydrodynamic waves (e.g. Emslie &
Sturrock, 1982; Fletcher & Hudson, 2008; Kerr et al., 2016; Reep et al.,
2018). Determining the relative contributions of each of these mechanisms is
an important goal in flare physics, which requires in part a solid
understanding of the observed radiative losses in the lower atmosphere.

Semi-empirical models have been used to estimate the radiative losses in the
lower solar atmosphere by various authors (Machado & Linsky, 1975; Machado et
al., 1980; Vernazza et al., 1981; Avrett et al., 1986; Gan et al., 1993;
Sobotka et al., 2016; Abbasvand et al., 2020). These models have been
performed using some form of non-local thermodynamical equilibrium data
inversion. Inversion codes adjust the stratification of physical parameters in
a model atmosphere in order to fit the observed spectra. The resulting model
can be assumed to be an approximation of the real thermodynamical state of the
solar atmosphere when the code includes all the physical ingredients that are
required to model the observations (del Toro Iniesta & Ruiz Cobo, 2016).

The radiative losses estimated from commonly used semi-empirical models of the
quiet-Sun chromosphere are 4.3 kW m-2, whereas in active regions the value is
20 kW m-2 (Withbroe & Noyes, 1977; Vernazza et al., 1981). In order to
maintain the energy balance in the chromosphere, some form of heating
mechanism is required, which is still under debate. The comparison of
observations with the semi-empirical models has revealed that the
chromospheric radiative losses are mainly dominated by strong lines of Ca ii,
Mg ii, Ly$\alpha$, and the H- continuum (Vernazza et al., 1981). Anderson &
Athay (1989) studied the role of Fe ii lines in the chromosphere and concluded
that Fe ii can be responsible for up to one half of the radiative losses if
they are included in the calculations from quiet Sun models. Despite this
potentially important contribution, Fe ii lines have not been traditionally
included in the calculation of the radiative losses because of the vast number
of transitions that would need to be calculated in NLTE, and the lack of
accurate atomic models. Semi-empirical models have also been used to
investigate various possible energy-transfer processes occurring during flares
(Machado et al., 1980).

Leenaarts et al. (2018) used the Ca ii K brightness as a proxy of the
radiative losses in the chromosphere in order to investigate the atmospheric
structure and heating of the solar chromosphere in an emerging flux region.
They reported that the frequency-integrated Ca ii K brightness correlates
strongly with the total linear polarization in the Ca ii 8542 Å, while the Ca
II K profile shapes indicate that the bulk of the radiative losses occurs in
the lower chromosphere. In a recent study, Díaz Baso et al. (2021) estimated
radiative losses from the inferred atmospheric model in a strong magnetic
reconnection event by integrating the net cooling rates in the chromosphere
from temperature minimum to the height where the temperature reaches 10 kK.
They reported integrated radiative losses up to 160 kW m-2, which is roughly
five times more than the typical value for the chromosphere in active regions.

Due to recent improvements in solar observations and in inversion techniques,
it is now possible to infer the stratification of physical parameters such as
temperature, line-of-sight (LOS) velocity, micro-turbulence velocity, magnetic
field, etc. in the solar atmosphere using multi spectral line observations (de
la Cruz Rodríguez et al., 2019; Vissers et al., 2019). Recently, Yadav et al.
(2021, hereafter Paper I) constructed the stratification of physical
parameters in a flaring atmosphere using multiple spectral lines forming in
different parts of the solar atmosphere. Such multi-line observations,
spanning different layers of the chromosphere and photosphere, allow us to
reconstruct the atmospheric parameters more accurately compared to a single
layer or spectral line.

The temporal analysis of the inferred parameters in Paper I revealed that the
mean temperature increased up to $\sim$11 kK in the chromosphere (between
$\log\tau_{500}\sim-3.5$ and $-5$), whereas no significant change was noticed
in the photosphere. As a consequence, intense brightening, which is a
signature of radiative losses, is observed in the chromospheric lines (Ca ii K
and Ca ii 8542 Å) at flare ribbons.

In this paper, we estimated the radiative losses from hydrogen, singly ionized
Ca and singly-ionised Mg using the semi-empirical flare model inferred from
multi-line observations. We report the temporal evolution of the integrated
radiative losses in the chromosphere during the flare. We have also compared
our results with radiative losses calculated from a flare simulation performed
with the RADYN code (Carlsson & Stein, 1997; Allred et al., 2015) in order to
demonstrate how our observationally derived results can attempt to set
constraints on the energy injection.

A brief overview of observations and inversions is given in Section 2. The
approach for estimating the cooling rates is described in Section 3. The
results are discussed in Section 4. Finally, the paper is summarized in
Section 5.

Figure 1: Overview of the flare observed at 08:46 UT on May 6, 2019. The
continuum intensity map (left panel) and the Chromospheric intensity map
(middle panel) near line core of Ca ii K observed at the SST. Right panel: the
integrated radiative losses estimated from the inferred semi-empirical model
atmosphere. Solar north, solar east, and the direction of the disk center are
indicated by “N,” “E,” and “DC,” respectively. The red asterisk symbol in the
middle panel shows the location of pixels analyzed to investigate the temporal
evolution of the integrated radiative losses. The mean values of the total
integrated losses at the A, B, C, and D locations (5$\times$5 pixels) are
listed in Table 1.

##  2 Observation and inversion

A C2-class solar flare (SOL2019-05-06T08:47) was observed between 08:34 and
9:33 UT on May 6, 2019, with the CRisp Imaging SpectroPolarimeter (CRISP;
Scharmer et al., 2008) and the CHROMospheric Imaging Spectrometer (CHROMIS;
Scharmer, 2017) instruments at the Swedish 1-m Solar Telescope (SST; Scharmer
et al., 2003). This flare occurred in the active region (AR) NOAA 12740, which
was located at N08E48 ($\mu$=0.62). The CRISP simultaneously recorded full
spectropolarimetric data in the Ca ii 8542 Å and Fe i 6173 Å spectral lines.
The Ca ii 8542 Å line scans consisted of 17 wavelength positions spanning a
range of 1.4 Å around line center, with steps of 75 mÅ in the inner wings and
125 mÅ at four outer wing positions, whereas the Fe i 6173 Å spectral line was
scanned at 15 wavelength positions spanning a range of 0.5 Å around line
center, with steps of 35 mÅ in the inner wings and 50 mÅ in the outer wings at
four different positions. The CHROMIS recorded Ca ii K intensity profiles at
28 wavelength positions spanning a range of 3 Å around line center, with steps
of 65 mÅ in the inner wings and an additional sampling at wavelength positions
$\pm{1.5}$ Å and $\pm{1.18}$ Å relative to the line center. In addition to
this, one point in the continuum at 4000 Å was also observed with the CHROMIS
instrument. The CRISP and CHROMIS data were obtained quasi-simultaneously with
a cadence of 21 sec and 15 sec, respectively.

The data were reduced and post-processed using the CRISPRED/SSTRED pipelines
(de la Cruz Rodríguez et al., 2015; Löfdahl et al., 2021). Further
observational details, data reduction, and post-processing procedures are
described in Paper I. See Figure 2 in Paper I for an overview of the observed
flare at the SST along with aligned ultraviolet (UV) and extreme ultraviolet
(EUV) images observed by the Atmosphere Imaging Assembly (AIA; Lemen et al.
2012) aboard the Solar Dynamic Observatory (SDO; Scherrer et al. 2012).

We employed the STiC inversion code (de la Cruz Rodríguez et al., 2016, 2019)
to infer the atmospheric stratification of temperature, velocity, and magnetic
field from the Fe i 6173 Å, Ca ii K, and Ca ii 8542 Å spectral lines. The STiC
code is built around a modified version of the RH code (Uitenbroek, 2001) in
order to derive the atomic populations by assuming statistical equilibrium and
a plane-parallel geometry. The equation of state is borrowed from the
Spectroscopy Made Easy (SME) computer code described in Piskunov & Valenti
(2017). The radiative transport equation is solved using cubic Bezier solvers
(de la Cruz Rodríguez & Piskunov, 2013).

For the inversion of the Stokes profiles, we considered the Ca II 8542 Å line
in non-LTE conditions, under the assumption of complete frequency
redistribution, while the Ca ii K line was synthesized in non-LTE conditions
with partial redistribution effects of scattered photons. We used the fast
approximation proposed by Leenaarts et al. (2012). The Fe I 6173 Å line was
treated under the assumption of LTE conditions. We inverted all the Stokes
parameters in the Fe I 6173 Å and Ca II 8542 Å lines, but only Stokes $I$ in
the Ca ii K line. The inversion setup and retrieved atmospheric parameters are
discussed in detail in Paper I.

##  3 Radiative cooling rate

The energy released by radiation is characterized by net radiative cooling
rates (radiative losses), which is evaluated as the frequency-integrated
radiative flux divergence (Rutten, 2003)

$\nabla\cdot
F=\int_{0}^{\infty}\int_{\Omega}\alpha_{\nu}(S_{\nu}-I_{\nu})\,\mathrm{d}\Omega\,\mathrm{d}\nu,$
(1)

where $\alpha$, $S$, and $I$ are the opacity, source function, and intensity
at the frequency $\nu$ and direction $\Omega$.

Equation 1 describes the total radiative loss rate. In this study we focus on
losses in individual spectral lines. The contribution of a single bound-bound
transition to Eq. 1 can be expressed as

$\displaystyle Q_{\mathrm{bb}}$ $\displaystyle=$
$\displaystyle\oint\int\frac{h\nu}{4\pi}\Bigg{[}n_{j}A_{ji}\psi_{ij}(\nu,\@vec{n})$
(2) $\displaystyle-
B_{ij}\bigg{(}n_{i}-\frac{g_{i}}{g_{j}}\rho_{ij}(\nu,\@vec{n})n_{j}\bigg{)}\phi_{ij}(\nu,\@vec{n})I(\nu,\@vec{n})\Bigg{]}\,\mathrm{d}\nu\,\mathrm{d}\Omega,$

where $h$ is the Planck constant, $n_{j}$ and $n_{i}$ are the population
densities of the upper levels $j$ and lower level $i$, $A_{ji}$ and $B_{ij}$
are Einstein coefficients, $g_{i}$ and $g_{j}$ are the statistical weights of
level $i$ and $j$, and $I(\nu,\boldsymbol{n})$ is the intensity for a given
direction $\@vec{n}$ at frequency $\nu$. This formulation takes PRD effects
into account through the frequency and angle-dependent ratio of the emissivity
to the absorption profile $\rho_{ij}=\psi_{ij}/\phi_{ij}$ (Uitenbroek, 2002).

We note that Eq. 1 includes many continuum opacity sources that are being
neglected when we make a selection of transitions and apply Eq. 2. In flares,
many of those continua can cool or heat even in the chromosphere. A limitation
is that most of them are computed under the assumption of LTE in the RH module
of STiC, and this assumption cannot be used to properly estimate atom
population densities in the chromosphere. Therefore, we have opted for a
conservative approach, and we are only accounting for the contributions that
we can model reasonably well. Therefore our results pose a lower-limit
estimate of the radiative losses.

The net cooling rates are evaluated from the inferred flare model by summing
the losses in the Ly$\alpha$, Ca II H & K, Ca IR triplet (Ca II 8542, 8498,
8662 Å), and Mg II h & k spectral lines as computed using Eq. 2. The radiative
rates and populations were computed using the STiC code after taking into
account the PRD effect and non-LTE conditions (see Section 2 for the STiC code
and inversion).

We decided to ignore the contribution from bound-free transitions in the
integrated losses. Because the transition region is not well constrained in
our inferred models as the observed spectral lines form in the lower/middle
chromosphere. Additionally, the node description of the atmosphere sets limits
on the steepness of the transition region by imposing a relatively smooth
stratification of the physical parameters with height (de la Cruz Rodríguez et
al., 2016, 2019). Because the code did not need to place a steep transition
region in order to fit the Ca II lines, the model predicts, close to the
transition region, contributions from the Lyman and Balmer continua that are
unrealistically high. Therefore our calculated radiative losses are only a
lower-limit to the real ones, which in turn are a lower-limit approximation to
the heating terms in the chromosphere. This is not a problem for our analysis
as we compute the radiative losses in the same way with the inversion models
and the simulations.

Although our flare was not observed in the Mg II h & k and Ca II H lines, we
assumed that the inferred model atmosphere retrieved from Ca ii 8542 Å, Ca II
K and Fe I 6173 Å observations represent the upper chromosphere correctly near
the formation height of Mg II h & k (Bjørgen et al., 2019). Our observed flare
is a relatively weak event and does not influence the photosphere like white-
light flares. Therefore, we neglected the contribution from weak spectral
lines in the photosphere and the H- continuum as their contributions in the
lower and middle chromosphere is small compared to the Ca II and Mg II lines
(Vernazza et al., 1981). We note that in our inferred model the mean
temperature in the lower/middle chromosphere is normally below 11 kK and the
uncertainty in the parameters is larger in the higher layers (¡
log${\tau_{500}}=-5$) where the observed lines are less sensitive.

Figure 2: Stratification of the net radiative cooling rates per unit volume
for a quiet-Sun pixel (left panel) and a pixel located at flare footpoint
(right panel) indicated by asterisk symbol in Figure 1. The colored solid
lines refer to different spectral lines. The black solid line represents the
total cooling rate, whereas the dashed curve refers to the temperature
stratification. The gray shaded area indicates the integrated heights.

##  4 Results

###  4.1 Stratification and distribution of radiative losses

We evaluated the bound-bound radiative losses per unit volume using the
approach outlined in Section 3. In order to estimate the total height-
integrated radiative losses in the chromosphere, we followed the approach
adopted by Díaz Baso et al. (2021). We integrated the contribution of selected
lines over geometrical heights, derived by assuming hydrostatic equilibrium
using the STiC code, between temperature minimum and the height where the
temperature reaches 10 kK. As an example, the total integrated radiative
losses for the field-of-view (FOV) comprising flare-ribbons (near flare peak
time) is shown in Figure 1. The figure illustrates that, during the flare peak
time, the bright flare ribbon regions display strong integrated radiative
losses.

The distribution of radiative losses across the FOV is greatly correlated with
the integrated intensity of the Ca ii H&K lines, similarly reported by
Leenaarts et al. (2018) and Díaz Baso et al. (2021). We note that the
integrated radiative losses are strongly dependent on the heights used to
integrate the net cooling rates. The contribution changes significantly in the
deeper layers compared to higher layers as density drops sharply outwards in
the solar atmosphere. The stratification of the net radiative cooling rate per
unit volume for a pixel located at the quiet-Sun and the flare footpoint is
shown in Fig. 2. The analysis of the quiet-Sun pixel exhibits that the
contribution of Ca ii IR lines is dominating compared to the Ly$\alpha$, Ca II
H & K and Mg ii h & k lines. On the other hand, during the flare peak time,
the net cooling rate at flare footpoint shows that the contribution from Ca ii
H & K and Mg ii h & k lines is also strong and comparable to the Ca IR triplet
($\sim$32 kW m-2). The stratification of net cooling rate of the Ly$\alpha$
shows that it increases around 1 Mm but remains below from that of the Ca ii
IR lines.

To compare the integrated radiative losses at different parts of the FOV (near
flare peak time), we selected small patches of pixels (5$\times$5) indicated
by A (Quiet-Sun region), B (bright fibrils), C (near polarity inversion line),
and D (flare footpoint) in the right panel of Fig. 1. The obtained integrated
radiative losses and contribution of all selected lines at different locations
are listed in Table 1. The total integrated value (4.9 kW m-2) in the Quiet-
Sun region is compatible with previous findings (Withbroe & Noyes, 1977;
Vernazza et al., 1981). At the bright fibrils, the total integrated value is
23.9 kW m-2. Near the polarity inversion line, the radiative losses are
$\sim$37 kW m-2, which is roughly twice the value reported for active regions
(Vernazza et al., 1981). At the flare footpoint, the mean integrated radiative
loss value is $\sim$133 kW m-2. As listed in the table, in all cases, the
majority of contribution comes from the Ca II IR triplet lines, with the
exception of the flare footpoint location where the contribution from the Ca
ii H & K and Mg ii h & k lines is similar to that from the Ca ii IR lines.

Figure 3 depicts the integrated radiative losses, at the flare peak time,
across a line passing through two footpoints. It shows that the maximum value
can go up to $\sim$175 kW m-2 for pixels located at the flare footpoint.
Similar high values are also reported by Díaz Baso et al. (2021) using a
similar approach and spectral lines (except Ly$\alpha$) in a strong magnetic
reconnection event. They reported that the chromospheric radiative losses at
the reconnection site are as high as 160 kW m-2.

Figure 3: The total integrated radiative losses evaluated along a blue solid line passing through two footpoints shown in Fig. 1. Figure 4: Temporal evolution of the total integrated radiative losses in the chromosphere at asterisk symbol highlighted in Fig. 1. The vertical solid, dashed, and dotted lines refer to the start, peak, and end times of the flare, respectively. Table 1: Mean contribution of Ly$\alpha$, Ca ii and Mg ii lines to the total integrated radiative losses (kW m-2) for the regions highlighted in Fig. 1. The Ca IR refers to the triplet: Ca II 8542, 8498, 8662 Å. A, B, C and D regions correspond to quiet-Sun, bright fibrils, pixels near PIL, and footpoint pixels at flare peak time, respectively. Spectral lines | Regions
---|---
| A | B | C | D
Ca II H | 0.75 | 3.9 | 6.1 | 17.6
Ca II K | 0.72 | 5.4 | 8.3 | 22.8
Ca IR | 3.2 | 9.1 | 14.9 | 32.5
Mg II h | 0.27 | 2.4 | 3.4 | 25.2
Mg II k | 0.21 | 2.8 | 3.7 | 28.2
Ly$\alpha$ | -0.17 | 0.3 | 0.76 | 7.05
Total | 4.93 | 23.9 | 37.16 | 133.35

Figure 5: Left and middle panels: The temporal evolution of temperature and
velocity in the RADYN simulation. The simulation were performed with
$\delta=7$, E${}_{c}=20$ keV and energy flux of $5\times 10^{10}$ erg s-1
cm-2. Right panel: Temporal evolution of the integrated radiative losses in
the chromosphere estimated using the simulated atmosphere. Solid black lines
refer to the total integrated radiative losses.

###  4.2 Evolution of the total integrated radiative losses

In order to investigate the temporal evolution of the integrated radiative
losses at flare footpoints, we selected a small patch of pixels (5$\times$5)
indicated by the asterisk symbol in the chromospheric intensity map of Fig. 1.
Figure 4 displays the obtained temporal evolution of the total integrated
radiative losses at the flare footpoint. The vertical lines in Fig. 4 indicate
the flare start, peak, and end times, which are determined using the X-ray
flux measured by the Geostationary Operational Environmental Satellite (GOES).
The curve showing the total contribution sharply increases after the flare
start time until the flare peak time with an average increase rate of
$\sim$0.4 kW m-2 s-1. We also note that the total integrated losses rise is
delayed by a few minutes from the flare start time. The flare start time
refers to the global flare, but the small patch we study evidently becomes
activated sometime after the initial onset of reconnection, as flare
reconnection propagates. The total integrated losses from our localized patch
then decay at a slower rate following the impulsive peak, ($\sim$0.23 kW m-2
s-1) in comparison to the sharp rise rate. During the flaring period, the main
cooling contribution comes from Ca ii and Mg ii lines. We also noticed that at
the flare peak time Ca ii H & K and Mg ii h & k lines also contribute
significantly, although their contribution before the start and after the end
time of the flare is less dominating compared to the Ca ii IR lines.

Moreover, during the flaring time, the radiative loss due to the Ly$\alpha$
line increases until the flare peak time, but its value stays below the Ca ii
and Mg ii lines. At the flare footpoint, during flare peak time, the maximum
value of integrated radiative losses from Ly$\alpha$ is around $\sim$7 kW m-2.
We note, however, that there might be a bias in this result. Our chromospheric
inversion is only based on spectra from Ca ii lines. For $T\gtrsim 10$ kK,
these lines rapidly lose sensitivity because of Ca ii ionization into Ca iii
(Pietarila et al., 2007). Therefore, our models typically predict a very cold
transition region in comparison with an inversion that includes lines that
sample the transition region. While the effect of a colder transition region
is not necessarily visible in the Ca ii lines, it can have a profound effect
in the Ly$\alpha$ line, which has a strong contribution from regions with
temperatures ranging from 20 to 40 kK (Machado et al., 1980; Vernazza et al.,
1981). This could explain why the response in this line is so subdued compared
to the expected stronger response (Milligan et al., 2014; Brown et al., 2018).
So it is not only a problem of restricting the integration limits for the
radiative losses to $T\leq 10$ kK, but of the intrinsic structure of the
reconstructed transition regions.

###  4.3 Comparison of radiative losses retrieved from observations and
simulations

In order to compare the values retrieved from observations with the
simulations, we also estimated the integrated radiative losses from the
atmospheric model synthesized using the radiative hydrodynamic (RADYN) code
(Carlsson & Stein, 1997; Allred et al., 2015). The RADYN simulation represents
a 1D strand or field-line in the flare that attempts to catch most of the
microphysics that are at work in a flare. RADYN was recently merged with the
code FP (Allred et al., 2020), offering an improved treatment of the non-
thermal particle propagation and dissipation. The RADYN code was employed to
simulate various flaring atmospheres, which are characterized by the electron
beam parameters: low-energy cutoff (Ec), spectral index ($\delta$), and non-
thermal energy flux. The $\delta$ and Ec values determine the spectral
‘hardness’ of the non-thermal electron distribution. A higher $\delta$ and
lower $E_{c}$ represent a softer distribution, that is a relatively larger
proportion of lower energy electrons. Lower energy particles thermalize more
easily in the upper chromosphere/lower transition region compared to higher
energy electrons that can penetrate deeper. We performed a small parameter
study, with a fixed $\delta=7$, but variable E${}_{c}=[10,15,20]$ keV. For
each value of Ec and $\delta$ we injected a beam with energy flux values F =
[0.75$\times 10^{10}$, 1.0$\times 10^{10}$, 2.5$\times 10^{10}$, and
5.0$\times 10^{10}$] erg s-1 cm-2, for a period of $t=30$ s. Following
cessation of the beam, each simulation was allowed to cool.

While generating STiC input model atmosphere we set a threshold temperature of
$40$ kK to the RADYN synthetic model atmospheres, which removes the
contribution of the corona. This is required as the STiC code is not well
constrained to model the coronal dynamics. We also investigated the influence
of temperature threshold on the radiative losses. We find that there is no
significant change in the radiative losses values if the temperature threshold
changed from 40 kK up to 60 kK. This could be related to our integrated
heights which are located in the deeper chromosphere.

As an example, the temporal evolution of the temperature and velocity
generated from one of the RADYN simulations is shown in Fig. 5. The obtained
stratification of temperature and velocity are qualitatively similar to the
inverted model parameters shown in Paper I, but there is a steeper rise in
temperature, and the velocity field has larger amplitudes in the chromosphere
compared to the inverted model retrieved from observations. Although this may
be a consequence of the starting atmosphere as well as the choice of electron
beam parameters, we note that node-based inversion tends to predict smooth
model atmospheres as sharp gradients and very complex stratifications require
a large number of nodes.

We have analyzed the radiative losses in the mid-chromosphere for the
different electron beam parameters and we find that the total radiative losses
in that part of the atmosphere are higher for larger values of Ec and energy
flux. In our case, the total radiative losses obtained with E${}_{c}=20$ keV
and energy flux of $5\times 10^{10}$ erg s-1 cm-2 are the most comparable to
the value obtained from the inverted model atmosphere, though differences do
exist. We stress that we only used a small sample of RADYN models as an
illustration of how such a model-data comparison may be performed. A more
comprehensive model-data comparison can follow in a future study.

The right panel of Fig. 5 displays the total integrated radiative losses
obtained from different lines. The maximum value reaches around 125 kW m-2,
which is comparable to the value ($\sim$135 kW m-2) obtained from the semi-
empirical model constructed from observations. Compared to the inverted
models, the RADYN curves have a much sharper rise and decay phase. While the
RADYN model is representative of one field line, the behavior of the time
evolution in Fig. 4 can be understood as a superposition of events like the
one seen in the rightmost panel of Fig. 5 in space and time. It is also known
that the modelled flare gradual phase is too short compared to observations
(see Kerr et al., 2020; Allred et al., 2022, and references therein).

Moreover, during the flare peak time in the RADYN simulations, in addition to
the Ca ii and Mg ii lines, we also noticed a sharp rise in the net cooling
rate of the L$\alpha$ line for a few seconds. In contrast to the simulation,
we did not observe the sharp rise in Ly$\alpha$ losses estimated from the
inferred model atmosphere, probably due to the colder and relatively poorly
constrained transition region in the inversions that we mentioned earlier.

As illustrated in Fig. 5, at the beginning of the simulation ($<$ 30 s), there
is a sharp rise in temperature which falls close to the temperature minimum
location. The right panel in Fig. 5 shows that the integrated radiative losses
show no response before $t<$ 30 s, even though the atmosphere is heated
strongly. This behavior is due to our integrated height selection criteria
(described in Section 3), which results in a very small range to integrate the
net radiative losses. Consequently, the integrated radiative losses do not
show high values before 30 s in the lower-limited chromosphere. In this model,
the radiative losses within the height (temperature) range selected actually
peak after the cessation of flare energy injection. The rapidly cooling upper
chromosphere/transition region results in a heat flux that propagates to a
greater depth. A small bubble forms, raising the temperature by up to a few
thousand degrees (but below $T=10$ kK), and elevating the electron density.
This produces a strong peak in radiative losses from the lower chromosphere.
Comparing the radiative losses derived from the inversion to those from the
forward model, the observed losses we report may reflect heating of the lower
chromosphere following cessation of the electron beam. A similar analysis with
lines that are sensitive to a larger temperature range, and that thus
constrain the mid-upper flaring chromosphere, are essential and should be a
priority in the present solar cycle.

##  5 Summary and Discussion

A C2-class solar flare (SOL2019- 05-06T08:47) that occurred on May 6, 2019,
was observed simultaneously in the Ca II K, Ca ii 8542 Å and Fe i 6173 Å lines
with the CRISP and CHROMIS instruments at the SST. We have inferred an
empirical flare model from the observational data by inverting all these lines
simultaneously with the STiC code (see Paper I for details). In this paper we
have utilized this model to estimate the net radiative losses, which
ultimately are a lower-limit approximation of the chromospheric heating terms.
Our analysis includes the main bound-bound chromospheric coolers: Ly$\alpha$,
Ca II H & K, Ca IR triplet and the Mg II h & k spectral lines. Due to the lack
of an accurate and realistically large atomic model, our calculations do not
include the contribution from Fe ii lines, which could be responsible for up
to one half of the radiative cooling in the chromosphere (Anderson & Athay,
1989). In order to estimate the radiative losses in the chromosphere, we
integrated the contribution of all selected lines over the geometrical heights
between the temperature minimum and the height below 10 kK temperature (e.g.,
indicated by gray area in Figure 2). Normally above 10 kK the observed
chromospheric spectral lines loose sensitivity as at such high temperatures
the Ca ii atom is ionized into Ca iii.

A novel aspect of our study is that we have calculated a very high spatial-
resolution map of net integrated losses around the peak of the flare and the
entire time-series for selected regions in the field-of-view of an
observational dataset. For the quiet-Sun pixels in the surroundings, the
estimated net radiative losses have values around 5 kW m-2, which is
comparable to previously reported values (Vernazza et al., 1981). We find that
the maximum value for a pixel located at footpoint can go up to 175 kW m-2. A
similarly high value is also reported by Díaz Baso et al. (2021), although in
a magnetic reconnection event that took place deeper down in the photosphere
and chromosphere. In the quiet-Sun the main contributor to the net radiative
losses are the Ca ii IR lines. Our analysis of the stratification of net
cooling rate shows that for a pixel located at a flare footpoint, the
contribution from Ca ii H & K and Mg ii h & k lines is also significant and
comparable to the Ca IR triplet (see Table 1). The obtained temporal evolution
of integrated radiative losses exhibit that the total integrated losses rise
sharply ($\sim$0.4 kW m-2 s-1) before the flare peak time and decreases
relatively slowly ($\sim$0.23 kW m-2 s-1) after the flare peak.

Due to sensitivity limitations of our setup in the transition region, the
inverted models do not catch the sharp and hot transition region that is
usually associated with simulations of flare atmospheres (see e.g., Rubio da
Costa et al., 2016; Simões et al., 2017; Kerr et al., 2021). Therefore, the
estimated integrated radiative losses in the Ly$\alpha$ line are very low
($\sim$7 kW m-2) compared to the contributions from the Ca ii and Mg ii lines.
We find that the contribution of Ly$\alpha$ is not significant in the
middle/lower chromosphere or below 10 kK temperature. This low contribution of
Ly$\alpha$, in the selected integrated heights, can also be related to our
relatively weak flare. A significant contribution can be noticed from
Ly$\alpha$ when the chromosphere is heated above 10 kK, which is possible in
strong flares. A recent study of the Ly$\alpha$ contrasts by Milligan (2021)
found that C-class flares, with some notable exceptions, produced a far
smaller contrast compared to M- or X-class flares. Note, however, that those
observations lacked spatial resolution and therefore don’t represent the
response of a single footpoint but of the global flare.

We have compared our inversion results with a flaring atmosphere generated by
the RADYN code. To facilitate the comparison of our semi-empirical model
retrieved from observations with the simulated flare models, we estimated the
radiative losses using the same approach described in Section 3 with the STiC
code. We find the best agreement with the RADYN simulations that have electron
beam parameters E${}_{c}=20$ keV and energy flux of $5\times 10^{10}$ erg s-1
cm-2. The RADYN atmosphere has a sharper rise in temperature compared to the
semi-empirical atmosphere inferred from observations. Given the 1D nature of
the RADYN simulation, the much wider peak in the net radiative losses derived
from the observations can be understood as a superposition in space and time
of events like the RADYN one. However, in both cases, the peak values of the
total integrated radiative losses in the chromosphere are comparable.

In this study, the radiative losses are calculated assuming statistical
equilibrium. This approximation works relatively well for Ca ii and Mg ii
transitions (Wedemeyer-Böhm & Carlsson, 2011; Kerr et al., 2019). Due to large
uncertainty in the parameters in the higher layers (¡ log${\tau_{500}}=-5$),
we did not consider the Balmer continuum in our study. However, the Balmer
continuum could also contribute significantly to the radiative cooling in
flares (Hong et al., 2022). To further investigate the role of Ly$\alpha$ and
other spectral lines that forms in the upper chromosphere we need high spatial
and temporal resolution observations spanning different layers, especially in
the upper chromosphere, which is now possible with ground-based (e.g., DKIST;
Tritschler et al. 2015 and SST) and space-based (e.g., IRIS; De Pontieu et al.
2014) observations. Additionally, the SPICE (SPICE Consortium et al. 2020)
instrument onboard the Solar Orbiter mission (Müller et al., 2020) has the
capability of sampling several strong emission lines formed at both upper
chromosphere/TR and coronal temperatures at a high spatial and temporal
resolution. Future solar missions like Solar-C/EUVST (Shimizu et al., 2019)
will also provide continuous spectral coverage throughout the solar
atmosphere. The coordinated observations among different instruments could
provide new insights to better constrain our flare models.

###### Acknowledgements.

We would like to thank the anonymous referee for the comments and suggestions.
The Swedish 1-m Solar Telescope is operated on the island of La Palma by the
Institute for Solar Physics of Stockholm University in the Spanish
Observatorio del Roque de los Muchachos of the Instituto de Astrofísica de
Canarias. The Institute for Solar Physics is supported by a grant for research
infrastructures of national importance from the Swedish Research Council
(registration number 2017-00625). RY is supported through the CHROMATIC
project (2016.0019) funded by the Knut och Alice Wallenberg foundation. This
project has received funding from the European Research Council (ERC) under
the European Union’s Horizon 2020 research and innovation program (SUNMAG,
grant agreement 759548). G.S.K. acknowledges support from NASA’s Early Career
Investigator Program (Grant 80NSSC21K0460). The inversions were performed on
resources provided by the Swedish National Infrastructure for Computing (SNIC)
at the High Performance Computing Center at Linköping University. This
research has made use of NASA’s Astrophysics Data System. We acknowledge the
community effort devoted to the development of the following open-source
packages that were used in this work: numpy (numpy.org), matplotlib
(matplotlib.org) and sunpy (sunpy.org).

## References

  * Abbasvand et al. (2020) Abbasvand, V., Sobotka, M., Heinzel, P., et al. 2020, ApJ, 890, 22
  * Allred et al. (2020) Allred, J. C., Alaoui, M., Kowalski, A. F., & Kerr, G. S. 2020, ApJ, 902, 16
  * Allred et al. (2022) Allred, J. C., Kerr, G. S., & Emslie, A. G. 2022, arXiv e-prints, arXiv:2204.11684
  * Allred et al. (2015) Allred, J. C., Kowalski, A. F., & Carlsson, M. 2015, ApJ, 809
  * Anderson & Athay (1989) Anderson, L. S. & Athay, R. G. 1989, ApJ, 346, 1010
  * Avrett et al. (1986) Avrett, E. H., Machado, M. E., & Kurucz, R. L. 1986, in Chromospheric flare models., ed. D. F. Neidig & M. E. Machado, 216–281
  * Benz (2017) Benz, A. O. 2017, Living Reviews in Solar Physics, 14, 2
  * Bjørgen et al. (2019) Bjørgen, J. P., Leenaarts, J., Rempel, M., et al. 2019, A&A, 631, A33
  * Brown (1971) Brown, J. C. 1971, Sol. Phys., 18, 489
  * Brown et al. (2018) Brown, S. A., Fletcher, L., Kerr, G. S., et al. 2018, ApJ, 862, 59
  * Carlsson & Stein (1997) Carlsson, M. & Stein, R. F. 1997, ApJ, 481, 500
  * Carmichael (1964) Carmichael, H. 1964, NASA Special Publication, Vol. 50, 451
  * de la Cruz Rodríguez et al. (2016) de la Cruz Rodríguez, J., Leenaarts, J., & Asensio Ramos, A. 2016, ApJ, 830, L30
  * de la Cruz Rodríguez et al. (2019) de la Cruz Rodríguez, J., Leenaarts, J., Danilovic, S., & Uitenbroek, H. 2019, A&A, 623, A74
  * de la Cruz Rodríguez et al. (2015) de la Cruz Rodríguez, J., Löfdahl, M. G., Sütterlin, P., Hillberg, T., & Rouppe van der Voort, L. 2015, A&A, 573, A40
  * de la Cruz Rodríguez & Piskunov (2013) de la Cruz Rodríguez, J. & Piskunov, N. 2013, ApJ, 764, 33
  * De Pontieu et al. (2014) De Pontieu, B., Title, A. M., Lemen, J. R., et al. 2014, Sol. Phys., 289, 2733
  * del Toro Iniesta & Ruiz Cobo (2016) del Toro Iniesta, J. C. & Ruiz Cobo, B. 2016, 13
  * Díaz Baso et al. (2021) Díaz Baso, C. J., de la Cruz Rodríguez, J., & Leenaarts, J. 2021, A&A, 647, A188
  * Emslie & Sturrock (1982) Emslie, A. G. & Sturrock, P. A. 1982, Sol. Phys., 80, 99
  * Fletcher et al. (2011) Fletcher, L., Dennis, B. R., Hudson, H. S., et al. 2011, Space Sci. Rev., 159, 19
  * Fletcher & Hudson (2008) Fletcher, L. & Hudson, H. S. 2008, ApJ, 675, 1645
  * Gan et al. (1993) Gan, W. Q., Rieger, E., & Fang, C. 1993, ApJ, 416, 886
  * Hirayama (1974) Hirayama, T. 1974, Sol. Phys., 34, 323
  * Holman et al. (2011) Holman, G. D., Aschwanden, M. J., Aurass, H., et al. 2011, Space Sci. Rev., 159, 107
  * Hong et al. (2022) Hong, J., Carlsson, M., & Ding, M. D. 2022, A&A, 661, A77
  * Hudson (1972) Hudson, H. S. 1972, Sol. Phys., 24, 414
  * Kerr et al. (2020) Kerr, G. S., Alaoui, M., Allred, J. C., et al. 2020, arXiv e-prints, arXiv:2009.08407
  * Kerr et al. (2019) Kerr, G. S., Carlsson, M., & Allred, J. C. 2019, ApJ, 885, 119
  * Kerr et al. (2016) Kerr, G. S., Fletcher, L., Russell, A. e. J. B., & Allred, J. C. 2016, ApJ, 827, 101
  * Kerr et al. (2021) Kerr, G. S., Xu, Y., Allred, J. C., et al. 2021, ApJ, 912, 153
  * Kontar et al. (2011) Kontar, E. P., Brown, J. C., Emslie, A. G., et al. 2011, Space Sci. Rev., 159, 301
  * Kopp & Pneuman (1976) Kopp, R. A. & Pneuman, G. W. 1976, Sol. Phys., 50, 85
  * Leenaarts et al. (2018) Leenaarts, J., de la Cruz Rodríguez, J., Danilovic, S., Scharmer, G., & Carlsson, M. 2018, A&A, 612, A28
  * Leenaarts et al. (2012) Leenaarts, J., Pereira, T., & Uitenbroek, H. 2012, A&A, 543, A109
  * Lemen et al. (2012) Lemen, J. R., Title, A. M., Akin, D. J., et al. 2012, Sol. Phys., 275, 17
  * Löfdahl et al. (2021) Löfdahl, M. G., Hillberg, T., de la Cruz Rodríguez, J., et al. 2021, A&A, 653, A68
  * Machado et al. (1980) Machado, M. E., Avrett, E. H., Vernazza, J. E., & Noyes, R. W. 1980, ApJ, 242, 336
  * Machado & Linsky (1975) Machado, M. E. & Linsky, J. L. 1975, Sol. Phys., 42, 395
  * Milligan (2021) Milligan, R. O. 2021, Sol. Phys., 296, 51
  * Milligan et al. (2014) Milligan, R. O., Kerr, G. S., Dennis, B. R., et al. 2014, ApJ, 793, 70
  * Müller et al. (2020) Müller, D., St. Cyr, O. C., Zouganelis, I., et al. 2020, A&A, 642, A1
  * Pietarila et al. (2007) Pietarila, A., Socas-Navarro, H., & Bogdan, T. 2007, ApJ, 670, 885
  * Piskunov & Valenti (2017) Piskunov, N. & Valenti, J. A. 2017, A&A, 597, A16
  * Reep et al. (2018) Reep, J. W., Russell, A. J. B., Tarr, L. A., & Leake, J. E. 2018, ApJ, 853, 101
  * Rubio da Costa et al. (2016) Rubio da Costa, F., Kleint, L., Petrosian, V., Liu, W., & Allred, J. C. 2016, ApJ, 827, 38
  * Rutten (2003) Rutten, R. J. 2003, Radiative Transfer in Stellar Atmospheres
  * Scharmer (2017) Scharmer, G. 2017, in SOLARNET IV: The Physics of the Sun from the Interior to the Outer Atmosphere, 85
  * Scharmer et al. (2003) Scharmer, G. B., Bjelksjo, K., Korhonen, T. K., Lindberg, B., & Petterson, B. 2003, in SPIE Conference Series, Vol. 4853, Proc. SPIE, 341–350
  * Scharmer et al. (2008) Scharmer, G. B., Narayan, G., Hillberg, T., et al. 2008, ApJ, 689, L69
  * Scherrer et al. (2012) Scherrer, P. H., Schou, J., Bush, R. I., et al. 2012, Sol. Phys., 275, 207
  * Shimizu et al. (2019) Shimizu, T., Imada, S., Kawate, T., et al. 2019, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 11118, UV, X-Ray, and Gamma-Ray Space Instrumentation for Astronomy XXI, 1111807
  * Simões et al. (2017) Simões, P. J. A., Kerr, G. S., Fletcher, L., et al. 2017, A&A, 605, A125
  * Sobotka et al. (2016) Sobotka, M., Heinzel, P., Švanda, M., et al. 2016, ApJ, 826, 49
  * SPICE Consortium et al. (2020) SPICE Consortium, Anderson, M., Appourchaux, T., et al. 2020, A&A, 642, A14
  * Sturrock (1966) Sturrock, P. A. 1966, Nature, 211, 695
  * Tritschler et al. (2015) Tritschler, A., Rimmele, T. R., Berukoff, S., et al. 2015, in Cambridge Workshop on Cool Stars, Stellar Systems, and the Sun, Vol. 18, 18th Cambridge Workshop on Cool Stars, Stellar Systems, and the Sun, 933–944
  * Uitenbroek (2001) Uitenbroek, H. 2001, ApJ, 557, 389
  * Uitenbroek (2002) Uitenbroek, H. 2002, ApJ, 565, 1312
  * Vernazza et al. (1981) Vernazza, J. E., Avrett, E. H., & Loeser, R. 1981, ApJS, 45, 635
  * Vissers et al. (2019) Vissers, G. J. M., de la Cruz Rodríguez, J., Libbrecht, T., et al. 2019, A&A, 627, A101
  * Wedemeyer-Böhm & Carlsson (2011) Wedemeyer-Böhm, S. & Carlsson, M. 2011, A&A, 528, A1
  * Withbroe & Noyes (1977) Withbroe, G. L. & Noyes, R. W. 1977, ARA&A, 15, 363
  * Yadav et al. (2021) Yadav, R., Díaz Baso, C. J., de la Cruz Rodríguez, J., Calvo, F., & Morosin, R. 2021, A&A, 649, A106


# Benchmarking MESA isochrones against the Hyades single star sequence

Wolfgang Brandner,1 Per Calissendorff2 and Taisiya Kopytova3,4,1
1Max-Planck-Institut für Astronomie, Königstuhl 17, 69117 Heidelberg, Germany
2Department of Astronomy, University of Michigan, Ann Arbor, MI 4810, USA
3Division of Medical Image Computing, German Cancer Research Center (DKFZ),
69120 Heidelberg, Germany
4Ural Federal University, Yekaterinburg, 620002, Russia E-mail:
<EMAIL_ADDRESS>
(Accepted 2022 August 5. Received 2022 July 29; in original form 2022 June 11)

###### Abstract

Based on GAIA EDR3, we revisit and update our sample of bonafide single stars
in the Hyades open cluster. The small observational uncertainties in parallax
and photometry of EDR3 result in a tightly defined stellar sequence, which is
ideal for the testing and calibration of theoretical stellar evolutionary
tracks and isochrones. We benchmark the solar-scaled MESA evolutionary models
against the single star sequence. We find that the non-rotating MESA models
for [Fe/H] = +0.25 provide a good fit for stars with masses above 0.85, and
very low mass stars below 0.25 M⊙. For stars with masses between 0.25 and 0.85
M⊙ the models systematically under predict the observed stellar luminosity.
One potential limitation of the models for partially convective stars more
massive than 0.35 M⊙ is the prescription of (superadiabatic) convection with
the mixing-length theory parameter $\alpha_{\rm ML}$ tuned to match the Solar
model. Below 0.35 M⊙, the increased scatter in the stellar sequence might be a
manifestation of the convective kissing instability, which is driven by
variations in the 3He nuclear energy production rate due to instabilities at
the convective core to envelope boundary. For a Hyades-like stellar
population, the application of solar-scaled models to subsolar mass stars
could result in a significant underestimate of the age, or an overestimate of
the metallicity. We suggest that future grids of solar-scaled evolutionary
stellar models could be complemented by Hyades-scaled models in the mass range
0.25 to 0.85 M⊙.

###### keywords:

open clusters and associations: individual: Hyades – convection – stars:
evolution – stars: fundamental parameters – stars: interiors – Hertzsprung-
Russell and colour-magnitude diagrams

††pubyear: 2020††pagerange: Benchmarking MESA isochrones against the Hyades
single star sequence–Benchmarking MESA isochrones against the Hyades single
star sequence

##  1 Introduction

Stellar cluster and individual binary stars constitute the most important
astrophysical calibration sources. At an average distance of 45 pc, the Hyades
open cluster is the closest (populous) stellar cluster to the Sun. The Hyades
have super-solar metallicity. Kopytova et al. (2016) derived [Fe/H]=+0.14 for
the best-fitting BT-Settl2010+PISA and DARTMOUTH isochrones (Allard et al.,
2013; Dotter et al., 2008; Tognelli et al., 2011; Da Rio & Robberto, 2012;
Tognelli et al., 2012), while Gossage et al. (2018) derived [Fe/H]=+0.10 to
+0.12 using MESA isochrones (Paxton et al., 2011; Paxton et al., 2013, 2015;
Dotter, 2016; Choi et al., 2016; Paxton et al., 2018) in near infrared (NIR).

At an age of $\approx 635\pm 135$ Myr, the Hyades open cluster comprises main
sequence and post-main sequence stars with initial masses in the range
$\approx$0.1 to 3.6 M⊙ (e. g. Perryman et al., 1998; de Bruijne et al., 2001;
Krumholz et al., 2019), which serve as benchmarks for models of stellar
evolution. Colour-absolute magnitude diagrams (CMD) (Castellani et al., 2001;
Röser et al., 2011) revealed discrepancies between stellar models and
observations for sub-solar mass stars. Castellani et al. (2001) attributed
this to limitations in the description of the efficiency of superadiabatic
convection in the outer layers of partially convective stars. Kopytova et al.
(2016) showed that the incorporation of updated input physics (equation of
state, opacities, etc., Degl’Innocenti et al., 2008) results in an improved
match between theoretical isochrones and 2MASS photometric measurements in the
mass range 0.6 to 0.8 M⊙. Below 0.6 M⊙, the close-to-vertical (i.e. constant
colour) stellar sequence in NIR CMDs makes isochronal fitting rather
insensitive to stellar luminosity, and hence metallicity or age.

In Brandner et al. (2022) we used MESA and BHAC2015 (Baraffe et al., 2015)
isochrones in the GAIA photometric system to determine the age of the nearby
exoplanet host star GJ 367. Both sets of isochrones suggested a young age in
the range of $\approx$30 to 60 Myr for the star, which is considerably younger
than its age suggested by gyro-chronology, and by its space motion and
galactic dynamics models. This and the unprecedented photometric and parallax
accuracy of GAIA EDR3 observations prompted us to benchmark the solar-scaled
MESA isochrones against the single star sequence of the Hyades open cluster.

The structure of the paper is as follows. In section 2 we present the update
sequence of bonafide single stars in the Hyades open cluster. In section 3 we
summarize literature age and metallicity estimates of the Hyades based on
photometric data. In section 4 we benchmark the MESA isochrones against the
single star sequence in the GAIA photometric system. In section 5 we discuss
potential short-comings of grids of stellar models, and suggest ways forward.

##  2 The Hyades single star sequence

Table 1: Single and multiple candidate members of the Hyades cluster,
classified according to GAIA EDR3, sorted by RA. The full table is available
online.

GAIA EDR3 ID | dpgeo | lo$\\_$dpgeo | hi$\\_$dpgeo | G | $\sigma_{\rm G}$ | BP | $\sigma_{\rm BP}$ | RP | $\sigma_{\rm RP}$ | flag
---|---|---|---|---|---|---|---|---|---|---
| [pc] | [pc] | [pc] | [mag] | [mag] | [mag] | [mag] | [mag] | [mag] |
395696646953688448 | 59.612 | 59.565 | 59.655 | 10.7791 | 0.0009 | 11.4353 | 0.0004 | 10.0031 | 0.0002 | 1
386277165192421120 | 30.935 | 30.908 | 30.962 | 14.7670 | 0.0089 | 16.6299 | 0.0030 | 13.4805 | 0.0009 | 2
393017579491591168 | 64.815 | 64.423 | 65.115 | 17.0372 | 0.0073 | 19.2764 | 0.0380 | 15.6597 | 0.0031 | 1
420637590762193792 | 83.610 | 82.524 | 84.646 | 18.5022 | 0.0007 | 20.9101 | 0.0994 | 17.0266 | 0.0068 | 1
385502112574538624 | 42.821 | 42.780 | 42.859 | 14.3729 | 0.0005 | 16.1299 | 0.0041 | 13.1075 | 0.0015 | 1
2860677398591440768 | 39.751 | 39.395 | 40.054 | 11.7972 | 0.0017 | 12.8881 | 0.0017 | 10.6002 | 0.0022 | 4

> flag: 1 - bonafide single, 2 - likely binary or multiple, 3 - white dwarf, 4
> - peculiar GAIA EDR3 BP-G vs. G-RP colours; Median , low (lo, 16th quantile)
> and high (hi, 84th quantile) of the photogeometric distance posterior dpgeo
> are from Bailer-Jones et al. (2021)

Figure 1: Colour-absolute magnitude diagram of members of the Hyades open
cluster based on GAIA EDR3 parallaxes and photometry. Blue dots indicate all
stars in the GAIA GCNS sample with good photometry. Red crosses indicate
members of our bonafide single star sequence. Overplotted are MESA isochrones.
While the isochrone matching the age of the Hyades ($\approx$710 Myr, dashed
line) provides a good fit to the upper and lower main sequence, it
underpredicts the luminosity of stars in the mass range $\approx$0.25 to 0.85
M⊙. The 1$\sigma$ error bars for the single stars are based on the
uncertainties in photometry and parallax.

Kopytova et al. (2016) defined a fiducial observational sequence of single
stars suitable for testing of stellar evolutionary and atmospheric models. The
stars were selected from a sample of 724 probable members of the Hyades open
cluster established by Röser et al. (2011) based on their proper motion
according to the PPMXL catalog (Röser et al., 2010). Using literature data and
high-angular resolution Lucky Imaging observations with AstraLux Norte
(Hormuth et al., 2008), the stars were screened for stellar binarity and
photometric blends to derive a sample of single stars. This single star sample
was selected quite conservatively. Considering the intrinsic 2′′ angular
resolution of the 2MASS Point Source Catalog (Cutri et al., 2003), and in
order to minimize the effect of photometric blends, Kopytova et al. (2016)
flagged all occurrences of another source within $\approx$4′′ as potential
binary companions, and excluded them from the single star sample.

GAIA EDR3 facilitates a refinement of the single star sequence from Kopytova
et al. (2016). Gaia Collaboration et al. (2021) published a GAIA Catalogue of
Nearby Stars (GCNS) listing 920 candidate members of the Hyades. In order to
reject photometric outliers caused, e.g., by blends in the GAIA BP and RP
bands, we first applied a colour cut-off: $-0.2$ mag $\leq$ BP - G $\leq$ 3.2
mag and $-0.3$ mag $\leq$ G - RP $\leq$ 1.7 mag. This rejected 30 sources,
resulting in a sample of 890 candidate members listed in Table 1. As a second
step we fitted a 4th order polynomial to the data in a G - RP vs. BP - G two-
colour diagram. Application of an iterative sigma-clipping resulted in a
sample of 783 candidate members with good photometric quality data. As a third
step, we used the Renormalized Unit Weight Error (RUWE, see GAIA-Collaboration
(2020); Lindegren et al. (2021)) to distinguish between bonafide single stars
and likely unresolved binary and multiple systems. RUWE values around 1.0
indicate that the GAIA astrometric observations are well fitted by the single-
star model. A significantly larger RUWE value indicates that the single-star
model does not provide a good fit to the astrometric solution due to, e.g.,
the non-single nature of the source.

Next we computed absolute Gabs magnitudes based on the apparent G magnitudes
and the EDR3 photogeometric distances according to Bailer-Jones et al.
(2021)111For 18 stars, including five of the Hyades white dwarfs (GAIA EDR3 ID
45980377978968064, 3313606340183243136, 3313714023603261568,
3306722607119077120, 3294248609046258048), we substituted the missing
photogeometric distance by the geometric distance.. In order to identify and
flag photometric binaries (i.e. sources falling on the binary sequence in the
CMD, but with separations too close to be identified by the RUWE selection),
we fitted an 8th order polynomial to the single star main-sequence, and
applied an iterative sigma clipping. This resulted in 616 sources classified
as bonafide single stars, 156 as likely binary or multiple systems, and 11 as
white dwarfs.

Figure 2: Same as Figure 1, but zoomed in on very low-mass stars with masses
$\lessapprox$0.35 M⊙.

Figure 1 shows the colour-absolute magnitude diagram of the Hyades open
cluster, covering the main sequence and part of the post-main sequence. Blue
dots mark the candidate members of the Hyades from the GCNS sample. Red
crosses mark bonafide single stars, with observational uncertainties derived
from uncertainties in GAIA photometry and parallax indicated. The median
uncertainty in BP-RP colour amounts to 4.1 mmag, and to 3.3 mmag in Gabs. In
particular for BP-RP $\leq 2.3$ mag, the single stars form a very tight
sequence, which is clearly distinct from the scatter of apparently
overluminous Hyades members located on the binary sequence. For redder (and
intrinsically fainter) stars, there is a larger scatter in the bonafide single
star sequence.

As Figure 2 highlights, the single star sequence becomes successively
incomplete for BP-RP $\geq 3.2$ mag. The sample of GCNS Hyades candidate
members extends to BP-RP = 4.7 mag. The reddest and lowest mass object
included in the single star sample is LSPM J0354+2316 (GAIA EDR3
65638443294980224), which is of spectral type M8 (Bardalez Gagliuffi et al.,
2014), and has a mass of $\approx$0.1 M⊙ (Goldman et al., 2013).

##  3 Age and metallicity of the Hyades

There is a vast literature on abundance estimates and the calibration of
astrophysical parameters for stars in the Hyades cluster (see, e.g., Perryman
et al. (1998); Tognelli et al. (2021), and references therein). Age estimates
for the Hyades are in general derived from the main sequence turn-off and
isochrone fitting. Abundance estimates rely both on spectral analysis and
isochrone fitting.

Table 2: Compilation of abundance, $\alpha_{\rm ML}$, and age estimates for
the Hyades based on isochrone fitting

[Fe/H] | X | Y | Z |  $\Delta$Y/$\Delta$Z | $\alpha_{\rm ML}$ | mass range | age | PD1 | reference
---|---|---|---|---|---|---|---|---|---
|  |  |  |  |  | [M⊙] | [Myr] |  |
$+0.14$${}^{+0.05}_{-0.05}$ | 0.716 |  $0.260$${}^{+0.020}_{-0.020}$ |  $0.024$${}^{+0.003}_{-0.003}$ |  | 1.64 | [0.8,1.6] | $625\pm 50$ | BD | Perryman et al. (1998)
+0.14 | 0.691 | 0.285 | 0.024 |  | 1.68 | [0.5,0.9] | $638\pm 13$ | TY | de Bruijne et al. (2001)
+0.14 | 0.716 | 0.260 | 0.024 |  | 1.64 | [0.9,1.6] | $638\pm 13$ | TY | de Bruijne et al. (2001)
+0.14 | 0.708 | 0.273 | 0.019 |  | 1.68 | [1.6,2.4] | $631$ | TY | de Bruijne et al. (2001)
$+0.14$ | 0.700 | 0.283 | 0.0175 | 2 | 1.74 | [0.13,2.30] | $726\pm 50$ | 2M | Kopytova et al. (2016)
$+0.24$${}^{+0.02}_{-0.02}$ |  |  |  |  | 1.82 | [0.5,2.4] | $726\pm 50$ | TY | Gossage et al. (2018)
$+0.10$${}^{+0.02}_{-0.02}$ |  |  |  |  | 1.82 | [0.5,2.4] | 741${}^{+36}_{-14}$ | 2M | Gossage et al. (2018)
$+0.169$${}^{+0.025}_{-0.025}$ | 0.6947 | 0.2867 | 0.01863 | $2.03\pm 0.33$ | $2.01\pm 0.05$ | [0.83,1.35] | 500 | G2 | Tognelli et al. (2021)

> 1 key to photometric data set (PD): 2M - based on 2MASS photometry (Cutri et
> al., 2003); BD - based on BDA (Mermilliod, 1995); G2 - based on GAIA DR2
> (Gaia Collaboration et al., 2016, 2018); TY - based on TYCHO photometry (Høg
> et al., 2000)

Table 2 summarizes some of the canonical estimates, including the helium-to-
metal enrichment ratio $\Delta$Y/$\Delta$Z, based on isochrone fitting. The
differences in the parameter estimates can in part be explained by variations
in the observational methods and data sets and their intrinsic uncertainties,
in part by differences and advances in the modelling of stellar interiors and
atmospheres, and in part by advances in the analysis of the solar elemental
abundances (see, e.g., Asplund et al. (2009)). The majority of the estimates
focused on post-main sequence stars and main sequence stars of spectral type K
and earlier (more massive than 0.5 M⊙). The sole exception is the study by
Kopytova et al. (2016), which includes stars with masses as low as 0.13 M⊙.
Some of the studies also consider variations of the mixing length (ML)
parameter $\alpha_{\rm ML}$ or stellar rotation (Gossage et al., 2018;
Tognelli et al., 2021). The latter effect appears to be most noticeable in the
colours and luminosity of post-main sequence stars.

Common to all studies is the derived (or assumed) super-solar metallicity of
the Hyades, with [Fe/H] estimates in the range $+0.10$ to $+0.24$. Age
estimates for the Hyades cover the range 500 to 770 Myr. The majority of the
studies also indicate a higher than solar He abundance for the Hyades (see
Table 2).

##  4 Benchmarking isochrones

Figure 3: Same as Figure 1, but varying the stellar metallicity of the 710 Myr
isochrones. Table 3: Model stellar parameters at inflection points between
MESA isochrone (age = 710 Myr, [Fe/H]=+0.25, v/vcrit = 0) and Hyades single
star sequence

Gabs1 | BP - RP1 | B - V | log Teff | log g | log L | Mass | note
---|---|---|---|---|---|---|---
[mag] | [mag] | [mag] | [K] | [cm/s2] | [L⊙] | [M⊙] |
$<$6.7 |  $<$1.2 |  $<$1.03 |  $>$3.676 |  $<$4.61 |  $>$-0.594 |  $>$0.85 | good fit
6.7 to 11.2 | 1.2 to 2.6 | 1.03 to 1.30 | 3.676 to 3.519 | 4.61 to 4.90 | -0.594 to -1.864 | 0.85 to 0.35 | model underluminous
11.2 to 13.5 | 2.6 to 3.1 | 1.30 to 1.36 | 3.519 to 3.485 | 4.90 to 4.97 | -1.864 to -2.229 | 0.35 to 0.25 | transition region
$>$13.5 |  $>$3.1 |  $>$1.36 |  $<$3.485 |  $>$4.97 |  $<$-2.229 |  $\leq$0.25 | good fit

> 1 Gabs and BP-RP refer to the observed stellar sequence, while the other
> quantities are according to the MESA isochrone for the corresponding BP-RP
> colour.

The solar scaled MESA isochrones and stellar tracks222We use the
MIST˙v1.2˙vvcrit0.0˙UBVRIplus packaged model grid, dated 2020-12-04, which
includes updated synthetic photometry for GAIA EDR3 based on Riello et al.
(2021). The grid steps are 0.25 dex in the range -2.00 $\leq$ [Fe/H]
$\leq$+0.50, and 0.05 dex in the range 5.0 $\leq\log_{10}$(age [yr]) $\leq$
10.3 use a Ledoux plus mixing length theory prescription of convection, with
$\alpha_{\rm ML}=1.82$ tuned to fit the Sun (Choi et al., 2016).

In Figure 1 we overlay two MESA isochrones for [Fe/H]=+0.25, and no rotation
(v/vcrit = 0) on the colour-absolute magnitude diagram of the Hyades. As
presented in Gossage et al. (2018), the best fitting MESA isochrones in the
optical yield systematically higher metallicities than the best fitting MESA
isochrones in the NIR for the Hyades, Praesepe, and Pleiades
clusters.333Gossage et al. (2018) attribute this to the small number of stars
in their optical samples ($<$40 stars for the Hyades, according to their
Figure 7), and suggest that their optical samples do not provide meaningful
constraints on the metallicity of either of these clusters. In the pre-
computed grid of solar-scaled MESA isochrones, [Fe/H]=+0.25 is closest to
[Fe/H]=+0.24$\pm$0.01 as deduced by Gossage et al. (2018) from the best
fitting MESA isochrone in the TYCHO BT ,VT photometric system. The choice of
non-rotating stellar models is based on the dearth of rapid rotators among
single stars with masses $\geq$0.3 M⊙ in the Hyades (Douglas et al., 2016).
They find that stars with masses of $\approx$0.4 M⊙ from the sample defined by
Kopytova et al. (2016) have typical rotational periods of $\approx$20 days.

We find that the isochrone for $\log_{10}$(age [yr]) = 8.85 ($\approx$710 Myr)
provides a better fit for stars with masses $>$1.35 M⊙ than the next younger
(630 Myr) or older isochrones (795 Myr). This age is in good agreement with
isochronal age determinations by Kopytova et al. (2016) and Gossage et al.
(2018). Tognelli et al. (2021) only considered stars with masses $<$1.5 M⊙ for
the age determination, and were thus less sensitive to the rapid evolution of
stellar luminosity near the upper end of the main sequence. For stars between
0.25 and 0.85 M⊙, the 710 Myr isochrone tends to underpredict the stellar
luminosity. In the colour range BP-RP = 1.2 to 2.6 mag, the 55 Myr isochrone
(dash-dotted line) provides a good fit to the observed sequence, but it
overpredicts the stellar luminosity for BP-RP $\geq$2.6 mag.

In Figure 3 we overlay three isochrones for an age of 710 Myr, and for [Fe/H]
= 0.00, +0.25, and +0.50. The highest metallicity isochrone (dash-dotted line)
provides a good fit the colour range BP-RP = 1.2 to 2.6 mag, but overpredicts
the stellar luminosity for bluer and redder stars.

Table 3 lists the corresponding B-V colour, effective temperature, $\log$ g,
$\log$ L, and stellar mass according to the MESA isochrone (age = 710 Myr
([Fe/H]=+0.25, and v/vcrit = 0) at the boundaries of the four BP-RP colour
regions marked by the vertical dotted lines in Figures 1 and 3. None of the
single-age, single-metallicity isochrones is capable of fitting the entire
single star sequence.

##  5 Discussion

In Figures 1 and 3 we have marked the BP-RP colour regions where the 710 Myr,
[Fe/H] = +0.25 isochrone provides a good fit to the observed sequence, and
where it significantly deviates from the observed sequence by predicting
fainter (underluminous) stars. The good fit of the MESA isochrones for solar-
type stars with masses above 0.85 M⊙, and for very low-mass stars with masses
between 0.1 and 0.25 M⊙ is highly encouraging, and speaks for the maturity of
stellar model grids.

Potential problem areas in modelling CMDs of the Hyades have been noticed in
the literature. In general, the challenges considered for the Hyades stellar
sequence are related to stellar rotation, elemental abundance ([Fe/H] and
$\Delta{\rm Y}/\Delta{\rm Z}$) and nuclear energy production rates, and the
description of convection. Gossage et al. (2018) discusses the effect of
rotation and variations in $\alpha_{\rm ML}$ for stars above 1.2 M⊙ in the
Hyades. They conclude that a larger $\alpha_{\rm ML}=2.0$ would provide a
better fit to the giant stars in the Hyades. For Hyades members with masses
less than 0.85 M⊙, Castellani et al. (2001) noted a discrepancy between
observed and theoretical optical colour and brightness, in particular in the
mass range where superadiabatic convection dominates the outer convective zone
of a star. They suggest to consider $\alpha_{\rm ML}$ as a free parameter in
this mass range which could be tuned to a (lower) value to better describe the
efficiency of superadiabatic convection. Based in part on the ideas initially
explored by Stevenson (1979) and MacDonald & Mullan (2014), Ireland & Browning
(2018) study the effects of rotation and magnetic fields on convection, and
investigate how they could be parameterized by a depth dependent $\alpha_{\rm
ML}$.

The good match of the synthetic photometry of the MESA isochrone to the
observed GAIA data for BP-RP $>$3.1 mag suggests that there is no generic
issue in the ATLAS12 and SYNTHE conversion (Choi et al., 2016) of luminosity
and temperature to synthetic photometry for cool stellar photospheres. Choi et
al. (2016) use different sets of opacity tables for the mass ranges 0.1 to 0.3
M⊙, 0.3 to 0.6 M⊙, and $>$0.6 M⊙ as boundary conditions for the atmospheres.
The discrepancy for 1.2 mag$<$BP-RP$<$3.1 mag (Teff = 3050 to 4750 K, m = 0.25
to 0.85 M⊙), thus might warrant a review of the opacity tables and their
transitions in this parameter range. Choi et al. (2016) also point out
systematic differences in the evolutionary track of a 0.3 M⊙ star between the
MIST and Lyon (Baraffe et al., 1998, 2003, 2015) models on the one side, and
the PARSEC (Girardi et al., 2002; Marigo et al., 2008; Bressan et al., 2012)
models on the other side. They attribute this difference to the modified
temperature-Rosseland mean optical depth (T-$\tau$) relation (Chen et al.,
2014) employed for low-mass stars by the PARSEC models.

The colour range around BP-RP $\approx$2.80 to 2.95 mag (G${}_{\rm abs}\approx
10.5$ to 11 mag) stands out in the CMD as the observed stellar sequence shows
a larger scatter than for stars with bluer colours or BP-RP $>$ 3.1 mag
(G${}_{\rm abs}>11.0$ mag, Figures 1 and 2). According to the MESA tracks,
this correspond to stellar masses just below 0.35 M⊙, which roughly coincides
with the fully convective boundary. As discussed by Baraffe & Chabrier (2018),
for stars in this mass range the energy productions rate of the proton-proton
I branch (3He + 3He $\longrightarrow$ 4He + 2 p) is crucial for a proper
description of the stellar luminosity. Depending on the precise stellar
properties and the initial He abundance, 3He in Hyades members in this mass
range might not yet have reached its equilibrium abundance. An overabundance
in 3He resulting in an enhanced energy productions rate could explain the
observed overluminosity of stars compared to the 710 Myr isochrone. The
increased scatter in the stellar magnitude-colour sequence could suggest the
presence of an instability in the stellar luminosity for this particular mass
and age range. van Saders & Pinsonneault (2012) identified a 3He-driven
instability for stars near the fully convective boundary, which they referred
to as convective kissing instability. Baraffe & Chabrier (2018) confirmed the
existence of this instability using a different evolutionary code. Caveats are
that the observed increase in the scatter of Gabs by $\approx$0.2 mag (20%) is
about a factor 2 to 4 larger than the variations in luminosity according to
the models by van Saders & Pinsonneault (2012). For stellar interior models of
solar metallicity, the convective kissing instability seems to be restricted
to a relatively narrow mass range of 0.34 to 0.37 M⊙ (van Saders &
Pinsonneault, 2012; Baraffe & Chabrier, 2018).

In the Hyades, a stellar mass of $\approx$0.30 M⊙ marks the boundary between
very low mass stars with rotational periods ranging from a fraction of a day
to 5 days, and more massive single stars with rotational periods in the range
of 10 to 20 days (Douglas et al., 2016; Douglas et al., 2019). Fast rotation
and the associated stellar activity has been associated with radius inflation
(Somers & Pinsonneault, 2015; Somers & Stassun, 2017). As discussed by Feiden
& Chaboyer (2014) and Feiden (2016) a significant radius inflation requires
strong interior magnetic fields in the range of 10 MG, which might be
difficult to maintain over the age of the Hyades. Douglas et al. (2016)
discuss that poloidal fields in Hyades members in the mass range 0.3 to 0.6 M⊙
resulted in effective magnetic braking, and strongly reduced stellar activity,
while fully convective stars of lower mass primarily rely on their (weak)
stellar winds for shedding angular momentum. We suggest that a future study
could focus on the rotation periods and activity levels of the stars in the
mass range 0.25 and 0.35 M⊙, and look for, e.g., correlations with their
luminosity, or photometric variability.

As of mid 2022, Choi et al. (2016) had more than 1200 citations, with the MESA
models being used to assess ages, metallicity, radii, effective temperatures,
luminosity, surface gravity, etc. for individual stars and (complex) stellar
populations. The primary science of more than 130 of these articles is on
exoplanets, where in general planetary properties are derived relative to the
astrophysical properties of the host star. Biased stellar properties could
thus directly bias the deduced properties of exoplanets.

The example of the Hyades single star sequence highlights the potential perils
in the analysis of low-mass (late-type) stellar populations. Even in the
presence of GAIA high-precision parallax and photometric information, missing
supplemental information could result in biased conclusions on absolute
stellar ages, metallicity, or the intrinsic spread of these properties. A
‘blind’ isochronal analysis of the Hyades sample in the BP-RP colour range 1.2
to 2.6 mag might underestimate the true age by more than a factor of 10 (55 vs
710 Myr, see Figure 1), or result in a significant overestimate of its
metallicity ([Fe/H] = +0.50 vs. +0.25, see Figure 3).

The updated single star sequence of the Hyades cluster with its accurate GAIA
EDR3 distance and photometric measurements could serve as a reference to tune
modelling parameters like, e.g., $\alpha_{\rm ML}$ to the efficiency of super-
adiabatic convection, or to tune astrophysical parameters like, e.g.,
$\Delta{\rm Y}/\Delta{\rm Z}$ (and in particular 3He abundances) to reflect
actual energy productions rates. We suggest that future grids of solar-scaled
evolutionary models should be tested against the Hyades single star sequence
presented in Table 1, or against comparable data sets for the Pleiades or
Praesepe open clusters.

## Acknowledgements

We thank H.-W. Rix for the initial discussion, which prompted this research.
We thank the anonymous referee for constructive comments, which helped to
improve the paper.

This work has made use of data from the European Space Agency (ESA) mission
Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing
and Analysis Consortium (DPAC,
https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has
been provided by national institutions, in particular the institutions
participating in the Gaia Multilateral Agreement.

## Data availability

The data underlying this article are available in the article and in its
online supplementary material. The online version of Table 1 includes
coordinates, which makes objects discoverable in VizieR.

## References

  * Allard et al. (2013) Allard F., Homeier D., Freytag B., 2013, Mem. Soc. Astron. Italiana, 84, 1053
  * Asplund et al. (2009) Asplund M., Grevesse N., Sauval A. J., Scott P., 2009, ARA&A, 47, 481
  * Bailer-Jones et al. (2021) Bailer-Jones C. A. L., Rybizki J., Fouesneau M., Demleitner M., Andrae R., 2021, AJ, 161, 147
  * Baraffe & Chabrier (2018) Baraffe I., Chabrier G., 2018, A&A, 619, A177
  * Baraffe et al. (1998) Baraffe I., Chabrier G., Allard F., Hauschildt P. H., 1998, A&A, 337, 403
  * Baraffe et al. (2003) Baraffe I., Chabrier G., Barman T. S., Allard F., Hauschildt P. H., 2003, A&A, 402, 701
  * Baraffe et al. (2015) Baraffe I., Homeier D., Allard F., Chabrier G., 2015, A&A, 577, A42
  * Bardalez Gagliuffi et al. (2014) Bardalez Gagliuffi D. C., et al., 2014, ApJ, 794, 143
  * Brandner et al. (2022) Brandner W., Calissendorff P., Frankel N., Cantalloube F., 2022, MNRAS, 513, 661
  * Bressan et al. (2012) Bressan A., Marigo P., Girardi L., Salasnich B., Dal Cero C., Rubele S., Nanni A., 2012, MNRAS, 427, 127
  * Castellani et al. (2001) Castellani V., Degl’Innocenti S., Prada Moroni P. G., 2001, MNRAS, 320, 66
  * Chen et al. (2014) Chen Y., Girardi L., Bressan A., Marigo P., Barbieri M., Kong X., 2014, MNRAS, 444, 2525
  * Choi et al. (2016) Choi J., Dotter A., Conroy C., Cantiello M., Paxton B., Johnson B. D., 2016, ApJ, 823, 102
  * Cutri et al. (2003) Cutri R. M., et al., 2003, VizieR Online Data Catalog, p. II/246
  * Da Rio & Robberto (2012) Da Rio N., Robberto M., 2012, AJ, 144, 176
  * Degl’Innocenti et al. (2008) Degl’Innocenti S., Prada Moroni P. G., Marconi M., Ruoppo A., 2008, Ap&SS, 316, 25
  * Dotter (2016) Dotter A., 2016, ApJS, 222, 8
  * Dotter et al. (2008) Dotter A., Chaboyer B., Jevremović D., Kostov V., Baron E., Ferguson J. W., 2008, ApJS, 178, 89
  * Douglas et al. (2016) Douglas S. T., Agüeros M. A., Covey K. R., Cargile P. A., Barclay T., Cody A., Howell S. B., Kopytova T., 2016, ApJ, 822, 47
  * Douglas et al. (2019) Douglas S. T., Curtis J. L., Agüeros M. A., Cargile P. A., Brewer J. M., Meibom S., Jansen T., 2019, ApJ, 879, 100
  * Feiden (2016) Feiden G. A., 2016, A&A, 593, A99
  * Feiden & Chaboyer (2014) Feiden G. A., Chaboyer B., 2014, ApJ, 789, 53
  * GAIA-Collaboration (2020) GAIA-Collaboration 2020, Gaia eDR3 source catalogue "light", VO resource provided by the GAVO Data Center, https://dc.zah.uni-heidelberg.de/tableinfo/gaia.edr3lite
  * Gaia Collaboration et al. (2016) Gaia Collaboration et al., 2016, A&A, 595, A1
  * Gaia Collaboration et al. (2018) Gaia Collaboration et al., 2018, A&A, 616, A1
  * Gaia Collaboration et al. (2021) Gaia Collaboration et al., 2021, A&A, 649, A6
  * Girardi et al. (2002) Girardi L., Bertelli G., Bressan A., Chiosi C., Groenewegen M. A. T., Marigo P., Salasnich B., Weiss A., 2002, A&A, 391, 195
  * Goldman et al. (2013) Goldman B., et al., 2013, A&A, 559, A43
  * Gossage et al. (2018) Gossage S., Conroy C., Dotter A., Choi J., Rosenfield P., Cargile P., Dolphin A., 2018, ApJ, 863, 67
  * Høg et al. (2000) Høg E., et al., 2000, A&A, 355, L27
  * Hormuth et al. (2008) Hormuth F., Hippler S., Brandner W., Wagner K., Henning T., 2008, in McLean I. S., Casali M. M., eds, Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series Vol. 7014, Ground-based and Airborne Instrumentation for Astronomy II. p. 701448 (arXiv:0807.0497), doi:10.1117/12.787384
  * Ireland & Browning (2018) Ireland L. G., Browning M. K., 2018, ApJ, 856, 132
  * Kopytova et al. (2016) Kopytova T. G., Brandner W., Tognelli E., Prada Moroni P. G., Da Rio N., Röser S., Schilbach E., 2016, A&A, 585, A7
  * Krumholz et al. (2019) Krumholz M. R., McKee C. F., Bland-Hawthorn J., 2019, ARA&A, 57, 227
  * Lindegren et al. (2021) Lindegren L., et al., 2021, A&A, 649, A2
  * MacDonald & Mullan (2014) MacDonald J., Mullan D. J., 2014, ApJ, 787, 70
  * Marigo et al. (2008) Marigo P., Girardi L., Bressan A., Groenewegen M. A. T., Silva L., Granato G. L., 2008, A&A, 482, 883
  * Mermilliod (1995) Mermilliod J.-C., 1995, in Egret D., Albrecht M. A., eds, ASSL Vol. 203, Information & On-Line Data in Astronomy. p. 127, doi:10.1007/978-94-011-0397-8_12
  * Paxton et al. (2011) Paxton B., Bildsten L., Dotter A., Herwig F., Lesaffre P., Timmes F., 2011, ApJS, 192, 3
  * Paxton et al. (2013) Paxton B., et al., 2013, ApJS, 208, 4
  * Paxton et al. (2015) Paxton B., et al., 2015, ApJS, 220, 15
  * Paxton et al. (2018) Paxton B., et al., 2018, ApJS, 234, 34
  * Perryman et al. (1998) Perryman M. A. C., et al., 1998, A&A, 331, 81
  * Riello et al. (2021) Riello M., et al., 2021, A&A, 649, A3
  * Röser et al. (2010) Röser S., Demleitner M., Schilbach E., 2010, AJ, 139, 2440
  * Röser et al. (2011) Röser S., Schilbach E., Piskunov A. E., Kharchenko N. V., Scholz R.-D., 2011, A&A, 531, A92
  * Somers & Pinsonneault (2015) Somers G., Pinsonneault M. H., 2015, MNRAS, 449, 4131
  * Somers & Stassun (2017) Somers G., Stassun K. G., 2017, AJ, 153, 101
  * Stevenson (1979) Stevenson D. J., 1979, Geophysical and Astrophysical Fluid Dynamics, 12, 139
  * Tognelli et al. (2011) Tognelli E., Prada Moroni P. G., Degl’Innocenti S., 2011, A&A, 533, A109
  * Tognelli et al. (2012) Tognelli E., Degl’Innocenti S., Prada Moroni P. G., 2012, A&A, 548, A41
  * Tognelli et al. (2021) Tognelli E., Dell’Omodarme M., Valle G., Prada Moroni P. G., Degl’Innocenti S., 2021, MNRAS, 501, 383
  * de Bruijne et al. (2001) de Bruijne J. H. J., Hoogerwerf R., de Zeeuw P. T., 2001, A&A, 367, 111
  * van Saders & Pinsonneault (2012) van Saders J. L., Pinsonneault M. H., 2012, ApJ, 751, 98