© ESO, 2016

1. Introduction

In a previous paper (Buldgen et al. 2016), we studied the binaries 16CygA and 16CygB using the full Kepler dataset from Davies et al. (2015). The system is in fact more complex since a red dwarf orbits the A component and a Jovian planet orbits the B component (Cochran et al. 1997; Holman et al. 1997; Hauser & Marcy 1999). We carried out a forward modelling process of both stars without taking into account binarity as a constraint and used our inversion techniques to further constrain their fundamental parameters, and demonstrated the importance of microscopic diffusion. The inversion technique provided strong constraints on the chemical composition and mixing of 16CygA, the brightest of the two components. However, when carrying out the same inversion for 16CygB, we faced the problem of the amplification of the observational error bars. The problem is well-known in the context of inversions, since the results are always a trade-off between amplifying the errors and fitting the target function of the inversion (Pijpers & Thompson 1994). In the context of asteroseismology, since more weight has to be given to the fit of the target function due to the small number of observed frequencies compared to the solar case, we are always limited in terms of error amplification. Trying to reduce the error bars by amplifying the trade-off parameters can result in a significant reduction of the quality of the fit, thus implying that what is gained by reducing the propagation of observational error bars is lost due to the poor quality of the averaging kernel.

In the following sections, we re-analyse the trade-off problem of 16CygB and show that the seismic information is sufficient to analyse this star independently with the t_u indicator. To explain the trend seen with the inversion, we try unsuccessfully to restore the agreement by modifying the surface chemical composition of this star. Since this leads to inconsistencies with the 16CygA results of our previous paper, we analyse the potential necessity of an additional mixing process, which has already been mentioned to explain the lithium depletion in this star (Deal et al. 2015). We emphasize that the solution we propose for consistency with the inversion result is hypothetical and is subject to the same limitations and model-dependencies as our previous study on 16CygA. We compute models using a parametrised approach of the extra mixing which should not be considered as a physical solution but rather a hint that a certain amount of mixing is required in deep regions of the B component in order to reconcile the modelling of both components.

The paper is structured as follows, we start by briefly presenting additional reference models in Sect. 2. We then present our inversion results as well as the regularisation in Sect. 3. These results are further analysed and discussed in Sect. 4 in light of the possible necessity for extra mixing in 16CygB. We then conclude with the implications and perspectives of this study in Sect. 5.

2. Reference models

In this section, we will describe the forward modelling process that has been carried out to obtain the reference models for the inversion. The process has been already described in Buldgen et al. (2016), but we recall it here for the sake of clarity. Nevertheless, the number of models computed has been increased to improve the diagnostic process of the inversion and to ensure unbiased results.

In practice, we computed these models independently from the modelling of 16CygA presented in our previous paper. We used the frequency spectrum from Davies et al. (2015), which was based on 928 days of Kepler data¹. A Levenberg-Marquardt algorithm was used to determine the optimal set of free parameters for our models. We used the CLES stellar evolution code and the LOSC oscillation code (Scuflaire et al. 2008b,a). The stellar models used the CEFF equation of state (Christensen-Dalsgaard & Daeppen 1992), the OPAL opacities from Iglesias & Rogers (1996) supplemented at low temperatures by the opacities of Ferguson et al. (2005). The nuclear reaction rates are those from the NACRE project (Angulo et al. 1999), including the updated reaction rate for the ¹⁴N(p,γ)¹⁵O reaction from Formicola et al. (2004) and convection was implemented using the classical, local mixing-length theory (Böhm-Vitense 1958). We also used the implementation of microscopic diffusion from Thoul et al. (1994), for which three groups of elements are considered and treated separately: hydrogen, helium and the metals (all considered to have the diffusion speed of ⁵⁶Fe). No additional transport mechanism, beside microscopic diffusion, was included in the models. No surface correction of the individual frequencies was used in this study since we used quantities that are naturally less sensitive to these effects.

Moreover, since the inversion results for 16CygA implied that microscopic diffusion had to be included in the stellar models and since both stars are very similar, we considered that we had to include atomic diffusion in the models of 16CygB. We also emphasize that obtaining consistent results in age for both components is impossible if one considers that one component of the binary system is subject to microscopic diffusion effects while the other is not. Yet, we also want to stress that the implementation of microscopic diffusion has its own uncertainties. First, we consider here the implementation from Thoul et al. (1994) which considers only three components to the mixing; secondly, in their own paper, Thoul et al. (1994) consider the diffusion velocities obtained to be accurate within approximately 15%; thirdly, it may be possible that radiative accelerations play a role in competing with gravitational settling effects. Thus, the use of microscopic diffusion as a solution to be consistent with the inversion results for 16CygA is a first hypothesis of this study. It does not mean that another combination of mixing processes could not successfully reproduce the trends previously seen with the inversion technique for this star.

In this study, we substantially increased the number of reference models used to carry out the inversions for 16CygB but did not use any hypothesis on the chemical composition of this star. In fact, surface chemical composition differences between the A and B components have been claimed by Tucci Maia et al. (2014) when carrying out a differential spectroscopy analysis between both stars. Moreover, although the centroid of the present surface helium abundance, Y_f, interval found by Verma et al. (2014) is the same, the scatter is larger for the B component, and if microscopic diffusion is included in the stellar models, it is also clear that surface chemical composition differences will be seen since this mixing will not have the same efficiency for stars of different masses².

Nevertheless, it should be noted that chemical composition differences between 16CygA and 16CygB are still under some debate since their existence has been claimed by Ramírez et al. (2009) and Tucci Maia et al. (2014) as well as by previous studies (see Deliyannis et al. 2000) but could not be confirmed by Schuler et al. (2011). In Tucci Maia et al. (2014), one finds [Fe/H]_A = 0.101 ± 0.008 and [Fe/H]_B = 0.054 ± 0.008 whereas Schuler et al. (2011) finds [Fe/H]_A = 0.07 ± 0.05 and [Fe/H]_B = 0.05 ± 0.05. These results are not totally incompatible, and what is more striking is the difference in error bars between various studies.

