Probing core overshooting using subgiant asteroseismology: the case of KIC10273246

The size of convective cores remains uncertain, despite its substantial influence on stellar evolution, and thus on stellar ages. The seismic modeling of young subgiants can be used to obtain indirect constraints on the core structure during main sequence, thanks to the high probing potential of mixed modes. We selected the young subgiant KIC10273246, observed by Kepler, based on its mixed-mode properties. We thoroughly modeled this star, with the aim of placing constraints on the size of its main sequence convective core. We first extracted the parameters of the oscillation modes of the star using the full Kepler data set. To overcome the challenges posed by the seismic modeling of subgiants, we proposed a method which is specifically tailored for subgiants with mixed modes and consists in a nested optimization. We then applied this method to perform a detailed seismic modeling of KIC10273246. We obtained models that show good statistical agreements with the observations, both seismic and non-seismic. We showed that including core overshooting in the models significantly improves the quality of the seismic fit, optimal models being found for $\alpha_{\mathrm{ov}} = 0.15$. Higher amounts of core overshooting strongly worsen the agreement with the observations and are thus firmly ruled out. We also found that having access to two g-dominated mixed modes in young subgiants allows us to place stronger constraints on the gradient of molecular weight in the core and on the central density. This study confirms the high potential of young subgiants with mixed modes to investigate the size of main-sequence convective cores. It paves the way for a more general study including the subgiants observed with Kepler, TESS, and eventually PLATO.


