Massive Interacting Binaries as an Enrichment Source for Multiple Populations in Star Clusters

We present a suite of binary evolution models with massive primaries (10 ≤ M 1 ≤ 40 M ⊙) and periods and mass ratios chosen such that the systems undergo nonconservative mass transfer while the primaries have helium cores. We track the total mass and chemical composition of the ejecta from these systems. This material shows the abundance signatures of hot hydrogen burning that are needed to explain the abundance patterns seen in multiple populations in massive star clusters. We then calculate the total yield of a population of binary stars with masses, mass ratios, and periods consistent with their distribution in a field population. We show that the overall abundance of this material is enriched in helium, nitrogen, sodium, and aluminum, and depleted in carbon, oxygen, and magnesium, by amounts that are consistent with observations. We also show that such a population of binaries will return approximately 25% of its mass in this ejecta (compared to 4% if all the stars were single), over a characteristic timescale of about 12 Myr. We argue that massive binaries must be seriously considered as a contributor to the source of enriched material needed to explain the multiple populations in massive clusters, since essentially all massive stars are formed in binaries or higher-order multiples, massive binaries are primarily formed in clusters, and massive binaries naturally produce material of the right composition.


INTRODUCTION
Massive star clusters were for many years thought to be simple stellar populations, consisting of stars formed at the same time, and out of material with a uniform composition.This assumption is still reasonable for the purposes of determining their ages or metallicities (e.g.Dickson et al. 2023;Adamo et al. 2023;Ying et al. 2023).As simple populations, clusters can encode the galactic conditions at the time and place when the cluster formed, and therefore can trace galaxy formation and assembly (Krause et al. 2020).Clusters have also been invaluable for our understanding of stellar evolution and stellar dynamics, since they were thought to provide constraints such as a common age, composition, and distance which are usually unavailable for individual stars.
However, the advent of exquisite Hubble Space Telescope photometry of large number of globular clusters in the Milky Way, and accompanying high resolution spectroscopy, confirmed the hints from many decades earlier that in fact cluster stars do not have uniform surface abundances.While many stars show the same chemical patterns as field stars of the same metallicity, approximately half of the stars show variations.In particular, we see correlations or anti-correlations in light elements but (generally) no variation in iron-peak or s-process elements.Specifically, the enriched population shows enhancement in helium, nitrogen, sodium; depletion in carbon and oxygen; and in a few clusters there is an increase in aluminum for either constant or decreased magnesium.The variations are consistent with material that has been processed by the nucleosynthetic cycles associated with hot hydrogen burning beyond the CNO cycle (at about 70 MK) but not subsequent burning.Essentially all massive clusters show evidence of these so-called "multiple populations" when studied carefully, regardless of cluster age or galactic environment.For recent reviews of the state of the observational evidence for multiple populations in globular clusters, see Gratton et al. (2019); Bastian & Lardo (2018); Milone & Marino (2022).
The origin of these multiple populations is still not clear.Any complete model must identify a nucleosynthetic site capable of processing hydrogen-rich material at 70 MK.That material must then be made available to either form new stars or be added to the surfaces of pre-existing stars.And finally there must be enough material available so that half or more of the presentday stars in a cluster have the modified abundances.The nucleosynthetic source must also primarily act in clusters, as only a small percentage of field stars show these same chemical patterns (Martell et al. 2011).It must act on a timescale that is quite short, as there is no demonstrable age difference between the various populations.In addition, the abundance anomalies are found in stars throughout the HR diagram, suggesting that the stars were formed with this composition rather than being polluted later in their lives.
There are handful of viable sources for the material with the appropriate composition (usually called 'enriched' material even though some elements are depleted).These include asymptotic giant branch stars (Ventura & D'Antona 2009), rapidly-rotating massive stars (Decressin et al. 2007), very massive (∼ 300 M ⊙ ) or supermassive (≳ 5000 M ⊙ Gieles et al. 2018) stars, stellar collisions (Sills & Glebbeek 2010), and interacting massive binary stars (de Mink et al. 2009).While many of these sources have been explored in detail (particularly asymptotic giant branch stars and rotating massive stars), the massive binary stars have been given less attention.The original paper that proposed these objects as a source of material de Mink et al. (2009) modelled only one binary system, essentially as a proof of concept.
However, massive binaries are an attractive solution to the multiple population problem.Multiple populations are found in essentially every massive cluster studied, and so the mechanism that creates them must be a normal outcome of star formation.Massive stars are preferentially found in binary or higher order systems, with a multiplicity fraction of 100% for stars more massive than ∼ 10 M ⊙ (Offner et al. 2023).Binary interactions (e.g.mass transfer, mergers, and common envelope evolution) dominate the evolution of massive stars (Sana et al. 2012).In dense regions, dynamical evolution of binaries will primarily drive them to shorter orbital periods, potentially increasing the number of interacting massive binaries in clusters.In a previous study, we have demonstrated the viability of a pollution source which is proportional to the number of massive stars and which acts during the cluster assembly process inside the cluster's natal giant molecular cloud (Howard et al. 2019).
The binary system investigated by de Mink et al. ( 2009) consisted of two massive stars with an orbital period such that mass transfer would occur after core hydrogen exhaustion but before core helium burning had begun.They demonstrated that, with an appropriate treatment of angular momentum which spins up the secondary and drives non-conservative mass transfer, significant amounts of processed material could be ejected from the system.In this paper, we revisit massive interacting binary stars as a pollution source for multiple populations in massive clusters.We present a suite of binary evolution models covering a larger range of masses and orbit properties, and quantify the yields from these systems.We then look at the ensemble yield from a population of such binaries.Finally, we comment on the viability of massive interacting binaries as a source of enriched material which should be considered in a solution to the multiple population problem.

Binary evolution calculations
We simultaneously evolve both components of our binary systems using Modules for Experiments in Stellar Astrophysics (MESA) version r15140 (Paxton et al. 2011(Paxton et al. , 2013(Paxton et al. , 2015(Paxton et al. , 2018(Paxton et al. , 2019;;Jermyn et al. 2023).We generally used the suggested parameters from the MESA Isochrones and Stellar Tracks (MIST) suite for single star models (Choi et al. 2016) and followed the example of Paxton et al. (2015) for binary evolution.A full description of the parameters used can be found in Nguyen (2023), and sample inlists, the custom reaction rate network file, and our run stars extras.f90&run binary extras.f90files are available at 10.5281/zenodo.10998057.
We constructed zero-age main sequence (ZAMS) models from the pre-main sequence with an initial alphaenhanced metallicity of [Fe/H] = −1.44,similar to that used by de Mink et al. (2009) and in the middle of the range of metallicities of Milky Way globular clusters.We follow the prescription of Reddy et al. (2006) for the alpha-enhancement, which results in [O/Fe]=0.44,Z=1.2 × 10 −3 and Y=0.25.These ZAMS models were evolved as single stars to allow us to determine appropriate periods for our binary systems, and as a comparison population for the nucleosynthetic yields.The ZAMS models also serve as the starting point for each of our binary evolution calculations.Mixing mechanisms appropriate for massive stars are included (convection, semi-convection, overshoot, thermohaline mix-ing, rotational mixing).Wind mass loss is treated using the 'Dutch' scheme, as outlined in Glebbeek et al. (2009).Because we wish to evaluate the nucleosynthetic yields relevant for the multiple population problem, we use a custom reaction network which is an expansion of MESA's 'h burn.net' and includes a more complete coverage of H, He, and C burning processes including the high-temperature hydrogen burning cycles.At each time step we compute and output the ejection of 50 isotopes of elements up to the iron peak.
Our choices for binary evolution generally follow those given in Paxton et al. (2015).We use the default Ritter (1988) mass transfer scheme, which was also used in de Mink et al. (2009).Our binaries are initialized on circular orbits and given an initial rotation rate as set by tidal synchronization.As the binary evolves, angular momentum loss from the orbit is computed based on mass loss unaccounted for by the stellar angular momentum loss.We also allow angular momentum transfer through accretion, both of which are necessary for following the non-conservative mass transfer which takes place in these systems.The maximum rotation rate that a star can have is set to 90% of its critical rotation rate, as MESA is not numerically stable beyond that rate (Paxton et al. 2019).

Initial Conditions
For this study, we aim to cover the full parameter space which may produce ejecta with appropriate abundances for the multiple populations problem, while reducing the time spent on simulations which do not add much material to the population.We chose our primary masses M 1 to be in the range of 10 − 40 M ⊙ since these stars are massive enough to reach high central temperatures.Higher mass stars also reach appropriate temperatures, but as will be described in more detail below, binaries with more massive primaries are unlikely to contribute significant amounts of material of the appropriate abundances for the multiple populations problem.Mass ratios q = M 2 /M 1 were evenly sampled over the range of 0.15 -0.9.This excludes twins, q ≳ 0.95, which if close enough to initiate mass transfer in our simulations, would become contact binaries.It also excludes very unequal mass systems, which behave more and more like single stars or conservatively-transferring systems.
Our binary periods were chosen such that the systems would undergo mass transfer between core hydrogen exhaustion and core helium ignition (i.e.case B mass transfer).Using our single-star models, we determined the radii of each of our stars at these two evolutionary stages, which turns to correspond to two peaks in the radius evolution, as shown in Figure 1.We then set these radii equal to the Roche lobe radius as given by Eggleton (1983), assuming a mass ratio of 0.5, and solved for the period given this semi-major axis.We also added one longer period for each primary mass, just to confirm that we had in fact bounded the period appropriately.The evolutionary endpoints correspond to periods ranging from approximately 2-700 days, and the longest period we simulate is ∼ 2800 days.The lower limit falls into a narrow range of around 2-3 days for the mass range we cover.The upper period limit varies sharply with mass.With a narrower range of periods for the high-mass primaries expected to act as enrichment sources, and with the relative rarity of high-mass primaries, the high-mass binaries in our sample are not expected to contribute significant enrichment.Figure 2 shows the parameter space sampled by the binaries, spanning period P , primary mass M 1 , and mass ratio q.
A total of 204 binaries were simulated.

Stopping Conditions
Periods, primary masses, and mass ratios of the binaries simulated in this work.Binaries of the same initial primary mass and initial orbital period are triangular segments of a hexagon centred on the primary mass and period.The relative orientation of the segment represents the mass ratio, q, shown in the legend.The grey dashed line represents the expected upper limit on the period for binaries of interest based on single-star models at the sampled primary masses.The orange triangle indicates the system modelled by de Mink et al. (2009).
Stars within a binary were simultaneously evolved from zero-age main sequence (ZAMS) models until they reached one of our four possible stopping conditions.The endpoint of all 204 of our systems are shown in Figure 3. First, we stopped our simulation when we saw a sudden and very large increase in the system mass loss rate, which signals the onset of dynamically unstable mass transfer and, potentially, the creation of a common envelope system.Second, we ended a run if the stars came into contact, i.e. when both stars are simultaneously overflowing their Roche lobes.These systems usually occur at the shortest periods.MESA has been used to model contact systems for primary masses below 20 M ⊙ (Henneco et al. 2024).In our higher mass models, we encountered both significant numerical instability and uncertainty in the physical validity of different approaches.Therefore, we made the choice to stop our simulations, knowing that we may be under-estimating the amount of material that could be ejected from these systems.
Third, some systems evolved to core carbon depletion, which we defined as the time when the central 12 C mass fraction drops to 5 × 10 −3 along with a low central He mass fraction (< 1 × 10 −4 ).These are typically the intermediate mass primaries (12 ≲ M 1 ≲ 25 M ⊙ ).The final stopping condition is called 'none' in the figure, because in these systems, mass transfer ceased and the stars' final masses were too low to reach core carbon burning.These primaries will likely eventually evolve to be white dwarfs, but we stopped their evolution when it was apparent that they had reached this state.
For a binary to be included in our completed sample, it must eject at least 0.1 solar masses of material and have the primary reach Roche lobe overflow during simulation time.This ensures that the binaries we examine undergo non-conservative mass transfer, the mechanism we are interested in for acting as an enrichment source.Many of the systems at the longest periods did not go through Roche lobe overflow, although some of the higher mass systems beyond our nominal period limit did contribute, confirming that the evolution of a primary star in a binary system is indeed slightly different than that of a single star, even before mass transfer is initiated.

Abundances of binary ejecta
The binary systems in our simulation suite follow the same general evolution as the system described in detail in de Mink et al. (2009): the primary exceeds its Roche lobe radius after core hydrogen exhaustion and begins to transfer material to the secondary.At first this mass transfer is conservative, but the material transfers significant angular momentum to the secondary and it spins up.When it nears its critical rotation rate, it can no longer accept any more material and so any additional material from the primary leaves the system.Typically this happens after relatively small amount of material (≲ 1 M ⊙ ) has been added to the secondary.We note that this is the default behaviour in MESA and results in a low accretion efficiency.Measurements of accretion efficiency in post-interaction systems, and also different treatments of mass transfer modelling, suggest higher efficiencies could be possible Marchant & Bodensteiner (2023).Since there is considerable uncertainty, we continue with the default MESA approach but note that our yields are dependent on the amount of material that can remain on the secondary during mass transfer compared to the amount that flows out of the system.
We track how much material is lost from the system, and also the abundance of each element in that material over time.We expect that much of this material will show signatures of hot hydrogen burning, as it will have originated fairly deep in the primary star.The outer, unmodified layers of the primary may have been shed through a normal stellar wind, but mostly will be transferred and accepted by the secondary during the first, conservative, phase of mass transfer.By the time the material starts leaving the system, some or all of that material will be removed from the primary, and the enriched layers will be available.
For each binary in our simulation suite, we show the abundance of the elements observed to vary in multiple populations (C, O, He, N, Na) in Figure 4.Here we show the mass fraction of that element in the ejecta relative to the initial surface abundance of that element, X/X i .Points that are darker have ejecta with stronger enrichment signals (associated with more extremely enriched populations) than lighter points.We have a clear signature of hot hydrogen burning, with depleted carbon and oxygen (top panels) and enhanced helium, nitrogen and sodium (bottom panels).Direct comparisons to the observed enhancement values seen in globular clusters are not straightforward since there is a wide variety of enhancement levels between clusters, and some abundances such as helium and often carbon are not measured directly but inferred from photometry or molecular bands.Nevertheless, some clusters have been studied in detail (e.g. the extreme NGC 2808, and the more typical NGC 6752, as described in Bastian & Lardo 2018;Gratton et al. 2019).The maximum observed enrichment levels are approximately X/X i =1.5 for helium, up to 100 for nitrogen, and 5 for sodium; oxygen depletion can reach X/X i =0.5, and carbon may be as low as 0.3.
Systems with high initial primary masses M 1 = 40 M ⊙ consistently show the strongest enrichment signals.Enriching binaries (darker points) tend to have higher primary masses M 1 ≳ 12, higher q ≳ 0.45 (bottoms and upper right of hexagons), and are at longer periods.These are the systems where most of the pristine outer layers can be accreted onto the secondaries so that it will not be ejected.These periods have also given the primary enough time to allow for significant hydrogen shell burning to proceed, adding to the total amount of hydrogen-processed material that is available.Binaries with shorter periods, particularly the high primary mass systems, do not contribute as much because they reach contact relatively quickly and stop contributing.We also note that there are systems that show essentially no signs of processing with X/X i ∼ 1.This suggests that non-conservative mass transfer can not only provide enriched gas but can also provide pristine gas for dilution.
The total amount of mass ejected by each system varies considerably, as shown in Figure 5.Here we show the enrichment in five representative elements as a function of the fraction of the initial system mass that has been ejected from the system.In general, less enriching binaries with X/X i ∼ 1 eject a smaller fraction of their mass.A transition appears when around 20% of the system mass has been ejected from the system.We also see that the maximum fraction of system mass lost is under 50% of the total mass of each binary system.
Our simulations followed many more elements than just those shown in the previous figures.The main ones are listed in Table 1, where we give the minimum and maximum enrichment or depletion X/X i in our suite of 204 systems, and the median X/X i of all systems.The iron peak elements (Ti, Fe, and Ni) do not change, and neither do the heavier alpha elements (Ca, Si).
The fragile light elements Li, Be, and B are all significantly depleted.In some systems, hydrogen is depleted which indicates the ejection of H-burned, He-enhanced material, but the overall depletion of hydrogen is very small.The other light elements show the expected signatures of multiple populations (and hot hydrogen burning): depletion in C, O, and Mg and enhancement in He, N, Na, and Al.We have also confirmed that the expected correlations (e.g.C-O, N-Na) and anticorrelations (e.g.C-N, Na-O) are found in the ejecta  from massive interacting binaries.Since we are essentially probing the abundances of hot-hydrogen-processed material, this is the expected outcome of our choices of binary parameters.
We make special note of the Mg-Al anti-correlation, since that is seen in only the most massive globular clusters, and is often used as an important diagnostic to test models of multiple population yields.Most of our binaries do not show significant change in their magnesium abundance.If magnesium does change, we see  M ⊙ are responsible for almost all the change in Mg, as the Mg-Al cycle requires very high temperatures to be significantly activated.

Timescales for binary ejecta
The timescale on which these systems eject mass is an important factor in their ability to act as an enrichment source.In all cases, the duration of the mass transfer event is short (≲ 0.5 Myr, and the non-conservative portion of the mass transfer event is half or less of that time), compared to the evolution of the system before the onset of mass transfer.Therefore, we assign a single time of mass ejection for all of the binaries in our samples, shown in Figure 6.We define this time as the age of the system when the primary's radius first exceeds the Roche lobe radius.The lowest mass primaries eject on timescales ≳ 23 Myr, and the highest mass primaries we simulate can eject their material as early as ∼ 4.3 Myr.The time of mass transfer onset increases with increasing period, but this results in less spread (∼ 5 Myr for our lowest mass systems than the difference caused by the evolutionary timescale of the primary.

Expected yields from a population of binaries
The total amount of available material, and its composition, to create the enriched population in star clusters will not come from a single binary, but instead from the ensemble of the binary population.Now that we have the nucleosynthetic yields for individual systems, we can combine them by weighting each binary by its likelihood of being part of the population.To do this,  2. Total enrichment or depletion relative to the initial mass fraction of each element, X/Xi, and also the new value of Y or [X/Fe], from a population of binary stars weighted by their likelihood of occurring in the field.
we need to make some assumptions about the properties of the binary population of a forming cluster.For this work, we make the simplifying assumption that the binary population in a cluster is the same as the field population.Binary populations can be modified dynamically through interactions with other stars and the natal gas in forming clusters (Cournoyer-Cloutier et al. 2021), but the initial distribution cannot be very different from the field population.We use the same binary sampling routine as in Cournoyer-Cloutier et al. ( 2021), which is based on the field binary properties from Moe & Di Stefano (2017).The binary fraction and period distribution are functions of the primary mass, and the mass ratio distribution depends on both primary mass and period.
We determine the probability of having a binary in a three-dimensional grid covering primary mass, mass ratio, and period.We do not consider the eccentricity distribution of the binaries as we only simulate circular binaries.We interpolate the yields of each element and the total amount of mass ejected from our binary system simulations onto the gridpoints of primary mass, mass ratio and period.We do not extrapolate any yields beyond the boundaries of our simulations (1.87 < P < 661.65 days, 10 < M 1 < 40 M ⊙ , and 0.15 < q < 0.9 in steps of ∆(log P ) = 0.01 [days], ∆M 1 = 0.1 M ⊙ , and ∆q = 0.01) but we do calculate the probabilities over the full binary population explored in Cournoyer-Cloutier et al. ( 2021) which covers masses from 0.08 to 150 M ⊙ , periods between 10 0.5 and 10 7.5 days, and all mass ratios.The total yield of each element is then the sum, over all grid points, of the individual binary yield weighted by the probability of having such a binary.We are not considering any material outside the boundaries of our simulations, so the results presented here should be considered as the yields from hot hydrogen burning only.
The results for the elements relevant for multiple populations are given in Table 2 in two different ways: the total enrichment/depletion relative to the initial mass fraction of the element, X/X i , and also the resultant [X/Fe] value (except for helium which is reported as the new value of Y).This weighted average of the individual binary yields is precisely what is expected to explain the multiple populations problem: enhancement in helium up to above Y=0.3,enhancement in oxygen, sodium, and aluminum accompanied by depletion in carbon, nitrogen, and magnesium.The levels of enhancement or depletion (e.g.large values of nitrogen, small values of magnesium) are also what is seen in detailed spectroscopic observations in globular clusters (e.g.Gratton et al. 2019).
Using the same weighting of the individual binaries, we calculate that the population of binary stars will release 23% of its total mass in ejecta with this composition.Earlier work investigating the binary ejection hypothesis (Howard et al. 2019) assumed 3 possible values for this quantity (10, 25, or 50%) and found that 25% was the best match to observations, so we can confirm this was a reasonable assumption.Finally, we calculate the binary-weighted average timescale for this ejecta to be released to be just under 12 Myr.For comparison, the same population of stars as our binaries, but evolved as two single stars, would eject approximately 3.8% of their mass on a timescale of 14 Myr, a reduction of about 1/6 in total mass compared to the binary population.

SUMMARY AND DISCUSSION
Using the stellar evolution code MESA, we simulated the simultaneous evolution of both stars in 204 binary systems with parameters that covered the range over which we expected to see non-conservative case B mass transfer: 10 ≤ M 1 ≤ 40 M ⊙ , 0.15 ≤ q ≤ 0.9, and 2 ≲ P ≲ 700 days.We tracked the composition of the material which is lost from these systems and compared the abundances to the observed and inferred abundances of multiple populations in globular clusters.The yields of our binaries showed the expected signature of hot hydrogen burning: they are enhanced in helium, nitrogen, sodium, and aluminium, and depleted in carbon, oxygen, and magnesium.Fragile light elements Li, Be, and B are destroyed, and the iron peak elements are unchanged.The material is ejected between when the stars are between 4 and 23 Myr old, depending on the mass of the primary.
We then combined these yields with our expectation for the distribution of these binaries in mass, mass ratio, and period space to calculate the total yield of a normal, field-like population of binaries.We find that the total population will release approximately 25% of its total mass in ejecta on a typical timescale of 12 Myr, and with overall abundances which match the most extreme populations in massive globular clusters.Dilution of this material with pristine gas will provide the range of abundances seen in clusters.This confirms the suggestion of de Mink et al. (2009) that a population of massive binaries are a viable source of enriched material out of which to form multiple populations.
The total mass ejected is consistent with the amount needed to produce multiple populations in the context of pollution during cluster assembly (e.g.Howard et al. 2019).If instead we assume that clusters self-enrich some time after all the un-enriched stars have formed, then binaries do not produce enough material to solve the mass budget problem.However, most other stellar populations, individually or even when considering multiple enrichment sources (Sills & Glebbeek 2010), are also not able to create enough enriched material.This problem, combined with the increasing evidence of no age difference between the different populations, points strongly toward an enrichment source which must be acting during cluster assembly, and enriching only some of the material out of which the cluster will assemble during that process.
If indeed enrichment must happen during cluster assembly, then the timescale to produce the enriched material must be shorter than the timescale on which cluster stars are forming.We find that binaries will produce their enriched material within about 12 Myr.Nearby young clusters have observed age spreads of at most a few Myr (see Massey & Hunter 1998;Da Rio et al. 2010, for studies of R136 and Orion respectively).Such a short formation time would pose a problem for the massive binaries to release their enriched material in time to be relevant.On the other hand, simulations of massive star cluster formation suggest that star formation durations up to 20 Myr (Andersson et al. 2024) are possible, with many simulations suggesting that about 10 Myr is needed to assemble massive clusters (e.g.Howard et al. 2018;Rizzuti et al. 2023).Dating very young stars is notoriously difficult, and simulations include different prescriptions for important physical processes.Efforts will be needed on both fronts to resolve these discrepancies.We also note that the most enriching systems are those with the highest primary masses, which evolve on shorter timescales and so may do most of the necessary enrichment on timescales which are as short as about 4 Myr.
In this paper, we have assumed that the binaries are evolving in isolation, and that the overall population is the same as that seen in the field.However, binaries in clusters are dynamically modified (Heggie 1975), even during the very earliest stages of cluster forma-tion (Cournoyer-Cloutier et al. 2021).Dynamical evolution, in general, reduces binary periods and increases binary mass ratios.All our binaries with mass ratios closer to 1 produce more enriched material than those with lower q.In general, a longer binary period is more conducive to producing more enriched material, mainly because a longer period allowed the primary time to process more of its material though hot hydrogen burning in its hydrogen-burning shell.How dynamics would affect the abundances depends on when the system's period is reduced.Some systems with periods beyond our current upper period boundary would be pushed into a region where they would now start contributing; others may have their period reduced early and may produce less enriched ejecta; but some are likely to have their periods reduced after most of the burning has happened, and then they might eject even more material than if they had stayed at their original period.If we had the complete dynamical evolution of a population of binaries inside a forming cluster, it may be possible to extrapolate along the binary tracks calculate in this paper to determine the final properties of a modified, non-field binary population.
Massive binaries provide a very attractive solution to the source of multiple populations in globular clusters, simply because they are so common and so commonplace.Massive stars are formed in clusters.Massive stars are formed in binaries and higher order systems such as triples or quadruples.Multiple star systems transfer mass and undergo non-conservative mass transfer.They are dynamically modified in a cluster environment.We have shown that a population of massive binaries will produce a significant amount of material with the right abundances, up to and including the Al-Mg anti-correlation.They eject this material on a timescale which is comparable to the assembly timescale of massive star clusters.Just as any models of star cluster formation should include binary stars, any model of multiple populations should include binary stars.It remains to be seen whether binaries alone can produce all the signatures of multiple populations, but since we know that binaries exist in significant numbers, we must include them in our models.

Figure 1 .
Figure 1.Evolution of an example single model (20 M⊙),showing the radius evolution in the upper panels, the stellar track in the lower left, and the central abundances in the lower right.The x axis on the right-hand side shows the model number, which increases with time but shows the interesting phases of evolution more clearly.Orange crosses and grey vertical lines correspond to the identified peaks in radius which determine the minimum and maximum periods we allow.

Figure 3 .
Figure 3. Stopping conditions reached by each binary in the sample of 204 binaries.Filled-in triangles are included in our main sample and empty triangles are omitted.The layout of this figure follows Figure 2. The colour of each point indicates the stopping condition reached by that binary.Out of the 190 binaries that we include in our main sample, 29 reach core-C exhaustion, 34 initiate contact, and 67 have sharp changes in their mass loss rates, and 14 out of 204 binaries are excluded for not ejecting at least 0.1 M⊙.

Figure 4 .Figure 5 .
Figure 4. Enrichment of the total ejected material from each simulated binary.The enrichment is the relative change in the mass fraction of the element compared to the initial surface abundances, X/Xi.Values less than 1 indicate that element is depleted in the ejecta from that binary; values above 1 indicate enhancement.The colour indicates the strength of enrichment;darker values correspond to stronger depletion (top panels) or stronger enhancement (bottom panels).The strongest signals of enhancement are present at higher q, higher P , and higher M1.

Figure 6 .
Figure 6.The time of onset of Roche lobe overflow as a function of initial primary mass.Systems with longer initial periods begin mass transfer at later times for a given mass.[Na/Fe] is shown on the colour bar as an indication of the enrichment level of the binary ejecta.