Moreover, these values depend on the reference solar metallicity assumed in the study since the observational constraint provided is the [Fe/H] value which must be translated in a value using the Sun as a reference. In our previous paper, we used the most recent abundance tables given by AGSS09 (Asplund et al. 2009) and found that they led to a better agreement with the inversion results for 16CygA. In this study, we computed most models with the AGSS09 abundances but also used some models with the older GN93 abundances (Grevesse & Noels 1993). We explain our motivations for using such models in Sect. 4.

We summarise the observational constraints used for 16CygB in Table 1 and the fundamental parameters obtained for some of the reference models in Table 2. In this table, we also recall the intervals from the forward modelling process of 16CygA obtained previously. The forward modelling was carried out starting from various initial conditions with the Levenberg-Marquardt algorithm. The set-up of the minimization process was the following:

• Constraints: individual small frequency separationsd_0,2d_1,3, inverted mean density for which conservative error bars of 0.005 g/cm³ were considered, acoustic radius (τ) for which conservative error bars of 30 s were considered, present surface metallicity (Z_f/X_f) from Ramírez et al. (2009), present surface helium abundance (Y_f) from Verma et al. (2014) and the effective temperature from Tucci Maia et al. (2014), for which we considered error bars of 30 K.

• Free parameters: Mass, age, initial hydrogen abundance (X₀), initial abundance of heavy elements (Z₀), mixing-length parameter (α_MLT).

In total, we had 5 free parameters for 31 constraints. In addition to these constraints, we checked the values of the luminosity L, surface gravity log g and radius R after the forward modelling to see if they were consistent with the constraints from the literature. Models which were completely inconsistent with these additional constraints were disregarded. An additional comment should be made on some error bars used in the forward modelling. Firstly, we considered the errors from Tucci Maia et al. (2014) to be unrealistic and assumed a conservative 30 K error bar which is already very accurate but more consistent with other studies. Secondly, both the inverted mean density and acoustic radius are known to have underestimated error bars with the SOLA method, from the multiple hare and hounds we performed to calibrate the inversion techniques, we noticed that error bars of 0.5% were to be expected as a conservative error bar for the inverted values of the mean density. For the acoustic radius, the precision has to be assessed from the dispersion of the inverted values, in this particular case this lead to a precision of around 0.7%. Consequently, we used these conservative error bars in the Levenberg-Marquardt algorithm rather than the error bars derived directly from the SOLA method.

We can see that the scatter of fundamental parameters is very similar to that obtained for 16CygA. However, we only give the results for models including diffusion in Table 2, as can be seen by looking at the values of the D parameter. This parameter is related to the implementation of diffusion we use, it is a multiplicative factor of the microscopic diffusion velocities such that if D = 1.0, one uses the diffusion velocities of standard solar models. We can see that some models have radii and luminosities that are below the observed values. Thus, these models can already be rejected or at least questioned in terms of quality. The age and chemical composition intervals are completely consistent with the values obtained for the reference models of 16CygA recalled in the third column of Table 2. We recall here that the models associated with ages above 7.4 Gy were rejected for 16CygA, based on the t_u inversion results and their implications on microscopic diffusion and chemical composition. A successful modelling of the binary system implies finding similar ages and initial chemical composition for both stars as well as being consistent with the seismic, spectroscopic and interferometric constraints at hand. Ultimately, the models shall also be compatible with the inversion results. This is not an easy task and requires a careful analysis and a good trade-off between all of the constraints.

3. Inversion results

In this section, we present updated inversion results for 16CygB. In our initial work, we faced the problem of large error bars for the t_u inversion. These error bars implied that we could not derive any additional constraints on the structure of 16CygB. In fact the inversion results showed that all models should be accepted, regardless of whether they included diffusion or not. However, we will show in the following sections that a more careful look at the frequency data can lead to an independent diagnostic with the inversion and provide additional interesting insights on the structure of this star.

The inversion technique we present is based on the linear integral equations presented in Gough & Thompson (1991) derived for the squared isothermal sound speed and the helium mass fraction, Y. The basic equation of the inversion is then written: (1)where the notation stands for the relative difference between observed quantities and quantities of the reference model, defined as follows: (2)The most striking difference between inversions in asteroseismology and inversions in helioseismology is the number of observed frequencies, leading to the fact that the classical linear kernel based inversion methods cannot be used to derive full structural profiles of observed stars. In previous studies, we have adapted the SOLA inversion techniques from Pijpers & Thompson (1994) to carry out inversions of structural integrated quantities (See Reese et al. 2012; Buldgen et al. 2015a,b, for various examples.). Amongst the indicators derived, we defined a core condition indicator in Buldgen et al. (2015a) as follows: (3)with f(r) = r(r−R)²exp(−7(r/R)²), the weight function used for this inversion with R the stellar radius and r the radial coordinate associated with each layer inside the model, u is the squared isothermal sound-speed previously defined.

First, we recall a few basic equations of seismic inversion techniques. It is important to remember that seismic diagnostics using classical inversion techniques involve individual relative frequency differences (defined as in Eq. (2)). In that sense, any inverted result is generated from a recombination of these frequency differences. When we use the linear SOLA technique (Pijpers & Thompson 1994), we build a linear combination of frequency differences. In the case of the t_u inversion, for example, we have: (4)with the c_i being the inversion coefficients, which are determined by finding the optimal value of the SOLA cost function for given trade-off parameters values. The SOLA cost function is defined as follows for the t_u indicator and denoted : (5)where is the target function associated with the indicator, K_Avg is the averaging kernel, and K_Cross the cross-term kernel, defined with respect to the fractional radius position . η is a Lagrange multiplier, k is a regularization factor related to the non-linear generalization of indicator inversions (see Buldgen et al. 2015a, for details.), σ_i are the errors associated with each individual frequency and β and θ are the free parameters of the SOLA method, related to the trade-off with the cross-term and the amplification of observational errors and the accuracy of the fit of the target function. Nevertheless, for this particular inversion, no additional terms used to deal with surface effects have been added since they often bias the results and reduce the quality of the fit of the target function. This is also justified by the fact that the t_u indicator probes core regions and that its target function has low amplitude in the surface.