Introduction
One of the most important current open questions in stellar physics is the extent of convective cores. Several physical processes are known to extend the convective core boundaries beyond the standard Schwarzschild limit. The most frequently quoted are overshooting of ascending blobs of fluids due to their inertia, rotational mixing or semi-convection. All these processes remain poorly described by theory and the way they interact is even less understood. They are therefore generally modeled, in stellar evolution codes, as an extension of the mixed core over a distance d ov , which is often referred to as the distance of overshoot, even though other processes can contribute as well. In practice, this can either be achieved by simply extending the fully mixed central region, or by treating overshooting as a diffusive process, following a behavior found in numerical simulations (Freytag et al. 1996). In both cases, a free parameter controlling the extent of the additional mixing is required. Observations are therefore necessary to better constrain those processes.
First constraints have been obtained thanks to the study of the HR diagram of clusters (see e.g. Maeder & Mermilliod 1981, VandenBerg et al. 2006, and the modeling of eclipsing binaries (Claret & Torres 2016). Most of those studies favor adding overshooting, with various values for d ov , generally around 0.2 H p , H p being the pressure scale height. Claret & Torres (2016) found Send offprint requests to: A. Noll e-mail: anoll@irap.omp.eu that α ov , the ratio between d ov and H p , increases with mass for stars under 2 M before reaching a plateau. However, this result is still debated (Constantino & Baraffe 2018).
Over the last decade, asteroseismology allowed us to probe the structure of stellar cores. Thanks to the data of CoRoT (Baglin et al. 2006), Kepler (Borucki et al. 2010) and now TESS (Ricker et al. 2014) missions, we have been able to precisely measure the oscillation frequencies of numerous pulsators. The study of pressure (p) modes in low-mass main sequence (MS) stars, showed the need for core overshooting to correctly reproduce the observed frequencies , Silva Aguirre et al. 2013. Deheuvels et al. (2016), modeling several MS stars, found that α ov increases with the mass. Moreover, gravity (g) mode pulsators, like slowly-pulsating B (SPB) stars, are interesting targets to constrain the additional mixing around convective cores. Indeed, gravity modes probe the inner chemical structure of the star and allow fine investigation of the convective core extensions. Moravveji et al. (2015Moravveji et al. ( , 2016, modeling SPB stars, found that overshoot was necessary and favored diffusive overshooting over a simple extension of the central mixed region. Post main sequence stars are another way to put constraints on the amount of overshooting. Once the central hydrogen is exhausted, nuclear energy production stops, leaving an inert radiative helium core. This core then contracts, heating the surrounding hydrogen layers of the star until shell burning starts. The star begins at that moment its subgiant phase, and evolves on a nu-Article number, page 1 of 15 arXiv:2101.11025v1 [astro-ph.SR] 26 Jan 2021 A&A proofs: manuscript no. main clear timesale for masses below about 1.5 solar masses (M ). For stars that are close to the Terminal-Age Main Sequence (TAMS), the structure and the evolution remain highly influenced by the properties of the MS convective core. Interestingly, the star begins at that moment to exhibit mixed modes. These modes behave like gravity modes in the internal propagation cavity and pressure modes in the outer one. Thus, they allow us to finely probe the deepest layers of the star, all the while being observable. This, and the closeness of the subgiant to the TAMS, make the mixed modes of young subgiants valuable data to study the extension of convective cores.
Another particularity of mixed modes is their very fast evolution, compared to the nuclear evolution timescale of the subgiant phase. Indeed, mixed modes frequencies change dramatically over the course of a few million years. This makes their seismic modeling challenging. Recently, increasing efforts have been made to model subgiants Stokholm et al. 2019;Metcalfe et al. 2020;Deheuvels et al. 2020;Li et al. 2019Li et al. , 2020, driven by both their great physical interest and the sudden increase of seismic data for these stars. Most of those works focused on finding the optimal stellar parameters for one or several subgiants. So far, few studies have used subgiants as tools to test stellar physics, mainly due to the challenges of their modeling, as mentioned above. Deheuvels & Michel (2011) successfully constrained α ov from a subgiant observed by CoRoT, HD 49385, which exhibits only one g-dominated mode and is therefore very close to the TAMS. They found that either no overshooting, or a model with α ov = 0.19 were giving equally good results. In this work, we model a young subgiant, KIC10273246, observed by Kepler during almost 1000 days. That star exhibits two g-dominated modes, which allows us to better constrain its inner structure. We perform a thorough seismic modeling of the star, in order to precisely estimate its stellar parameters and to place constraints on the extension of its MS convective core.
In Section 2, we show the interest of having access to two g-dominated mixed modes in young subgiants. In Sect. 3, we present the surface observables of KIC10273246, and perform a fresh analysis of its oscillation spectrum using the full Kepler data set. We then describe in Sect. 4 the modeling technique that we adopted, which is an improved version of the method developed by Deheuvels & Michel (2011). Sect. 5 presents our optimal stellar models and the constraints that were obtained on the extent of the MS convective core for KIC10273246. We discuss these results in Sect. 6 and Sect. 7 is dedicated to conclusions.

Probing potential of mixed modes
Just after the main sequence, the oscillation spectra of solarlike pulsators show the presence of mixed modes, which are due to the coupling between the observed p-modes and low radialorder g-modes (n g = 1, 2, 3, n g being the number of nodes in the g-mode cavity). The frequency spacing between low-order g-modes is large (several times the large separation of p modes), so that only a few are in the observable frequency window during the subgiant phase. Moreoever, n g being low, the pure g modes that couple to p modes do not follow an asymptotic behavior (as described in Shibahashi 1979, Tassoul 1980. The oscillation spectra of subgiants therefore constrast with those of more evolved stars, which typically have more g-dominated modes than p-dominated modes, and for which n g is of the order of several tens (e.g. Mosser et al. 2012).
The frequency of mixed modes are mostly determined by two structural features of the star: gφ Hp ∇ µ Propagation zone Evanescent zone Fig. 1: Typical propagation diagram of a low-mass subgiant star. The Brunt-Väisälä frequency is represented in blue, delimiting the g-mode cavity (light blue). The Lamb frequency, in orange, delimits the p-mode cavity (light orange). Two gdominated mixed modes angular frequencies, with n g = 1, 2, are represented (solid lines in propagation zones, dotted lines in evanescent zones). The G cavity turning points are noted r i1 , r i2 and r o1 , r o2 . Finally, the thermal and chemical contributions to the Brunt-Väisälä frequency are represented in green dashed and red dot-dashed, respectively.
1. The g-mode (G) cavity, delimited by the Brunt-Väisälä frequency N. 2. The evanescent zone between the g-mode and p-mode (P) cavities, the latter being delimited by the Lamb frequency S l .
The G cavity affects the frequency of the g-mode that is involved in the mixed mode frequency. G-modes frequencies, in the asymptotic theory, can be approximated by l being the degree of the mode, r 1 and r 2 the turning points in the G cavity, and r the local radius of the star. In our case, n g is low for the observed modes, so that the asymptotic expression given in Eq. 1 should not apply. However, it has been shown that it can provide qualitative information about the behavior of the mixed modes frequencies . It tells us that g-dominated modes should give strong constraints on the Brunt-Väisälä frequency in the G cavity. One can write it in the following form (e.g. Kippenhahn et al. 2012): where g is the local gravity, δ = −(∂ ln ρ/∂ ln T ) P,µ , φ = (∂ ln ρ/∂ ln µ) P,T , ∇ ad = (∂ ln T/∂ ln P) ad , ∇ = ∂ ln T/∂ ln P and ∇ µ = ∂ ln µ/∂ ln P. The Brunt-Väisälä frequency consequently carries information about both the thermal (two first terms in the parenthesis) and compositional structure (last term) of the star. The evanescent zone affects the coupling between the two cavities, whose strength is linked to the size of this region and the value of N inside it (e.g. Unno et al. 1989). Using a toy model, Deheuvels & Michel (2011) showed that a strong coupling induces a shift of the l ≥ 1 p-dominated frequencies that are close to a g-mode. The p-dominated frequencies therefore provide complementary information about the internal structure of the subgiant.
In this work, we investigated whether having several gdominated modes in the observed oscillation spectrum can bring more information about the extension of the MS core. From above, we know that the frequencies of the g-dominated mixed modes are related to the N/r integral between the turning points of the G cavity. Fig. 1 shows the propagation diagram of a subgiant close to the TAMS, highlighting the frequencies of the two first g-dominated l = 1 mixed modes, that is those that arise due to the coupling of p modes with g modes or radial orders n g = 1 and 2. We denote their turning points in the G cavity as r i1 , r o1 for the mode with n g = 1, and r i2 , r o2 for the mode with n g = 2. The difference between the two frequencies is thus mainly related to the Brunt-Väisälä frequency value between r o1 and r o2 (as the one in the [r i1 , r i2 ] region is negligible). This region, as it can be seen in Fig. 1, is dominated by the µ-gradient contribution. This gradient is related to the characteristics of the hydrogen-burning shell, especially the nuclear energy generation and thus its temperature and composition. It has been shown that H-burning shell structure depends on the MS core features, and especially on the amount of core overshooting. One can see this in Fig. 5 of , which exhibits two Brunt-Väisälä profiles of stars having the same evolutionary stage and position in the HR diagram, but computed with different α ov . The two profiles differ mainly in the peak caused by the µ-gradient, and these structural differences are large enough to have a significant impact on the eigenfrequencies of the star.
For all those reasons, stars with two visible g-dominated modes are therefore expected to be interesting targets to place constraints on the efficiency of core overshooting.
That criterion led us to the choice of KIC10273246, a subgiant which has two g-dominated modes and a mass of 1.26 ± 0.10 M (Creevey et al. 2012) which places it safely in the mass range where stars are expected to have a convective core during the MS.

Observational properties of KIC 10273246
3.1. Surface constraints 3.1.1. Constraints from spectroscopy Surface constraints have been obtained for KIC10273246 by Creevey et al. (2012). The authors used different algorithms on the same spectra obtained with the FIES spectrograph. For our target, they found effective temperatures (T eff ) ranging from 5933 to 6380 K. A weighted mean gives us a value of 6150±100 K, which we have used to constrain our seismic modeling. The star was also found to have to have a sub-solar metallicity, with [Fe/H] = −0.13 ± 0.1 dex.

Constraints from broadband photometry
To obtain a reliable value of the luminosity of the star, we performed a spectral energy distribution (SED) fit, following the procedure of Stassun & Torres (2016).
We extracted photometric values using the VizieR photometry tool. Those data come from the NUV filter from GALEX (Martin et al. 2005), the B T and V T filters from Tycho-2 (Høg et al. 2000), the J,H and K s filters from 2MASS (Skrutskie et al. 2006), the gri filters from SDSS (Skrutskie et al. 2006), the W1-W4 filters from WISE (Wright et al. 2010) and the G magnitude from Gaia (Evans et al. 2018). The atmosphere model comes from the Kurucz atmosphere grid (Kurucz 2005), with the surface gravity (log g) derived from the seismic relations (see 3.2.4), and the metallicity coming from spectroscopic measurements. We then fitted the photometry points to the spectrum, with the T eff and extinction A v as free parameters. We also used the spectroscopic data from Creevey et al. (2012) and the extinction from Green et al. (2019) as priors. With a reduced χ 2 of 0.7, we found T eff = 6000 ± 33 K, and A v = 0.00 +0.037 −0.000 mag. The fitted spectrum and the photometric data are represented in Fig. 2. Finally, we integrated the flux over all the wavelengths and used the distance coming from Gaia to obtain the luminosity of the star. According to Zinn et al. (2019), a parallax bias exists in the Kepler field, which depends on the G-band magnitude and the pseudo-color ν eff (effective wavenumber of the photon flux distribution in the Gaia band) of the star. We found − Gaia = 39.15 ± 9.46 µas, which gives L = 5.74 ± 0.17 L . This result is as expected lower than the Gaia archive value (5.92 ± 0.13 L ) due to the parallax offset.

Preparation of Kepler light curve
KIC10273246 was observed with Kepler between quarters Q0 and Q11 (total duration of 978 days) in short cadence (58.85 s). An early seismic analysis of the target was performed by Campante et al. (2011) using the first four quarters of Kepler observations (325 days of data). They estimated the frequencies of oscillation modes of degrees l = 0, 1, and 2 over eight overtones. We revisited this analysis using the complete Kepler data set. The light curve of the star was processed using the Kepler pipeline developed by Jenkins et al. (2010). Corrections from outliers, occasional drifts and jumps were performed following the method of García et al. (2011). The power density spectrum (PSD) was then obtained by applying the Lomb-Scargle periodogram (Lomb 1976, Scargle 1982. The PSD is shown in the shape of an échelle diagram in Fig. 3. We recall that the échelle diagram is built by dividing Fig. 3: Échelle diagram of KIC10273246, folded with ∆ν = 48.2 µHz. For clarity, the spectrum has been smoothed over a 0.2-µHz boxcar. The white symbols indicate the frequencies that have been extracted for modes of degree l = 0 (crosses), l = 1 (diamonds), and l = 2 (triangles) in Sect. 3.2.2. The two l = 1 g-dominated mixed modes have been circled in red.
the PSD in consecutive chunks with a length corresponding to the large separation of acoustic modes ∆ν, and piling them up. Here, we used the estimate of ∆ν = 48.2 µHz obtained by Campante et al. (2011). The main interest of échelle diagrams is that acoustic modes of same degree align in nearly straight ridges, which eases mode identification. The neighboring l = 0 and l = 2 ridges are readily identified on the left part of the échelle diagram (modes indicated by crosses and triangles, respectively, in Fig. 3).
The ridge of l = 1 modes (indicated by diamonds in Fig. 3) deviates from the nearly-vertical line that is expected for purely acoustic modes. This behavior is known to arise for dipolar modes in the presence of avoided crossings between low-order g modes and p modes in subgiants. Each avoided crossing is characterized by the presence of an additional mode, which lies away from the ridge of theoretical purely acoustic l = 1 modes (which would be a nearly-vertical line at an abscissa of about 35 µHz in Fig. 3). This mode is most strongly trapped in the core and is thus g-dominated. The neighboring l = 1 modes are p-dominated, but their frequencies are nevertheless affected by the presence of the g-dominated mode. The modes with frequencies larger than the g-dominated mode are shifted to higher frequencies (to the right in the échelle diagram) and those with frequencies below the gdominated mode are shifted to lower frequencies (to the left in the échelle diagram). These features are clearly visible in Fig.  3, corresponding to two l = 1 avoided crossings. The l = 1 g-dominated modes associated to these avoided crossings have been circled in red in Fig. 3.

Extraction of oscillation mode parameters
To extract the parameters of the oscillation modes, we followed the method of Appourchaux et al. (2008). We here briefly recall the main steps of the procedure and refer the reader to the latter paper for more details.
Prior to fitting the individual oscillation modes, we modeled the background of the PSD. The contribution from granulation was modeled by two Harvey-like profiles, following the pre-scription of Karoff et al. (2013), and we added a white noise component to account for photon noise. The overall contribution from the oscillations was modeled as a Gaussian function. We fitted this model to the PSD using maximum-likelihood estimation (MLE). The central frequency of the Gaussian function gives an estimate of the frequency of maximum power of the oscillations ν max . To determine the error on this quantity, we subdivided the Kepler light curve in ten chunks of equal duration, and fitted the background model on the PSD calculated with these time series. The error on ν max was taken as the standard deviation of the measurements of this quantity for each chunk. We thus obtained ν max = 843 ± 20 µHz. The PSD was then divided by the optimal background model.
We then performed a fit of the oscillation modes, which were modeled as Lorentzian functions to account for their finite lifetimes. Each mode profile of degree l, radial order n, and azimuthal order m was characterized by its central frequency ν n,l,m , its height H n,l,m and its line width Γ n,l,m . Since dipolar modes have a mixed character, it cannot be assumed that they share similar heights and line widths with the neighboring radial modes, as is often done for main sequence solar-like pulsators. Most quadrupolar modes are expected to be p-dominated, owing to the weak coupling between the P and G cavities for these modes. We therefore assumed that the l = 2 modes have the same heights and widths as their closest l = 0 modes, with the exception of one g-dominated l = 2 mode, which is discussed below and in Sect. 3.2.3.
Non-radial modes are split into multiplets by rotation. Owing to the slow rotation of subgiants, the effects of rotation on the mode frequencies can be found by applying a first-order perturbation analysis. The components of a rotational multiplet are thus expected to be equally spaced by the rotational splittings δν n,l . We also assumed that they share similar line widths, and that their height ratios depend only on the inclination angle i of the star following the expressions given by Gizon & Solanki (2003). In principle, mixed modes can have different rotational splittings because they probe the rotation at different depths in the star. This has been used to probe the internal rotation of subgiants (e.g., Deheuvels et al. 2014).
To test whether individual rotational splittings can be measured in KIC10273246, we first performed local fits of the nonradial modes. Following the method described by Deheuvels et al. (2015), we fitted each mode using two different models: one under the H 0 hypothesis (no rotation, so that each mode is modeled as a single Lorentzian), and one under the H 1 hypothesis (rotation is considered and each mode is modeled as a set of 2l+1 Lorentzians separated by the rotational splitting). It is clear that hypothesis H 1 necessarily provides better fits to the data than hypothesis H 0 since it involves two additional free parameters (inclination angle and rotational splitting). The significance of hypothesis H 1 can be tested using the likelihoods 0 and 1 of the best fits obtained under the H 0 and H 1 hypotheses, respectively. As shown by Wilks (1938), the quantity ∆Λ ≡ 2(ln 1 −ln 0 ) follows the distribution of a χ 2 with ∆n degrees of freedom, where ∆n is the difference between the number of free parameters involved in hypotheses H 1 and H 0 (here, ∆n = 2) 1 . For each multiplet, we thus obtained a value of ∆Λ. The false-alarm probability was then given by the p-value p = P(χ 2 (2 dof) ∆Λ), which corresponds to the probability that a mode under the null hypothesis can produce such a high value of ∆Λ. Fig. 4: Oscillation spectrum of KIC102773246 in the vicinity of a quadrupolar mode that was found to be significantly split by rotation (see Sect. 3.2.2). The thick red curve corresponds to our best-fit model of the spectrum. Two quadrupolar mixed modes are visible (around 779.4 µHz and 783.9 µHz) and one radial mode (around 785.6 µHz).
For dipolar modes, the lowest p-value that we found is 0.08, which is too high to consider the measurement as significant. This means that we cannot reliably extract individual rotational splittings for dipolar modes in this star. The most likely explanation is that the modes have large line widths compared to the rotational splitting. For quadrupolar modes, only one mode (the one with a frequency around 779.4 µHz) was found to have a low p-value, of about 4 × 10 −5 , which shows a very high significance level. A rotational splitting of 0.53 ± 0.03 µHz was obtained for this mode (see Fig. 4). This mode is in fact a g-dominated mixed mode, as we show in Sect. 3.2.3.
We then performed a global fit of the modes (all the modes are fitted simultaneously). Since individual splittings cannot be measured, we assumed a common rotational splitting for all l = 1 and l = 2 modes (except for the aforementioned l = 2 mode around 779.4 µHz). Since most non-radial modes are pdominated, we expect the common rotational splitting to essentially measure the rotation in the envelope. The best fit corresponds to a rotational splitting of δν = 0.45 ± 0.02 µHz for nonradial modes and an inclination angle of i = 55±6 • . As was done for local fits, we also performed an additional fit of the modes without including the effects of rotation (null hypothesis). We could therefore estimate the p-value corresponding to the measurement of a mean rotational splitting. We found p ∼ 10 −4 , which indicates a high level of significance. Our results are compatible with the estimates of Campante et al. (2011), who had found i 20 • for this star, and optimal values of the rotational splitting slightly below 0.5 µHz.
The best-fit parameters for the oscillation modes (frequencies, heights, and line widths) are given in Table A.1. The uncertainties of the fitted dipolar mode frequencies range from 0.08 to 0.50 µHz. The measured mode frequencies are in quite good agreement with the ones found by Campante et al. (2011). Discrepancies at the level of 3 σ were found for only two modes (the dipole mode around 1055 µHz and the quadrupole mode around 880 µHz). Using the complete Kepler data set enabled us to detect l = 0 and l = 2 modes over three additional radial overtones compared to Campante et al. (2011). Our results are also in very good agreement with the recent measurements of mode frequencies for KIC10273246 by Li et al. (2020) using the com-plete Kepler data set (agreement at the level of 2 σ or better for all oscillation modes).

Detection of an l = 2 mixed mode
We mentioned above that the l = 2 mode with a frequency of about 779.4 µHz is the only mode for which an individual rotational splitting could be measured. This mode also has other distinctive features. It is separated from the closest radial mode by 6.1 ± 0.2 µHz. By comparison, for the other radial orders, the average separation between the l = 2 mode and the neighboring l = 0 mode is 4.4 µHz, with a standard deviation of 0.4 µHz. This suggests that this mode might be an l = 2 mixed mode, whose frequency is modified by the coupling between the p-and g-mode cavities. This hypothesis is strengthened by the fact that it has a short line width (0.26 ± 0.08 µHz) compared to the width of the neighboring l = 2 modes (between 1.7 and 2.4 µHz). Indeed, if the mode under study is a g-dominated mixed mode, it should have a higher inertia than p-dominated l = 2 modes, and therefore a shorter line width. Figure 4 shows the profile of the radial mode that is closest to the l = 2 mode under study. There appears to be an additional mode in the left wing of the radial mode, at a frequency of about 783.9 µHz. To determine the significance of this mode, we performed local fits assuming either its presence (H 1 hypothesis) or absence (H 0 hypothesis). We found a p-value of 0.01, indicating a reasonable significance level. This also supports the identification of the l = 2 mode at 779.4 µHz as a mixed mode. In this case, the additional mode at 783.9 µHz would also be an l = 2 mixed mode undergoing an avoided crossing with its close neighbor. As is shown in Sect. 5, the best-fit models for KIC10273246 do show a pair of mixed modes in the vicinity of these two modes, which confirms our identification.

First estimates of stellar parameters using seismic scaling relations
To obtain first estimates of the global stellar parameters of the star, we used seismic scaling relations, which relate the global seismic parameters ∆ν and ν max to stellar properties such as the mass, radius and surface gravity (Brown et al. 1991). These relations could be derived because ν max scales to good approximation as the acoustic cut-off frequency (Brown et al. 1991;Stello et al. 2008;Belkacem et al. 2011).
To estimate the asymptotic large separation of acoustic modes, we followed the prescription of Mosser et al. (2013). We fitted an expression of the type to the observed radial modes, where ∆ν obs is the observed large separation around ν max , α measures the curvature corresponding the to the second-order term in the asymptotic development, ε p is an offset, and n max = ν max /∆ν obs . We thus obtained ∆ν obs = 48.47 ± 0.02 µHz, which translates into an asymptotic large separation of ∆ν as = 50.63 ± 0.02 µHz, following Mosser et al. (2013).
Using our estimates of ∆ν as , ν max from Sect. 3.2.2, and T eff from Sect. 3.1, we could apply seismic scaling relations to derive preliminary estimates of the star's mass, radius and surface gravity. We obtained M = 1.24 ± 0.12 M , R = 2.10 ± 0.07 R , and log g = 3.88 ± 0.03.

Physics of the models
We used MESA v10108 (Paxton et al. 2015) evolution models, with OPAL equation of states and opacity tables , with the solar mixture from Asplund et al. (2009). The models have been computed with an Eddington-gray atmosphere. The convection regions are treated using the standard mixing length theory (MLT) as prescribed in Cox & Giuli (1968), with a free parameter α conv corresponding to the ratio between the mixing length and the pressure scale height. Microscopic diffusion has been taken into account, unless otherwise specified, using the Burgers formalism (Burgers 1969) and diffusion coefficients from Stanton & Murillo (2016). However, radiative accelerations have not been included in the computed models, as the increase in computational time could not be afforded in this study. The impact of those processes are discussed in Sect. 6.2.
As raised by Gabriel et al. (2014), for stars that have a growing convective core, it is necessary to use the Ledoux criterion to determine the radius of the convective core R cc . This way, we avoid the creation of unphysical convective zones outside the core in strongly chemically stratified regions, which may have an impact on the composition profile of the star and thus on its evolution. Moreover, we used the predictive mixing scheme (Paxton et al. 2018).
Core overshooting was modeled as a step extension of the convective core, over a distance where d ov is the distance of instant mixing overshooting, H p the pressure scale height, and α ov a free parameter quantifying the phenomenon. Eq. 4 replaces the traditional expression d ov = α ov H p in order to prevent d ov from becoming unphysically large when the core is small (H p → ∞ when r → 0). It is important to note that this definition varies from one evolution code to another (see, e.g., Eq. 1 of Deheuvels & Michel 2011 for Ce-sam2K). Low-mass stars have small convective cores, therefore those differences must be kept in mind when comparing models coming from different codes. Additionally, the impact on our results of using a diffusive overshooting, as proposed by Freytag et al. (1996), is discussed in Sect. 6.1. The adiabatic oscillations of the models have been computed using ADIPLS (Christensen-Dalsgaard 2008), and the surface effects have been corrected for using the cubic term of the prescription of Ball & Gizon (2014).

Why modeling subgiants is challenging
The frequencies of g-dominated mixed modes evolve on a very short time, with a non-linear change of several µHz per million years, which is much larger than the usual uncertainties coming from Kepler data. As this timescale is several orders of magnitude shorter than the typical nuclear evolution time of low-mass subgiants, reproducing the mixed modes with a traditional grid technique requires extremely fine steps in mass and age. This makes this method prohibitive when the number of free parameters is large, as is required to test the model physics. Interpolation in age is possible (Li et al. 2020), but somehow difficult for l = 2 g-dominated modes, which we observed in KIC10273246. Interpolation across tracks (as used e.g. in AIMS, Rendle et al. 2019) could mitigate the need for extremely fine steps in mass, iso-∆ν iso-ν g Fig. 5: HR-diagram representing the evolution tracks of stellar models with masses varying from 1.2 M (lightest gray) to 1.3 M (darkest gray) and otherwise identical physics. Each evolution is stopped when ν g is correctly reproduced. but needs to be tested for subgiants, especially regarding the extreme sensitivity of the g-dominated frequencies to the masses of the models. Additionally, an "on-the-fly" optimization technique may perform badly due to the highly non-linear behavior of the mixed modes, especially during the computation of the derivatives in the case of a gradient-descent kind of algorithm.
To overcome those difficulties, a new approach is necessary. We thus developed a nested optimization, where we optimize the physical parameters of models (e.g. metallicity, initial helium abundance...) that have been beforehand optimized in mass and age. This way, we can handle those two sensitive parameters using a dedicated and specific procedure, separately from the other ones for which a more traditional technique is possible. This modeling method originates from Deheuvels & Michel (2011) and has been adapted to make it more robust. We recall in the following the basic principles of this method and highlight the differences with the one used in the present study.

Optimization in mass and age
In that part of the optimization process, we compute models with only two free parameters, the mass and the age of the star, the rest being fixed. The optimization of those two parameters can be made easier thanks to the fact that, if all the other physical parameters (such as metallicity, mixing-length parameter...) are fixed, reproducing only ∆ν and the frequency ν g of a g-mode is enough to precisely constrain the mass and the age.
A physical justification of that approach can be found in Deheuvels & Michel (2011). We remind it here using an HRdiagram represented in Fig. 5. It shows the iso-∆ν line, as L ∝ T 5 eff for models with the same large separation, and the iso-ν g line, computed from stellar evolution models. The two lines meet at a unique point, that can be reached by tuning only the mass (i.e., choosing the "right" evolution path) and age (i.e., stopping at the right moment on that path).
In concrete terms, our first step is, at a given mass, to find the age that correctly reproduces the ν g frequency.
As we only see mixed modes and not pure g-modes, we cannot directly measure ν g . A possible solution would be to Here, like in our modeling of KIC10273246, δν is defined as ν 1,11 − ν 0,13 , and the observed value is represented by the dotted line. The plot has been strongly magnified in order to see the 1-σ uncertainties from the data.
choose a g-dominated mode (i.e., a non-radial mode far from its ridge) frequency. Unfortunately, such a frequency does not evolve monotonously with age, as it can be seen on the upper panel of Fig. 6. We then preferred to look at the distance between that g-dominated mode and its closest radial mode, which we call δν.
As we can see in the top panel of Fig. 6, this value always decreases with age, but also keeps the interesting properties of the mixed modes as it evolves very quickly during an avoided crossing, allowing us to tightly constrain the age.
This step would be equivalent, on Fig. 5, to follow an unique evolution path and stop it when it crosses the iso-ν g line.
We then optimize on those "good-age" models in order to correctly reproduce the large separation. In practice, to do this we minimize the χ 2 of only the radial modes, defined as We do not take into account the non-radial modes at this stage to get rid of the complexity of behavior of the mixed modes. This approach differs from the one followed by Deheuvels & Michel (2011), who at this stage searched for models that minimized the difference between the observed average large separation and the one coming from the models. By using here all the radial modes instead, we found that the optimization process is more robust to the correction of near-surface effects. It may be observed that the behavior of ∆ν (and so of the radial frequencies) is close to linear when varying the mass. Then, a simple Newton-type algorithm (such as a Levenberg-Marquard algorithm) is enough to quickly find the optimal mass. This step would then be equivalent, on Fig. 5, to the right evolution path that leads to the meeting points of the two iso-lines. Figure 7 shows the échelle diagram of a model that we can obtain after that first step, with arbitrary physical parameters: metallicity [Fe/H] = −0.2 dex, mixing-length parameter α conv = 1.5, initial helium abundance Y 0 = 0.28. We can see that the radial modes and δν = ν 1,11 − ν 0,13 (the proxy for ν g ) are, by construction, correctly reproduced. However, the other frequencies are far from the observed ones. Especially, the g-dominated mode ν 1,18 is about 10 µHz away from the observations. Thus, to find a better matching model, we adjust the other parameters.

Optimizing the other parameters
Now that we have a method to correctly find the mass and the age of a star at a given physics, we must find the other parameters of the stars, taking into account this time all the observational constraints. Thus, we define a new χ 2 as where N obs is the total number of observational constraints, both seismic and non-seismic, x obs i , x mod i the values of those observed constraints or their computed equivalent, and σ i their respective uncertainties. We also introduced the quantities ∆ i ≡ (x obs i − x mod i ) 2 /σ 2 i , which indicate the contributions of every observable to the χ 2 , for later purposes.
As those parameters have a lower impact on the frequencies than the mass and age, it is now possible to use more traditional approaches. One possibility is to compute grids of models, where each model of the grid is optimized in mass and age. Another option is to perform an optimization using an iterative method, where again each iteration consists in an optimization of the mass and age. To model KIC10273246, we opted for a hybrid method, which is described in the following section.

Fitting procedure adopted for KIC10273246
For the modeling of KIC10273246, we left the initial metallicity [Z/X] 0 , the mixing-length parameter α conv , the initial helium abundance Y 0 and, of course, the overshoot parameter α ov as free parameters. At first, to explore the global behavior of the χ 2 , we computed a very loose grid ([Fe/H] between -0.2 and 0.2, step 0.1; Y 0 between 0.24 and 0.28, step 0.02; α conv between 1.5 and 2.1, step 0.2 and α ov between 0.0 and 0.2, step 0.05). We recall that each model of this grid is optimized in mass and age as explained in Sec. 4.3. The purpose of this loose grid was to investigate whether there exist double solutions or local minima. No such features have been found. Moreover, those grids allowed to identify the region of the best parameters.
We thereafter refined those parameters. As mentioned in Sec. 4.4, optimizing [Z/X], α conv , Y 0 , α ov can be performed either through a grid approach or an iterative procedure. We therefore made several tests, using stellar models as mock observations, to determine which method is preferable. We found the best robustness when following a hybrid approach: for given values of Y 0 (0.26 through 0.31, step 0.01) and α ov (0 through 0.25, step 0.05), we conducted iterative optimizations with the Levenberg-Marquardt algorithm to find optimal values of [Fe/H] and α conv . This method differs from the one used in Deheuvels & Michel (2011) where a single grid was used for all the free parameters.
Among those models, we considered only the models that were compatible with the observational estimates of the chemical enrichment law ∆Y 0 /∆Z 0 . Typical values quoted for ∆Y 0 /∆Z 0 range from 1.4 to 4 (e.g. Chiappini et al. 2002, Balser 2006, Casagrande et al. 2006. We consequently had a conservative approach and took into account all models having ∆Y 0 /∆Z 0 < 5.

Results
In this section, we describe the general characteristics of the best models, before commenting the constraints on α ov . We finally investigate at the internal structures of the best models and the constraints brought by the mixed modes.

General characteristics of the best models
Following the method described in the Sec. 4.4.1, we obtained optimal models for each value of α ov , whose caracteristics are in Table 1. The best model, with α ov = 0.15, has a reduced χ 2 of 3.2. The échelle diagram of this model is represented in Fig. 8. Also, surface observables are consistent with the ones found in the literature, or with the SED fitting formerly described: we find a less than 1-σ difference for the effective temperature T eff , the metallicity [Fe/H] and the luminosity L. We found a radius and a mass that are consistent with the seismic scaling relations as well. We can note that, as expected from the mass-based prediction, all good models had a convective core during the MS. This supports our choice of using this star to constrain α ov .
We note that our best-fit models are significantly less massive and older than those of Li et al. (2020), who also performed a seismic modeling of KIC20173246 and found M = 1.49 ± 0.08 M and an age of 2.84 ± 0.60 Gyr. These discrepancies could be partially explained by the different assumptions made on the input physics. For instance, Li et al. (2020) considered a solar-calibrated mixing length parameter, while we left this parameter free in our study. Also, contrary to us, Li et al. (2020) neglected element diffusion and adopted the mixture of Grevesse & Sauval (1998). Finally, we stress that the agreement with the observed dipole mode frequencies, in particular for the g-dominated mode ν 1,18 , is significantly better in the present study than it is for the best-fit models of Li et al. (2020) (compare Fig. 10 of Li et al. (2020) to Fig. 8 of the present paper). These mismatches between models and observations for dipole modes are acknowledged by Li et al. (2020), and the authors attribute them to an imprecise modeling of the core structure. For each combination of (α ov , Y 0 ), our minimization using the LM algorithm can be used to derive uncertainties in the stellar parameters. The error bars in the free parameters of the fit ([Fe/H] and α conv ) are obtained as the diagonal coefficients of the inverse of the Hessian matrix. The uncertainties on the other parameters can then be obtained using Eq. 10 of Deheuvels et al. (2016). We thus obtained very small error bars, of the order of 0.007 for [Fe/H] and 0.002 for α conv , which translates into uncertainties of approximately 0.004 M for the stellar mass and 0.04 Gyr for the age. This means that for a given combination of (α ov , Y 0 ), the available observables provide very strong constraints on the stellar parameters. By comparison, we found that optimal models with different Y 0 can yield similar agreement with the observations (statistically equivalent χ 2 ) but have quite different stellar parameters. This degeneracy of stellar models with respect to Y 0 is addressed in more detail in Sect. 6.3. It thus seems that the uncertainties in the stellar parameters are dominated by the model degeneracy in Y 0 . We thus used the optimal models with different Y 0 to estimate uncertainties in the stellar parameters. For a given α ov , we fitted a second order polynomial to the χ 2 curve and retrieved the interval of values corresponding to χ 2 min + 1. This gave us the 1-σ uncertainty, which are reported in Table 1. Figure 9 shows the variation in the χ 2 of the optimal models with α ov . We can see that adding overshooting allows us to reproduce significantly better the observed frequencies, the χ 2 of the models without overshoot and with α ov = 0.15 being 315 and 127, respectively.

Constraints on core overshooting
To better understand what frequencies allow us to favor models with overshoot, we investigated the contributions of the different observables to the χ 2 (∆ i in Eq. 6). We denote as ∆ nov i and ∆ ov i the χ 2 contributions of the observables coming from the  optimal models without overshoot and with α ov = 0.15, respectively. Figure 10 represents the differences ∆ nov i − ∆ ov i , positive values meaning that the observable is better reproduced by the α ov = 0.15 model. As expected, we can see that the main χ 2 difference is due to the dipolar modes, which have a mixed behavior. However, we observe that the g-dominated modes (indicated by the dotted vertical lines) hardly contribute to distinguishing the models with and without overshooting. Either type of models fit well the g-dominated frequencies. The main contributors to the χ 2 differences are in fact the l = 1 p-dominated modes in the neighborhood of the g-dominated modes. As explained in Sect. 2, the frequencies of these modes are mainly influenced by the coupling between the P and the G cavities. The intensity of that coupling thus accounts for the main part of the differences between the models. We note that all those models correctly reproduce ν 1,18 , as the high sensitivity of this g-dominated mode strongly constrains the region of parameters of the models with the smallest χ 2 . The major role played by the dipolar modes in the constraints on α ov is also illustrated in Fig. 9, where the colored regions indicate the contributions to the χ 2 of the surface observables and modes depending on their degree.
Moreover, Fig. 9 indicates that the contribution to the χ 2 of the l = 2 modes hardly changes with α ov . This was partly expected, because their evanescent zone is larger than that of the dipole modes, making the coupling between the G and P cavities weaker. Most of the detectable modes are therefore very strongly p-dominated and do not constrain the deep structure of the star, hence α ov . Yet, one g-dominated l = 2 mode has been detected (ν 2,10 , see Sec. 3.2.3). It is interesting to see that, in a similar way to the previous paragraph, its frequency is equally well reproduced by models with and without overshooting. One can see this in Fig. 10, where ∆ nov i − ∆ ov i of that mode is less than 1-σ. On the other hand, the (2,11) mode, whose frequency is close enough to ν 2,10 to be influenced by the coupling, varies substantially with α ov . Figure 10 shows a 3-σ difference, despite the high 0.65 µHz observational uncertainty, confirming the key role of the coupling in the constraint on α ov . Interestingly though, whereas the α ov = 0.15 model reproduces better the dipolar modes, the (2,11) mode is better fitted in the model without overshooting. However, its large observational uncertainty prevents it from being too constraining.
Finally, we notice that adding a larger amount of overshooting (α ov > 0.15) strongly worsens the quality of the fit, placing a strong maximum limit to the value of α ov . To better understand this behavior, we investigate in the next section the seismic constraints on the internal structure of the models.

Contraints on central density
A tight link is expected between the g-dominated frequencies and the central density ρ c . This comes from the relation between the g-mode frequencies and the Brunt-Väisälä frequency, which approximately scales as ρ c (see eq. 15 from Deheuvels & Michel 2011). To verify this, we investigated the constraints placed on ρ c by the frequency of a g-dominated dipolar mode. For this purpose, we considered the values of ρ c in the models computed in the loose grids defined in Sec. 4.4.1, in which models are optimized in mass and age in order to reproduce ∆ν and the frequency of a g-dominated mode (here ν 1,11 , as described in Sec. 4.3). We observed that, despite the wide range of parameters considered, the standard deviation of ρ c among those models is as low as 32.4 g.cm −3 , which represents around 1% of the mean central density ρ c = 2100 g.cm −3 . This standard deviation even decreases to 10 g.cm −3 if we only keep the 200 best models, illustrating the tight relation between ρ c and the frequency of a g-dominated mixed mode.
This plays a role in the increase of χ 2 for α ov > 0.15. To illustrate this point, we computed models with the same values of [Z/X], Y 0 and α conv , but with different values of α ov . Each of these models was optimized in mass and age, as described above. Fig. 11 shows the evolution of ρ c with age for those models. One can see that they all reach approximately the same value of central density, ρ c , in accordance with the previous paragraph. Moreover, the intensity of the ρ c jump that is due to the post-MS core contraction increases with α ov . We can explain this by the fact that for bigger cores, the layers where hydrogen remains are more distant from the center. They are colder and require a stronger core contraction, hence a stronger jump in ρ c , to reach the fusion temperature. Therefore, when α ov increases, models that reach ρ c are closer to the TAMS. When the model gets too close to the TAMS, this impacts the whole stellar structure. Internally, the µ-gradient in the core is very much affected, because the nuclear reactions in the H-burning shell did not have enough time to smooth its shape. In the outer layers, the convective envelope, which expands during the subgiant phase, has a different size. Those processes alter the frequencies (the g-and p-dominated ones, respectively), which are thereby not compatible with the observations.

Constraints on the Brunt-Väisälä profile.
Based on Eq. 1, we expect the frequency of mixed modes to place constraints on the integral of the Brunt-Väisälä frequency in the G cavity. To investigate this, we plotted in Fig. 12 the N 2 profiles for the models of the loose grid defined in Sect. 4.4.1, which all reproduce the low frequency g-dominated mode (ν 1,11 ) and ∆ν. We observe that reproducing both already strongly constrains the part of the Brunt-Väisälä frequency dominated by the thermal term (∇ ad − ∇), which corresponds to the most central layers (r < 0.05 R ). This was expected because of the 1/r factor in the Eq. 1 integral. On the contrary, the part of N 2 that is dominated by the µ-gradient changes significantly within the grid.
We expect that part to be strongly determined by the dipolar modes. We therefore investigated the constraints brought by the two most determining seismic features (see Sect. 2): the coupling intensity and the frequency of pure g-modes. As those two are not directly measurable, we used observational proxies. The intensity of the coupling can be quantified by δ ≡ ν 1,12 − ν 1,11 , coming from Deheuvels & Michel (2011). That value is the difference between the low frequency g-dominated mode (ν 1,11 ) and the following dipolar p-dominated mode (ν 1,12 ). Thus, δ increases with the coupling. The frequency of a pure g-mode is measured through ν 1,18 , which is the high-frequency g-dominated mode. For those two values we color-coded, in Fig. 12, the profiles of the Brunt-Väisälä and Lamb frequencies based on their agreement with the observations (left panel for δ and right panel for ν 1,18 ).
One can see on the right panel that models correctly reproducing the coupling (i.e. dark profiles) have very similar Hburning shell position (N 2 peak around r = 0.05 R ). However, the Brunt-Väisälä profiles become more degenerate for higher r: several different profiles can reproduce δ within 1-σ. This degeneracy is lifted thanks to the high-frequency g-dominated mode: on the left panel, models closely reproducing ν 1,18 all have similar Brunt-Väisälä profiles. This corroborates our theoretical approach in Sect. 2, the high-frequency g-dominated mode adding tight constraints on the shape of the µ-gradient. Important gains in structural constraints are therefore obtained from having a second detectable g-dominated mode.

Diffusive overshooting
During this study, overshooting has been modeled as a step extension of the convective core. However, Freytag et al. (1996) proposed another prescription, based on results coming from 2D simulations of A-star and white dwarfs. Overshoot is then modeled as a diffusive process, with a coefficient D ov exponentially decaying from the boundary of core, following the expression with D conv the MLT derived coefficient taken just below R cc and f ov a free parameter that tunes the length scale of the overshooting region. In order to compare the results coming from the two types of overshooting, we first found the f ov that is equivalent to a given value of α ov . In order to do this, we searched for the value of f ov that gives models reaching the TAMS at the same age as models computed with a step overshooting and α ov = 0.15. We found f ov = 0.01. After that, we modeled it using a method similar to the one used in Sect. 5 and compared the best model with the one computed with a step overshoot.
As we can see in Fig. 13, the differences between the frequencies and observables of the best models with step and diffusive overshoot are mainly less than 1-σ. We note the exception of the g-dominated ν 2,10 mode, which is better reproduced by the diffusive overshoot model. However, its impact on the global χ 2 is counter-balanced by the generally better reproduced dipolar frequencies of the step overshoot model. Moreover, the difference between the characteristics of the two models are within the uncertainties of Table 1. Therefore, we cannot discriminate between the two kinds of modelings with the current set of data. Fig. 14: Brunt-Väisälä profiles of the best models with gravitational settling and chemical diffusion (blue) and without those two processes (orange).

Effect of microscopic diffusion
The models presented in Sect. 5 of our study include both gravitational settling and chemical diffusion. Such processes, which happen to be necessary in order to correctly reproduce the helioseismic observations (Christensen-Dalsgaard et al. 1993), are expected to have an impact on our results for two main reasons. The first is the sinking during the main sequence of heavier elements because of the gravitational settling. This reduces the hydrogen abundance in the core and shortens the main sequence, which eventually decreases the age of models with same mean density. The second is the smoothing of the structure of the sub-  Table 2: Characteristics of the two models of the Fig. 15 giant. High µ-gradient peaks, like the one produced by the withdrawal of the convective core, are strongly smoothed out by the chemical diffusion, which impacts the mixed mode frequencies.
Thus, it is interesting to see how gravitational settling and chemical diffusion change the characteristics and the quality of the fit. We therefore modeled the star following a similar methodology, but without including those two processes. We found that the best model also has α ov = 0.15, but provides a significantly worse agreement with the observations than the best diffusion model, with χ 2 nodiff − χ 2 diff = 71. It is more massive (1.29 M ) and older (4.85 Gyr), as expected from the fact that heavy elements do not sink during the MS. We note that the surface observables are less well reproduced, with a best model that is too cold (5874 K, which represents 1.84 times the observational uncertainty) and that has a too low luminosity (5.0 L , which represents 4.3 σ).
Moreover, similarly to what we found in Sect. 5, the quality of the fit improves as α ov increases for α ov ≤ 0.15. However, this is less significant (χ 2 α ov =0 − χ 2 α ov =0.15 = 24). For higher values, the quality of the fit strongly worsens, in a comparable way to what has been found with gravitational settling and chemical diffusion. Figure 14 illustrates the differences in the Brunt-Väisälä profiles between the two best models, with and without both diffusive processes. One can see the high ∇ µ peak (at r = 0.06 R ) is as expected smoothed out by the chemical diffusion. Otherwise, the two profiles are remarkably similar, despite the different physics of the models, which highlights the robustness of the constraints coming from the mixed modes.
Finally, the effects of gravitational settling are expected to be somewhat counter-balanced in the envelope by radiative accelerations, which can have a significant impact on stars with masses greater than 1.1 M (Deal et al. 2018). However, including the effects of radiative accelerations in the models increases the computational time of stellar evolution calculation by orders of magnitude, and it could not be afforded in the present study. To test the impact of this process on our modeling, we computed a model which takes into account radiative accelerations, with the same parameters as the best model for α ov = 0.15. We obtained slightly different large separations (∆ν rad − ∆ν norad = 0.12 µHz) but very similar frequencies, once normalized by the large separation. Radiative accelerations are therefore not expected to change the conclusions of this study.

Helium degeneracy
To model KIC10273246, we initially performed optimizations with fixed α ov and considering Y 0 , α conv and [Fe/H] as free parameters. In this case, we observed an unwanted sensitivity of the Y 0 parameter to the guess value of our optimization process. This We found that optimal models with different values of Y 0 have indeed surprisingly close frequencies, despite their wide range of mass. This is illustrated in Fig. 15, which shows the échelle diagrams of the best model with Y 0 = 0.28 (blue) and the best model with Y 0 = 0.24 (orange), both having α ov = 0.15. Those models have quite different characteristics, as reported in Table 2. However, their frequencies are almost indistinguishable, despite the very small uncertainties on the mode frequencies coming from Kepler data. Only the g-dominated l = 2 mode allows us to slighlty favor the Y 0 = 0.28 model. Such degeneracy is related to the anti-correlation between mass and Y 0 , that has been observed in MS stars (see e.g. Lebreton & Goupil 2014) as well as subgiant stars (Li et al. 2020). Additionally, we note that no monotonic behavior has been found between the age and Y 0 . We therefore conclude that the seismic modeling of subgiants, despite bringing strong constraints on the deep structure of the star, does not lift the degeneracy between Y 0 and the mass.

Internal rotation
We mentioned in Sect. 3.2.2 that a rotational splitting of 0.53 ± 0.03 µHz could be measured with a high level of significance for the l = 2 mode at 779.4 µHz. This is obviously not enough to probe the internal rotation in detail. However, since this mode is g-dominated, it can be used to place approximate constraints on the rotation in the core of the star. Using our best-fit model from Sect. 5, we could compute the rotational kernel K(r), which relates the splitting δν s of this mode to the rotation profile Ω(r) This can be re-written as δν s = K g Ω g + K p Ω p , where Ω g and Ω p are average rotation rates in the g-and p-mode cavities, respectively, and K g (resp. K p ) correspond to the integral of K(r) in the g-mode (resp. p-mode) cavities. For the l = 2 mode under study, we found using our best-fit stellar model that 84% of the kernel energy is enclosed in the g-mode cavity, which confirms that the mode is indeed g-dominated.
Campante et al. (2011) found a clear rotational modulation in the Kepler light curve of KIC10273246. They thus estimated the surface rotation rate of the star to about 0.5 µHz (rotation period of about 23 days). This value is comparable to the average rotational splitting of 0.45 ± 0.02 µHz that we obtained in this study. As mentioned in Sect. 3.2.2, this average splitting is dominated by the contribution of p-dominated modes, so that it essentially measures the envelope rotation rate. Taken together, these two measurements suggest a low rotation contrast within the p-mode cavity, in line with the conclusions of Benomar et al. (2015), Nielsen et al. (2015) for main-sequence solar-like pulsators.
The splitting measured for the l = 2 g-dominated mode is close to the rotation rate inferred for the envelope, which suggests a low amount of radial differential rotation in the star. If we take Ω p /(2π) ≈ 0.45 µHz, we obtain a core rotation rate of about 0.65 µHz (rotation period of about 18 days). Clearly, more rotational splittings would be required to precisely measure the core rotation rate. However, our results indicate that KIC10273246 could be rotating nearly as a solid-body, like the two Kepler subgiants whose internal rotation profiles have been recently measured (Deheuvels et al. 2020).

Conclusion
In this study, we performed a seismic modeling of KIC10273246, a subgiant observed by Kepler, and obtained strong constraints on the size of its main sequence (MS) convective core. We chose this target because it exhibits two dipolar g-dominated modes, which we expected to bring stronger constraints on the internal structure. We extracted the mode parameters from the oscillation spectrum of the star using the full Kepler data set and thus updated the mode frequencies that were previously obtained by Campante et al. (2011).
The seismic modeling of subgiants is notoriously complex. We here elaborated on the algorithm proposed by Deheuvels & Michel (2011). This method consists in a two-step approach. The first step consists in finding the mass and age that match the large separation of p modes and the frequency of one g-dominated mixed mode. The second step optimizes the other free parameters ([Fe/H], Y 0 , α conv and α ov ) to reproduce the other frequencies as closely as possible. In this work, we improved this algorithm to make it more robust. This enabled us to perform a detailed seismic modeling of KIC10273246, with a particular emphasis on the investigation of the size of the MS convective core.
We found models in good agreement with the observations, with a reduced χ 2 of 3.2 for the best model, and with surface observables that are reproduced to within less than 1 σ. One key result of this study is that models with core overshooting during the MS reproduce significantly better the observations, with an optimal value of α ov = 0.15. For higher values of α ov , the quality of the fit strongly worsens. We found that such models are very close to the terminal-age main sequence (TAMS). Their internal structure thus differs from that of the lower-α ov solutions, and their seismic properties show strong mismatch with the observations. We tested the robustness of our conclusions by considering other choices for the input physics. No significant difference was found when modeling core overshooting as a diffusive process. Models computed without microscopic diffusion also favor models with α ov = 0.15, albeit less significantly, and show strong mismatch with the observations for higher values of α ov . However, they yield poorer agreement with the seismic and surface observables compared to the models computed with microscopic diffusion. This study thus confirms the high potential of young subgiants with mixed modes to measure the extent of the MS convective cores.
We also investigated the information conveyed by the mixed modes about the core structure. We showed that the combined knowledge of the large separation ∆ν and the frequency of one g-dominated mixed mode is enough to estimate the central density ρ c to a precision of about 1%. This helps us understand why models with a larger amount of core overshooting (α ov > 0.15) are not compatible with the observations. Because of their larger MS convective core, they have a higher ρ c just after the end of the MS and they thus reach the optimal central density closer to the TAMS. We then studied the roles of the different mixed mode frequencies in determining the profile of the Brunt-Väisälä frequency inside the star. While the first g-dominated mixed mode strongly constrains the thermal part, the second one helps constrain the part of the Brunt-Väisälä frequency that is dominated by the µ-gradient. We therefore confirm that having access to two g-dominated mixed modes helps better characterize the Brunt-Väisälä profile.
Also, despite the strong constraints that were obtained on the internal structure, we noted the existence of a degeneracy between the stellar mass and the initial helium abundance Y 0 . This degeneracy, which is already well known for MS stars (e.g., Lebreton & Goupil 2014), is not lifted by the mixed modes. We found that it is in fact the main source of uncertainties in the determination of the stellar parameters. This should be kept in mind when modeling subgiants. Current modeling techniques, such as traditional grid-based methods, tend to miss a significant fraction of the best-fit models because of the size of the mesh. This could lead to explore only partially the degeneracy between Y 0 and mass, and thus to underestimate the uncertainties on the stellar parameters.
As a byproduct of this work, we obtained partial constraints on the internal rotation of KIC10273246. We could not measure individual rotational splittings for the dipolar mixed modes, but obtained a splitting of 0.53 ± 0.03µHz for the only g-dominated l = 2 mixed mode in the spectrum of the star. Interestingly, this value is close to the surface rotation rate of 0.5µHz that was found for this star by Campante et al. (2011) using photometric data from Kepler. This suggests that this star might be rotating nearly as a solid-body, similarly to the two young subgiants recently studied by Deheuvels et al. (2020).
This work highlights the large potential of the seismic modeling of young subgiants to indirectly obtain constraints on the core structure of the star during the MS. The next step will be to use this method on a larger sample of stars drawn from the targets observed with Kepler and TESS, and therefore place quantitative constraints on the overshooting process in low-mass stars. The data coming from the upcoming PLATO mission (Rauer et al. 2014) will add a large amount of potential targets for this type of analysis.
Moreover, we have shown in this study that detecting several g-dominated dipole modes places stronger constraints on the shape of the Brunt-Väisälä profile, and therefore on the µgradient in the stellar core. It could thus be anticipated that more evolved subgiants, which show a larger number of g-dominated mixed modes, would be more favorable targets for our purpose. However, these stars are also further from the end of the MS and it is to be feared that the chemical composition in the core might be less dependent on the properties of the MS core. We plan therefore to study this effect as well in a subsequent study.