Bridging theory and observations in stellar pulsations: The impact of convection and metallicity on the instability strips of Classical and Type-II Cepheids

The effect of metallicity on the theoretical and empirical period-luminosity (PL) relations of Cepheid variables is not well understood and remains a highly debated issue. Here, we examine empirical colour-magnitude diagrams (CMDs) of Classical and Type-II Cepheids in the Magellanic Clouds and compare those with the theoretically predicted instability strip (IS) edges. We explore the effects of incorporating turbulent flux, turbulent pressure, and radiative cooling into the convection theory on the predicted IS at various metallicities using MESA-RSP. We find that the edges become redder with the increasing complexity of convection physics incorporated in the fiducial convection sets, and are similarly shifted to the red with increasing metallicity. The inclusion of turbulent flux and pressure improves the agreement of the red edge of the IS, while their exclusion leads to better agreement with observations of the blue edge. About 90% of observed stars are found to fall within the predicted bluest and reddest edges across the considered variations of turbulent convection parameters. Furthermore, we identify and discuss discrepancies between theoretical and observed CMDs in the low effective temperature and high luminosity regions for stars with periods greater than ~ 20 days. These findings highlight the potential for calibrating the turbulent convection parameters in stellar pulsation models or the prediction of a new class of rare, long-period, 'red Cepheids', thereby improving our understanding of Cepheids and their role in cosmological studies.