The averaging and cross-term kernels are defined as follows for the (u,Y) structural pair, with Y the helium mass fraction and , the squared isothermal sound speed and the functions and the structural kernels associated with u and Y respectively: The fact that we have two free parameters in the SOLA cost function is due to the ill-posed nature of the problem and leads to the well-known trade-off problem when using inversion techniques. In this particular case, the question of the trade-off is particularly important since we have three oscillation modes in particular that have larger error bars than the all the others and two of these could sometimes see their individual frequencies fitted within their error bars.

Another specificity of asteroseismic inversions is that they are performed with little or no knowledge of the radius of the observed target, noted R_tar. In Sect. 2.1 of Buldgen et al. (2015a), we analysed the impact of this problem on equations of the type of Eq. (1). It was then shown that the inversion implicitely scaled the observed target to the same radius as the reference model used to perform the inversion while keeping its mean density constant. This meant that the target studied by the inversion was not defined by a mass M_tar and a radius R_tar but was a scaled target defined by a mass and a radius R_ref.

This does not restrict the diagnostic potential of the inversion technique but means that if we want to compare results from various reference models, we need to compare values of to get rid of the implicit scaling process introducing a dependency in R_ref in the inversion process.

3.1. Analysis of the error contributions

In Fig. 1, we illustrate in orange the initial inversion results of with their quite large error bars, R_Tar being the target photospheric radius. They seemed disappointing since the kernel fits were excellent and implied that there were enough kernels to fit the target function of the t_u inversions.

Fig. 1 In orange, inversion results for the tu indicator and with the full set of modes for 16CygB. In green, inversion results for the same models excluding the modes with the largest error bars in the frequency set. The blue squares are associated with inversion results for which the trade-off parameter θ has been slightly enhanced. In red, blue and magenta, tu and values in the reference models (See text for the explanation of the colour code).

This implied that the problem was simply stemming from the observational errors propagation term in the cost function of the SOLA method. The classical way to deal with this problem is to increase the θ parameter in the cost function thus reducing the propagation of observational errors. While this may be a solution, changing the θ parameter can lead lead to a much less accurate fit of the target function and thus reduces the quality of the inversion. This implies larger errors on the inverted result coming from the kernel fit as shown in Buldgen et al. (2015a). From our previous test cases, we also know that around 50 frequencies is sufficient to obtain an inverted value for t_u, especially if octupole modes are available. Consequently, we looked at the observed frequencies for which there were large uncertainties and found that the ℓ = 3, n = 14 mode, the ℓ = 3, n = 15 mode and the ℓ = 3, n = 16 had much larger uncertainties than the other modes of similar radial order. The error bars on the individual frequencies were sometimes even larger than the frequency differences between 16CygB and the computed reference models with the Levenberg-Marquardt algorithm. This is of course somewhat inefficient since it implies that we are using frequency differences that cannot be exploited by the inversion techniques.

In fact, frequency differences with large error bars can dominate the error contribution in the inversion results, especially if the inversion coefficient associated with the particular mode is important. This is in fact simply due to the form of the term associated with the error propagation in the SOLA cost-function which is written: (8)It is thus clear that modes with high inversion coefficients and large uncertainties contribute the most to the error propagation. Although the SOLA method tends to mitigate the impact of the modes with large uncertainties, the result is always a compromise between precision and accuracy. This trade-off is realized through the change of the free parameters of the inversion. In the context of asteroseismic inversions, the fact that each oscillation spectra has its own error bars, that each star is fitted individually within a given accuracy that can be variable and that each star occupies a different position in the HR diagram for which the linear approximation might be irrelevant to a certain degree, makes each inversion process unique. Therefore, from the mathematical point of view, each inversion has to be analysed differently, although trends in terms of inversion parameters can be seen and are understandable since they are linked to the data and model quality which can be objectively assessed.

The trade-off problem of inversion techniques is illustrated by the so-called trade-off curves that can be seen in the original paper on the OLA method by Backus & Gilbert (1967) or Pijpers & Thompson (1994) for the SOLA method. Typically, each frequency set defines the number of coefficient available, thus the resolution of the inversion. However, this resolution is mitigated by the error bars of these individual modes which limit the amplitude of the coefficient that can be built to fit the target function. The trade-off curve materializes this competition with respect to the parameter θ of the inversion. We describe a little bit more in depth the trade-off problem and the effect of eliminating modes in the frequency spectrum in Sect. A.

Important error bars can indeed be seen for the ℓ = 3, n = 14 mode, which is the octupole mode of lowest radial order. We know indeed from our previous test cases (see Buldgen et al. 2015a) that the t_u inversion uses preferentially the low order modes and tends to benefit from the presence of octupole modes and use them as much as possible. Since this particular mode has the highest error bar, we wanted to see how eliminating it from the frequency set used for the inversion could help us obtain a smaller error propagation. As previously explained, inversion techniques use individual frequencies to extract information. However, this is only possible if the frequencies used by the inversions are not fitted within their observational error bars. Typically, if one eliminates a mode with large error bars, one reduces the amplification of the errors but also the resolution of the inversion. Ultimately, eliminating a mode from the frequency set is only justified if its detection is arguable or if it is already fitted within the error bars. Otherwise, reducing the error bars is more efficiently done by increasing slightly the value of the θ parameter.

In the particular case of 16CygB, some individual modes could be fitted within their error bars and thus could not bring any additional seismic constraints if used in an inversion process. Finally, eliminating the worst offenders in terms of error bars is a process that has also been described in helioseismic inversions (see Basu et al. 2009), since they can have strong impact on SOLA inversions when adjusting the trade-off parameters for the inversion.

Fig. 2 Left panel: averaging kernels for the core conditions indicator (tu) for various θ values and reduced frequency spectrum. Right panel: same figure for the cross-term kernels of the tu inversion. We used the (u,Y) structural pair in both plots.

