Open Access
Issue
A&A
Volume 676, August 2023
Article Number A87
Number of page(s) 23
Section Planets and planetary systems
DOI https://doi.org/10.1051/0004-6361/202346604
Published online 11 August 2023

© The Authors 2023

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.

1 Introduction

The formation of exoplanets takes place in protoplanetary disks around young host stars, consisting of mainly hydrogen and helium gas, but also heavier elements in both solid and gaseous form. Their presence is a natural outcome of star formation (for a review, see Williams & Cieza 2011). In these disks, planet cores grow by accreting material from the disk. As this process takes place around the young host star, it is apparent that stellar evolution cannot be treated as taking place in an isolated system. While stellar irradiation is a common aspect considered in planet formation models as a form of stellar influence on the disk (e.g., Chiang & Goldreich 1997; Dullemond & Dominik 2004; Bitsch et al. 2015; Savvidou et al. 2020); the reverse impact of the surroundings on the star is not to be neglected.

The large fraction of stars hosting at least one planet naturally leads to the conclusion that planet formation is a ubiquitous phenomenon, further arguing for its consideration in the study of young stars. The protostar and its disk initially share the same chemical composition, having formed from collapsing molecular cloud material. With the evolution of the protoplanetary disk, possibly including planet formation, and the rapid dispersal of the disk at the end of its lifetime due to photoevaporation (for a review, see Alexander et al. 2014), the distribution of material in the disk may change. In turn, this may lead to the star changing its chemical composition as it accretes disk material. For example, neglecting the effects of planet formation, material accreted onto the star is initially enriched in millimeter- to centimeter-sized pebbles, because they experience radial drift and have a radial velocity that exceeds that of the gas (Weidenschilling 1977; Brauer et al. 2008). The accreted enriched material cannot be diluted by the depleted gas remaining in the disk, as it is dispersed before the entire disk mass can be accreted. Planet formation can then influence the accretion of these solids, either through taking away material from the disk that can then not be accreted onto the star, or by exerting pressure bumps in the disk that block inward drifting pebbles. The diversity of the exoplanet sample, which in turn raises ideas of different formation scenarios and different formation locations of planets, could mean different abundance changes of the host star depending on the planetary system they harbor.

Unlike for exoplanet atmospheres, where observations of atmospheric chemical abundances prove to be challenging at this time, stellar spectra and therefore elemental abundances can be readily observed. With the above concepts in mind, observations of young stars could provide indirect information about the protoplanetary disks they host, and the formed or forming planets therein. Previous studies have employed such an approach to study the depletion of refractory elements in the Sun compared to solar twins (Meléndez et al. 2009) by considering relative differences in the accretion rate of gas and dust (e.g., Hoppe et al. 2020). Such a difference was previously considered as a result of locking up solids in rocky planets (Chambers 2010), but it may more likely be explained by pressure bumps created as a result of giant planet formation (Booth & Owen 2020). The latter idea is supported by observations of transition disks, which is a group spanning about 10–20% of protoplanetary disks characterized by large gaps (e.g., Kenyon & Hartmann 1995; Koepferl et al. 2013; van der Marel et al. 2018), which could be linked to planets (e.g., Owen & Clarke 2012). Additionally, the observed gaps and rings in millimeter observations with the Atacama Large Millimeter/submillimeter Array (ALMA; e.g., Andrews et al. 2018) are also caused by pressure perturbations preventing a smooth inward drift of pebbles (e.g., Pinilla et al. 2012). The link between gaps and stellar abundance differences was previously explored for young A stars, which offer direct insight into the recent accretion history due to their small convective zones (Jermyn & Kama 2018). This also follows from the correlation between stellar mass and convective zone size during premain-sequence (PMS) evolution, as shown in Fig. 2. It was shown that there is a correlation between large disk cavities and depletion of refractory elements (Kama et al. 2015), which can be observed as a prominent feature due to the absence of convective mixing at the investigated stellar masses. For solar-mass and subsolar-mass stars, convective mixing reduces the impact of accreted material, but studies investigating the composition of gas accreting onto T Tauri stars and TW Hya directly have shown a similar link (McClure 2019; McClure et al. 2020).

The same concept can also be applied to stellar binary systems. Stellar binaries with wide separations (from ~ 100 AU to ~ 1 pc) form from the same molecular cloud at the same time (Kouwenhoven et al. 2010; Reipurth & Mikkola 2012) and are therefore expected to exhibit the same chemical abundances. Observations show, however, that this is not always the case (Teske et al. 2016b; Saffe et al. 2019; Liu et al. 2021). Liu et al. (2021) used high quality spectra from VLT/UVES and Keck/HIRES to derive detailed abundance differences of a wide variety of chemical elements for several wide binary systems containing a different number of planets. The abundance differences could be linked to the presence of a planet around one star, but not the other, in particular for HD106515, where only one constituent harbors exactly one planet (Li et al. 2021). Since close binaries would influence each other’s disk evolution, for example, by stellar radiation, such studies are limited to wide binaries where these effects are negligible. Previously, differences in chemical abundances in wide binaries were used to investigate the formation location of giant planets with respect to the water ice line, under the premise that locking up material in those planets is the dominant effect (e.g., Tucci Maia et al. 2014; Ramírez et al. 2015; Teske et al. 2016a; Bitsch et al. 2018a). This approach was also applied to the formation location of super-Earths with respect to the CO ice line (Bitsch et al. 2018a).

In this work, the convective envelope elemental abundance change after the accretion of disk material is obtained and compared to the initial condition. In a parameter study, we aim to give a perspective on initial conditions commonly applied for planet formation models. Furthermore, the investigation of the influence of planets on stellar abundances is continued, with a particular focus on the wide binary system HD106515. However, instead of considering locking up material in the planets as the main mechanism, pressure bumps acting as traps for dust, but not for gas are considered (Rice et al. 2006; Pinilla et al. 2012; Zhu et al. 2014). These bumps can be created by giant planets reaching a limiting mass to carve a partial gap into the disk’s surface density, the so-called pebble-isolation mass (Morbidelli & Nesvorny 2012; Lambrechts et al. 2014; Bitsch et al. 2018b; Ataiee et al. 2018; Paardekooper & Mellema 2006), typically also opening a deep gap into the disk in subsequent evolution (Lin & Papaloizou 1986; Crida et al. 2006). The approach is in principle similar to previous works on the elemental abundances of the Sun (Booth & Owen 2020; Hoppe et al. 2020). Key differences lie in the differentiation of solid material accreted onto the star into specific chemical species, as well as a more detailed treatment of planetary growth by pebble accretion. The basic processes that lay the foundation for this work are shown conceptually in Fig. 1. The applied model employs a chemical partitioning model to distribute the various chemical elements into molecules and track their movement through the disk as part of the overall viscous evolution in both gaseous and solid phases.

In Sect. 2, we give a brief overview over some core models implemented in the chemcomp code designed to represent the processes described above. Section 3 is used to present the results of the study, where the change to convective zone abundances due to disk accretion are considered. First, we perform a parameter study for cases without a growing giant planet, discussing viscosity and fragmentation velocity (Sect. 3.1), initial stellar mass and metallicity (Sect. 3.2) and the impact of planetesimal formation (Sect. 3.3). After that, results from simulations including a giant planet are presented (Sect. 4), where the influence of the planet’s formation location (Sect 4.1), as well as implications for Jupiter in the solar disk (Sect. 4.2) are discussed. In Sect. 5, we present the application to the HD106515 binary system. Finally, Sect. 6 discusses the results and Sect. 7 gives a brief summary.

2 Methods

We used the chemcomp code (Schneider & Bitsch 2021) to simulate the diffusive evolution of the disk, including radial drift of the dust. To be able to find the individual elemental abundances in the convective zone of the host star, the chemical composition in both gas and dust in the disk is tracked. In addition, the code includes a model for the growth of a giant planet opening a gap in the disk. In this section, the core aspects of these models are highlighted. More details can be found in Schneider & Bitsch (2021).

2.1 Dust drift

As the main ingredient for the time evolution of the disk’s solid component, the two-component model by Birnstiel et al. (2012) is used. Based on numerical and analytical considerations, this model simplifies the treatment of the grain growth and size distribution to save computation time by only treating a population of the smallest grains and one of grains close to size limit set by fragmentation and drift. A more detailed description of the applied dust drift model is given in Appendix A.

The dust velocity is set by representative sizes of the two populations and is evolved with a single advection-diffusion equation, ΣZt+1rr[ r(ΣZu¯ZvΣgasr(ΣZΣgas)) ]=SSaccSpla,${{\partial {{\rm{\Sigma }}_Z}} \over {\partial t}} + {1 \over r}{\partial \over {\partial r}}\left[ {r\,\left( {{{\rm{\Sigma }}_Z}{{\bar u}_Z} - v{{\rm{\Sigma }}_{{\rm{gas}}}}{\partial \over {\partial r}}\left( {{{{{\rm{\Sigma }}_Z}} \over {{{\rm{\Sigma }}_{{\rm{gas}}}}}}} \right)} \right)} \right] = - S - {S_{{\rm{acc}}}} - {S_{{\rm{pla}}}},$(1)

where v is the kinematic viscosity and ūZ is the mass-weighted velocity, u¯Z=(1fm)usmall+fmularge,${\bar u_Z} = \left( {1 - {f_m}} \right){u_{{\rm{small}}}} + {f_m}{u_{{\rm{large}}}},$(2)

with usmall and ularge the radial drift velocities for the representative sizes of the smaller and large population, respectively, as described by Eq. (A.2). The model parameter fm describes the mass distribution of the populations, that is Σlarge = ΣZfm, implemented as fm = 0.97 in the drift limited case and fm = 0.75 in the fragmentation limited case. The right-hand side of Eq. (1) contains the source term describing evaporation of dust, S (see Sect. 2.3), and the sink term Sacc related to pebble accretion of a potential giant planet. In addition, it contains the sink term related to planetesimal formation Spla, discussed in the following subsection.

It is imperative to note here that, despite the original division of the dust into two populations with distinct size and resulting velocity, the dust component as a whole is advected with a singular, mass-averaged, velocity, as given by Eqs. (2) and (1). For the dynamics, that is a reasonable assumption, as the large grain population accounts for 95% or 75% of the dust mass in the drift and fragmentation limited regime, respectively. However, this approximation prevents a treatment of small grains being accreted through the planetary pressure bump (Stammler et al. 2023; Ataiee et al. 2018; Bitsch et al. 2018b; see also Sect. 6.1).

thumbnail Fig. 1

Qualitative illustration of the impact of dust drift and planet formation on the composition of material accreted by the host star. The top part shows the drift of pebbles consisting of three different chemical species, a volatile one (blue), a moderately refractory one (orange) and a highly refractory one (green). The ice lines of the blue and orange species are shown with vertical lines of the corresponding color, while the background color of the disk indicates the gas components. The green species does not evaporate in the disk. More refractory species drift further before evaporating at the ice line, enriching the accreted material more strongly than volatile species. The bottom illustration shows how a giant planet core at pebble isolation mass changes the accreted material’s composition. The volatile component evaporates before the pebbles reach the planet’s pressure bump, allowing diffusion past the planet’s orbit and onto the star. The orange and green species are solid at the pressure bump location and blocked from further accretion. Therefore, the accreted material is no longer enriched in the moderately (orange) and highly refractory (green) species.

2.2 Planetesimal formation

In the chemcomp code, a pebble flux-regulated planetesimal formation model is used (Lenz et al. 2019; Voelkel et al. 2020), which assumes a continuous mechanism that converts a fraction of the local pebble flux into planetesimals, preventing that material from partaking in the enrichment of the central star’s convective envelope. The sink term Spla in Eq. (1) is given by Spla=Σ˙pla=ϵd(r)M˙peb2πr.${S_{{\rm{pla}}}} = {{\rm{\dot \Sigma }}_{{\rm{pla}}}} = { \over {d\left( r \right)}}{{{{\dot M}_{{\rm{peb}}}}} \over {2\pi r}}.$(3)

Hence, planetesimal formation can prevent heavy elements from being accreted onto the central star by locking up solids as they are created, because they are decoupled from the protoplanetary disk and do not drift inward. The sink term in Eq. (3) includes the pebble flux peb, the radial separation of pebble traps, d(r)=dplaHgas,$d\left( r \right) = {d_{{\rm{pla}}}}{H_{{\rm{gas}}}},$(4)

and the efficiency parameter ϵ, which describes the fraction of the pebble flux used for the formation of planetesimals after drifting the distance d. Both the separation in terms of the gas pressure scale height, dpla, and efficiency ϵ are model parameters. The former is kept constant at dpla = 5, the latter at ϵ = 0.05. The efficiency is varied in Sect. 5. Different from the original model by Lenz et al. (2019), only the large grain population is considered for the pebble flux, but there is no specific Stokes number cutoff. Furthermore, a threshold pebble flux must be reached, which is set by the condition that the mass fraction which can be converted at a given location is high enough such that at least one planetesimal can be formed, M˙crit=Mplaϵτt,${\dot M_{{\rm{crit}}}} = {{{M_{{\rm{pla}}}}} \over {{\tau _t}}},$(5)

where Mpla=43πρplaRpla3${M_{{\rm{pla}}}} = {4 \over 3}\pi {\rho _{{\rm{pla}}}}R_{{\rm{pla}}}^3$ is the mass of one planetesimal, with ρpla = 1 g cm−3 and Rpla = 50 km kept constant, corresponding to the characteristic size of planetesimals formed by the streaming instability (Johansen et al. 2015; Simon et al. 2015; Klahr & Schreiber 2020), and the trap lifetime, τt=100×2πΩK.${\tau _t} = 100 \times {{2\pi } \over {{{\rm{\Omega }}_K}}}.$(6)

At disk locations where the critical flux is not reached, the sink term is identically zero.

Table 1