INTRODUCTION
Pulsating stars serve as excellent astrophysical laboratories for testing theories of stellar structure, evolution, and pulsation (Cox 1980).Radially pulsating stars such as Classical Cepheids (CCs), RR Lyraes (RRLs), and Type-II Cepheids (T2Cs) provide stringent constraints for these models by enabling a quantitative comparison of their observed and predicted pulsation properties.These stars are of additional interest as they obey well-defined Period-Luminosity (PL) relations, making them important distance indicators as well as proxies for studying galactic structures (e.g., Deka et al. 2022a;Bhuyan et al. 2023).
CCs are intermediate-mass evolved stars with ∼ 3 − 20 M ⊙ (Bono et al. 2000b;Anderson et al. 2016;Musella 2022;Espinoza-Arancibia et al. 2022).They are Population I stars located in the instability strip ★ E-mail:mamideka8@gmail.com (IS) of the Hertzsprung-Russell diagram (HRD).For stars within the IS, small departures from hydrostatic equilibrium are amplified because of the increase in the opacity of the gas located in the ionization zones of the outer envelope on compression (Cox 1974).The pulsations are mainly driven by the kappa mechanism and reinforced by the gamma mechanism (Catelan & Smith 2015).
CCs are used as the primary calibrators for the first rung of the cosmic distance ladder (Riess et al. 2022, and references therein) and play a crucial role in the determination of the local value of the Hubble constant ( 0 ).However, the ongoing discord between  0 values derived from cosmic microwave background data (Planck Collaboration et al. 2020) and measurements based on CCs and Type Ia supernovae (SNe) (Riess et al. 2022) warrants a detailed investigation of possible systematics in the CC PL relations.For example, the addition of metallicity term to the CC PL relations is found to improve the accuracy of the distance ladder (Breuval et al. 2022;Riess et al. 2022;Bhardwaj et al. 2023).The effect of metallicity on the PL relation has been investigated using both observed (Ripepi et al. 2020;Riess et al. 2021;Romaniello et al. 2022;Breuval et al. 2022;Trentin et al. 2023;Bhardwaj et al. 2024) and theoretical data (Bono et al. 1999(Bono et al. , 2008;;Caputo et al. 2000;Marconi et al. 2005;Anderson et al. 2016;De Somma et al. 2022).The value of the coefficient in the metallicity term (also known as the  term in literature) in the empirical period-luminosity-metallicity (PLZ) relation is found to vary over a wide range from −0.5 to 0.0 mag dex −1 depending on wavelength from near-infrared to optical band (Udalski et al. 2001;Riess et al. 2021;Ripepi et al. 2022;Breuval et al. 2022;Freedman & Madore 2023;Trentin et al. 2023).On the other hand, most of the aforementioned studies based on non-linear convection models have found  > 0, except Anderson et al. (2016) who obtained  < 0 using linear models that incorporate the effect of rotation.The debate surrounding the sign of  term from the empirical and theoretical points of view thus remains unresolved and necessitates further investigations.
The slope and intercept of PL relation are strongly coupled with the slope and intercept of IS edges (Sandage & Tammann 2008).A steeper/shallower instability strip implies a steeper/shallower PL slope (Simon & Young 1997).To improve the accuracy of the PL relation, a precise estimation of IS edges is essential.Besides, the study of the IS can offer valuable insights into the physical mechanisms of the stars residing within this region.Also, the effect of metallicity on the PL relation is strongly correlated with the effect of metallicity on the IS edges.The IS boundaries of CCs have been studied extensively using theoretical models (Bono et al. 2000b;Fiorentino et al. 2002;Anderson et al. 2016;De Somma et al. 2022).The works of Bono et al. (2000a), Fiorentino et al. (2002) and De Somma et al. (2022) have found a strong dependence of metallicity on the topology of the IS based on the theoretical formulation as outlined by Bono et al. (2000a,c).However, Anderson et al. (2016) found a weak dependence of metallicity on the IS topology based on linear non-adiabatic (LNA) pulsation models when including the effect of rotation in their models.Therefore, in conjunction with different chemical compositions, varying input physics applied in computing the models also influences the location of the IS edges, thereby impacting the corresponding PL relationship.In this work, we carry out a detailed study of the effect of turbulent convection theory (more specifically turbulent convection parameters) and metallicity on the IS morphology of CCs.
In addition to CCs, we also include T2Cs in this work to study the effect of metallicity and turbulent convection parameters on the topology of their IS.As the Population-II counterparts of CCs, they present us with a unique opportunity to improve our understanding of these effects in the old, low-mass, metal-poor stellar populations.Similar to CCs, they also obey a tight PL relation, enabling their use as distance indicators (Smolec 2016;Das et al. 2021;Bhardwaj et al. 2017;Bhardwaj 2022;Das et al. 2024), particularly in regions where CCs are not abundant and RRLs are very faint, thereby underscoring their crucial role in the cosmic distance scale.Nevertheless, the evolutionary and pulsation characteristics of T2Cs have been relatively understudied compared to other radial pulsators, despite their useful role in the distance determinations within the Milky Way and other nearby/satellite galaxies (e.g., Bhardwaj et al. 2017;Wielgórski et al. 2022;Sicignano et al. 2024).
T2Cs are low-mass stars (∼ 0.5 − 0.8 M ⊙ ) with periods between 1 − 50 days, although there remains some uncertainty regarding the precise upper and lower period limits (Catelan & Smith 2015;Smolec 2016).These stars are categorized into three distinct classes based on their periods: BL Herculis (BL Her), with periods between 1−4 days; W Virginis (W Vir), with periods between 4 − 20 days; and RV Turbulent pressure,   0 0 Tauris (RV Tau), with periods above 20 days and correspond to different evolutionary stages.We exclude RV Tauris stars based on their periods due to the ambiguity surrounding their upper mass limit (Bódi & Kiss 2019) and diversity in their progenitors-ranging from young-massive to old stars in binary systems (Manick et al. 2018).
Henceforth, when we refer to T2Cs in this work, it encompasses only BL Her and W Vir stars.The present study makes use of both observed and theoretical data of CCs and T2Cs.This combined study will not only give us observationally consistent IS, but also put constraints on the input physics involved.The observed data are taken from OGLE-IV (Soszyński et al. 2015(Soszyński et al. , 2018)).We generate the theoretical data using Modules for Experiments in Stellar Astrophysics (MESA, Paxton et al. 2010Paxton et al. , 2013Paxton et al. , 2015Paxton et al. , 2018Paxton et al. , 2019;;Jermyn et al. 2023).
We generate the theoretical data using the new Radial Stellar Pulsations (RSP; Paxton et al. 2019) module in MESA for modelling linear/non-linear radial pulsations exhibited by variable stars like CCs, T2Cs, and RRLs (Smolec & Moskalik 2008;Paxton et al. 2019).We carry out an LNA radial pulsation stability analysis to obtain the growth rates of radial pulsation modes, their linear periods, and the mean luminosity.From the growth rates, we then determine the boundaries of the instability strip (IS).
MESA-RSP enables us to investigate the impact of various turbulent convection parameters on the properties of the models.It implements the time-dependent turbulent convection model by Kuhfuss (1986) following the stellar pulsation treatment from Smolec & Moskalik (2008) Four sets of turbulent convection, labeled A, B, C, and D, are defined based on different combinations of these parameters (Paxton et al. 2019).Set A corresponds to the simplest convection model, Set B adds radiative cooling; Set C adds turbulent pressure and turbulent flux; and set D includes all these effects simultaneously (Paxton et al. 2019).The free parameters that enter the convective models are kept the same as provided in Table 3 and 4 of Paxton et al. (2019) (Smolec & Moskalik 2008).Radiative cooling   is set to 2/3 following Wuchterl & Feuchtinger (1998).For further details on these parameters, we refer the interested reader to Kuhfuss (1986); Wuchterl & Feuchtinger (1998); Yecko et al. (1998); Smolec & Moskalik (2008).
The location of red and blue edges as well as the width of the IS are highly influenced by the value of ,   and   .Generally,   is attributed to the excitation, and   and   to the damping of stellar oscillations.For example, Yecko et al. (1998) found the turbulent eddy viscosity as the primary agent for damping pulsations.In case of red giant variables, Xiong (2021) have found that their contributions to the excitation and damping effect change depending on the stellar parameters ,  and  eff and the frequency of the mode.However, it is important to note that the turbulent convection model (Kuhfuss 1986) integrated in MESA-RSP is not equivalent to that of Yecko et al. (1998) and Xiong (2021), and that the assumptions and simplifications made in these models are different as well.Furthermore, the driving and damping of stellar pulsation cannot be exclusively attributed to the free parameters of the convection model.These might vary from model to model and a detailed study on a specific model is necessary to find out the main agents for the damping and driving mechanism (for details see Houdek & Dupret 2015).
In addition to the IS, the resulting pulsation growth rate, radial velocity and light curve shapes are also highly sensitive to these parameters (Paxton et al. 2019).It must be noted that the fiducial turbulent convection parameters given by Paxton et al. (2019) are not calibrated to observations.In a recent study, Kovács et al. (2023) and Kovács et al. (2024) thoroughly examined the impact of these turbulent convection parameters on the radial velocity and light curves of several RRLs and proposed revised values for certain parameters.The present study instead focuses on the impact of the fiducial turbulent convection parameters as outlined in Paxton et al. (2019) on the ensemble-based properties of the observed stars.We plan to extend this work in a future paper by considering further variations in these parameters individually for each star, aiming for a more rigorous comparison with observations.The remainder of the paper is organized as follows.In Section 2, we discuss the theoretical and observed data.Section 3 presents the methodology to carry out the analysis.We discuss the results of the analysis in Section 4. Finally, a summary and discussion of the implications of our study are presented in Section 5.

Theoretical Data
We have computed static models using the following input stellar parameters: mass (), luminosity (), effective temperature ( eff ), metal fraction () and hydrogen fraction ().We have constructed two fine grids of pulsation models for CCs and T2Cs after a thorough study of the previous literature (Bono et al. 2000a,c;Anderson et al. 2016;De Somma et al. 2022;Bhardwaj 2020;Das et al. 2021).
Although the typical lower mass limit for CCs is reported to be ∼ 3M ⊙ in most literature (e.g., Pilecki et al. 2022;Kurbah et al. 2023), our investigation reveals that a grid spanning a mass range of (3 − 11)M ⊙ fails to reproduce observations corresponding to lower magnitudes/periods.Through a rough estimation, we find that the mass of a CC with a period of 0.8 days (the lowest period observed for an FU CC from OGLE-IV data) is approximately ∼ 2.3M ⊙ , using the period-mass relation from Turner (1996).Following this, we adopt a wide range of stellar masses (2−11 M ⊙ ) to closely match the observed properties.We also adopt an effective temperature range of (4000 − 8000 ) in the CC grid.
The lower/upper limit of luminosity for a particular mass is set by decreasing/increasing the mass-luminosity relation by 0.25 dex (Anderson et al. 2016).Three different metallicities [Fe/H] = (0.00, −0.34, −0.75) suitable for the Galaxy, Large Magellanic Cloud (LMC) and Small Magellanic Cloud (SMC), respectively have been adopted for the input grid following Romaniello et al. (2008); Gieren et al. (2018).
Similarly, the masses of T2Cs are quoted to be in the range (0.5−0.8) M ⊙ in literature (Smolec 2016;Bhardwaj 2020;Das et al. 2021).However, we have considered a wider range of (0.3−1.0) M ⊙ in order to match with observations (also see Das et al. 2024).Additionally, we have considered a range of luminosities and effective temperatures spanning (10 − 1000) L ⊙ and (4000 − 8000) K, respectively.The upper limit for luminosity is slightly extended to account for a few very high-luminosity T2Cs observed in the data.T2Cs from different regions are found to have different metallicities (Welch 2012).A few of them as obtained in previous literature are listed as follows.T2Cs from all the Galactic globular clusters are found to have metallicities (in dex): [Fe/H] < −1.0 (Clement et al. 2001); Galactic field region: −1.0 < [Fe/H] < 0.0 (Schmidt et al. 2011) and for the Bulge: −1.4 < [Fe/H] < 0.6 (Harris & Wallerstein 1984;Wallerstein 2002) and the majority of field T2Cs from the Galactic Halo are found to have [Fe/H] > −0.9.Hence, we adopt a wide range of metallicities spanning −2.0 < [Fe/H] < 0.0 for T2C input grid.Since our model-to-observation comparisons are exclusively focused on data from the LMC and SMC, this adopted range of [Fe/H] will be sufficient for the scope of this study.
Using the aforementioned ,  and  eff values, we have constructed two fine grids for both the classes of variables (shown in Fig. A1).The grid is constructed using Sobol numbers, i.e., quasirandom low discrepancy numbers (Sobol 1967).The advantage of using such a grid over a linear grid or a random grid is that Sobol numbers are uniformly distributed and fill the input parameter space more rapidly.For more details on the use of Sobol numbers to construct a grid, the readers are referred to Bellinger et al. (2016).Final grids are constructed by mapping the compositions  and  for each class with the Sobol grid (, ,  eff ).Table 3 presents a summary of the input parameters of the constructed grids.
RSP solves the pulsation equations on a Lagrangian mesh.The mesh structure is divided into inner and outer zones with respect to a specified anchor temperature (see Fig. 2 of Paxton et al. 2019).Initially, we have kept the total number of Lagrangian mass cells  total = 200 with  outer = 60 cells above the anchor.However, if this configuration failed to create the initial model, then we have considered  total = 150 and  outer = 40.The choice of numerical parameters: inner boundary temperature ( inner = 2 × 10 6 K) and anchor temperature ( anchor = 11, 000 K) are kept as outlined in Paxton et al. (2019) (for more details, see Appendix B).These parameters are defined in such a way that we achieve a good resolution of the pulsation driving region.An example inlist used for the LNA computation is provided in https://github.com/mami-deka/LNA_CMD.

Observational Data
In order to compare the theoretical findings with the observed data, the observed optical (, )-band light curves of fundamental (FU), first-overtone (FO) mode CCs and T2Cs are taken from the OGLE-IV database (Soszyński et al. 2015(Soszyński et al. , 2018)).The light curves are available for both LMC and SMC.In the present analysis, well-sampled light curves with more than 30 data points are chosen both for CCs and T2Cs.The stars marked as uncertain in the database are also excluded.A concise overview of the observed data is outlined in Table 4.

Theoretical Data Processing
The theoretical stellar grids (discussed in Section 2) are used as input pulsation models for the LNA computation.From the computed models, we choose those models with a positive growth rate of the radial fundamental mode and linear periods between 0.8 − 200 days for the FU CCs.Similarly for the FO CCs, a positive growth rate of the radial first overtone mode and inear periods less than 10 days are chosen for further analysis.Likewise, T2Cs models are chosen from the models which have a positive growth rate of the radial fundamental mode and linear periods between 1−20 days.The period range for theoretical models is chosen based on the observed periods.We use the pre-computed bolometric corrections (BCs) tables from Lejeune et al. (1998)

Observational Data Processing
To compare theoretical prediction with observation in CMD, we require absolute magnitudes and extinction corrected colour.First, we phase the light curves using the following equation: where  0 ,  and  represent the epoch of maximum light for the -band, the pulsation period of a star in days and the times of observations, respectively.These values are taken from the OGLE-IV database.
Then to obtain the mean magnitude from each light curve, we have fitted the phased light curves with a Fourier sine series (Deb & Singh 2009;Deka et al. 2022b) of the following form: where  0 and  = 2   denote the mean magnitude and the angular frequency, respectively.The th order Fourier coefficients are represented by   and   .The order of the fit  is obtained using Baart's criteria (Deb & Singh 2009, and references therein) by varying it between (6 − 10), (6 − 10), (4 − 8) for the -band light curves of FU CCs, FO CCs and T2Cs respectively, while it is varied from 4 − 8 for the corresponding -band light curves.
The optical reddening map of Skowron et al. (2021) for the Magellanic Clouds is used to obtain the extinction-corrected observed apparent mean magnitudes.The  ( − ) reddening values obtained from Skowron et al. (2021) map are converted into the extinctions   and   using the relations   = 1.5  ( − ) and   = 2.5  ( − ) (Skowron et al. 2021).The true apparent magnitudes are then converted into the corresponding absolute magnitudes assuming a mean distance modulus of  = 18.48 mag (Pietrzyński et al. 2019) and  = 18.977 mag (Graczyk et al. 2020) for the LMC and SMC, respectively.
Furthermore, some of the observed data are discarded based on the approach following Madore et al. (2017) due to their higher horizontal deviations in the CMD.The CMDs in Fig. A2 highlight these deviated points in black colour.These deviated flares are presumably due to uncertainties associated with the magnitude and colour determinations and the reddening values adopted for each star.We employ a 2.5 clipping to the linear fit on the magnitude residuals of the period-luminosity (PL) relation versus the corresponding residuals of the period-Wesenheit (PW) relation to remove the stars with high vertical dispersion.

Instability strip edges
We now aim to determine the IS blue and red edges using the theoretical data generated by MESA-RSP.We use machine learning to train a linear classifier using the Support Vector Machine algorithm that we optimize using Stochastic Gradient Descent (SGD).
We note that while a non-linear classifier would be appropriate for some of the models considered, it is challenging to establish a simple relation coming out of such a classifier.Hence, we have enforced linear IS edges to provide simple relations that can be adopted for a direct comparison with observations.We train the classifier on the theoretical data to distinguish models with negative and positive growth rates, and hence to find the boundaries of the IS.We use the implementation from the scikit-learn python package (Pedregosa  Romaniello et al. (2008), the metallicity range of the CCs span from −0.18 dex to +0.25 dex for the Galaxy, −0.62 dex to −0.10 dex for the LMC, and −0.87 dex to −0.63 dex for the SMC.This variation may further expand with the inclusion of larger spectroscopic data.In the recent study by Bhardwaj et al. (2024), the metallicity range for the Galactic sample has already extended from −0.18 dex to −1.1 dex for metal-poor CCs and from +0.25 dex to +0.6 dex for the metal-rich ones.Hence, we have considered the entire range of chemical compositions from the input grid while estimating IS edges.The fraction of observed CCs and T2C variables (with and without outliers removal) within predicted IS boundaries are shown in Fig. 4. It is evident from this figure that the IS boundaries predicted from our model grid of CCs are able to include the observed stars within the IS very well.However, more complex convection physics captures more FU CCs while these numbers decrease for FO CCs and T2Cs with the increasing complexity of the convection physics.We emphasize that the "increasing complexity of convection physics" here refers to the increasing number of turbulent convection parameters integrated into the fiducial convection sets of MESA-RSP.For instance, set A comprises 5 turbulent convection parameters, set B includes 6, set C incorporates 7, and set D comprises 8 parameters.In the calculation of the fraction of observed stars, we have taken into account the uncertainties in intrinsic colour and absolute magnitude.We have propagated these uncertainties into our calculation by averaging over 10 5 random realizations.Specifically, we took into account uncertainties in apparent mean magnitude (derived from the Fourier fit), distance modulus (obtained from Pietrzyński et al. (2019); Graczyk et al. (2020)), and reddening (obtained from Skowron et al. (2021)).We note that the dominant source of uncertainty in absolute magnitude and intrinsic colour is coming from the uncertainties in reddening.
Similarly, the CMD for LMC and SMC FO CCs are displayed in Fig. 2 and Fig. C2, respectively.The IS edges exhibit a nearly parallel alignment for FO CCs.
We further present the comparative CMD for LMC/SMC T2Cs along with their predicted edges in Fig. 3 and Fig. C3, respectively.It is evident that the edges are becoming redder with the increasing complexity of convection sets in this case as well.It should be noted that the number of unstable T2C theoretical models generated in the low magnitude regime is also small.This may hint at constraints within MESA-RSP or possibly highlight gaps in our understanding and selection of turbulent convection parameters.
The slopes and intercepts of the predicted IS edges in colourmagnitude and  −  eff -plane are listed in Tables 5 and D1, respectively.It is noteworthy that the proportion of observed stars within each set is fairly consistent, except for the FO CCs set C, which contains nearly 70% of the total observed stars.

Effect of turbulent convection parameters on the IS
Here, we discuss the effect of turbulent convection parameters on the predicted IS edges.We have found that the IS edges become redder with the increasing complexity of the convection sets for both CCs and T2Cs.In other words, the set A blue edge is the bluest one and the set D red edge is the reddest one among the IS edges.
Within the set of eight turbulent convection parameters, the values of four parameters (,   ,   and  d ) are kept the same for all the four sets.The remaining four parameters (  ,   ,   , and   ) assume distinct values for each set.The effect of set A and set B on the IS edges looks similar, while the same is true for the set C and set D. Furthermore, the difference in the topology of IS edges with sets (A and B) and sets (C and D) is very distinct, attributed mainly to the parameters   and   .Turbulent pressure   contributes to the driving of pulsation, while convective flux   and eddy viscosity   act as damping mechanisms.We find the inclusion of   and   shifts the IS towards a relatively lower temperature as compared to set A and B, reducing the  eff obtained after LNA calculation by ∼ 100 − 200 K.Additionally, they also damp pulsations in some models near the blue edge while driving the same in some models near the red edge.Furthermore, Kovács et al. (2023) have found a correlation between   (in their case, it is denoted as   ) and  eff , in which increasing value of   is suitable for decreasing  eff in case of RR Lyrae stars.In Table 1, it can be seen that   increases from set A to set D except set C. Since we have also found that the IS edges become redder from set A to set D, it may indicate that for CCs also,   is dependent on  eff .Nevertheless, we still do not have a clear picture of the exact contributions of these turbulent convection parameters at low/high temperatures, which we plan to follow up in a future study.
To quantify the shift towards lower temperature with increasing complexity of the convection theory, we have determined the mean effective temperature from each distribution (for each convection set) by fitting a Gaussian.The mean  eff of set A is relatively similar to that of set B, while the mean  eff of set C is similar to that of set D. The difference in mean  eff from set (A,B) to set (C,D) is (150 − 200) K.The variations in the mean effective temperature of the models across each convection sets and as a function of metal fractions  are shown in Fig. 5.

Comparison with literature
We perform a comparative analysis between our predicted IS edges and those derived from recent investigations by De Somma et al. (2022) and Espinoza-Arancibia et al. ( 2023) in the  eff -plane to better constrain our findings.The study by De Somma et al. (2022) uses non-linear convective theoretical pulsation models, incorporating a mass-luminosity relation based on Bono et al. (2000a) that excludes the effects of mass loss, rotation, and core-overshooting (referred to as 'A').Additionally, they consider two alternative luminosity levels for the same mass by increasing the luminosity level by 0.2 dex ('B') and 0.4 dex ('C').The corresponding ISs are computed, accounting for the influences of superadiabatic convection efficiency with  ML = (1.5, 1.7, 1.9), and varying metallicities (,  ) = (0.004, 0.25), (0.008, 0.25), and (0.03, 0.28).The slopes and intercepts of the IS edges overplotted in Fig. 6 are taken from Table .10 of De Somma et al. ( 2022) with ML relation 'C' and 'B' and  ML = 1.5 (as referred in their paper) for FU and FO CCs, respectively, as these IS edges show the closest agreement with our predictions.The slopes show a close match for FU CCs with that of our work, while disparities in intercepts arise from variations in the adopted chemical compositions.Interestingly, for FO CCs, the slope of their red edge aligns more closely with ours, in contrast to their blue edge.For both FU and FO CCs, the blue edge is observed to be significantly cooler, while the width of the IS appears to be wider as compared to that of De Somma et al. (2022).Comparison with observed OGLE-IV LMC T2Cs with our estimated set A/B blue edges, while the lower empirical blue edge corresponds closely to our set C/D blue edges.However, we can see some distinct differences in our predicated red edges with their empirical lower/upper red edges.This also highlights the importance of optimising the current sets of turbulent convection parameters in MESA-RSP as the position of the red edge is highly dependent on the complex interaction between convection and pulsation.

Effect of metallicity
We have also investigated the effect of metallicity on the IS boundaries for all convection sets.The slopes and intercepts of IS edges corresponding to different metallicities are listed in Table D2.It has been observed that the IS edges generally exhibit a redder shift with increasing metallicities for both CCs and T2Cs, with a few exceptions in specific convection set scenarios.This is found to be consistent with De Somma et al. (2022).This is a consequence of the enhanced opacity in stellar envelopes with higher metal content.Specifically, an increase in metal fraction () shifts the blue edge towards lower temperature, while the increase in helium fraction ( ) shifts the blue edge towards higher temperature (Cox 1980).In our model grids, the  values increase with an increase in the  values.Therefore, we observe a very subtle effect of metallicity on the computed IS edges.The mean effective temperature ( eff ) is quantified for the model distribution for the individual metallicity values and is shown in Fig. 5.

Flare in the theoretical colour-magnitude diagram
Unstable models show "flaring" in the low  eff and high  region beyond the predicted IS red edge in the present work, as shown in Fig. 7.These flares are seen in models with periods greater than 20 days with all four convection sets, A-D.The flare is relatively larger in sets C and D as compared to sets A and B. However, the observed data do not show such flares.The effects of convection become very significant near the red edge of IS (Smolec & Moskalik 2008).Some preliminary findings from the flared models: (i) The flared models consist of ∼ 6% of the total pulsating models generated in this work.
(ii) The coolest of the flared models extend nearly to the Hayashi limit.(iii) With increasing metallicity, the number of flared models tends to increase.
(iv) Inclusion of radiative cooling   reduces the number of flared models.
(v) The pulsation can be damped in the flared models by increasing the value of eddy viscosity   and convective flux   from the fiducial values.Notably, adjusting   has minimal impact on the linear pulsation period, while changes in   significantly influence the period.This also suggests the necessity of establishing   as a function of  eff for CCs similar to RRLs (Kovács et al. 2023(Kovács et al. , 2024)).

Cepheid evolutionary models
To show the consistency of the predicted IS edges with the evolutionary models, we have computed stellar evolutionary tracks for Cepheids with mass range (3 − 11) M ⊙ in steps of 1 M ⊙ using MESA r15140.We have adopted the LMC composition for the calculation (see Section 2.1 and Table 2).We treat convection using the mixing-length theory (MLT) from Henyey et al. (1965) and the boundary of the convective region is determined using the Ledoux criterion.We use the step overshooting prescription integrated in MESA.We use a value of 0.2 pressure scale heights for incorporating overshooting of core and envelope convection zones.We have implemented the overshooting during the core hydrogen burning and shell non-burning stages.To include the effects of semi-convection (Langer et al. 1983), we use an efficiency parameter of 0.1 for models with mass less than 9 M ⊙ and 100 for masses 9 − 11 M ⊙ .The evolution is computed from the pre-main sequence stage till the completion of the blue loop.The computed tracks are shown in Fig. 8 along with the predicted IS edges for CCs.For this computation, we have used the 5M_cepheid_blue_loop inlist from the mesa-r15140 test_suite with a few modifications.For details, see the inlists available at https://github.com/mami-deka/LNA_CMD.We note that the evolutionary tracks computed for this project are indicative of the evolutionary path followed by CCs.The calibration of  MLT and overshooting parameters with observation are beyond the scope of the present study.Additionally, to explain the dearth of long-period observed FU CCs and the flared theoretical models, we determine the speed of the evolutionary models in the HRD using the formulation as given in Bellinger et al. (2023): Here, we set the normalization factor  to be the zero-age main sequence (ZAMS) velocity of a 3 M ⊙ star.The "speed" here is not a physical velocity but rather indicates how quickly the star moves through the HRD throughout its evolution, and hence indicates the relative longevity of stellar models in different phases of evolution.
For instance, a log  HR value of 0 within a particular region of the HRD implies an expectation of the same stellar population density as that of a 3M ⊙ star at the ZAMS, assuming the birthrate is the same.Conversely, a log  HR value of 3 indicates an expectation of a thousand times fewer stars in that region (Bellinger et al. 2023).We can thus explain the existence of CCs in the blue loops from Fig. 8 as the models in this region have the lowest speeds.However, for high-mass (long-period) CCs, the speed of the models increases relatively even during the blue-loop phase, which accounts for the fewer observed CCs in that region in addition to the observed initial mass function which has few massive stars.This could also offer a plausible explanation for the absence of observed stars in the flared region.

SUMMARY AND CONCLUSIONS
We have carried out a comparative study between the observed and theoretical CMDs of FU CCs, FO CCs and T2Cs.We used the observed data from the OGLE-IV database and generated the theoretical data using MESA-RSP.We constructed dense theoretical grids of pulsation models with stellar parameters , ,  eff , ,  and different sets of turbulent convection parameters for both CCs and T2Cs.The grids are constructed in such a way that they produce approximately similar values for the minimum and maximum magnitude/colour limits to those of the observed target stars.Then, using an ML algorithm, we estimated the IS edges from the theoretical data.The results obtained from the present analysis are summarized as follows: (i) We find that ∼ 90% of observed stars (CCs and T2Cs) fall within the predicted Set A blue edge and Set D red edge.We find the following theoretical IS boundaries to include the majority of the observed stars: (ii) We find discrepancies between theoretical and observed CMDs of FU CCs in the low-effective temperature and high luminosity regions for stars with periods greater than ∼ 20 days.We find unstable theoretical models in that region where no observed stars are found yet.These flared models might be useful for optimizing the turbulent convection parameters for CCs, or for potentially predicting a class of rare, long-period, red Cepheids.However, these predictions may only be a consequence of the chosen turbulent convection parameters, and may not be realized in nature.It is worth noting that if these red Cepheids do exist, then these results suggest they may arise through an evolutionary path that is distinct from that of CCs, namely, outside of blue loops while high up on the red giant branch.They may represent a new class of variable stars, or potentially belong to a previously recognized class of pulsating stars, such as the one found by Clark et al. (2015).A detailed investigation will provide a deeper insight into their nature which we intend to carry out in a follow-up study.However, regardless of their classification, they are expected to exhibit a PL relation that is noticeably less steep than that of FU CCs (see Fig. D1).
(iii) The IS edges tend to become redder as the complexity of the convection sets increases.The locations of the IS edges (more specifically the red edge) are highly influenced by the turbulent pressure and turbulent flux parameters.The present analysis also indicates a correlation between   and  eff .Hence, this indicates that the appropriate turbulent convection parameters may vary as a function of stellar parameters.
(iv) The IS edges also become redder with increasing metallicities, except in a few cases involving specific convection sets.
(v) MESA-RSP is not able to generate T2C theoretical pulsating models similar to the observations in lower magnitude-lower temperature region with the fiducial turbulent convection and numerical parameters.Detailed investigation and calibration of the parameters involved are required in this regard.
The present work is based on a linear analysis of CC and T2C models obtained using MESA-RSP.We have found that the present two grids of models are able to reproduce the average properties of the observed stars very well.Additionally, it is evident that the combination of different sets of turbulent convection parameters is necessary to reproduce the observed data.We plan to investigate the effect of turbulent convection parameters on the light curves using non-linear computation in a future work.We also plan to extend this analysis from the optical to the near-infrared band in the future in order to better constrain the physics of these stars.Comparison with observed OGLE-IV SMC FU CCs Comparison with observed OGLE-IV SMC FO CCs to obtain the absolute magnitude of a model in the  and -bands.Detailed conversion steps are given in Paxton et al. (2018) and briefly summarized in Das et al. (2021).

Figure 1 .
Figure 1.Comparison of the color-magnitude diagram between observed fundamental-mode Cepheids in the LMC (purple dots) and theoretical models (orange dots) for four different prescriptions of convection.The theoretical blue and red edges from each convection set are shown with lines.The black cross in the top-left panel corresponds to the median of the typical uncertainties in the observed data.Set A refers to the simplest convection model.Sets B, C and D additionally include radiative cooling; turbulent pressure and turbulent flux; and their combination, respectively.

Figure 2 .
Figure 2. Same as Fig. 1, here with first-overtone Classical Cepheids in the LMC.

Figure 3 .
Figure 3. Same as Fig. 2, here with Type-II Cepheids in the LMC.

Figure 4 .Figure 5 .
Figure 4. Percentage of observed Cepheids within the predicted instability strip across different assumptions in input physics in the convection theory.The upper/lower panels show the LMC/SMC.Fundamental and first-overtone Cepheids are in the left and middle panel, and Type-II Cepheids are in the right panel.

Figure 6 .Figure 7 .
Figure 6.Hertzsprung-Russell diagram of the instability strip (IS) edges for Cepheids across different convection sets (A, B, C, D) from this work compared with recently determined empirical and theoretical IS edges from the literature.Red and blue colours represent the theoretical red and blue edges obtained in this work.The left panel shows fundamental-mode Classical Cepheids and the right panel shows first-overtone Classical Cepheids.The black lines show the theoretical IS edges from De Somma et al. (2022).The orange and magenta lines represent the empirical IS edges from Espinoza-Arancibia et al. (2023) with and without a break, respectively.

Figure 8 .
Figure 8. Cepheid evolutionary tracks are shown from core H exhaustion to core He exhaustion.The tracks are computed with masses (3 − 11) M ⊙ in steps of 1 M ⊙ .The FU (solid lines) and FO (dashed lines) IS edges as given in Section 5 are overplotted, where blue and red colours represent blue and red edges respectively.

Figure A1 .Figure A2 .
Figure A1.Scatterplot matrix of the input stellar parameters: mass (), luminosity () and temperatures ( eff ) for CCs (left panel) and T2Cs (right panel).The observed underdensity of low-mass stars in the CC grid (left panel) can be attributed to the grid's construction in two segments.Specifically, the grid was divided into two parts: the first part encompassing stellar masses ranging from (3 − 11)M ⊙ and the second part covering (2 − 3)M ⊙ .

Figure B1 .
Figure B1.The difference in growth rates and periods is plotted against the different input parameters taken for set A. Models (1 − 3) for fundamental and first-overtone mode Classical Cepheids are chosen randomly from the Classical Cepheids input grid, while models (1 − 3) for Type II Cepheids are selected from Type II Cepheids input grid.

Figure D1 .
Figure D1.The plot shows the theoretical pulsating models of FU classical Cepheids in period-luminosity plane.The orange dots are the models within our predicted IS edges while the red dots represent the flared models.

Table 1 .
Paxton et al. (2019) parameter sets fromPaxton et al. (2019).Set A is the simplest one; set B includes radiative cooling; set C includes turbulent pressure and turbulent flux, and set D includes all of them.
. It is mainly governed by eight free dimensionless parameters: : Mixing length of the convection   : Eddy viscosity   : Source of turbulence   : Convective energy transport   : Dissipation of turbulence   : Turbulent pressure   : Turbulent flux   : Radiative losses of convective eddy energy.
. For convenience, the parameters are listed in Table 1.When turbulent pressure, turbulent flux, and radiative cooling are neglected (i.e.,   =   =   = 0) and values of   ,   and   are kept the same as in Table 1, the time-independent version of the Kuhfuss (1986) model reduces to the standard mixing-length theory (MLT; Böhm-Vitense 1958).The turbulent pressure (  ) and convective flux (  ) parameters were introduced by Yecko et al. (1998) which were not present in the original Kuhfuss (1986) model.These values are set to   = 2/3 and   ≡

Table 2 .
Chemical compositions of the adopted models.

Table 3 .
Stellar parameter range for CCs and T2Cs used in constructing the pulsation model grids.

Table 4 .
Stars used in the present study.Besides the determination of IS edge, we have also compared the theoretical data with the observed data in CMD to see how well our theoretical models are able to reproduce the observed properties.The CMD for LMC and SMC FU CCs are as shown in Figs. 1 and C1.In

Table 5 .
Slopes and intercepts of IS edges for FU CCs, FO CCs and T2Cs in the CMD.These IS boundaries are estimated including the entire range of chemical compositions considered while constructing the two input grids in this work.Class of Variables Convection Sets IS edges ( =  ( −  ) + )

Table B1 .
The stellar parameters used for investigating the sensitivity of the numerical parameters ( total ,  outer ,  inner and outer  total ) on the growth rates and linear periods.These models are selected randomly from the original input grids.

Table D1 .
Slopes and intercepts of IS edges for FU CCs, FO CCs, and T2Cs in the HRD.
⊙ ) = c log T eff + d) Width of the IS in  eff (in K) *

Table D2 .
Slopes and intercepts of IS edges for FU CCs, FO CCs and T2Cs at different metalicities in the CMD.Class of Variable Metal Fraction () Convection SetIS edges ( =  ( −  ) + ) FU CC theoretical models in PL plane