In Fig. 1, we show the impact of the modes with the largest error bars on the final inversion error propagation and values of t_u and for reference models. The new result with reduced error bars are illustrated in green and light-blue. The green results are obtained by eliminating problematic modes and the light-blue results are those obtained by slightly increasing the value of the θ parameter. We associated the following colour code for the values inferred from the reference models: a blue cross means that the t_u value agrees will all inversion results, a magenta cross that it agrees with some inversion results and a red cross that it does not agree with any result. We can see that eliminating the modes with large uncertainties reduces significantly the error bars on the inverted result, without changing much the actual inverted value. A change in the inverted value would have meant that the mode had a significant impact on the inversion result. In practice, this could be seen by a change in the fit of the target function by the averaging kernel. This could be the case if one had fewer individual frequencies and that the problematic oscillation mode was used by the inversion despite its large error bars. In Fig. 2, we illustrate the change in the averaging and cross-term kernel fit that is induced by the elimination of the most problematic modes in terms of observational error bars and an increase of the trade-off parameter θ. As was the case for the inverted t_u values, the differences on the averaging kernels are minimal. Hence, an independent study of 16CygB in terms of t_u can be performed. In the next section, we present new inversions results using a greater number of models for different surface chemical compositions, yet within the observational constraints, and for different diffusion coefficients, in much the same way as what was done in our previous study, more focused on 16CygA.

3.2. t_u inversion for 16CygB

In this section, we present the results for the t_u inversions for 16CygB. Using the reference models computed with our Levenberg-Marquardt algorithm and the more regularized inversions, we were able to obtain a value of t_u for 16CygB with lower error bars. However, the uncertainties are still non-negligible. Thus, we have to combine our analysis with other diagnostics and carefully discuss our final results, as was done in our previous study of 16CygA. In Fig 3, we present our results for various models with various surface chemical composition and changes of the factor D associated with atomic diffusion in CLES. The results we obtain are slightly model-dependent, which is very similar to what was obtained for 16CygA, but the trend is in this particular case opposite to what was seen before. Indeed, in Buldgen et al. (2016), we saw that including microscopic diffusion provided much more consistent values of the t_u indicator when compared to the inverted values. For 16CygB, models with lower helium surface abundances, higher surface metallicities and less diffusion are favoured. In fact, reducing the t_u value is directly related to a reduction of the gradient of , with T the temperature and μ the mean molecular weight. Consequently, reducing t_u implies reducing the mean molecular weight gradient within the star or changing the temperature gradient in the regions where the t_u indicator is sensitive. Reducing the mean molecular weight gradient can first be done by eliminating microscopic diffusion in the models. Indeed, this process tends to accumulate heavy elements in the deeper regions since for stars around 1.0 M_⊙, gravitational settling dominate the transport mechanism in the deep radiative regions. However, as stated before, not including this process leads to inconsistent ages and chemical compositions for both stars. Therefore, the reason for this discrepancy has to be explained using a more subtle effect.

We show in Fig. 4 the differences in chemical composition and in the weight function involved in the integral expression for the t_u indicator for two of our reference models in the chemical composition box. Model₁ is a model with a higher helium content (Y_S = 0.26), lower metallicity ((Z/X)_S = 0.0208) and microscopic diffusion (D = 1.0), thus following the prescription derived from our previous study. Consequently, it is also less massive (M = 0.9 M_⊙) and within the “young” range of our reference models (Age = 7.32 Gy). Due to the higher helium content and efficient microscopic diffusion, this model is rejected by the t_u inversion. Model₂ has a low helium content (Y_S = 0.22) and a higher metallicity ((Z/X)_S = 0.0214) and a less efficient microscopic diffusion (D = 0.5). This model is significantly more massive than Model₁(M = 1.01 M_⊙) but has a quite similar age of 7.54 Gy. The strong difference in mass is due to the well known degeneracy associated with the helium abundance. It should be noted that this model is validated by the t_u inversion.

This illustrates the fact that simply changing the surface chemical composition or microscopic diffusion has a strong impact on the fundamental parameters of the star and implies strong changes in the internal structure even if the model fits all the observational constraints (although Model₁ should be rejected due to its lower radius). Both models were chosen because they were extreme cases and illustrated well the strong degeneracy due to helium abundance.

3.3. Comparison with 16CygA

If we consider again 16CygA, the models with masses around 1.01 M_⊙, high helium content and ages 7.2Gy around were considered to be the best models of this star since they reproduced the t_u trend seen in our previous paper. This would mean that we would chose a model closer to Model₁ to be consistent in terms of the initial chemical composition of both components. However, since in this case we have to reduce the t_u values, and thus apply the opposite changes to the chemical composition and microscopic diffusion, the 16CygB have higher masses and ages (like Model₂ mentionned above), going up to 1.03 M_⊙ and 8.0 Gy.

The fact that the inversion is able to distinguish between Model₁ and Model₂ proves again the diagnostic potential of this approach. In this particular case, due to the fact that both stars are within a binary system, we are even able to see whether our selected result will be consistent with the previously determined parameters for 16CygA. Due to the very similar chemical composition derived spectroscopically and seismically, due to the results of independent forward modelling of both components leading to similar ages and initial chemical composition, we rather consider that the differences seen with the inversion technique is to be explained by inaccuracies in the models rather than considering the binary system to have merged from two isolated stars.

Fig. 3 Left panel: mean density vs. core conditions indicator (tu) plot. The inversion results are plotted in orange with their respective error bars. The crosses are values for the reference models computed with the Levenberg-Marquardt algorithm with AGSS09, the black × shows one example of a model computed with GN93 and Yf = 0.25. Right panel: surface chemical composition box for 16CygB. The colour code used allows direct trend comparisons between the surface chemical composition and the tu values. The size of the symbols is related to the intensity of microscopic diffusion, the smaller the symbol, the smaller the D coefficient. The + and the ◇ illustrate the impact of the metallicity on the tu value.

We also illustrate in Fig. 3 the results for one model using the GN93 abundances and models which were computed using the AGSS09 abundances and assuming a similar initial chemical composition to what was derived for 16CygA in our previous study. These models show values of around 3.7 g²/ cm⁶ whereas the model with GN93 is more consistent with the inversion results of 3.0 ± 0.5 g²/ cm⁶. It is clear that models computed assuming the same ingredients as 16CygA are incompatible from the point of view of the inversion. However, since these stars form a binary system and thus are thought to have formed together, we should be able to derive similar values of the initial chemical abundance and similar ages for both components of the system. This problem is also reflected in the effective temperature and radii determination. The well-known helium mass degeneracy leads to smaller radii for models with higher helium abundances, for example. We also tried using larger error bars on the effective temperature and looked at models with T_eff between 5600 K and 5900 K to see if this could affect the results. Ultimately, no trend was found since they are ultimately related to the chemical abundances and the way the elements are mixed within the star. These effects are well-known to affect the position of the models in the HR diagram at the end of its evolution. Thus, in what follows, we will focus on these aspects to try to reconcile our models of 16CygB with those of 16CygA and the inverted results.