Volume mixing ratios for the considered molecular species, given by the applied chemical partitioning model, along with the condensation temperatures of the respective species.

2.3 Disk chemical composition

Initially, the chemical composition of the protoplanetary disk is assumed to be identical to the composition of the star, which is motivated by the disk formation process, where both are made of the same molecular cloud material. The initial composition is defined by the initial [X/H] values defined relative to the solar measurements (Asplund et al. 2009). Together with the solid-to-gas ratio ϵd, taken at a temperature where all available solids are condensed (T < 20 K), the absolute mass of every chemical element is set.

The chemical elements are distributed into molecular species Y according to the chemical partitioning model described in Schneider & Bitsch (2021). This simple model does not treat changes in the solid composition caused, for instance, by chemical reactions, but tracks the relevant main carriers. This approach is justified, because normally the chemical reaction time scales are longer than the drift timescales, leaving inward drifting pebbles unaltered (Booth & Ilee 2019; Eistrup & Henning 2022). Table 1 provides the corresponding volume mixing ratios based on the disk elemental abundances, along with the condensation temperatures of the molecules (Lodders 2003). Initially, the entire mass of molecule Y is in the gas phase in the interior disk region where the temperature is larger than the condensation temperature, and completely in the solid phase exterior to that. During the simulation, transitions between solid and gas phases of the dust species in both directions are possible, accounting for evaporation as molecules cross their ice lines from colder to warmer areas in the disk, as well as for condensation in the opposite case1. This is reflected in the source term S (Eq. (1)) and a corresponding term for the gas diffusion equation, given for species Y as SY=Σ˙Y={ Σ˙Yevapr<rice,Y,Σ˙Ycondrrice,Y, ${S_Y} = {{\rm{\dot \Sigma }}_Y} = \left\{ {\matrix{ {{\rm{\dot \Sigma }}_Y^{{\rm{evap}}}} \hfill &amp; {r lt; {r_{{\rm{ice,Y}}}},} \hfill \cr {{\rm{\dot \Sigma }}_Y^{{\rm{cond}}}} \hfill &amp; {r \ge {r_{{\rm{ice,Y}}}},} \hfill \cr } } \right.$(7)

where rice,Y describes the radial coordinate where T(rice,Y) = Tcond,Y in the disk. Under the assumption that dust can only condensate onto preexisting grains, the condensation term reads Σ˙Ycond=3ϵp2πρsΣgas,Y(ΣZasmall+Σlargealarge)ΩKμμY,${\rm{\dot \Sigma }}_Y^{{\rm{cond}}} = {{3{_p}} \over {2\pi {\rho _s}}}{{\rm{\Sigma }}_{{\rm{gas}},{\rm{Y}}}}\left( {{{{{\rm{\Sigma }}_Z}} \over {{a_{{\rm{small}}}}}} + {{{{\rm{\Sigma }}_{{\rm{large}}}}} \over {{a_{{\rm{large}}}}}}} \right)\,{{\rm{\Omega }}_K}\sqrt {{\mu \over {{\mu _Y}}}} ,$(8)

with µ the mean molecular weight of the disk, µY the mass of the molecular species Y, a the representative size of the small and large dust component, respectively, ΩK the Keplerian angular frequency, and ϵp = 0.5 the sticking efficiency. The evaporation term reads Σ˙Yevap=ΣZ,Yu¯Z1×103AU,${\rm{\dot \Sigma }}_Y^{{\rm{evap}}} = {{{\Sigma _{Z,Y}}{{\bar u}_Z}} \over {1 \times {{10}^{ - 3}}{\rm{AU}}}},$(9)

where it is assumed that solids passing their ice line evaporate within 1 × 10−3 AU.

The integration of Eq. (1) and the analogous equation for the gas is done using a modified version of the flux-conserving donor-cell scheme (Birnstiel et al. 2010). The particular implementation utilized in chemcomp is adapted from the unpublished DISKLAB code (Dullemond & Birnstiel 2018, priv. comm.). The disk radial grid extends from 0.1 AU to 1000 AU and is divided into 500 cells whose sizes are distributed logarithmically.

thumbnail Fig. 2

Time evolution of the stellar convective envelope mass obtained from models by Hoppe et al. (2020) for various [Fe/H] (left panel) and M (right panel), indicated by a change in line color. The ordinate shows 1 – MCE/M, so that if the star has a more massive convective zone, a lower value is shown, to indicate that the convective zone abundances adapt to the composition of the accretion flux more slowly.

2.4 Convective zone abundances

Solid and gaseous disk material crossing the inner edge is accreted onto the convective zone of the star and mixed with the material therein. To investigate this, we added routines to chemcomp that calculate the change of the convective zone abundances based on the flux onto the star and its mass. To obtain the mass of material that is accreted, the flux of gas and solids at the inner edge is integrated. The gas flux is given by M˙gas=2πr3rr(vΣgasr).${{\dot M}_{{\rm{gas}}}} = - 2\pi r{3 \over {\sqrt r }}{\partial \over {\partial r}}\left( {v{{\rm{\Sigma }}_{{\rm{gas}}}}\sqrt r } \right).$(10)

As seen in Eq. (1), the mass flux for solids reads M˙dust=2πrΣZu¯Z2πrΣgasvr(ΣZΣgas).${{\dot M}_{{\rm{dust}}}} = 2\pi r{{\rm{\Sigma }}_Z}{{\bar u}_Z} - 2\pi r{{\rm{\Sigma }}_{{\rm{gas}}}}v{\partial \over {\partial r}}\left( {{{{{\rm{\Sigma }}_Z}} \over {{{\rm{\Sigma }}_{{\rm{gas}}}}}}} \right).$(11)

The abundance of the convective envelope was advanced in time as follows. First, the mass that was accreted during a given time step from tk to tk+1 was added to the envelope. In detail, the elemental abundances and total mass of the envelope at time tk were used to calculate the total mass per element at time tk, Mj,k, Mi,k=MCE,k(NiNH)mimH×10[ i/H ]1+j{ X }\{ H }(NjNH)mjmH×10[ j/H ],${M_{i,k}} = {{{M_{{\rm{CE}},k}}{{\left( {{{{N_i}} \over {{N_H}}}} \right)}_ \odot }{{{m_i}} \over {{m_H}}} \times {{10}^{\left[ {{i \mathord{\left/ {\vphantom {i H}} \right. \kern-\nulldelimiterspace} H}} \right]}}} \over {1 + \sum\limits_{j \in \left\{ {\rm{X}} \right\}\backslash \left\{ {\rm{H}} \right\}} {{{\left( {{{{N_j}} \over {{N_H}}}} \right)}_ \odot }{{{m_j}} \over {{m_H}}}} \times {{10}^{\left[ {{j \mathord{\left/ {\vphantom {j H}} \right. \kern-\nulldelimiterspace} H}} \right]}}}},$(12)

where the sum is taken over all modeled elements except hydrogen, and mi indicates the atomic mass of element i. This approach holds true under the assumption that the envelope preserves its composition during its evolution. After adding the accreted mass, Mi,k+1=Mi,k+tktk+1(M˙gas,k+M˙dust,k) dt,${M_{i,k + 1}} = {M_{i,k}} + \int_{{t_k}}^{{t_{k + 1}}} {\left( {{{\dot M}_{{\rm{gas,k}}}} + {{\dot M}_{{\rm{dust,k}}}}} \right)} \,{\rm{d}}t,$(13)

the elemental abundances were recalculated. This implies that the mixing is instantaneous. During the next time step, the procedure was repeated for a new convective envelope mass MCE,k+1 as given by the convective zone evolution. Overall, material accreted at a time when the convective zone is massive has a smaller impact on its abundance. Therefore, as material with a particular [X/H] is accreted, it takes longer for the abundances in the convective zone to approach that value if the convective zone is massive than if it is light.

For the mass evolution of the convective envelope, models from Hoppe et al. (2020) were used. In essence, the mass fraction of the envelope decreases over time, allowing late accretion (or absence thereof) to have a larger impact on [X/H] than earlier on. The particular time evolution depends on two main parameters considered in this work, [Fe/H] and M, which are shown in Fig. 2. It is illustrated how for constant M, higher stellar metallicity results in the envelope mass decreasing slower and at later times, such that the mass is higher at all times. For constant metallicity, a higher stellar mass leads to the opposite effect, where the decline of mass in the convective zone is faster and occurs earlier, so that the mass in the convective zone is lower at all times. We note that the accretion of material has an impact on the stellar evolution itself (e.g., Serenelli et al. 2011), also modifying the convective zone evolution (Baraffe & Chabrier 2010), which is neglected for this study. Furthermore, the stellar models from Hoppe et al. (2020) were used to compute the average luminosity in the first 10 Myr, which in turn was used to calculate the temperature profile of the disk.

thumbnail Fig. 3

Surface density evolution of carbon, oxygen, iron and nitrogen, for α = 1 × 10−4. The sum of the gas and dust component is shown. The color indicates the surface density of the element with respect to the hydrogen surface density, compared to the initial fraction at a given location on a logarithmic scale, such that a blue color indicates that the surface density of the element is enriched, and vice-versa for a red color. Due to this choice of norm, a blue color at the inner edge corresponds to an increase in [X/H] for that element up to the magnitude of the enrichment. All evaporation front locations of species where the corresponding element is a component in are indicated with black dashed lines in the relevant panels.

Table 2

Model parameters varied in the simulations for the parameter study, and their default value.

3 Results

In this section, we present simulations of the evolution of the protoplanetary disk’s dust and gas components, focusing on the impact of disk accretion on the chemical abundances of the host star’s convective envelope. We explore the effect of the variation of parameters described in Table 2. If not otherwise specified, the default values indicated in the table were used, and planetesimal formation was not considered. Furthermore, we assumed a disk radial surface density profile with an exponential cutoff at R0 = 75 AU and a total disk mass of M0 = Mdisk/M = 0.1 unless specified otherwise, which results in a gravitationally stable disk.

3.1 Abundance changes due to dust drift

Based on the growth and radial drift of material in the disk, the surface density of the various disk components naturally changes over time. This surface density evolution can be seen in Fig. 3, for α = 1 × 10−4, where the dust grains can grow to large sizes and thus drift inward very fast. The chemical model we used contains numerous elements whose evolution is very similar. Therefore, the figure presents the evolution of the surface density of carbon, oxygen and nitrogen, while iron is chosen to represent all refractory species. The ice lines of all molecular species an element is part in are shown as well. The surface densities are presented relative to the surface density of hydrogen and to the initial fraction of the elements with respect to hydrogen. This achieves a view that is consistent with the analysis of [X/H] in the star’s convective zone in the following sections.

As dust grains grow to large sizes, their radial drift speed surpasses the gas radial speed significantly. This leads to a pile-up of material at the ice lines, interior to which a particular species moves together with the slow gas, because it is now in a gaseous form after its evaporation. Large pebbles drift from the outer disk to the location of the ice lines quickly, and therefore create a large enrichment at those locations. One can observe that for elements mainly represented by volatile species (C, N, O), the region outside the ice line, and consequently the mass that can drift toward the ice line, is much smaller than for those present in refractory species (e.g., Fe). Altogether, there is a period of significant enrichment for refractory species close to the inner edge, which is very brief, after which the species is mostly depleted throughout the disk. In contrast, the enrichment occurs farther away from the inner edge for volatile species. Hence, the material needs to be viscously accreted in gaseous form to reach the inner edge and be removed from the disk, which means it is present also at later times in the disk. As oxygen and carbon are contained in both volatile and refractory species, they are an in-between case. This gives rise to both an early, albeit weaker, enrichment near the inner edge and a presence in the disk at later times.

Material that crosses the inner edge is assumed to be accreted onto the star’s convective zone. The difference in the surface density evolution between the elements, and its dependence on particle size, is therefore reflected in the composition of the accreted material over time, (see Eq. (11)). This is displayed in Fig. 4. For α = 1 × 10−3, the slow but constant enrichment of every element but nitrogen corresponds to a similar behavior of the accretion rate, which remains largely constant, close to the initial fraction of material. For α = 1 × 10−4, though, the enrichment of refractories corresponds to a strongly increased accretion rate at early times and a strongly decreased accretion rate at later times, where the refractories are depleted. In the case of nitrogen, the accretion rate does not increase as strongly and remains almost constant after an initial jump, owning to the slow depletion of the mostly gaseous nitrogen in a disk with low viscosity. The initial jump stems from the time needed for piled-up material at the far-out ice lines of N2 and NH3 to diffuse to the inner edge after evaporating.

To start, we discuss the time evolution of the convective zone abundances caused by the composition of the accreted material while focusing on the impact of particle size and disk viscosity. Figure 5 shows the time evolution of [X/H] – [X/H]0 for all four combinations of α and ufrag. In addition, the evolution of the mass of the convective envelope is shown. Generally, the evolution of the chemical abundances can be explained by the envelope’s [X/H] steadily approaching the value of the material accreted from the disk, with the speed of the approach being set by the mass of the envelope. For α = 1 × 10−4 and both values of ufrag, refractory abundance initially rises sharply, as expected in a disk that allows dust to grow to larger sizes.

Carbon abundance also rises steeply at the beginning, but falls short of reaching the same maximum abundance that the refractories and sulfur do. While the carbon grains also evaporate close to the star and are accreted quickly, only 60% of the mass of carbon is in the form of grains, with the rest in volatile molecules accreted with the gas phase. The oxygen abundance rises less steeply than that of carbon or the refractories, coming close to the refractory abundance at the end of the simulation. The difference to carbon stems from the fact that less oxygen mass is distributed toward refractory species, which causes the less steep rise. Furthermore, the ice line of the most massive volatile component, water, is closer to the inner edge of the disk than the ice lines of the carbon-containing volatiles are. With water being present in the disk for a more extended period of time than the refractories, the lighter convective envelope at later times allows the abundance to reach a level close to the refractories toward the end of the simulation.

The difference between ufrag = 1 m s−1 and ufrag = 10 m s−1 in this case is the steepness of the initial slope of the abundance increase, which is greater for uſrag = 10 m s−1; an effect that can be attributed to the difference in maximum grain size and the corresponding maximum radial drift speed. For the same reason, the maximal refractory abundance is also reached earlier for ufrag = 10 m s−1. The time axis of Fig. 5 can also be viewed as the disk lifetime, so that the abscissa would indicate the final abundance difference of the central star after having accreted a disk with a lifetime as indicated on the ordinate. In turn, this interpretation leads to the conclusion that a disk that is short-lived leaves behind a star whose convective zone is more enriched in refractories and less enriched in volatiles than a long-lived disk would.

As discussed above, the case of α = 1 × 10−3, ufrag = 1 m s−1 exhibits small particles drifting slowly through the disk. Therefore, the behavior of Δ[X/H] is substantially different for this case, with it increasing monotonically throughout the simulated time frame. While the dust grains still drift at speeds exceeding the gas, buildup of refractory grains is hindered, preventing the quick early accretion of grains. Instead, grains are accreted steadily as they reach their ice lines and get evaporated. As fully refractory elements are evaporated closer to the star, they can be accreted more quickly, thereby allowing the abundance difference to rise more quickly as well. Equally, oxygen and carbon rise slightly slower due to being part of both refractory and volatile species, and nitrogen rises the slowest of all considered elements as a fully volatile element.

Increasing the fragmentation velocity to ufrag = 10 m s−1, particles can grow to larger sizes, and previous considerations can be applied. When comparing this case to the one of lower viscosity, an even steeper initial slope is shown, and the peak of the abundance difference is reached earlier on. This is because the small distance between the refractory species’ ice lines and the inner edge of the disk can be traversed faster by evaporated material at this viscosity. However, there is a strong difference in evolution after the initial peak of fast refractory accretion. It relates to the state of the disk after the initial, fast growth and subsequent accretion of refractory dust has concluded. At this time, it is depleted of the elements contained in those grains (see Fig. 3). Moreover, gas is accreted at a high rate provided by the high viscosity. Combined, this leads to gas depleted in refractory elements being accreted rapidly, reducing the overall abundance of these species in the convective envelope. The same effect is responsible for the strong rise of nitrogen at later times, too. Nitrogen is enriched in the outer parts of the disk, and fast accretion of gas at high viscosity allows the enriched parts to be accreted onto the star during the later stages. The nitrogen abundance is higher at the end of the simulation for higher uſrag as well, due to a stronger enrichment at the nitrogen-related ice lines.

The relevance of late accretion in the simulation is amplified by the change in mass of the convective zone. The composition of material accreted onto to the star at later times plays a bigger role at late than at earlier times, as the shrinking convective envelope amplifies the rate of the abundance change. While this is relevant for all scenarios, it is best seen in the case of α = 1 × 10−3, ufrag = 10, where the abundance of nitrogen is almost the same as the other elements at the end of the simulation.

We note that the simulation time was picked as the typical upper limit for the lifetime of protoplanetary disks (Ribas et al. 2015). If the disk were to be dispersed earlier than at 10 Myr, the abundance change would stop at that time. From this perspective, Fig. 5 indicates the final envelope abundances as a function of disk lifetime. Furthermore, Fig. 5 also presents a comparison to simulations with M0 = 0.05. Naturally, a smaller disk mass diminishes the effect that disk accretion can have on the stellar abundances. For an even lighter disk with M0 = 0.01, the resulting abundance differences are of the order Δ[X/H] ~ 10−3. However, in this study, we focus on the maximal influence of disk accretion on the stellar composition and use M0 = 0.1 for the remainder of the study.

thumbnail Fig. 4

Accretion flux of disk material onto the host star for small particles (α = 1 × 10−3, left panel) and for large particles (α = 1 × 10−4, right panel) as a function of time. On the ordinate, the accretion rate with respect to the hydrogen accretion rate, i/MH, for various elements i is shown, indicated by different colors. The accretion rate is shown normed to the mass fraction with respect to hydrogen of the element as defined by the initial condition. Therefore, if the ordinate value is larger than 1, the material that is accreted results in a growth of the abundance of element i in the convective zone up to the magnitude of the enrichment. The accretion rate is directly related to the surface density evolution (see Fig. 3).

thumbnail Fig. 5

Stellar convective envelope elemental abundances as a function of time relative to the initial abundances, [X/H] – [X/H]0. Four different simulations are shown, varying α and ufrag to consider different maximum grain sizes in combination with different gas accretion speeds. Differently colored solid lines represent the abundances of the elements carbon (yellow), oxygen (red), iron (green), and nitrogen (blue), while the black dashed line shows the convective zone mass evolution as 1MCEM$1 - {{{M_{{\rm{CE}}}}} \over {{M_ \star }}}$, as in Fig. 2, corresponding to M = 1 M and [Fe/H] = 0. Solid lines indicate results for M0 = 0.1, while dashed lines show a comparison to simulations with M0 = 0.05.

3.2 Influence of different stellar metallicity and mass

In order to discuss the influence of different initial stellar metallicities, and therefore also different initial disk elemental abundances (see Table 2), it is necessary to consider how a change in metallicity affects the abundances of the chemical species. While typically the iron abundance [Fe/H] is used as a proxy for the overall metallicity, the elemental abundances scale with [Fe/H] separately for each element (e.g., Burbidge et al. 1957; Bitsch & Battistini 2020). Analogously to the methods employed in Bitsch & Battistini (2020, hereafter BB20), the scaling of abundances of the chemical model’s elements was obtained from the GALAH+ survey’s third data release (Buder et al. 2021). Due to changes in the reduction pipeline and analysis workflow between the second and third data release, the selection criteria were adjusted compared to BB20. We required a star to have Teff > 5000 K and log g > 4, as well as quality flags of flag_sp = 0, flag_fe_h = 0, flag_x_fe = 0 and snr_c3_iraf > 30. As in BB20, stars were grouped in [Fe/H] bins with size 0.1. The corresponding abundance [X/H] was then obtained as the mean [X/H] per bin, shown in Fig. 6. For sulfur, the same scaling as for silicon was assumed (Chen et al. 2002). While the uncertainty of the elemental abundances is large, especially for elements not considered by BB20, the resulting scaling of the elements was employed whenever [Fe/H] was adjusted in a simulation. We note that there are slight differences in the error bars for the individual elements between our analysis using GALAH DR3 compared to BB20, who used GALAH DR2. In addition, the error bars and trends are also recovered with the Apogee survey, see Cabral et al. (2023) for a comparison.

The time evolution of the convective zone abundances for low metallicity ([Fe/H] = −0.4) and high metallicity ([Fe/H] = 0.4) are shown in Fig. 7 for two values of ufrag. When considering the results of the simulations with uſrag = 1 m s−1, where grains cannot grow to large sizes, the final values for [X/H] are smaller for higher initial metallicity, while the overall trend remains similar. In contrast, simulations with large dust grains show larger final abundances for a higher metallicity. In addition, the species that exhibits the strongest enrichment at the end of the simulation changes between refractories having the highest abundance for [Fe/H] = +0.4 and nitrogen, carbon and oxygen having the highest abundance for [Fe/H] = −0.4.

Because the initial metallicity was changed for both the disk and the star, both being made out of material from the same molecular cloud, there are several changes to the disk evolution when compared to solar metallicity. First, the total surface density of dust in the disk increases for a higher metallicity, which decreases the growth timescale of dust grains. Second, a higher metallicity star has a lower luminosity averaged over 10 Myr, as computed from the Hoppe et al. (2020) stellar models. As the outer region’s temperature profile is dominated by stellar irradiation, this leads to a colder outer disk, shifting ice lines of volatile molecular species closer to the star. At the same time, the inner disk becomes hotter, as the dust optical depth increases with the dust-to-gas ratio. Third, as the elements scale differently with [Fe/H] (see Fig. 6) and the relative distribution of elements determines the partitioning in the chemical model as per Table 1, mass contributions to elements by chemical species can shift to more volatile or refractory species, for example, decreasing the water fraction by a factor of ~2 when increasing [Fe/H] from −0.4 to 0.4.

However, the difference in the evolution of the convective zone abundances is governed by the altered convective zone evolution. Material accreted late during the disk’s lifetime has a particularly strong impact on the changes caused by the metallicity variation. Figure 2 shows that the convective zone mass fraction is larger for all times when increasing the metallicity. Therefore, the rate of abundance change is lowered during late accretion phases for high metallicities. As disks that do not allow significant growth of dust grains do not exhibit a drop in abundance due to the accretion of depleted material, this results in smaller final abundances in such a case. If grains can grow larger and depletion becomes relevant, the drop of abundances caused by the accretion of depleted material is diminished by the more massive convective envelope. In such cases, high metallicities create large final abundances.

If the central star is more massive, the mean luminosity is higher over the disk lifetime, while the convective envelope mass is reduced (see Fig. 2). Both of these aspects are also modified by a change in initial [Fe/H], making the discussion of the M parameter analogous to the discussion of metallicity changes. With the disk mass having been defined relative to the stellar mass (see Table 2), the surface density increases for increasing stellar mass, while keeping the same dust-to-gas ratio, unlike for a change in metallicity. Similar to the case of metallicity variation, changes to the abundance evolution are dominated by the change in convective zone size.

Figure 8 shows the abundance changes for stellar masses of 0.6 M, 1M and 1.5 M. For the highest stellar mass, the disk lifetime was shortened to 5 Myr to account for observations of protoplanetary disks around more massive stars (Ribas et al. 2015), as well as to allow the star to still have a convective envelope at the end of the disk’s life. For refractory elements, there is only a small difference between the cases of a low-mass and a solar-mass star. This is because the refractory dust grains are accreted very early, at a time when the star is still fully convective in both cases. On the other hand, for the high-mass star, the early accretion of refractories is still commencing as the mass of the convective zone shrinks, allowing for a higher abundance peak. For the volatile nitrogen and the large water component of oxygen, the accretion is prolonged over the disk’s lifetime. Therefore, the presence of a more massive convective zone at later times prevents those elements from achieving an abundance difference as high as for the solar-mass star. As the most massive star has the smallest envelope, the opposite is true here, where oxygen can reach higher values of Δ[X/H].

thumbnail Fig. 6

Scaling of initial [X/H] with [Fe/H] for the elements considered in the disk model based on GALAH+ DR3. The scaling of the abundances is indicated by solid colored lines. In addition, a reference line of 1:1 scaling is shown as a black dashed line, and a dummy data point indicating the mean error of the measurements is shown in the top left corner of each panel. In the left panel, elements that were previously considered by Bitsch & Battistini (2020) are presented, while the remaining elements not previously considered, but used in the model, are shown in the left panel. Sulfur is not indicated; here, the same scaling as for silicon is assumed.

thumbnail Fig. 7

Same as Fig. 5, but the initial [Fe/H] is varied instead of α. It is kept constant at α = 1 × 10−3. On the left-hand side, the case for an initial metallicity of [Fe/H] = −0.4 is shown for both ufrag = 1 m s−1 and ufrag = 10 m s−1, whereas the right-hand size shows simulations with [Fe/H] = +0.4.

thumbnail Fig. 8

Same as Fig. 5, but the stellar mass is varied instead of α and ufrag. The three panels show simulations with M ∈ {0.6 M,1.0 M,1.5 M}. For the simulation with M = 1.5 M, the simulation run time is reduced to 5 Myr in accordance with observations of accretion disks around massive stars.

3.3 Influence of planetesimal formation

We conclude the discussion of planet-less accretion by considering the formation of planetesimals in the disk. Since planetes-imals form from a fraction of the pebble flux, dust locked up in planetesimals cannot be accreted onto the star, as planetesimals do not experience radial drift. Therefore, the final stellar abundances are lower if planetesimal formation is included. Two examples, for high and low viscosity, are shown in Fig. 9. With the direct dependence of the planetesimal formation sink term on the pebble flux, the absolute speed of the dust grains sets how strongly the enrichment is diminished, rather than the speed difference to the gas that is responsible for the enrichment itself. Hence, at high viscosity α = 1 × 10−3, even though the dust grains do not grow to large sizes, a large mass of heavy elements gets locked up and prevented from contributing to enrichment. As, in this case, the over-densities at the ice lines are already small, this results in a negative value of Δ[X/H] at the end of the simulation for refractory species. In disks with lower viscosity, however, the convective envelope is still enriched in refractories after 10 Myr, as strong enrichment is still possible even if some material is locked up. For all simulations presented in Fig. 9, the inclusion of planetesimal formation results in a maximal abundance difference for all elements below 10−2, an order of magnitude lower than differences obtained without it, which is lower than the precision that can be achieved with observations at this time.

Planetesimals can be formed anywhere in this model, and the total pebble flux is the largest for refractory species because they are in solid form even in the inner part of the disk. Therefore, they are subject to more removal of material. This is reflected in iron being the most depleted at the end of the simulation. It is also apparent when comparing the final abundances of oxygen and carbon. In our model, 60% of carbon is in the form of grains, whereas most oxygen is in water (see Table 1). Since the carbon grain evaporation front is closer to the star than the water ice line, more carbon is locked up in planetesimals than oxygen. Increasing the planetesimal formation efficiency ϵ would shift the location of planetesimal formation further outward, as pebbles in the outer disk would be converted into planetesimals before they can contribute to formation in regions of the disk closer to the star (Voelkel et al. 2020). Consequently, this means that planetesimals would include an even higher relative mass fraction of solids with lower condensation temperatures, as opposed to the less efficient case where the formation would be shifted more toward the inner disk regions, where they are present in gaseous form.

thumbnail Fig. 9

Same as Fig. 5, but only α is varied. For the shown simulations, planetesimals are formed in the disk.

4 Influence of a giant planet

In this section, simulations that include the formation of one giant planet, accreting pebbles onto a protoplanetary seed to from its core and accreting gas after reaching pebble isolation mass, if applicable, are presented. A brief summary about the employed planet formation model is given in Appendix B.

4.1 Planet location

The placement of the planet was varied to study how blocking different molecular species in the disk from accreting affects the final abundances. To achieve consistency between disks with different temperature profiles as model parameters were varied, the planet was always placed just outside one of three ice lines, at an orbital radius 110% of the ice line’s radial distance to the star. The considered ice lines were those of Fe3O4, to study a case where only highly refractory elements are blocked, H2O and CO2. For simplicity, the planet was prevented from migrating.

Once the giant planet reaches the pebble isolation mass, solid material weakly coupled to the gas is trapped in a pressure bump outside the planet’s orbit and prevented from moving inward, which blocks its accretion onto the star. Hence, the position of the planet relative to the ice lines impacts the composition of material that is added to the envelope. Figure 10 shows the pile-up of oxygen-containing solids outside the planet’s orbit. If the planet forms at the innermost location, only oxygen contained in refractory components is blocked, whereas a larger fraction of oxygen is prevented from accreting onto the star if the water and CO2 ice lines are blocked.

The impact of the planet’s positioning is reflected in the time evolution of [X/H] of the convective zone, shown in Fig. 11. The figure displays both the three different planet positions relative to the aforementioned ice lines, and cases for both ufrag = 1 m s−1 and ufrag = 10 m s−1, as a proxy for the dust grain size. Most notably, the position of the planet relative to the water ice line is reflected in the final oxygen abundance in the convective zone. For the innermost planet, this abundance is highest, as the water ice line is exterior to the planet’s orbit. For the other two farther outward planets, the opposite holds true, reducing the final abundance.

While the planet’s location relative to the ice lines defines what chemical species can potentially be blocked, some material can still be accreted. Material present interior to the planet before core formation is always accreted onto the star. If the core forms at a large distance to the star, more material is unaffected and can contribute to the abundance increase in the convective envelope. Furthermore, the growth time of the planet is prolonged in the inner regions of the disk due to a reduced pebble flux. Also, the magnitude of the isolation mass increases with radial distance (Lambrechts et al. 2014; Bitsch et al. 2018b). As a result, a planet that forms furthest from the star takes longest to reach pebble isolation, and more grains can be accreted before that happens, creating a larger abundance in the envelope. The middle planet takes the least time to reach the isolation mass, whereas the innermost planet’s time lies in between the other two cases. In disks where grains remain small, this difference is most noticeable. Overall, the magnitude of the pebble flux is indicative of both the growth time of the planet and the total grain mass that drifts past the planet before reaching the isolation mass.

Disks with small grains form less massive planets. The pressure bump that is created by the planet to reach pebble isolation is therefore weaker. If the planet is placed at the innermost location, the low pebble isolation mass gives rise to a particularly weak pressure bump. Combined with significant gas drag experienced by small particles, this low pressure bump does not successfully prevent the grains from drifting past the planet’s orbit, even after the pebble isolation mass is reached. This is reflected in Fig. 11. Comparing the top left and top center panels, the former shows a significantly stronger enrichment of, for example, iron. All iron related ice lines are blocked in both cases. Here, the stronger enrichment is caused by the failure of the weak pressure bump to prevent pebble drift past its orbit.

thumbnail Fig. 10

Same as Fig. 3, but the surface density evolution of oxygen in solid and gaseous form is shown for three different planet formation locations. Additionally, the H2O, CO2 and CO ice lines are labeled.

thumbnail Fig. 11

Same as Fig. 5, but now, a planet is included in the disk. The viscosity is kept constant, but ufrag is varied. Simulations where the planet was placed at three different locations are shown, relative to the Fe3O4, H2O and CO2 ice lines. Each column represents a different planet placement, with it being placed further out going from the left-most to the right-most column. The top row depicts disks with ufrag = 1 m s−1 and the bottom row those with ufrag = 10 m s−1.

4.2 Implications for Jupiter and the Solar System

Another interesting application of the methods presented here is to study how Jupiter’s growth and the accretion of the solar protoplanetary disk influenced the composition of the Sun. We do not aim for great detail in this section, but to give order-of-magnitude estimates using a simple approach.

As the most massive planet of the Solar System, Jupiter has the largest impact not only on the dynamics, but also on planet-disk interactions and potential influence on the accretion of solids onto the Sun. Since the chemcomp code can only consider a single planet in the disk, it was the only planet considered for determining the values of Δ[X/H] expected for the Sun. Several further simplifications are made. First, Jupiter was kept fixed at its present-day location of ~5.2 AU, neglecting migration scenarios where it has experienced both significant inward and outward migration in conjunction with Saturn (Walsh et al. 2011). Second, the gas accretion of the planetary core was stopped once it reached the mass of Jupiter. Changing the accretion rate of gas by a constant factor ≥ 0.1 was found to have no impact on the result. Furthermore, simplifications were made for the initial conditions of the Sun and its disk. In the previous sections, it was shown that the initial disk and star abundances are different from the present-day values after the disk has been accreted in part and dispersed. Despite this, rather than making assumptions about the initial abundances of the Sun and trying to reproduce present-day solar values, the simulation used the observed solar values (e.g., Asplund et al. 2009) as a starting point to find the abundance changes created by disk accretion. The main difference to the real situation created by this approach is the starting metallicity, which sets the evolution of the convective zone mass and the initial heavy element surface density in the disk. In Sects. 3.1 and 4, it was shown that the differences created by disk accretion are on the order of 10−2 in [X/H], which can be viewed as the expected difference in metallicity for the stellar model. The metallicity resolution of the stellar models used in Hoppe et al. (2020) was 0.05 dex, so that the change to the solar model by considering a modified initial metallicity of the Sun would in most cases be outside the scope of the model resolution.

The expected change of the convective zone abundances is presented in Fig. 12, where the final abundance difference after 4.5 Myr is shown. The choice of disk lifetime was based on information from the formation of chondrites (Johnson et al. 2016; Kleine et al. 2020). The abscissa shows the results for four values of α sampling the log space between 10−4 and 3 × 10−3 regularly. Again, cases for low (ufrag = 1 m s−1) and high (ufrag = 10 m s−1) fragmentation velocity are considered. In the applied model, Jupiter blocks the water ice pebbles, but not the CO2 pebbles, because they have already evaporated upon reaching Jupiter’s orbit. We note that Fig. 12 shows the abundance difference compared to the initial solar value, not to a case where Jupiter does not form in the disk.

The disk cutoff radius R0 and mass M0 affect the evolution of the disk, mainly due to a change in the surface density and total dust mass. For higher M0, the higher available mass of dust creates a stronger enrichment at the ice lines, which in turn is reflected by the convective zone abundances reaching a higher peak value. A smaller R0 increases the dust surface density Σd. While this increases disk temperature, the dominant effect is the reduced dust growth and drift timescale. As a result, the accretion flux onto the central star becomes enriched already at earlier time, when the convective zone is still massive. Depending on the subsequent convective zone evolution, given by stellar mass and metallicity, as well as the general accretion timescale, given by α and ufrag, this can increase or decrease the final and peak elemental abundances. Therefore, final solar abundance differences caused by disk accretion are computed for R0 ∈ {75 AU, 150 AU} and M0 ∈ {0.05,0.1}. The expected change of abundance is then taken as the average.

Because nitrogen-related ice lines are outside Jupiter’s orbit, [N/H] is not affected by the pressure bump caused by it. As such, the nitrogen abundance in the Sun’s convective envelope is governed by disk effects. At low viscosity, the pile-up at the N2 ice line is not accreted onto the star. Combined with the fact that only 10% of the total nitrogen mass is in the form of NH3, changing the size of the particles without changing the viscosity, that is due to variation of ufrag, has little effect. This changes once the viscosity is high enough for the enrichment at the N2 ice line to reach the inner edge of the disk during the simulation. In that case, a stronger enrichment corresponding to a larger fragmentation velocity results in a higher nitrogen abundance. However, for ufrag = 1 m s−1 and α = 3 × 10−3, the accretion of gas is more efficient than evaporation of material at the nitrogen ice lines, where no significant pile up occurs due to small grain sizes. In this scenario, the gas becomes depleted in nitrogen and the final value for [N/H] is lower than initially.

For the other elements, though, at least some molecules that contain them have their ice line blocked once the planet carves the shallow gap in the disk. The time it takes for Jupiter to grow to pebble isolation mass and the mass of material that drifts by its orbit before that is indicative for the final elemental abundances, as discussed above. The relative change in time scales when raising the fragmentation velocity allows more material to pass the planet’s orbit and accrete onto the star then for a low ufrag (see Fig. 11). In the high ufrag case, Jupiter always reaches the isolation mass. For high disk mass, it does so faster for higher viscosity until α = 1 × 10−3. For α = 3 × 10−3, the combined effect of a high pebble flux and smaller pebble sizes lead to Jupiter growing for longer before the isolation mass is reached. If the disk mass is lower, it takes Jupiter a much shorter amount of time to reach isolation for α = 3 × 10−3, but the difference between light and massive disks decreases as the viscosity decreases. For α = 1 × 10−3, the planetary core is able to reach the isolation mass earlier to block refractories if the disk is massive. The trends of isolation mass timing are directly reflected by a higher or lower final abundance of elements with refractory components, respectively. In addition, the iron abundance drops more significantly for cases where Jupiter is quick to form the pressure bump, because all ice lines of solids that contain those elements are inside the planets orbits, unlike for carbon or oxygen.

A value of ufrag = 1 m s−1 does not always allow Jupiter to grow to the pebble isolation mass during the disk’s lifetime. In fact, for α ≥ 1 × 10−3, regardless of disk mass or size, a partial gap to block solid grains is not opened. This is related to the fact that the gap opening is more difficult in disks with larger viscosity (e.g. Crida et al. 2006; Kanagawa et al. 2018; Bergez-Casalou et al. 2020). The final abundances of all elements present in refractory species for those high viscosities therefore show similar values, caused just by the growing core of Jupiter accreting part of the pebble flux. Otherwise, the final abundance is created by disk evolution effects independent of a planet (see Sect. 3.1). With dust growth being negligible for α = 3 × 10−3, no enrichment of material takes place in the disk and the abundance in the convective envelope remains at solar value. For low viscosity, Jupiter reaches pebble isolation, albeit late during the disk’s lifetime. At α = 3 × 10−4, Jupiter takes long to create a pressure bump, so the abundance of refractories in the envelope is higher in those cases.

thumbnail Fig. 12

Solar convective zone abundance changes created by accretion of a disk forming Jupiter as a function of disk viscosity α for ufrag = 1 m s−1 (blue) and ufrag = 10 m s−1 (orange). The data points show the mean value for simulations with R0 ∈ {75 AU, 150 AU} and M0 ∈ {0.05, 0.1}, and the error bars indicate the standard deviation.

5 Application to HD 106515

Finally, the aim of this section is to compare abundance differences created by planet formation and disk evolution to observations and study whether differences between wide binaries can be recreated. In particular, the wide binary system HD106515 was considered. It consists of two solar-like stars with a semi-major axis of 34547+95$345_{ - 47}^{ + 95}$ AU (Rica et al. 2017). The primary constituent HD 106515 A has a mass of 0.888 ± 0.018 M and a metallicity of [Fe/H] = +0.016 ± 0.009, while the companion HD 106515 B has a mass of 0.861 ± 0.015 M and a metallicity of [Fe/H] = +0.022 ± 0.010 (Saffe et al. 2019).

Observations of this system show that the two binary stars do not have the same elemental abundances. A clue to the origin of this difference could lie in the fact that the primary star harbors a confirmed giant planet, HD106515A b. Its mass is determined to be 18.91.4+1.5MJ$18.9_{ - 1.4}^{ + 1.5}{M_{\rm{J}}}$, and it has an observed semi-major axis of 4.48 ± 0.05 AU. In Sect. 4, we showed how the growth of a giant planet in the disk can have an impact on the accretion of solids onto the convective envelope of the disk. Therefore, we investigate whether the formation of a giant planet around one of the constituents of the HD106515 binary system is able to reproduce the observed differences in elemental abundances. For ease of comparison, similar to the considerations in Sect. 4.2, the present day abundances were applied as the initial conditions of the system. In addition to the reason presented in that section, two simulations were compared in this one, with and without the creation of the planet, such that only the difference between those two cases is relevant. Therefore, being off from the real initial condition could be neglected here for simplicity. We note that, while it seems possible to draw a straight line at Δ[X/H] = 0.001 to explain the abundance differences within the measurement uncertainties (see Fig. 13), it might be questionable why, in this case, the measurement of oxygen should be that far off, meaning why it is mostly a negative abundance. We thus think that these measurements point to an interesting scenario for planet formation and future observations will show if these features can be recovered in other binary systems.

Even though the two stars are not identical, chemical abundances of star B were used for the cases with and without a giant planet as initial conditions, and small effects created by the change of these properties between the two stars were neglected. In particular, the properties of star B were chosen as opposed to star A because the effects of the formation of a planet in the disk were considered, which has already influenced the abundances of star A. High precision measurements used to obtain the abundance differences do not give the individual stellar abundances (Liu et al. 2021), hence the measurements of the individual stars by Saffe et al. (2019) were used, with the abundance based on the Ti I line was applied for titanium. The mass of the central star in the simulations was chosen to be 0.88 M. While the stellar model that is in the closest agreement with these parameters employed that mass, it used [Fe/H] = 0. As for previous considerations, the planet’s location was kept fixed. The effect of planetary migration is shown in Appendix C. The employed approach assumes that all model parameters, apart from the presence of the planet, are identical for the accretion disks around both stars. We chose a disk radius of 75 AU and disk mass of 0.1 M. This was done for simplicity, but in reality, the disks around the stars might have been different. A scenario where the two disks are of unequal initial mass is presented in Appendix D. As discussed in Sect. 3.1, a difference in disk parameters between the two disks can have an impact on the difference of the final convective zone abundances. However, the focus here lies on the impact of planet formation in otherwise identical disks.

First, the planetary seed in system A was placed at the present-day observed location. For different sets of grain size-altering parameters, including only those that allow growth to pebble isolation, the top of Fig. 13 shows the resulting contrast between the convective zone abundances of the stars whose disks did and did not form a giant planet, respectively. It is apparent that in this case, the produced contrast is too large for the majority of elements, even in the case of the largest particles, where the impact of the forming giant planet is the smallest (see Sect. 4).

Furthermore, a key trend is found in the observations that is not reproduced with this setup: The abundance difference between stars B and A not only drops substantially compared to carbon and most refractories, but also likely negative. In fact, the present-day location of the planet corresponds to a location outside the water ice line, so that the planet blocks the accretion of water ice alongside the refractories. Therefore, the final stellar abundance of oxygen does not differ significantly from other elements in simulations where the planet is fixed at that location. The drop in abundance difference can, however, be reproduced in simulations if the planetary seed is placed such that the accretion of water ice is not affected by the pressure bump the planet creates. In that case, the difference of the oxygen abundance drops.

Reproducing Δ[O/H] < 0 is not possible if planet formation, without migration, in one of the disks is the only considered effect, and the disks have the same initial mass. The negative difference implies that the star whose disk did not form a planet experienced the accretion of less enriched material than when a planet is formed, contradicting results of previous sections. In addition, the simulated abundance differences are too large even if the planet forms inside the water ice line.

A solution to these issues for disks of the same mass could lie in the formation of plantesimals, since that also prevents the accretion of solids from the disk by locking them those objects. While a realistic model of planet formation involves forming planetesimals in the disks of both binary constituents, the abundance differences found when their formation is not considered suggests a higher efficiency of locking-up material in the disk where no planet forms. Therefore, in the interest of simplicity, planetesimal formation is only included in the simulation of the disk around star B.

Figure 14 presents the scenario of planetesimal formation in star B’s disk, where three different planetesimal formation efficiencies ϵ ∈ {0.01,0.05,0.1} were considered. Two different scenarios are shown. In the first scenario, shown on the left-hand side of Fig. 14, planetesimals could form anywhere in the disk given a high-enough pebble flux, while they could not form interior to the water ice line in the second scenario, shown on the right-hand side of Fig. 14. The second scenario is motivated by the planetesimal formation model by Drążkowska & Alibert (2017), where pebbles are too small to form planetesimals inside the water ice line. We note that this is caused by a radial change in fragmentation velocity in their model, whereas we employed a constant value. If planetesimals can form anywhere, but the pebble flux is low or planetesimal formation is inefficient, it predominantly occurs in the inner disk regions, where the dust composition shifts toward higher fractions of refractories. If planetesimal formation can only take place in the outer disk, where also volatiles are in solid form, all species are locked-up in equal fractions. Since planet formation also prohibits the accretion of refractories whose evaporation front locations are blocked equally, the second case shows a constant abundance difference per condensation temperature for these elements.

In simulations where star A’s planet forms at the present-day location, the observed oxygen abundance trend is not reproduced, regardless of the planetesimal formation scenario. If planetesimal formation can occur anywhere in the disk, the resulting abundance difference is negative or close to zero for all elements at medium and high formation efficiency (ϵ ≥ 0.05) and large particles (ufrag = 10 m s−1), not just for oxygen. In these simulations, locking-up material in planetesimals is more efficient at preventing enrichment of heavy elements than the planet’s pressure bump. In fact, the case where the largest pebbles are formed benefits locking up material due to planetesimal formation the most compared to blocking by the planet. Here, larger particles are less efficiently blocked by a pressure bump (see Sect. 4), while a large mass is locked up in planetesimals due to the high pebble flux. The pressure bump of the planet, located outside the water ice line, blocks the flux of icy pebbles, and the low condensation temperature of water does not result in significant inclusion in planetesimals. Hence, the relative importance of the two effects is reversed for large particles, so that those simulations show an upward bump for Δ[O/H], opposite to the observations. On the other hand, for the lowest efficiency ϵ = 0.01 and smaller particles, the influence of the planet is larger than that of planetesimal formation. If planetesimals do not form in the inner disk, their influence is also diminished compared to the pressure bump. Therefore, the abundance differences shift more toward positive final values. The difference in inclusion efficiency per species is reduced, as most species are solid in the outer disk. Here, only configurations of ϵ = 0.1 and ufrag = 10 m s−1 reach a negative final difference. Although the contrast between the convective zone abundances of the disk with a planet and without agrees for some configurations with the measured 1σ interval regardless of whether planetesimal formation is limited to the outer disk region or not, it fails to reproduce the observed negative difference between stars B and A for oxygen while also producing a positive difference for the other elements.

Results where such a trend is reproduced can be found if the planet is placed inside the water ice line, not hindering the accretion of oxygen in water. If planetesimal formation is not restricted, the most favorable results are produced for ϵ = 0.01 and for the largest grains. Here, both locking up material by planetesimal formation and pressure bump trapping is not very efficient, so that the contrast between the two star’s convective zone is close to zero, while the fact that the disk around star B forms planetesimals allows the contrast in [O/H] to be negative. Higher efficiencies lead to the other abundance difference values becoming negative as well. If planetesimal formation occurs only in the outer disk, a higher efficiency parameter is required to match observations. A medium efficiency of ϵ = 0.05 is found to be in good agreement with the measurements for the largest grains. The value for Δ[O/H] is more negative in this case, while the differences for all other elements is positive. This is because restricting the planetesimal formation to the outer disk means that fewer refractories are locked up in planetesimals overall, moving their abundance difference toward positive values, while the high formation efficiency and pebble flux in the outer disk create a negative contrast for oxygen, not affected by the planet. Disks with ϵ = 0.1 and the smallest grain size can agree with the measured differences, too.

The simulations presented in Fig. 14 show disks with M0 = 0.1 and R0 = 75 AU. A reduced disk mass moves the final abundance difference closer to zero overall, allowing more configurations, mainly those featuring lower formation efficiencies, to match the observations. A larger disk radius of 150 AU is not investigated here, because such large disks are close to overlapping, given the binary semi-major axis of 34547+95AU$345_{ - 47}^{ + 95}{\rm{AU}}$.

The results found in this section show that it is possible for planet formation to be the cause for the elemental abundance difference measured for the HD 106515 wide binary system by influencing the material that is accreted onto the convective envelope of the central star. However, to reach appropriate levels for the differences, especially Δ[O/H] < 0, additional requirements have to be met. One possibility is including the formation of planetesimals. In that case, the observations are reproduced only for cases where the disk that does not form a planet (around star B) forms planetesimals and the disk around star A does not, representing the physical case where planetesimal formation is more efficient in the disk that was unable to form a planet. This could possibly indicate that efficient planetesimal formation is not beneficial for giant planet formation. Alternatively, the disk that does not form a planet has to be lighter than the one that does (see Appendix D).

While the planet at its present-day location, outside the water ice line, can match the observations in some scenarios, the fact that the oxygen abundance difference drops significantly compared to the other elements cannot be reproduced this way. For that, the planet would need to be formed inside the water ice line. In this scenario, simulations with equal-mass disks and various disk and planetesimal formation parameters, as well as a simulation with disks of unequal mass, reproduce the observations for most elements.

thumbnail Fig. 13

Abundance differences between HD106515 B and A caused by planet formation taking place only around star A. The top shows the case of the planet forming where it is observed today at 4.5 AU, which lies outside the water ice line, while the bottom considers the planet forming inside the water ice line, just outside rFe3O4${r_{{\rm{F}}{{\rm{e}}_3}{{\rm{O}}_4}}}$. The gray area delimited by black lines corresponds to the measurements including the 1σ error (Liu et al. 2021). In blue, the disks around the constituents are parameterized by α = 1 × 10−3, ufrag = 10 m s−1. The red and green lines describe disks with a viscosity of α = 1 × 10−4 with ufrag = 1 m s−1 and ufrag = 10 m s−1, respectively.

thumbnail Fig. 14

Same as Fig. 13, but including planetesimal formation around star B. For simulations shown on the left-hand side, plantesimals can form in any disk region where the pebble flux is above a threshold value, in accordance with the model by Voelkel et al. (2020). For those on the right-hand side, planetesimal formation is turned off for r<rH2O$r < {r_{{{\rm{H}}_2}{\rm{O}}}}$. The planetesimal formation efficiency ϵ is varied between lines with different markers and styles. Sold lines with a circular marker show the results for ϵ = 0.05, dashed lines with a cross marker show ϵ = 0.01 and dotted lines with triangular markers show ϵ = 0.1.

6 Discussion

In this work, several simplifications were made to balance the required computational power and the level of physical detail the applied model can capture. Treatment of the evolution of the stellar convective zone was done by using stellar evolution models treating young stars in isolation. It was then assumed that the resulting evolution of the convective zone is applicable for calculating the change of abundances due to accretion of the disk. To compute that change, instantaneous mixing was assumed. Moreover, it was presumed that the convective envelope keeps its composition while transitioning away from a fully convective star. Regarding the composition, Kunitomo et al. (2018) show that accurate abundances are reproduced under such simplifications. However, this approach falls short of considering that the accretion has an impact on the stellar evolution itself, hence changing the mass evolution of the convective envelope (Baraffe & Chabrier 2010) and affecting how disk accretion is reflected on the elemental abundances. While Kunitomo et al. (2017, 2018) show that accretion complicates the PMS evolution of the star, the implication for this work is that the stellar interior is not necessarily fully convective even at early times, indicating that abundance differences presented here could be higher, emphasizing the importance of connecting disk and stellar evolution.

6.1 Grain size model

The model used to compute the Stokes number of the dust in the disk, in turn setting the radial drift speed and the efficiency of interactions with an accreting planet, uses an approach treating only two distinct populations with representative sizes, using fits to find good agreement with simulations that treat the growth and fragmentation of dust in detail (Birnstiel et al. 2010). A full treatment of dust is computationally involved (e.g., Stammler & Birnstiel 2022) and is not feasible to be included in this model, which also considers planet formation and the evaporation and condensation of material. However, the downside of using a mass-averaged velocity is that differently-sized grains cannot separate in the disk. In reality, the small dust-grain population is not affected by a planet at pebble isolation strongly enough to block it (Ataiee et al. 2018; Bitsch et al. 2018b; Picogna et al. 2018), allowing part of the solids outside the planet’s orbit to pass it even at later times. Stammler et al. (2023) find significant “leaking” of a planet-induced dust trap when grain fragmentation plays an important role, with 80% of the dust outside the planetary orbit being accreted through the gap after 10 Myr for a 200 M planet, depending on parameter choices. This is not reproduced in this model, where it is only possible for the entire dust population to avoid being blocked if the resulting mass-averaged outward-pointing velocity is smaller than the inward-pointed gas drag. Furthermore, the 1D-approach allows the trapping of all pebbles if the planet is massive enough, which does not match a more realistic 3D scenario, where small pebbles could be transported through the planetary orbit via meridional flows enhanced by vertical shear (e.g., Picogna et al. 2018). Altogether, the effects of pebble trapping by a massive planet would be diminished if such a more detailed treatment was included in the model.

6.2 Disk temperature and mass loss

The location of the ice lines in the disk play a crucial role for the strong enrichment of refractories early on in the simulations and the late accretion of material enriched in volatiles. For simplicity, the model used here employs a temperature profile constant in time. However, a more realistic scenario needs to treat changes of temperature by a change in the stellar luminosity, the surface density of the disk, the dust-to-gas-ratio change due to dust radial drift and radiative cooling (e.g., Bitsch et al. 2015; Savvidou et al. 2020; Savvidou & Bitsch 2021). Including such condensations could shift the ice lines in the disk based on the evolutionary stage of the system. Ultimately, the details of the temperature structure of the disk can be relevant to find a realistic placement of the ice lines over time, especially relative to a massive planet in the disk. If the ice line is, for example, outside a planet’s orbit at first but moves inside it as the disk cools, a stellar abundance on a level between the case without any planet and one with full trapping can be achieved.

Stellar irradiation does not only heat the disk, but photoevaporation creates steady mass-loss. In the model used here, the mass-loss is only considered in the form of rapid depletion once the disk lifetime is reached. However, the dust surface density is not strongly affected by the removal of disk mass, as large grains settle toward the mid-plane. Only smaller particles coupled to the gas could be carried away by photoevaporation winds, making the dust mass loss negligible (e.g., Owen et al. 2011). However, photoevaporation can carve a gap in the disk, thereby creating a pressure bump that can block pebbles (e.g., Ercolano & Picogna 2022). The difference to a gap formed by a giant planet is that a gap carved by photoevaporation also prevents further gas accretion, preventing a further change of convective envelope abundance. Disk winds, caused by photoevaporation, hydrodynamical or magneto-hydrodynamical effects have been found to influence disk evolution greatly by allowing for strong angular momentum transport (e.g., Chambers 2019; Lesur et al. 2022). While a full-on treatment of disk winds is computationally unfeasible for this work, a change in the accretion speed due to more efficient angular momentum transport without stronger turbulence would affect the timing of accretion, diminishing enrichment of the stellar envelope if it is still in a stage where it is massive. A disk-wind driven scenario would require little turbulence in the mid-plane, in turn allowing pebbles to grow to large sizes. If the disk winds are strong enough, outward migration is possible even for single giant planets (Kimmig et al. 2020), giving rise to various migration pathways crossing ice lines in both directions. Furthermore, the treatment of magnetic fields would constitute a more complex magneto-hydrodynamical model and can affect the accretion of the ferromagnetic iron.

6.3 Giant planet formation

The timing of the formation of giant planets that are massive enough to have an impact on the stellar envelope abundances is an uncertainty with this approach. The planet needs to form early enough to have an impact by trapping remaining material in the disk, for which it needs to reach the pebble isolation mass quickly before most material is accreted onto the star. While this is typically a requirement for giant planets to form their massive core, the question remains about how early the formation occurs relative to start of the PMS evolution, which in turn defines the impact disk accretion has by setting the mass of the convective zone. If planet formation generally occurs early on during the PMS evolution, the impact is hindered by a more massive convective envelope. On the other hand, for super-Earths, the efficiency of trapping material is likely small, as they might not form early enough and not grow massive enough to block a large pebble mass. As was shown in Sect. 3.1, the disk lifetime is indicative for the final elemental abundances, too. Refractory abundances typically reach their maximum value within the first ~2 Myr, while a longer lived disk allows for the enrichment of volatile material at the expense of diluting the refractory abundance.

A correlation between the occurrence rate of giant planets and the metallicity of stars around which they occur was observed previously (Johnson et al. 2010; Santos et al. 2004; Fischer & Valenti 2005), where more metal rich stars are more likely to harbor giant planets. Since giant planets trap refractory pebbles, stars where they form experience less enrichment in iron by disk accretion than stars without giant planets. Therefore, the measured correlation could be stronger if the initial metallicity was considered instead of present-day values.

6.4 Abundance difference trend with condensation temperature

Observations of the HD 106515 system show no clear trend of the elemental abundance difference with condensation temperature. However, such trends are found for other systems with no known planets (Liu et al. 2021) and commonly discussed in the context of the solar refractory depletion (e.g., Meléndez et al. 2009). In this work, most of the observed elements with high condensation temperature are not included in the chemical model, and a trend with condensation temperature was not investigated. This is because the evaporation fronts of the refractory molecules lie close together and close to the star given the temperature profile at hand, so that no significant difference can be created, neither by a difference in accretion timing nor by planet-induced trapping. Furthermore, the species with the highest condensation temperature in the model is TiO (see Table 1), and the aluminum-containing Al2O3 has a higher condensation temperature than VO, which is the only molecule containing vanadium here. Observations show Δ[Al/H] < Δ[Ti/H] < Δ[V/H]. This contradicts the general condensation temperature trend and the ordering based on molecular condensation temperature we use here. Using this approach, those difference therefore cannot be reproduced accurately. Even in cases where the trend matches the order of the partitioning model, differences between refractories would require a hotter disk so that the evaporation fronts are farther out and have greater spacing, allowing accretion timing or a pressure bump to influence refractories with higher condensation temperature differently from those with a lower one.

Generally, the technique we applied to find signatures of planet formation around a star and obtain information about the formation location can be applied to systems other than HD106515, too. However, it is not possible to find signatures of planet formation around field stars this way, because a comparison to a star who did not form a planet in its disk and had the same initial abundance is necessary. GALAH+ data show that Galactic variations of stellar initial conditions are considerably stronger than differences induced by planet formation, making statistical approaches not viable. This is also reflected by the error bars in Fig. 6. Wide binary systems are favorable, as close binaries exhibit interactions whose treatment would warrant for a more complex model of planet formation and disk evolution. Including the formation of multiple planets in the simulations would allow, for example, the investigation of the HD 133131 system, where the one constituent harbors one, and the other harbors two planets, or the abundance pattern of TW Hya, where a multigap model was proposed to explain observations (McClure et al. 2020).

7 Summary

In this work, we considered the growth and drift of millimeter-sized dust grains in protoplanetary disks, separated into chemical species using a partitioning model, to find the change of the elemental abundances in the shrinking convective envelope of the central star due to the accretion of disk material. In addition, we investigated the impact of a growing giant planet in the disk and the order of magnitude of the expected abundance change for the Sun. We now summarize our findings as follows:

  1. For large dust grains, the efficient radial drift creates a strong enrichment of refractories in the accretion steam at early times, in turn raising the respective elemental abundances. At that point in time, the convective zone of the star is still very massive, so that the accreted material’s composition is imprinted only slowly onto it. After the pebble reservoir of the outer disk is depleted, the material left in the disk is depleted in refractory species. Now, the disk is accreted onto a shrinking envelope which adapts to the disk accretion flux, only enriched in volatiles, more quickly. Due to the absence of refractories at that stage, the late accretion increases the volatile abundance while simultaneously diluting the refractory enrichment from the early accretion phase. The relative significance of this process is directly related to the disk lifetime, where a disk that is short-lived leaves behind a star whose convective zone is more enriched in refractories and less enriched in volatiles than if the disk was long-lived;

  2. The abundance changes of the stellar convective zone caused by disk accretion can reach up to ~5 × 10−2 for refractories and volatiles, where the highest level enrichment of refractories is reached if particles in the disk are large, so for low viscosity or high fragmentation velocity. For volatiles, the highest enrichment is reached if the disk viscosity is high. There are other disk and stellar parameters that affect the final abundance difference, too. Stellar mass and metallicity affect the mass evolution of the convective envelope, in turn changing the relative significance of early and late accretion. More massive disks can create larger enrichment when accreted. Smaller disks reduce the enrichment in most cases by reducing dust growth and drift timescales, shifting the period of enriched accretion toward times when the convective envelope is more massive and adaptation to the composition of the accreted material is slow. Overall, the elemental abundances of the stellar convective envelope after the disk has dissipated differ from the initial ones present in the natal molecular cloud, and in turn initially in the star and the disk;

  3. The formation of a massive planet changes the material accreted onto the convective envelope after reaching the pebble isolation mass by trapping solids outside its orbit. The effect of the growth of a planet in the disk is defined by the time it takes the initial protoplanetary seed to grow to this mass limit, and the amount of material that is still left outside the planet’s orbit at that point. Moreover, the composition of the blocked solids is related to the position of the planet relative to the evaporation fronts. Molecular species with ice lines outside the planet’s orbit are not affected by its presence, having already evaporated before they pass its location. Also, dust grains need to be of a certain size to be trapped, as the strong coupling to the gas and the resulting gas drag of small particles helps them overcome the pressure bump that is created outside the planet’s location;

  4. As water is the most massive heavy component in the disk, the oxygen abundance in the convective envelope is particularly strongly affected by the formation of a planet outside or inside the water ice line, creating changes in the stellar envelope [O/H] between those two locations of ~10−2. Given observations of this precision, the elemental abundances of stars can be used to constrain the formation location relative to ice lines of planets that formed in their nebula. Refractory material is affected by the formation of a planet in most cases, as their evaporation fronts are close to the star and likely to be inside the planets orbit in giant planet formation scenarios;

  5. The effect of disk accretion on the Sun’s elemental abundances was investigated with the simplification of only considering Jupiter as the most significant gravitational influence in the disk. After 4.5 Myr, the change in refractory abundances can reach ~2 × 10−2 for massive disks with high fragmentation velocity and low viscosity, but can be as low as ~10−3 at α = 1 × 10−3 with low fragmentation velocity and disk mass, and even lower at α > 1 × 10−3. In contrast, such highly viscous disks can enrich the volatile nitrogen up to [N/H] ~ 1 × 10−2. Depending on the parameters of the solar protoplanetary disk, the initial abundances of the Sun therefore differ in an observationally discernible matter from today’s values;

  6. Observations have revealed chemical abundance differences between stars in the HD106515 wide binary system. One of its constituents, HD106515A, harbors a massive giant planet at ~4.5 AU, while the other has no confirmed planets. A key aspect is that Δ[O/H] < 0, while all other Δ[X/H] > 0. This is best reproduced if the planet, despite its present-day location, formed inside the water ice line, and if either planetesimals are formed more efficiently around star B, which does not harbor a planet, or the disk around star B was lighter. For the former, the requirement of more efficient planetesimal formation in disk B compared to A, where no planet has been observed as of today, might provide hints about the role planetesimals can play in the broader picture of planet formation. This approach shows that stellar abundances can give information about the formation history of planetary systems similar to what measurements of the atmosphere of the planets could provide, such as constraining the position of the planet relative to the water ice line or the efficiency of planetesimal formation. A combination of the constraints obtained by a combination of both types of observations in a single system could allow for further insights in the future.

Our model presented here opens up new opportunities to constrain planet formation theories via stellar abundances and emphasizes the need for detailed stellar abundances of exoplanet host stars.

Acknowledgements

L.-A.H. acknowledges financial support from the European Research Council via the ERC Synergy Grant ECOGAL (grant 855130). B.B. thanks the European Research Council (ERC Starting Grant 757448-PAMDORA) for their financial support.

Appendix A Dust drift model

The chemcomp model, in addition, accounts for grain density (ρs) changes that arise as volatile and refractory elements condense and evaporate (see Sect. 2.3) by dynamically calculating it during each time step following Drążkowska & Alibert (2017), ρs=(mref+mvol)(mrefρref+mvolρvol)1,${\rho _s} = \left( {{m_{{\rm{ref}}}} + {m_{{\rm{vol}}}}} \right)\,{\left( {{{{m_{{\rm{ref}}}}} \over {{\rho _{{\rm{ref}}}}}} + {{{m_{{\rm{vol}}}}} \over {{\rho _{{\rm{vol}}}}}}} \right)^{ - 1}},$(A.1)

with ρref = 3 g cm−1 and ρvol = 1 g cm−1, as well as mref and mvol the mass of refractory and volatile material in solid form, respectively. The grains are subject to radial drift, gas drag and turbulent mixing. The relative significance of these effects depends on the coupling strength of the dust to the gas, described by the Stokes number St. Additionally, the grains grow in size as they are subject to sticking collisions, increasing their Stokes number. As the gas particles move with sub-Keplerian speed due to the pressure gradient in the disk, the dust grains experience the headwind effect so that their radial drift speed exceeds that of the gas. It is given by uZ=11+St2ugas2St1+StΔυ,${u_Z} = {1 \over {1 + {\rm{S}}{{\rm{t}}^2}}}{u_{{\rm{gas}}}} - {2 \over {{\rm{S}}{{\rm{t}}^{ - 1}} + {\rm{St}}}}{\rm{\Delta }}\upsilon ,$(A.2)

where ugas is the radial velocity of the gas and Δυ is the difference of the azimuthal gas velocity to the Keplerian value. As long as St < 1, the difference to the gas velocity becomes very large as the particles increase in size. Grain growth is subject to size limits, with the dominant growth-limiting factors being fragmentation and radial drift (Dominik & Dullemond 2008; Brauer et al. 2008).

Treating the former, the Stokes number at which relative velocities between similar-sized grains equal the fragmentation velocity ufrag sets an upper limit for the particle sizes (Birnstiel et al. 2009), Stfrag=ff13ufrag2αcs2,${\rm{S}}{{\rm{t}}_{{\rm{frag}}}} = {f_f}{1 \over 3}{{u_{{\rm{frag}}}^2} \over {\alpha c_s^2}},$(A.3)

where α is the turbulent viscosity parameter and cs is the sound speed. Lab experiments reveal ufrag values between ~1 m s−1 for silicates (Blum & Wurm 2008) and ~10 m s−1 for icy grains (Gundlach & Blum 2015). Both values were considered in this work, but a radial change of fragmentation velocity based on solid composition was not considered. Analytical models of the grain size distribution show that a significant fraction of the dust has sizes slightly below the limiting Stokes number in the fragmentation limited regime (Birnstiel et al. 2011). Therefore, an additional model parameter ff is introduced in Eq. A.3, which is set to ff = 0.37.

Furthermore, the size of particles at a given radial distance from the central star is limited by the drift of the particles. As they grow in size, their drift speed increases, making them drift away from their original location, thereby limiting the maximum particle size at that location. It is characterized by a maximum Stokes number, Stdrift=fdϵdυK2cs2γ1,${\rm{S}}{{\rm{t}}_{{\rm{drift}}}} = {f_d}{_d}{{\upsilon _K^2} \over {c_s^2}}{\gamma ^{ - 1}},$(A.4)

with ϵd=ΣzΣgas${_d} = {{{{\rm{\Sigma }}_z}} \over {{{\rm{\Sigma }}_{gas}}}}$ the vertically integrated solid-to-gas mass ratio, ΣZ the solid surface density, γ=| dlnPdlnr |$\gamma = \left| {{{d\,\ln \,P} \over {d\,\ln \,r}}} \right|$ the absolute value of the power-law index of the gas pressure profile, and υK the Keplerian velocity. Analogously to Eq. A.3, a model parameter fd = 0.55 is introduced to fit more detailed models.

The small population, being one of the two components of the model, is always characterized by the size of the smallest monomers (a0), being subject to gas drag while exhibiting significant coupling to the gas. Lastly, to find the time and location dependent representative size of the large population, delayed drift has to be accounted for, as grains first have to grow to large enough sizes to be influenced by drift, growing slower in outer regions of the disk. The size of the large population is hence given by the minimum size set by the fragmentation, drift, and growth limit.

Appendix B Planet growth and pebble isolation model

The growth of one giant planet can be modeled using chemcomp, starting from a planetary seed and growing through pebble accretion (for a review, see Johansen & Lambrechts 2017), which allows fast growth of massive planetary cores (Ormel & Klahr 2010; Johansen & Lacerda 2010; Lambrechts & Johansen 2012; Morbidelli & Nesvorny 2012), up to the pebble isolation mass (Morbidelli & Nesvorny 2012; Lambrechts et al. 2014; Ataiee et al. 2018; Bitsch et al. 2018a), and subsequent envelope gas accretion (e.g., Ikoma et al. 2000; Machida et al. 2010; Ndugu et al. 2021).

In the context of a planet’s influence on the composition of material accreted onto the stellar convective envelope, the timescale of core growth by pebble accretion is responsible for allowing it to reach the pebble isolation mass. Therefore, it describes at what point during the disk’s lifetime the planet creates a pressure bump strong enough to alter the dust evolution. Initially, the planetary seed is placed in the disk at the transition mass where pebble accretion becomes efficient (Lambrechts & Johansen 2012; Johansen & Lambrechts 2017). During the subsequent time evolution, accretion rates from Johansen & Lambrechts (2017) are applied. They strongly depend on the Stokes number of the pebbles, and are the largest close to St = 0.1. To account for the evaporation of pebbles during accretion, which contributes material to the planetary gaseous envelope, only 90% of the accreted pebble mass was used for the core growth, while the other 10% were added to the envelope.

As the planetary core grows in mass, gravitational interactions with the disk lead to a partial gap being opened, which has an impact on the movement of dust grains (Paardekooper & Mellema 2006; Rice et al. 2006). If the core grows massive enough, the accretion of pebbles is halted, preventing grains from drifting to further-in regions of the disk and reaching the central star, and stops the core from growing further. The limiting mass for this effect to occur is the pebble-isolation mass. For its prescription, findings from Bitsch et al. (2018b) are employed, Miso=25ffitM,${M_{{\rm{iso}}}} = 25{f_{{\rm{fit}}}}{{\rm{M}}_ \oplus },$(B.1)

with ffit=[ Hgas/r0.05 ]3[ 0.34(3log(α))4+0.66 ].${f_{{\rm{fit}}}} = {\left[ {{{{{{H_{{\rm{gas}}}}} \mathord{\left/ {\vphantom {{{H_{{\rm{gas}}}}} r}} \right. \kern-\nulldelimiterspace} r}} \over {0.05}}} \right]^3}\,\left[ {0.34\,{{\left( {{{ - 3} \over {\log \left( \alpha \right)}}} \right)}^4} + 0.66} \right].$(B.2)

Due to the pressure gradient changes caused by pile-up at the ice lines, the dependence of the pebble isolation mass on the pressure gradient was not included for simplicity. The formation of a pressure bump and subsequent opening of a gap due to gravitational interaction of the planet with the disk is modeled by a change in viscosity, where the gap profile is described as a Gaussian around the planet’s location.

Appendix C Planet migration for HD106515A b

When treating planet formation around HD106515A without considering planetary migration, the negative abundance difference of oxygen between star B and A is reproduced only if star B, without any confirmed planets, forms planetesimals efficiently. In that scenario, it is possible for star A to accrete material more enriched in oxygen than star B, despite forming a planet in the disk. However, a second scenario is possible if the planet migrates. Here, a negative oxygen abundance can be created if the planet crosses the water ice line well into the disk lifetime. Initially, the position of the planet was chosen such that its pressure bump prevents the accretion of water. As it migrates inward, it eventually crosses the water ice line, allowing the icy pebbles that follow the planets migration path in the pressure bump to evaporate and be accreted onto the star. At this time, the stellar convective envelope has already shrunk, resulting in a faster adaptation of its composition to that of the accreted material.

thumbnail Fig. C.1

Migration track of HD106515A b for ap,0 = 2.5 AU, α = 1 × 10−4 and ufrag = 1 m s−1, shown in blue. The black dashed lines indicate ice lines that are crossed during the inward movement.

thumbnail Fig. C.2

Same as Fig. 13, but including migration of the planet around star A, with an initial position of ap,0 = 2.5 AU.

In principle, due to the reduction in envelope mass, it is possible for the accretion of the previously blocked water to create an enrichment stronger than in a case where water was never blocked. However, this is very sensitive to the initial planetary location. If the ice line is crossed early, the blocked material does not give a significant contribution compared to a case where the initial position is inside the water ice line. In addition, the convective envelope is still close to its original mass at that time. On the other hand, if the ice line is crossed late, the oxygen-enriched gas is not accreted for a substantial amount of time before the disk dissipates. Therefore, even though the envelope has shrunk in mass at this time, a significant enrichment is not possible.

Due to this sensitivity to the initial position, only configurations with 2 AU ≲ ap ≲ 2.5 AU produce Δ[O/H] < 0. In Fig. C.1, the migration track of a planet with aP,0 = 2.5 AU is shown, for α = 1 × 10−4 and ufrag = 1 m s−1. It crosses the water ice line at t ~ 2 Myr, which is also the time when the convective envelope starts shrinking. At t ~ 4 Myr, the Fe3O4 evaporation front is crossed, which allows a fraction of iron to be accreted onto a stellar envelope which is lighter due to its evolution.

In Fig. C.2, the final abundance difference are shown, including planetary migration. The best fit to the observations is produced by the α = 1 × 10−4, ufrag = 1 m s−1. As in the case of a fixed planet and planetesimal formation (see Sect. 5), individual trends for the highly refractory elements are not reproduced. In this migrating scenario, no substantial abundance differences are produced for those elements, because they are released and accreted onto the star at the end of the simulation, when the planet has reached the inner edge. In the case of ap = 2 AU, the highly refractory abundance differences even reach Δ[X/H] ~ −0.01, which contradicts observations. As the planet ends up at the inner edge of the disk, significant outward migration would be required to commence to find the planet at its present day location of 4.5 AU. This problem is similar to the one of the Sect. 5 scenario, where the planet needs to have formed inside the water ice line, corresponding to a formation inside 4.5 AU. While in the former case, a solution could lie in the inward movement of the ice line due to disk cooling, outward movement of the ice line would be required to obtain results corresponding to the migration scenario while matching the observed planet location. Therefore, even though such outward movement can occur during the first few 100 kyr of disk evolution (e.g., Drążkowska & Dullemond 2018), the migration scenario might be less realistic.

Appendix D HD106515 model with different disk masses

While the investigation of every possible variation of initial disk parameters between the binary constituents of HD106515 is outside the scope of this work, an interesting scenario is a difference in initial disk mass. The case where M0,A < M0,B can be excluded, because the formation of a ~19 MJ planet is difficult to explain in less massive disks with current formation models (e.g., Savvidou & Bitsch 2023). Additionally, reproducing Δ[O/H] < 0 would be even more difficult in this case.

On the other hand, the case of M0,A > M0,B is similar to the scenario where the disk around star B formed planetesimals more efficiently than the disk around star A. Instead of planetesimal formation removing material from the disk (see Fig. 14), it has less mass initially. The main difference is that a change of M0 affects the abundance of all chemical species equally, while planetesimal formation only affects materials that are predominantly in a solid phase. Because most species are in a solid phase in the outer disk, the described scenario is therefore most comparable to the bottom right scenario of Fig. 14. The resulting abundance differences are shown in Fig. D.1. The best fit to the observations is produced when using the disk parameters α = 10−4, ufrag = 10 m s−1. In particular, Δ[O/H] < 0 is reproduced, and the abundance differences are within the measurement errors except for Na, Ti and V. In principle, further lowering M0,B could allow simulations with higher viscosities or lower fragmentation velocities to match the observations as well. However, further investigation would increase the parameters space, which is outside the scope of this work and would not change the conclusion that a lower disk mass around constituent B could help explain the observed abundance differences.

thumbnail Fig. D.1

Same as Fig. 13, showing only the case where the planet formed inside the water ice line. Here, M0,A = 0.1 and M0,B = 0.05

References

  1. Alexander, R., Pascucci, I., Andrews, S., Armitage, P., & Cieza, L. 2014, in Protostars and Planets VI, eds. H. Beuther, R. S. Klessen, C. P. Dullemond, & T. Henning, 475 [Google Scholar]
  2. Andrews, S. M., Huang, J., Pérez, L. M., et al. 2018, ApJ, 869, L41 [NASA ADS] [CrossRef] [Google Scholar]
  3. Asplund, M., Grevesse, N., Sauval, A. J., & Scott, P. 2009, ARA&A, 47, 481 [NASA ADS] [CrossRef] [Google Scholar]
  4. Ataiee, S., Baruteau, C., Alibert, Y., & Benz, W. 2018, A&A, 615, A110 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  5. Baraffe, I., & Chabrier, G. 2010, A&A, 521, A44 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  6. Bergez-Casalou, C., Bitsch, B., Pierens, A., Crida, A., & Raymond, S. N. 2020, A&A, 643, A133 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  7. Birnstiel, T., Dullemond, C. P., & Brauer, F. 2009, A&A, 503, L5 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  8. Birnstiel, T., Dullemond, C. P., & Brauer, F. 2010, A&A, 513, A79 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  9. Birnstiel, T., Ormel, C. W., & Dullemond, C. P. 2011, A&A, 525, A11 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  10. Birnstiel, T., Klahr, H., & Ercolano, B. 2012, A&A, 539, A148 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  11. Bitsch, B., & Battistini, C. 2020, A&A, 633, A10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  12. Bitsch, B., Johansen, A., Lambrechts, M., & Morbidelli, A. 2015, A&A, 575, A28 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  13. Bitsch, B., Forsberg, R., Liu, F., & Johansen, A. 2018a, MNRAS, 479, 3690 [NASA ADS] [CrossRef] [Google Scholar]
  14. Bitsch, B., Morbidelli, A., Johansen, A., et al. 2018b, A&A, 612, A30 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  15. Blum, J., & Wurm, G. 2008, ARA&A, 46, 21 [Google Scholar]
  16. Booth, R. A., & Ilee, J. D. 2019, MNRAS, 487, 3998 [CrossRef] [Google Scholar]
  17. Booth, R. A., & Owen, J. E. 2020, MNRAS, 493, 5079 [NASA ADS] [CrossRef] [Google Scholar]
  18. Brauer, F., Dullemond, C. P., & Henning, T. 2008, A&A, 480, 859 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  19. Buder, S., Sharma, S., Kos, J., et al. 2021, MNRAS, 506, 150 [NASA ADS] [CrossRef] [Google Scholar]
  20. Burbidge, E. M., Burbidge, G. R., Fowler, W. A., & Hoyle, F. 1957, Rev. Mod. Phys., 29, 547 [NASA ADS] [CrossRef] [Google Scholar]
  21. Cabral, N., Guilbert-Lepoutre, A., Bitsch, B., Lagarde, N., & Diakite, S. 2023, A&A, 673, A117 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  22. Chambers, J. E. 2010, ApJ, 724, 92 [NASA ADS] [CrossRef] [Google Scholar]
  23. Chambers, J. 2019, ApJ, 879, 98 [NASA ADS] [CrossRef] [Google Scholar]
  24. Chen, Y. Q., Nissen, P. E., Zhao, G., & Asplund, M. 2002, A&A, 390, 225 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  25. Chiang, E. I., & Goldreich, P. 1997, ApJ, 490, 368 [Google Scholar]
  26. Crida, A., Morbidelli, A., & Masset, F. 2006, Icarus, 181, 587 [Google Scholar]
  27. Dominik, C., & Dullemond, C. P. 2008, A&A, 491, 663 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  28. Drążkowska, J., & Alibert, Y. 2017, A&A, 608, A92 [Google Scholar]
  29. Drążkowska, J., & Dullemond, C. P. 2018, A&A, 614, A62 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  30. Dullemond, C. P., & Dominik, C. 2004, A&A, 417, 159 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  31. Eistrup, C., & Henning, T. 2022, A&A, 667, A160 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  32. Ercolano, B., & Picogna, G. 2022, Eur. Phys. J. Plus, 137, 1357 [NASA ADS] [CrossRef] [Google Scholar]
  33. Fischer, D. A., & Valenti, J. 2005, ApJ, 622, 1102 [NASA ADS] [CrossRef] [Google Scholar]
  34. Gundlach, B., & Blum, J. 2015, ApJ, 798, 34 [Google Scholar]
  35. Hoppe, R., Bergemann, M., Bitsch, B., & Serenelli, A. 2020, A&A, 641, A73 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  36. Ikoma, M., NakaZawa, K., & Emori, H. 2000, ApJ, 537, 1013 [NASA ADS] [CrossRef] [Google Scholar]
  37. Jermyn, A. S., & Kama, M. 2018, MNRAS, 476, 4418 [Google Scholar]
  38. Johansen, A., & Lacerda, P. 2010, MNRAS, 404, 475 [NASA ADS] [Google Scholar]
  39. Johansen, A., & Lambrechts, M. 2017, Annu. Rev. Earth Planet. Sci., 45, 359 [Google Scholar]
  40. Johansen, A., Mac Low, M.-M., Lacerda, P., & Bizzarro, M. 2015, Sci. Adv., 1, 1500109 [CrossRef] [Google Scholar]
  41. Johnson, J. A., Aller, K. M., Howard, A. W., & Crepp, J. R. 2010, PASP, 122, 905 [Google Scholar]
  42. Johnson, B. C., Walsh, K. J., Minton, D. A., Krot, A. N., & Levison, H. F. 2016, Sci. Adv., 2, e1601658 [NASA ADS] [CrossRef] [Google Scholar]
  43. Kama, M., Folsom, C. P., & Pinilla, P. 2015, A&A, 582, A10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  44. Kanagawa, K. D., Tanaka, H., & Szuszkiewicz, E. 2018, ApJ, 861, 140 [Google Scholar]
  45. Kenyon, S. J., & Hartmann, L. 1995, ApJ, 101, 117 [Google Scholar]
  46. Kimmig, C. N., Dullemond, C. P., & Kley, W. 2020, A&A, 633, A4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  47. Klahr, H., & Schreiber, A. 2020, ApJ, 901, 54 [NASA ADS] [CrossRef] [Google Scholar]
  48. Kleine, T., Budde, G., Burkhardt, C., et al. 2020, SSR, 216, 55 [NASA ADS] [Google Scholar]
  49. Koepferl, C. M., Ercolano, B., Dale, J., et al. 2013, MNRAS, 428, 3327 [Google Scholar]
  50. Kouwenhoven, M. B. N., Goodwin, S. P., Parker, R. J., et al. 2010, MNRAS, 404, 1835 [NASA ADS] [Google Scholar]
  51. Kunitomo, M., Guillot, T., Takeuchi, T., & Ida, S. 2017, A&A, 599, A49 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  52. Kunitomo, M., Guillot, T., Ida, S., & Takeuchi, T. 2018, A&A, 618, A132 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  53. Lambrechts, M., & Johansen, A. 2012, A&A, 544, A32 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  54. Lambrechts, M., Johansen, A., & Morbidelli, A. 2014, A&A, 572, A35 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  55. Lenz, C. T., Klahr, H., & Birnstiel, T. 2019, ApJ, 874, 36 [NASA ADS] [CrossRef] [Google Scholar]
  56. Lesur, G., Ercolano, B., Flock, M., et al. 2022, ArXiv e-prints [arXiv:2203.09821] [Google Scholar]
  57. Li, Y., Brandt, T. D., Brandt, G. M., et al. 2021, AJ, 162, 266 [NASA ADS] [CrossRef] [Google Scholar]
  58. Lin, D. N. C., & PapaloiZou, J. 1986, ApJ, 307, 395 [NASA ADS] [CrossRef] [Google Scholar]
  59. Liu, F., Bitsch, B., Asplund, M., et al. 2021, MNRAS, 508, 1227 [NASA ADS] [CrossRef] [Google Scholar]
  60. Lodders, K. 2003, ApJ, 591, 1220 [Google Scholar]
  61. Machida, M. N., Kokubo, E., Inutsuka, S.-I., & Matsumoto, T. 2010, MNRAS, 405, 1227 [Google Scholar]
  62. McClure, M. K. 2019, A&A, 632, A32 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  63. McClure, M. K., Dominik, C., & Kama, M. 2020, A&A, 642, A15 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  64. Meléndez, J., Asplund, M., Gustafsson, B., & Yong, D. 2009, ApJ, 704, L66 [Google Scholar]
  65. Morbidelli, A., & Nesvorny, D. 2012, A&A, 546, A18 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  66. Ndugu, N., Bitsch, B., Morbidelli, A., Crida, A., & Jurua, E. 2021, MNRAS, 501, 2017 [Google Scholar]
  67. Ormel, C. W., & Klahr, H. H. 2010, A&A, 520, A43 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  68. Owen, J. E., & Clarke, C. J. 2012, MNRAS, 426, L96 [NASA ADS] [CrossRef] [Google Scholar]
  69. Owen, J. E., Ercolano, B., & Clarke, C. J. 2011, MNRAS, 411, 1104 [NASA ADS] [CrossRef] [Google Scholar]
  70. Paardekooper, S. J., & Mellema, G. 2006, A&A, 453, 1129 [CrossRef] [EDP Sciences] [Google Scholar]
  71. Picogna, G., Stoll, M. H. R., & Kley, W. 2018, A&A, 616, A116 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  72. Pinilla, P., Benisty, M., & Birnstiel, T. 2012, A&A, 545, A81 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  73. Ramírez, I., Khanal, S., Aleo, P., et al. 2015, ApJ, 808, 13 [Google Scholar]
  74. Reipurth, B., & Mikkola, S. 2012, Nature, 492, 221 [NASA ADS] [CrossRef] [Google Scholar]
  75. Ribas, Á., Bouy, H., & Merín, B. 2015, A&A, 576, A52 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  76. Rica, F. M., Barrena, R., Henríquez, J. A., Pérez, F. M., & Vargas, P. 2017, PASA, 34, e004 [NASA ADS] [CrossRef] [Google Scholar]
  77. Rice, W. K. M., Armitage, P. J., Wood, K., & Lodato, G. 2006, MNRAS, 373, 1619 [Google Scholar]
  78. Saffe, C., Jofré, E., Miquelarena, P., et al. 2019, A&A, 625, A39 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  79. Santos, N. C., Israelian, G., & Mayor, M. 2004, A&A, 415, 1153 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  80. Savvidou, S., & Bitsch, B. 2021, A&A, 650, A132 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  81. Savvidou, S., & Bitsch, B. 2023, A&A, submitted [Google Scholar]
  82. Savvidou, S., Bitsch, B., & Lambrechts, M. 2020, A&A, 640, A63 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  83. Schneider, A. D., & Bitsch, B. 2021, A&A, 654, A71 [Google Scholar]
  84. Serenelli, A. M., Haxton, W. C., & Peña-Garay, C. 2011, ApJ, 743, 24 [NASA ADS] [CrossRef] [Google Scholar]
  85. Simon, A. A., Wong, M. H., & Orton, G. S. 2015, ApJ, 812, 55 [NASA ADS] [CrossRef] [Google Scholar]
  86. Stammler, S. M., & Birnstiel, T. 2022, ApJ, 935, 35 [NASA ADS] [CrossRef] [Google Scholar]
  87. Stammler, S. M., Lichtenberg, T., Drążkowska, J., & Birnstiel, T. 2023, A&A, 670, A5 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  88. Teske, J. K., Khanal, S., & Ramírez, I. 2016a, ApJ, 819, 19 [NASA ADS] [CrossRef] [Google Scholar]
  89. Teske, J. K., Shectman, S. A., Vogt, S. S., et al. 2016b, AJ, 152, 167 [NASA ADS] [CrossRef] [Google Scholar]
  90. Tucci Maia, M., Meléndez, J., & Ramírez, I. 2014, ApJ, 790, L25 [NASA ADS] [CrossRef] [Google Scholar]
  91. van der Marel, N., Williams, J. P., Ansdell, M., et al. 2018, ApJ, 854, 177 [Google Scholar]
  92. Voelkel, O., Klahr, H., Mordasini, C., Emsenhuber, A., & Lenz, C. 2020, A&A, 642, A75 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  93. Walsh, K. J., Morbidelli, A., Raymond, S. N., O’Brien, D. P., & Mandell, A. M. 2011, Nature, 475, 206 [Google Scholar]
  94. Weidenschilling, S. J. 1977, MNRAS, 180, 57 [Google Scholar]
  95. Williams, J. P., & Cieza, L. A. 2011, ARA&A, 49, 67 [Google Scholar]
  96. Zhu, Z., Stone, J. M., Rafikov, R. R., & Bai, X.-N. 2014, ApJ, 785, 122 [NASA ADS] [CrossRef] [Google Scholar]

1

To guarantee mass conservation, only a maximum of 90% of the mass of a given element and at a given location is allowed to change phase per time step Δt.

All Tables

Table 1

Volume mixing ratios for the considered molecular species, given by the applied chemical partitioning model, along with the condensation temperatures of the respective species.

Table 2

Model parameters varied in the simulations for the parameter study, and their default value.

All Figures

thumbnail Fig. 1

Qualitative illustration of the impact of dust drift and planet formation on the composition of material accreted by the host star. The top part shows the drift of pebbles consisting of three different chemical species, a volatile one (blue), a moderately refractory one (orange) and a highly refractory one (green). The ice lines of the blue and orange species are shown with vertical lines of the corresponding color, while the background color of the disk indicates the gas components. The green species does not evaporate in the disk. More refractory species drift further before evaporating at the ice line, enriching the accreted material more strongly than volatile species. The bottom illustration shows how a giant planet core at pebble isolation mass changes the accreted material’s composition. The volatile component evaporates before the pebbles reach the planet’s pressure bump, allowing diffusion past the planet’s orbit and onto the star. The orange and green species are solid at the pressure bump location and blocked from further accretion. Therefore, the accreted material is no longer enriched in the moderately (orange) and highly refractory (green) species.

In the text
thumbnail Fig. 2

Time evolution of the stellar convective envelope mass obtained from models by Hoppe et al. (2020) for various [Fe/H] (left panel) and M (right panel), indicated by a change in line color. The ordinate shows 1 – MCE/M, so that if the star has a more massive convective zone, a lower value is shown, to indicate that the convective zone abundances adapt to the composition of the accretion flux more slowly.

In the text
thumbnail Fig. 3

Surface density evolution of carbon, oxygen, iron and nitrogen, for α = 1 × 10−4. The sum of the gas and dust component is shown. The color indicates the surface density of the element with respect to the hydrogen surface density, compared to the initial fraction at a given location on a logarithmic scale, such that a blue color indicates that the surface density of the element is enriched, and vice-versa for a red color. Due to this choice of norm, a blue color at the inner edge corresponds to an increase in [X/H] for that element up to the magnitude of the enrichment. All evaporation front locations of species where the corresponding element is a component in are indicated with black dashed lines in the relevant panels.

In the text
thumbnail Fig. 4

Accretion flux of disk material onto the host star for small particles (α = 1 × 10−3, left panel) and for large particles (α = 1 × 10−4, right panel) as a function of time. On the ordinate, the accretion rate with respect to the hydrogen accretion rate, i/MH, for various elements i is shown, indicated by different colors. The accretion rate is shown normed to the mass fraction with respect to hydrogen of the element as defined by the initial condition. Therefore, if the ordinate value is larger than 1, the material that is accreted results in a growth of the abundance of element i in the convective zone up to the magnitude of the enrichment. The accretion rate is directly related to the surface density evolution (see Fig. 3).

In the text
thumbnail Fig. 5

Stellar convective envelope elemental abundances as a function of time relative to the initial abundances, [X/H] – [X/H]0. Four different simulations are shown, varying α and ufrag to consider different maximum grain sizes in combination with different gas accretion speeds. Differently colored solid lines represent the abundances of the elements carbon (yellow), oxygen (red), iron (green), and nitrogen (blue), while the black dashed line shows the convective zone mass evolution as 1MCEM$1 - {{{M_{{\rm{CE}}}}} \over {{M_ \star }}}$, as in Fig. 2, corresponding to M = 1 M and [Fe/H] = 0. Solid lines indicate results for M0 = 0.1, while dashed lines show a comparison to simulations with M0 = 0.05.

In the text
thumbnail Fig. 6

Scaling of initial [X/H] with [Fe/H] for the elements considered in the disk model based on GALAH+ DR3. The scaling of the abundances is indicated by solid colored lines. In addition, a reference line of 1:1 scaling is shown as a black dashed line, and a dummy data point indicating the mean error of the measurements is shown in the top left corner of each panel. In the left panel, elements that were previously considered by Bitsch & Battistini (2020) are presented, while the remaining elements not previously considered, but used in the model, are shown in the left panel. Sulfur is not indicated; here, the same scaling as for silicon is assumed.

In the text
thumbnail Fig. 7

Same as Fig. 5, but the initial [Fe/H] is varied instead of α. It is kept constant at α = 1 × 10−3. On the left-hand side, the case for an initial metallicity of [Fe/H] = −0.4 is shown for both ufrag = 1 m s−1 and ufrag = 10 m s−1, whereas the right-hand size shows simulations with [Fe/H] = +0.4.

In the text
thumbnail Fig. 8

Same as Fig. 5, but the stellar mass is varied instead of α and ufrag. The three panels show simulations with M ∈ {0.6 M,1.0 M,1.5 M}. For the simulation with M = 1.5 M, the simulation run time is reduced to 5 Myr in accordance with observations of accretion disks around massive stars.

In the text
thumbnail Fig. 9

Same as Fig. 5, but only α is varied. For the shown simulations, planetesimals are formed in the disk.

In the text
thumbnail Fig. 10

Same as Fig. 3, but the surface density evolution of oxygen in solid and gaseous form is shown for three different planet formation locations. Additionally, the H2O, CO2 and CO ice lines are labeled.

In the text
thumbnail Fig. 11

Same as Fig. 5, but now, a planet is included in the disk. The viscosity is kept constant, but ufrag is varied. Simulations where the planet was placed at three different locations are shown, relative to the Fe3O4, H2O and CO2 ice lines. Each column represents a different planet placement, with it being placed further out going from the left-most to the right-most column. The top row depicts disks with ufrag = 1 m s−1 and the bottom row those with ufrag = 10 m s−1.

In the text
thumbnail Fig. 12

Solar convective zone abundance changes created by accretion of a disk forming Jupiter as a function of disk viscosity α for ufrag = 1 m s−1 (blue) and ufrag = 10 m s−1 (orange). The data points show the mean value for simulations with R0 ∈ {75 AU, 150 AU} and M0 ∈ {0.05, 0.1}, and the error bars indicate the standard deviation.

In the text
thumbnail Fig. 13

Abundance differences between HD106515 B and A caused by planet formation taking place only around star A. The top shows the case of the planet forming where it is observed today at 4.5 AU, which lies outside the water ice line, while the bottom considers the planet forming inside the water ice line, just outside rFe3O4${r_{{\rm{F}}{{\rm{e}}_3}{{\rm{O}}_4}}}$. The gray area delimited by black lines corresponds to the measurements including the 1σ error (Liu et al. 2021). In blue, the disks around the constituents are parameterized by α = 1 × 10−3, ufrag = 10 m s−1. The red and green lines describe disks with a viscosity of α = 1 × 10−4 with ufrag = 1 m s−1 and ufrag = 10 m s−1, respectively.

In the text
thumbnail Fig. 14

Same as Fig. 13, but including planetesimal formation around star B. For simulations shown on the left-hand side, plantesimals can form in any disk region where the pebble flux is above a threshold value, in accordance with the model by Voelkel et al. (2020). For those on the right-hand side, planetesimal formation is turned off for r<rH2O$r < {r_{{{\rm{H}}_2}{\rm{O}}}}$. The planetesimal formation efficiency ϵ is varied between lines with different markers and styles. Sold lines with a circular marker show the results for ϵ = 0.05, dashed lines with a cross marker show ϵ = 0.01 and dotted lines with triangular markers show ϵ = 0.1.

In the text
thumbnail Fig. C.1

Migration track of HD106515A b for ap,0 = 2.5 AU, α = 1 × 10−4 and ufrag = 1 m s−1, shown in blue. The black dashed lines indicate ice lines that are crossed during the inward movement.

In the text
thumbnail Fig. C.2

Same as Fig. 13, but including migration of the planet around star A, with an initial position of ap,0 = 2.5 AU.

In the text
thumbnail Fig. D.1

Same as Fig. 13, showing only the case where the planet formed inside the water ice line. Here, M0,A = 0.1 and M0,B = 0.05

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.