In terms of precision and accuracy, it should be noted that neither the model-dependency, nor the regularisation can be held responsible for an inaccurate result. Hence, as shown in this section, in particular thanks to the large variety of reference models, we can see that surface chemical composition changes are not sufficient to explain the inverted values of t_u. In fact, taking θ = 10^-4 still implies very similar inversion results with reduced error bars and a slightly worse fit to the target function. Moreover, we know from our previous numerous test cases that the t_u inversion provided accurate seismic diagnostic of core regions (See Buldgen et al. 2015a).

Fig. 4 Left panel: in blue, helium abundance profile (Y) for one model with a lower surface helium abundance, around 0.22. In red, Y profile for a model with a higher surface abundance, around 0.26. Right panel: the profile of the target function of the core conditions indicator (tu) is plotted in corresponding colours for both models.

Fig. 5 Left panel: plot showing the difference of the gradient of the natural logarithm of temperature (T) and that of the mean molecular weight (μ) for models including different mixing processes: the green curve is for a model with microscopic diffusion, the blue curve is for a model with a constant turbulent diffusion coefficient and the red curve is for an exponentially decaying turbulent diffusion coefficient. Right panel: the gradient of the natural logarithm of the mean molecular weight and of the temperature for the same models as in the left panel, the colour code has been respected.

3.4. Influence of physical parameters on t_u

When analysing the effects of microscopic diffusion, the problem is even worse, since if we trust the values of Y_f ∈ [0.240.25] for the final surface helium abundance of 16CygA, we should obtain higher Y_f values for its less-massive counterpart due to the fact that its convective envelope goes slightly deeper and implies less-efficient microscopic diffusion. One should note that similar conclusions can be drawn for the surface heavy element abundance of this star. In fact, increasing the amount of heavy elements in the stars increases the opacity in the deep radiative regions where the t_u indicator is sensitive (see Fig. 2). Thus, it implies an increase of the temperature gradient, . Now, since (see Eq. (3)), it is worth looking more in depth at the behaviour of this indicator with changes in the stellar structure. Using the ideal gaz approximation, we have a straightforward relation between u, T and μ(9)This formula implies that the behaviour of the indicator depends on the values of the gradients themselves. As can be seen in the right panel of Fig. 5, it is not always straightforward to say whether an increase of the mean molecular weight gradient through diffusion will imply an increase of t_u. For example, below 0.1 R, with R the stellar radius, it will be the case because diffusion will increase the depth of the minimum just below 0.1 R. However, adding extra mixing around 0.2 R or 0.3 R will smooth the transition towards the chemically homogeneous convective envelope (around 0.7 for this model) thereby decreasing the value of t_u. Similarly, increasing the temperature gradient below 0.1 will reduce t_u, and reducing above 0.1 (thus sharpening the transition towards the convective envelope) will imply the same reduction for the indicator. One can see these effects in Fig. 5 where we illustrated the impact of different types of mixing on the temperature and μ gradients and thus on the t_u indicator. This gives us a clue as to what could be modified in the models to reconcile the inversion results with the other constraints. However, it does not mean that this is the only solution to the problem we presented previously. For the sake of illustration, we also illustrate the hydrogen gradient and the Brunt-Väisälä frequency of these models, showing the change of the slope of the hydrogen gradient at the bottom of the convective zone but also a significant deplacement of the base of the convective zone for these models due to the use of the new opacities from the OPAS project. The OPAS opacities are new opacity tables specifically designed for solar-like conditions, where great care has been given to the details of the absorption lines considered. These models also used the latest version of the OPAL equation of state (Rogers & Nayfonov 2002). These changes of course affect the stratification below the convective zone and thus the behaviour of the t_u indicator. Turbulent diffusion also implied a change of the Brunt-Väisälä frequency in the very deep regions (below 0.1 R), this is particulary seen for the model associated with constant turbulent diffusion.

Fig. 6 Left panel: plot showing the hydrogen gradient of the same models as in Fig. 5. Left panel: plot showing the Brunt-Väisälä frequency of these models.

In practice, all thermodynamic quantities are coupled through the equation of state. For example, adding a mixing process will affect the chemical composition, thus the mean molecular weight, but it will also affect the opacity and indirectly the temperature gradient. Consequently, the t_u inversion offers a new insight on some differences between the target and the reference model, but does not provide the physical cause of the observed differences in structure.

4. Impact of physical ingredients on the core conditions indicator

4.1. Adding extra mixing

Because of the t_u inversion results, we are faced with a very peculiar problem. We have two stars, in a binary system, with very similar surface chemical composition, similar masses and radii, that show significantly different seismic behaviours when carrying out inversions of their structure. The problem is that the models for both stars cannot be consistent with the inversion results and simultaneously present similar chemical composition and age. Small discrepancies in chemical composition between both stars have proven not to be sufficient to eliminate the discrepancy with the inverted t_u values. Therefore, we had to assume that something was neglected in the models for 16CygB, or 16CygA, or for both stars. In what follows, we study supplementary models including a parametrized approach for an additional mixing process. The physical nature of this mixing process is not discussed here, but we demonstrate that the t_u indicator is, as expected, able to discriminate between various processes inside the star. Figure 7 shows various t_u inversion results for different implementations of diffusion yielding different chemical compositions. At first, we still wish to see whether there is a way reconcile the chemical composition of 16CygB with that of 16CygA.

The parametrization of this additional mixing is based on an implementation of turbulent diffusion used in previous studies (see Miglio et al. 2007, for details). We tested different implementations of this mixing. First, we added a constant turbulent diffusion coefficient of around 20 cm² s^-1 acting in the entire stellar structure and computed a few models fitting the observational constraints for 16CygB. The impact of the constant turbulent diffusion coefficient is quite strong. Indeed, gradients are quickly attenuated and the t_u value decreases, as can be seen in Fig. 7 with the positions of the blue ◇ in the left panel. However, disagreement with other constraints is quickly found if this mixing is further increased. For example, it is impossible to fit the individual small frequency separations when the extra mixing is too important although the acoustic radius, the mean density and other constraints of the cost function of the forward modelling can be accurately fitted.

We also computed models with a diffusion coefficient implemented as an exponential decay starting either from the bottom of the convective envelope or from the surface. Two parameters are used for this formalism, one multiplicative constant and the rate of exponential decay. From previous studies (Miglio et al. 2007), we know that a multiplicative coefficient of around 100 cm² s^-1 is consistent with the effects of rotation expected in solar-like stars. This value was used as a benchmark for the order of magnitude of the mixing, but we did not limit ourselves to this value since we wanted to investigate the effects of this parametric implementation on the t_u indicator. We thus allowed changes of up to ±50 cm² s^-1 in the value of this diffusion coefficient. From Fig. 7, where the models with the implementation of turbulent diffusion as an exponential decay starting from the surface are represented by blue +, we can see that it can indeed help to reconcile the models with the t_u values for 16CygB, even if a higher present surface helium value is considered, as had to be done for 16CygA. The fundamental parameters of these models are presented in Table 3, we note that they have slightly higher masses and ages than the models without turbulent diffusion for the same chemical composition.

Fig. 7 Left panel: mean density vs core conditions indicator (tu) plot. The inversion results are plotted in orange with their respective error bars. The × are values for the reference models computed with the Levenberg-Marquardt algorithm without extra mixing. The □ are related to models with a decaying exponential turbulent diffusion coefficient starting at the bottom of the convective enveloppe. The + show models with a decaying exponential diffusion coefficient starting from the surface and the ◇ use a constant turbulent diffusion coefficient. The ∗ depict models using the new OPAS opacities and the decaying exponential coefficient starting from the surface. Right panel: surface chemical composition box for 16CygB. The colour code allows direct trend comparisons between the surface chemical composition and the tu and values as in Fig. 3.

As expected, additional mixing can indeed help to reconcile the chemical compositions of both stars, but does not reconcile them in age since some of the models computed with the extra mixing have ages up to 8Gy even if most are still around 7.4−7.7 Gy. It is also noticeable that the masses of the models in the present study tend to be slightly higher than those previously found in Buldgen et al. (2016) using constraints from the 16CygA modelling. However, all of these models are still consistent with the radius, luminosity and log g constraints from the litterature. A clear trend is also seen in the fact that increasing the mixing improves the agreement between the reference models and the inversion. However, as the models come closer to the inverted values for t_u, they tend to be less consistent with the small frequency separation values, meaning that the extra mixing should not be too intense. Indeed, reducing the rate of exponential decay (thereby extending the effects of extra mixing to lower regions) or increasing directly the turbulent diffusion coefficient leads to the same disagreement with the small frequency separations. To better understand the problems here, we plot the effects of the extra mixing on both the metallicity and helium profiles in Fig. 8. We see that the main effect is to reduce a metallicity peak right under the convective region. The more reduced the peak is, the closer the t_u values to the inverted ones. But in the meantime, we also degrade the agreement with the small frequency separations. Changes are also seen for the helium profile right under the convection zone. During the fitting process with the Levenberg-Marquardt algorithm, this affects the initial helium abundance required to be within the constraints from Verma et al. (2014) and thus indirectly the hydrogen profile and thus the metallicity constraint.

Fig. 8 Left panel: helium abundance plot for three models with different implementations of turbulent diffusion. For the red curve, a constant mixing coefficient was applied throughout the structure and the evolution of the model. For the blue curve, we used an exponential decay starting from the base of the convective envelope of the model. For the green curve, we used an exponential decay starting from the surface of the model. For the magenta curve, we used a model including only microscopic diffusion. Right panel: the metallicity profiles of the three models described above.

4.2. Changing the opacities

To further investigate the problem, we computed additional models with the new OPAS opacities (Mondet et al. 2015; Le Pennec et al. 2015) including atomic diffusion and the implementation of turbulent diffusion using an exponential decaying function starting from the surface. We used our forward modelling approach to compute these models and analysed whether they could agree better with both the small frequency separations and the inverted values of t_u. As stated before, adding extra mixing could reduce the agreement with the small frequency separations if its intensity was too high. However, when using the new opacities, we were able to further increase the intensity of the extra mixing, and thus the agreement with the t_u inversion, without degrading the agreement with the small frequency separations. As such, they provide a partial help to the problem of fitting all the constraints, as can be seen from the position of the blue ∗ in Fig. 7, but do not solve completely the problem. We can see the influence of these increased opacities, slightly moving deeper the base of the convective envelope and changing the efficiency of the mixing right below the base of the convective zone. We recall here that these models were still selected with the Levenberg-Marquardt algorithm using the observational constraints of 16CygB. The fundamental parameters of these models computed with were very similar to those obtained previously, but these models tended to show a slightly lower luminosity around 1.18 L_⊙ and mass around 0.98 M_⊙ and are thus responsible for the “lower” part of the intervals given in Table 3.

It can also be seen that when using a turbulent diffusion coefficient decaying from the lower boundary of the convective region, the effect on t_u is slightly more efficient, as illustrated by the position of the blue □ in the diagram of Fig. 7. Nevertheless, we did not seek here to fine-tune the parameters in this study since we are using a parametric approach to the problem without any physical background.

At this stage, we can already conclude that reconciling both models in terms of chemical composition and age will also probably need to remodel 16CygA to analyse whether effects other than diffusion could not be held responsible for the trend in t_u previously observed. In that sense, looking at constraints from the lithium abundance (King et al. 1997) and combining these constraints in the modelling of both stars might change the derived fundamental parameters by a few percents.

5. Conclusion

In this paper, we updated our study of the 16Cyg binary system by focusing our attention on 16CygB. From a re-analysis of the data, we were able to extract information from the t_u inversion and analyse the impact of extra mixing on the t_u values and other classical seismic indicators. First, we illustrated and solved the problem associated with the propagation of observational errors for inversions in 16CygB by analysing the impact of trade-off parameters and the presence of modes which in some cases were useless for the inversion technique. Ultimately, this approach could be used in similar situations for other observed targets.

From the t_u inversion, we were able to expose a problem in the surface chemical composition of 16CygB when compared to its companion. We computed a new set of models for this star, varying the surface chemical composition and restricting the effect of diffusion. We then observed that when the models were consistent with the inversion results, they were systematically inconsistent with the surface chemical composition we obtained for 16CygA. Since changing the chemical composition was not the solution, we sought to implement an extra mixing process in the models of 16CygB and tried to analyse its potential impact on the t_u values. As intuitively guessed, an extra mixing in the form of turbulent diffusion was found to be able to reconcile the models both with the surface chemical composition of 16CygA and the inversion results. Furthermore, using the new OPAS opacity tables further improved the agreement with the inversion. One could argue that other implementations could be tested, such as extra mixing in the form of undershooting using the prescription of Zahn (1991) as was found in HD 52265 by Lebreton & Goupil (2012). However, as was described in Lebreton & Goupil (2012), this extra mixing would leave an oscillatory pattern in the rr₀₁ and rr₁₀ seismic indicators. Due to the quality of the seismic data of 16CygB, we were able compute these indicators and found no evidence for an oscillatory pattern but rather a decreasing trend with frequency that is well reproduced by models without undershooting.

To conclude, we can state that various physical processes could improve the agreement. For example, a change in opacity would further change the results of the forward modelling process and thus the stellar parameters obtained with the Levenberg-Marquardt algorithm, these new models could potentially be in agreement with the inversion of the t_u indicator. In this study, extra mixing in the form of turbulent diffusion was invoked to reduce the disagreement between the inversion and the models. However, we did not seek to provide a physical explanation for this mixing and while it helps reducing the disagreement, further studies need to be performed to completely solve the problem. One first point would be to re-analyse 16CygA in the scope of the impact of extra mixing. Indeed, we have shown here that turbulent diffusion can change the t_u values. It is also well-known that rotation induces such type of extra mixing and it is believed to be responsible for the destruction of lithium in stars. Therefore, a first step would be to perform a thorough study of the impact of extra mixing on lithium abundances and inversion results for 16CygA. The case of 16CygB should be re-analysed afterwards, since it is well-known that this star shows even lower lithium abundances and is believed to have triggered thermohaline diffusion by accreting planetary matter (See Deal et al. 2015). As such, combining spectroscopic and seismic constraints in this binary system may provide new insights on stellar modelling of solar-like stars.

Moreover, additional indicators obtained through inversions seem to be a promising way to analyse the boundaries of convective envelopes. Consequently, from the sensitivity of seismic inversions and the quality of additional constraints, we are convinced that a re-analysis of the 16Cyg binary system with new stellar models should shed new lights on extra mixing processes in stellar interiors.

Acknowledgments

G.B. is supported by the FNRS (“Fonds National de la Recherche Scientifique”) through a FRIA (“Fonds pour la Formation à la Recherche dans l’Industrie et l’Agriculture”) doctoral fellowship. This article made use of an adapted version of InversionKit, a software developed in the context of the HELAS and SPACEINN networks, funded by the European Commissions’s Sixth and Seventh Framework Programmes.

References

Appendix A: The trade-off problem for t_u inversions

Appendix A.1: Origin of the trade-off problem

As described in Sect. 3, the SOLA inversion technique we use to obtain values of the t_u indicator computes a linear combination of individual frequency differences. These coefficients are obtained through the minimization of the cost-function defined as Eq. (5). We recall this definition here to better analyse the different contributions: (A.1)The first integral is associated with the averaging kernel, denoted K_Avg, this term defines the accuracy of the inversion technique, the better the fit of the target function, here , the more accurate the inversion is.

The second integral is associated with the cross-term kernel, denoted K_Cross. The cross-term stems from the presence of a second integral in Eq. (1), here for example associated with Y. Since the inversion only wants to extract information from the variable u, the contribution associated with Y as here to be damped. The trade-off between the reduction of the cross-term and the fit of the target function is calibrated by the free parameter β. In the case of the t_u inversion, the use of Y means that cross-term kernels have naturally smaller amplitudes due to the intrinsic small amplitude of helium kernels, as can be seen by comparing the right and left panels of Fig. 2. From previous hare and hounds, we know that the cross-term contribution is much smaller then the errors from the averaging kernels and have a negligible effect on the inversion results.

The third term of Eq. (A.1) is associated with the propagation of observational errors. This term regulates the precision of the inversion technique by damping the coefficients associated with large error bars. Since large coefficients are required for the t_u inversion, this term has an important impact on the final outcome of the inversion and is at the centre of the trade-off problem we will discuss. The importance given to the observational error bars of individual oscillation mode is materialized by the free parameter θ. Ultimately, the SOLA method comes down to a trade-off between precision and accuracy. In practice, a large value of θ will imply small error bars, but also potentially a very bad fit of the target function and the reduction of its accuracy. On the opposite, a very small value of θ means that the target function is well-fitted, but the result cannot be trusted due to its large error bars.

The fourth term is associated with an additional regularization based on homologous relations. The proper justification of the value of the coefficient k can be found in Sect. 3.2 of Buldgen et al. (2015a) and additional examples can be found in Reese et al. (2012) and Buldgen et al. (2015b) for other indicators. Eta is thus no free parameter but a Lagrange multiplier.

Appendix A.2: Effects of θ variations and mode suppression

As we stated in the previous section, the SOLA inversion is a compromise between precision and accuracy. This compromise is materialized by what is called a trade-off curve. It presents the accuracy of the result, in the form of the fit of the target function in abscissa plotted against the observational error amplification in ordinate for different values of the θ parameter. An example of a trade-off curve for the full-set of observation is plotted in Fig. A.1.

Fig. A.1 Trade-off curves for the tu inversion using SOLA for the full dataset and for the dataset not using the three octupole modes with lowest radial order.

As we can see, this trade-off curve is an L-shaped curve and it is also quickly understood that the optimal value of the free parameter θ is found nearby the edge of the trade-off curve. This position is associated with the best compromise achievable between precision and accuracy given a set of observational data. For the particular case of 16CygB, we started with values of θ = 10^-5 and found out that values around 5 × 10^-5 were better in terms of compromise between precision and accuracy. This is indeed seen in the plot of the trade-off curve were we zoomed on regions associated with θ = 10^-5. The green line is the trade-off curve obtained with the full set of data while the blue curve is associated with the trade-off curve when the ℓ = 3 and n = 14,15 and 16 modes have been suppressed from the dataset. The green vertical line indicates the position on the green trade-off curve (full dataset) associated with θ = 5 × 10^-5 while the blue vertical line indicates the position on the blue trade-off curve (restricted dataset) for θ = 10^-5. We can see that both positions are very close to each other in terms of error bars and fit of the target function values. However, the fact that the blue curve is always above and to the right of the green curve means that the compromise achieved with the restricted dataset will always be sub-optimal when compared to the compromise achieved with the full dataset. One can also see that the changes in error bars are quite quick when reducing θ to lower values. For example, if one considers the initial value of θ = 10^-5, the error bars are 30% larger then at θ = 5 × 10^-5 (which is even more striking then the example given in Fig. 1).

To illustrate the reason why we tried to eliminate the modes associated with ℓ = 3 and n = 14,15 and 16, we plot in Fig. A.2 the individual relative frequency differences with increasing frequencies. It can be seen that these three modes are well fitted and have larger error bars, this is why suppressing them helped us find a better compromise for the inversion technique. However, as stated above, this compromise is still sub-optimal in the strict mathematical sense due to the positions of the trade-off curves with respect to each other.

Fig. A.2 Relative frequency differences plotted with respect to the observed frequencies.

Fig. A.3 Effects of modes suppression on the averaging kernel of the tu inversion, illustrating the decrease in quality of the target function (in green) fit for various frequency sets. The degradation of the kernel fit is associated with instabilities in the tu values given in Table A.1.

Further illustrations are provided in Fig. A.3 and Table A.1. In Fig. A.3, we illustrate the variation of the fit of the averaging kernel for various sets of observed frequencies. Each fit is also associated with a result in Table A.1. For the set of 47 frequencies and 39 frequencies, we suppressed modes with small error bars and low n, that are known to be used by the inversion. It can be seen that the degradation of the kernel is correlated with a reduction in accuracy and some instability of the inversion results. This is basically due to the fact that each time we change the dataset, we are on a different position on a different trade-off curve. Also, this does not mean that change θ will always be a solution, because at some point the seismic information will simply be insufficient to infer some diagnostic using the t_u inversion.

All Tables

All Figures

Fig. 1 In orange, inversion results for the tu indicator and with the full set of modes for 16CygB. In green, inversion results for the same models excluding the modes with the largest error bars in the frequency set. The blue squares are associated with inversion results for which the trade-off parameter θ has been slightly enhanced. In red, blue and magenta, tu and values in the reference models (See text for the explanation of the colour code).

In the text

	Fig. 2 Left panel: averaging kernels for the core conditions indicator (tu) for various θ values and reduced frequency spectrum. Right panel: same figure for the cross-term kernels of the tu inversion. We used the (u,Y) structural pair in both plots.
In the text

Fig. 3 Left panel: mean density vs. core conditions indicator (tu) plot. The inversion results are plotted in orange with their respective error bars. The crosses are values for the reference models computed with the Levenberg-Marquardt algorithm with AGSS09, the black × shows one example of a model computed with GN93 and Yf = 0.25. Right panel: surface chemical composition box for 16CygB. The colour code used allows direct trend comparisons between the surface chemical composition and the tu values. The size of the symbols is related to the intensity of microscopic diffusion, the smaller the symbol, the smaller the D coefficient. The + and the ◇ illustrate the impact of the metallicity on the tu value.

In the text

	Fig. 4 Left panel: in blue, helium abundance profile (Y) for one model with a lower surface helium abundance, around 0.22. In red, Y profile for a model with a higher surface abundance, around 0.26. Right panel: the profile of the target function of the core conditions indicator (tu) is plotted in corresponding colours for both models.
In the text

Fig. 5 Left panel: plot showing the difference of the gradient of the natural logarithm of temperature (T) and that of the mean molecular weight (μ) for models including different mixing processes: the green curve is for a model with microscopic diffusion, the blue curve is for a model with a constant turbulent diffusion coefficient and the red curve is for an exponentially decaying turbulent diffusion coefficient. Right panel: the gradient of the natural logarithm of the mean molecular weight and of the temperature for the same models as in the left panel, the colour code has been respected.

In the text

	Fig. 6 Left panel: plot showing the hydrogen gradient of the same models as in Fig. 5. Left panel: plot showing the Brunt-Väisälä frequency of these models.
In the text

Fig. 7 Left panel: mean density vs core conditions indicator (tu) plot. The inversion results are plotted in orange with their respective error bars. The × are values for the reference models computed with the Levenberg-Marquardt algorithm without extra mixing. The □ are related to models with a decaying exponential turbulent diffusion coefficient starting at the bottom of the convective enveloppe. The + show models with a decaying exponential diffusion coefficient starting from the surface and the ◇ use a constant turbulent diffusion coefficient. The ∗ depict models using the new OPAS opacities and the decaying exponential coefficient starting from the surface. Right panel: surface chemical composition box for 16CygB. The colour code allows direct trend comparisons between the surface chemical composition and the tu and values as in Fig. 3.

In the text

Fig. 8 Left panel: helium abundance plot for three models with different implementations of turbulent diffusion. For the red curve, a constant mixing coefficient was applied throughout the structure and the evolution of the model. For the blue curve, we used an exponential decay starting from the base of the convective envelope of the model. For the green curve, we used an exponential decay starting from the surface of the model. For the magenta curve, we used a model including only microscopic diffusion. Right panel: the metallicity profiles of the three models described above.

In the text

	Fig. A.1 Trade-off curves for the tu inversion using SOLA for the full dataset and for the dataset not using the three octupole modes with lowest radial order.
In the text

	Fig. A.2 Relative frequency differences plotted with respect to the observed frequencies.
In the text

	Fig. A.3 Effects of modes suppression on the averaging kernel of the tu inversion, illustrating the decrease in quality of the target function (in green) fit for various frequency sets. The degradation of the kernel fit is associated with instabilities in the tu values given in Table A.1.
In the text