© ESO 2019

1. Introduction

Spectra of stellar populations provide us with powerful tools to study stars and their host galaxies across cosmic time. Existing surveys and those anticipated with future facilities, such as the James Webb Space Telescope (JWST, Gardner et al. 2006), are expected to deliver a wealth of observational data that can potentially revolutionize our understanding. Translating these data into measurements of the physical quantities of interest, such as star formation rates, requires the use of theoretical or semi-empirical models for the spectra of stellar populations. Accurate models for the spectra of stellar populations, and in particular the ionizing radiation, are therefore indispensable (Conroy 2013).

Ionizing photons are of primary interest for two reasons. Ionizing photons from stellar sources can be reprocessed by nearby gas and dust, giving rise to infrared excess and various prominent emission lines (Charlot & Longhetti 2001). These include lines that are used as diagnostics to infer star formation rates, metallicities, as well as possible variations in the initial stellar mass function and to infer the presence or absence of active galactic nuclei (AGN). The reliability of these measurements depends on the accuracy of the underlying models. Ionizing photons from stellar populations also play a crucial role as a source of stellar feedback. For example, it is thought to be important in regulating the efficiency of star formation (Krumholz et al. 2006; Dale et al. 2013). On a larger scale, ionizing radiation from stellar populations that escapes the host galaxies can ionize gas in the intergalactic medium (IGM). These galaxies are generally held responsible for the reionization of the Universe (Barkana & Loeb 2001; Robertson et al. 2010; Lapi et al. 2017). What type of sources that caused cosmic reionization is uncertain. However, recently, the role of interacting binary stars as sources of ionizing photons has been shown to potentially be important during the reionization (Ma et al. 2016; Wilkins et al. 2016).

Traditionally, massive single stars are considered to be the main producers of ionizing photons in stellar populations. They are rare and short-lived, but they have high luminosities of ∼10⁴ − 10⁶ L_⊙. Furthermore, with temperatures higher than ∼25 000 K, they emit most of their radiation at energies above the threshold for hydrogen ionization (e.g., Smith et al. 2002; Martins et al. 2005; Ekström et al. 2012; Köhler et al. 2015; Schneider et al. 2018). The most massive stars can eventually lose their hydrogen-rich envelopes as a result of strong stellar winds or eruptions creating Wolf-Rayet (WR) stars (Meynet & Maeder 2005), which can be so hot that they emit photons sufficiently energetic to ionize even helium (Crowther 2007).

Recent studies of nearby resolved stellar populations show that massive and intermediate mass stars are often accompanied by a companion star that orbits so close that interaction between the two stars is inevitable as the stars evolve and swell up (Sana et al. 2012; Moe & Di Stefano 2017). Interaction in binary systems can completely change the future evolution of the stars in the system, leading to mass accretion, rejuvenation, rapid rotation, mass loss, and possibly even coalescence (e.g., Vanbeveren et al. 1979; Podsiadlowski et al. 1992; Wellstein & Langer 1999; Eldridge et al. 2008; de Mink et al. 2014; Schneider et al. 2015; Stanway & Eldridge 2018). Including, or improving the treatment of, the effects of binary interaction in models for the spectra of stellar populations is therefore warranted.

In this study, we focus on stars that have been stripped of their envelope as a result of interaction with a binary companion (Kippenhahn & Weigert 1967; Paczyński 1967). This is expected to be a very common outcome of binary interaction (e.g., Sana et al. 2012), resulting in very hot and compact stars (van der Linden 1987; Yoon et al. 2010, 2017). Due to their high temperatures, they are expected to emit the majority of their radiation at ionizing wavelengths (Götberg et al. 2017, hereafter Paper I). This makes them potentially important, but still often neglected, contributors to the budget of ionizing photons produced by stellar populations, as argued, for example, in Van Bever & Vanbeveren (2003) and Vanbeveren et al. (2007, see also Belkus et al. 2003 and Stanway et al. 2016).

One of the challenges to properly account for the impact of these stripped stars on the ionizing spectra of stellar populations, was, until recently, the scarcity of appropriate atmosphere models. Direct observations of stars stripped in binaries are still scarce, likely because they are typically outshone by their companion (see, however, Gies et al. 1998; Groh et al. 2008; Peters et al. 2008, 2013; Wang et al. 2017, 2018; Chojnowski et al. 2018). As a result, there had been little motivation to calculate atmosphere models for stripped stars. The earliest estimates of the integrated spectra of stellar populations including stripped stars, therefore, relied on black-body approximations or rescaling of atmosphere models that were originally computed for more luminous WR stars. Modern simulations make use of spectral libraries, but these typically do not cover the full parameter space of interest for stripped stars nor consider the effect of metallicity.

This motivated us to undertake the effort of computing a library of evolutionary and spectral models custom-made for stripped stars that result from binary interaction for a range of masses and metallicities (Götberg et al. 2018, hereafter Paper II). In Paper I we showed how metallicity affects the stripping process. At lower metallicity, the progenitor star remains more compact leading to an incomplete stripping of the envelope. The resulting star is therefore partially stripped and more luminous, but also cooler. In Paper II we presented model grids that cover a range of masses and metallicities. We showed how these stars span a continuous range of spectral types, from WR-like spectra characterized by emission lines formed in the winds of the more massive and metal-rich stripped stars, to subdwarf-like spectra dominated by absorption features resulting from the photosphere of stripped stars with transparent outflows. We further predicted the existence of a hybrid intermediate class of spectra showing a combination of absorption and emission lines, very similar to those observed for the recently discovered new stellar class labeled WN3/O3 (Massey et al. 2014, 2015, 2017; Neugent et al. 2017, 2018; Smith et al. 2018).

The aim of this work is to estimate the radiative contribution from stripped stars to the spectral energy distribution of stellar populations and measure the additional emission rate of ionizing photons that stripped stars produce. For this purpose, we performed a population synthesis to estimate the number and type of stripped stars that are present in a population at any given time. We used the custom-made grid of spectral models for individual stripped stars that we published in Paper II. We use this to investigate the impact on the integrated spectra. We discuss the integrated spectra, the emission rate of H I-, He I-, and He II-ionizing photons, and further quantities that can be derived from observations, namely the production efficiency of ionizing photons, the ionization parameter, the UV luminosity and continuum slope, and stellar spectral features.

The first theoretical and semi-empirical studies of the integrated spectra of stellar populations primarily focussed on the effects of single stars. These include codes based on single star models, such as STARBURST99 (Leitherer et al. 1999, 2010, 2014 GALAXEV (Bruzual & Charlot 2003), and PEGASE (Fioc & Rocca-Volmerange 1997, 1999; Le Borgne et al. 2004). However, more recently, several groups have considered the effects of binary interaction. These include the Brussels code (Van Bever et al. 1999; Belkus et al. 2003; Vanbeveren et al. 2007), the Yunnan simulations (Zhang et al. 2004, 2012, 2015; Han et al. 2010; Chen & Han 2010; Li et al. 2012), and the BPASS code (Eldridge & Stanway 2009, 2012; Eldridge et al. 2017).

The diversity of codes has helped improve our understanding of the importance of binary interactions throughout the years. However, the outcome of binary interactions and their spectra is still uncertain (see e.g., De Marco & Izzard 2017) and detailed studies of each binary evolutionary channel is required. One important part for spectral synthesis models that has not been investigated in detail is how the spectra of stars stripped in binaries look like and what their contribution is to the spectrum of the full stellar population.

In this study, we use our spectral models custom-made for stars stripped in binaries to make predictions for their contributions to the spectrum of a full stellar population. We choose to focus on this particular binary evolutionary channel, both so that we can carefully examine the evolution and spectra of these objects and so that we can understand in detail the effects of stripped stars in particular. We compare our results with the results from STARBURST99 to understand the effect of envelope-stripping compared to a population containing only single stars. The code BPASS takes a sophisticated approach, combining very large grids of binary evolution models with spectral libraries with the aim to account for all binary evolutionary channels at once. The developers of BPASS have provided testable predictions for a large range of observable phenomena (Eldridge & Stanway 2012, 2016; Eldridge & Maund 2016; Eldridge et al. 2019; Stanway et al. 2014, 2016; Stanway & Eldridge 2018; Xiao et al. 2018, 2019). By comparing our results to those from BPASS, we can identify effects that come from stripped stars. This is not only valuable for understanding the effects of stripped stars and how various uncertain binary parameters would alter the results for stripped stars, but this approach will also enable the identification of parts that require careful modeling and which effects that are bluntly unreliable because of uncertainties in the models for stripped stars.

Our hope is that the predictions provided in this work will be of use for the interpretation of the data that is resulting from several recent and ongoing surveys that are probing the ionizing properties of stellar populations. This includes the anticipated James Webb Space Telescope (Gardner et al. 2006, see for example Feltre et al. 2016, Gutkin et al. 2016, Xiao et al. 2018), but our simulations are also interesting for various surveys currently conducted from the ground. For example, the MOSDEF survey (Kriek et al. 2015), the KBSS (Rudie et al. 2012; Steidel et al. 2014), the KLCS (Steidel et al. 2018), the GLASS (Treu et al. 2015), the VUDS (Le Fèvre et al. 2015), and the ZFIRE (Nanayakkara et al. 2016). We, therefore, make our predictions available at the CDS and they are also available through the STARBURST99 web-portal.

The structure of this paper is as follows. Section 2 describes the models for the evolution and spectra of stripped stars that we presented previously in Paper II and how we use these to construct a model for the integrated spectrum of stripped stars in a stellar population. In Sect. 3, we present the model predictions for how many stripped stars are present in stellar populations. In Sect. 4, we show how the presence of stripped stars affects the total spectral energy distribution of a stellar population. We quantify the contribution from stripped stars to the emission rate ionizing photons in Sect. 5. In Sect. 6, we discuss the impact of stripped stars on observable quantities. In Sect. 7 we summarize our findings and conclusions. The current paper is the third in a series with Paper I and Paper II, but it can be read independently.

2. Accounting for the radiative emission from stripped stars in stellar populations

We created an estimate of the radiative contribution from the stripped stars in stellar populations. We first describe the models for the evolution of individual stripped stars and their spectra in Sect. 2.1. We then describe the assumptions that we make to model stripped stars in stellar populations in Sect. 2.2 and how we represent the emission from a full stellar population in Sect. 2.3.

2.1. Models of individual stripped stars

2.1.1. Binary evolutionary models

We used the models presented in Paper II to describe the evolution of stars that lose their envelope through interaction with a binary companion (see also Paper I, for an in-depth discussion). These are models of binary systems in which stable mass transfer is initiated early during the Hertzsprung gap, after which the H-rich envelope is stripped off (commonly referred to as Case B type mass transfer, Kippenhahn & Weigert 1967). Stripped stars can also result from stable mass transfer initiated during the main-sequence evolution of the most massive star in the system (Case A type mass transfer, Kippenhahn & Weigert 1967), or from unstable mass transfer and a subsequent successful ejection of the common envelope (Paczynski 1976; Ivanova 2011). The contribution from the different formation channels vary depending on the progenitor mass, but we expect that Case B type mass transfer is responsible for the majority of the stripped stars. Despite the variety of evolutionary histories, the properties of the stripped stars are remarkably similar. These properties are primarily dependent on the mass of the stripped stars alone. This assumption works well for most systems, see however Claeys et al. (2011), Yoon et al. (2017), Siess & Lebreuilly (2018), Sravan et al. (2018) for evolution that leads to a larger fraction of the H-rich envelope is left, which is the case for systems with long initial periods or very low metallicity. We used our models of stripped stars created through Case B type stable mass transfer as an approximation for stripped stars formed via any evolutionary channel. This approximation is sufficient for our current purposes.

The evolutionary models were computed using the binary stellar evolution code MESA (Paxton et al. 2011, 2013, 2015, 2018, 2019). The models have the metallicities Z = 0.014, 0.006, 0.002, and 0.0002, which are representative of the metallicity of the Sun, the Large and Small Magellanic Clouds, and an environment with very low metallicity that may be representative of early stellar populations. Each metallicity grid consists of 23 models, which are different from each other by the initial mass of the donor star, M_1, init. The initial donor star masses were chosen between M_1, init = 2 M_⊙ and 20 M_⊙ with equal spacing in the logarithm of the mass. Since, in most cases, the properties of the stripped star are not sensitive to the formation channel or initial orbital period and initial mass ratio, one binary model was computed for each initial mass of the donor star. We chose the initial mass of the accretor star by setting the initial mass ratio to q ≡ M_2, init/M_1, init = 0.8. We chose the initial periods such that mass transfer was initiated early during the Hertzsprung gap. This corresponded to initial periods between P_init = 3 days and 35 days depending on the initial mass of the donor star (see Table 1 of Paper II). The resulting stripped stars have masses between 0.35 M_⊙ and 7.9 M_⊙. The wind mass-loss from stripped stars is an uncertain parameter because few stripped stars have been observed. The models, therefore, employ extrapolations of the wind mass-loss algorithm for hot and subluminous stars of Krtička et al. (2016) from the low mass end and of the empirical WR mass-loss recipe of Nugis & Lamers (2000) from the high mass end. The switch between the two wind regimes is chosen to occur for stripped stars with progenitor masses of 6 M_⊙. Low-mass stripped stars are affected by diffusion processes that impact the surface composition (Heber 2016). An algorithm accounting for the effect of gravitational settling (see Paxton et al. 2011, also Thoul et al. 1994; Paquette et al. 1986) is included when modeling the evolution of stripped stars. It has strong effects for stripped stars with masses below 0.7 M_⊙ (see Paper II).

Our evolutionary model grids cover the evolution of stripped stars from low to high mass. They stretch from stripped stars at the lower mass limit of central helium burning (∼0.35 M_⊙, Han et al. 2002) up to close to the mass limit where massive stars lose their envelope via their own wind (e.g., Chiosi & Maeder 1986). It is likely that stars of higher mass than what we consider experience envelope-stripping in binaries, for example, through mass-transfer initiated on the main sequence evolution of the donor star (see e.g., Yoon et al. 2010). However, these stars would primarily contribute at early stages because the progenitor stars are more massive and, therefore, live for a shorter time. Here, we focus on the contribution from stripped stars that can not have been created by strong wind mass-loss and, therefore, they have lower masses than most WR stars. For this mass range, we consider that our models are appropriate to use as a representation of the stellar evolution of stripped stars given the careful choices for the wind mass-loss rates and the treatment of diffusion processes on the stellar surfaces. We use the properties of the stripped stars that were computed in the evolutionary models as input for the spectral models described below. We also adopt the time of stripping and the duration of the stripped phase computed in the evolutionary models in our simulations of the integrated spectra of stellar populations, see Sect. 2.2.

2.1.2. Spectral models

The spectral models for individual stripped stars that we used in this work were computed with the non-LTE radiative transfer code CMFGEN (Hillier 1990; Hillier & Miller 1998) and were custom-made for stripped stars by employing the surface parameters given by the evolutionary models as input at the base of the atmosphere. The evolutionary models were used at the time when the stripped star had reached half-way through central helium burning (X_He, c = 0.5). We can use these models as a good approximation for the spectral properties throughout most of the stripped star phases. This is because the luminosity and the effective temperature do not change significantly during most of the time the stars are stripped. We show this in Appendix C. The spectral models for individual stripped stars are publicly available at the CDS database¹.

The shape of the spectral energy distribution and the emission rates of ionizing photons (see Table 1 of Paper II) depend on the assumed wind mass-loss rates, wind speeds, and wind clumping. These parameters are uncertain. Theoretical predictions are now available (e.g., Krtička et al. 2016; Vink 2017), but they have not yet been thoroughly tested against observations, because only very few stripped stars with sufficiently strong wind mass-loss have been identified and studied in detail so far (e.g., Groh et al. 2008). In Paper I, we showed that variations in wind mass-loss rate primarily affect the predicted emission rate of He II-ionizing photons, while the emission rates of H I- and He I-ionizing photons are not significantly affected. The mass-loss rates assumed in our models were chosen to smoothly connect the mass-loss rates of subdwarfs (Krtička et al. 2016) with the observed mass-loss rates of WR stars (Nugis & Lamers 2000). Our assumed mass-loss rates also match well with the observed mass-loss rate of the stripped star in the binary system HD 45166 (Groh et al. 2008). The recent theoretical predictions by Vink (2017) suggest that the mass-loss rates of stripped stars may be ten times lower than what we assume in this paper. The winds of stripped stars are likely not reaching close to the Eddington limit, in contrary to massive main-sequence and WR stars (cf. Bestenlehner et al. 2014). This suggests that the wind mass-loss rate from stripped stars is lower than that from WR stars and thus not well-described by the recipe for WR stars of Nugis & Lamers (2000). To establish which are the wind mass-loss rates from stripped stars, observations of a sample of stripped stars are necessary. If, as suggested by Vink (2017), the mass-loss rates from stripped stars indeed are lower than what the recipe from Nugis & Lamers (2000) predicts, it would likely imply an increase of the emission rates of He II-ionizing photons presented in this work. The emission rates of H I- and He I-ionizing photons are robust against wind uncertainties.

The wind parameters also affect the stellar spectral features. Higher wind mass-loss rate, slower winds, or more clumping result in stronger emission features as the stellar wind becomes denser. The wind speed could be slower than our assumptions. We assumed terminal wind speeds of 1.5 times the escape speed of the surface of the stripped star, which resulted in values of ∼1500 − 2500 km s⁻¹. The observed stripped star HD 45166 has an anisotropic wind that partially is slow, which gives rise to the strong emission lines the star exhibits (Groh et al. 2008). For a better understanding of the spectral features from stripped stars, an observed sample is necessary.

2.2. Modeling the contribution of stripped stars to a stellar population

We modeled the number and type of stripped stars that are present in a population as a function of time by taking a Monte Carlo approach. We first created a sample of stars by randomly drawing initial masses, M_init, from the initial mass function of Kroupa (2001), $\mathrm{d}N/\mathrm{d}{M}_{\mathrm{init}}\propto M_{\mathrm{init}}^{\alpha }$, where α = −1.3 for M_init <  0.5 M_⊙ and −2.3 for M_init >  0.5 M_⊙ until we reached a total stellar mass of 10⁶ M_⊙. We assumed mass limits of 0.1 M_⊙ and 100 M_⊙. We then chose which stars that have a companion star using the mass-dependent binary fraction of Moe & Di Stefano (2017) that follows closely the linear function f_bin = 0.09 + 0.63log₁₀(M_init/M_⊙) (we fit this polynomial to the data presented in Fig. 42 of Moe & Di Stefano 2017, for comparison purposes also see van Haaften et al. 2013). The mass of the companion stars were randomly drawn, such that the mass ratio, q = M_2, init/M_1, init, followed a flat distribution sampled between 0.1 and 1 (consistent with the observations of Kiminki & Kobulnicky 2012, Sana et al. 2012, and Moe & Di Stefano 2017, for early-type stars). After assigning companions and their masses, the total stellar mass increased and we therefore randomly removed single or binary stars until we reached a total stellar mass of 10⁶ M_⊙ again. The initial orbital periods of the binary systems were randomly drawn from a distribution that is flat in the logarithm of the period (e.g., Öpik 1924; Kouwenhoven et al. 2007; Moe & Di Stefano 2017). For systems where the most massive star of the system has a mass M_1, init ≥ 15 M_⊙, we used the distribution by Sana et al. (2012), which favors short-period systems. As a lower limit for the initial period, we chose the shortest period that allows both stars to fit inside their Roche-lobes at zero-age main-sequence. For the upper limit of the initial period, we followed Moe & Di Stefano (2017) and set 10^3.7 days. Moe & Di Stefano (2017) also presented period distributions for binaries, but with lower resolution in the period than Sana et al. (2012) do for the massive stars. We, therefore, combined the period distribution for massive stars presented by Sana et al. (2012) with the flat distribution that is consistent with the data presented by Moe & Di Stefano (2017) for lower mass stars.

We used evolutionary models of single stars to follow the radius evolution of the donor star and to determine when it will start to interact with its companion star. We created these models with MESA, using the same physical assumptions as we adopted for the models of binary stars (see Paper II). The moment mass transfer starts can then be determined by comparing the radius evolution of the most massive star with its Roche radius (Eggleton 1983). We predicted the further evolution from the initial mass ratio of each binary system and whether the donor star had developed a deep convective envelope at the time of interaction. We assumed that stable mass transfer occurs in systems with an initial mass ratio larger than a critical value (q_crit, MS = 0.65 and q_crit, HG = 0.4 for interaction initiated on the main-sequence and Hertzsprung gap following de Mink et al. 2007, and Hurley et al. 2002, respectively). For systems with a smaller initial mass ratio than the critical value and systems which have donor stars that have a deep convective envelope, we assumed that a common envelope developed. Since stars do not have well-defined cores during the main sequence, we assumed that stripped stars are not formed through common envelope evolution initiated during the main sequence evolution of the donor star. We also did not include stripped stars formed after or during central helium burning because the product is likely short-lived.

We used the classical α-prescription to determine whether the common envelope is successfully ejected or not (Webbink 1984). For this, we employed the standard value of one for the efficiency of the ejection of the common envelope, α_CE (see e.g., Hurley et al. 2002). We also assumed that the parameter λ_CE, which describes how strongly the envelope is bound to the core of the donor star (Dewi & Tauris 2000; Tauris & Dewi 2001), is 0.5, which is the average for Hertzsprung gap stars (see Appendix E of Izzard 2004). Dewi & Tauris (2000) showed that λ_CE is dependent on the stellar mass and evolutionary stage. They found that λ_CE decreases with increasing mass and that λ_CE typically is between 0.1 and 1 for Hertzsprung gap stars of initial masses between 3 M_⊙ and 10 M_⊙. The mass dependence of λ_CE suggests that lower mass donors easier eject a common envelope. However, since α_CE is considered to be a very uncertain parameter, we assumed the simple, mass-independent value for λ_CE.

We assumed that the stars that shared a common envelope coalesced if either the core of the donor star or the companion star filled their Roche-lobe during the in-spiral inside the common envelope. To control whether this occurs, we assumed that the core of the donor star had the same mass and radius as a stripped star that originated from a progenitor star with the same initial mass as the donor (see Table 1 of Paper II). For the companions, we assumed that they have the same mass and radius as on the zero-age main sequence. We used our models for stripped stars and for single stars and interpolated over the mass to find the radii of the stars in each system. In cases when the companion star had lower mass than our lowest mass model, we extrapolated to smaller radii, however, this never led to negative radii. Once mass transfer initiates within a common envelope, it likely leads to coalescence since the stars spiral closer together due to friction from the surrounding material. If coalescence occurred, we did not consider the further evolution of the star. The α-prescription likely oversimplifies the evolution of stars in a common envelope (see e.g., Ivanova et al. 2013). We therefore consider our approach approximate. To better constrain the formation of stripped stars through common envelope evolution, a large sample of observed stripped stars is needed.

Summarizing the above descriptions, we considered that stripped stars are formed through three main evolutionary channels: stable mass transfer initiated on the main-sequence evolution of the donor star, stable mass transfer initiated on the Hertzsprung gap evolution of the donor star, and un-stable mass transfer initiated on the Hertzsprung gap evolution of the donor star followed by a successful ejection of the common envelope. We assigned the masses of the stripped stars by interpolating between the masses of the progenitor stars and assuming that the stripped star masses are constant throughout the stripped phases, taken at X_He, c = 0.5 from the evolutionary models (see Sect. 2.1). We also interpolated over the progenitor masses in the evolutionary grid to determine the time of envelope-stripping and the duration of the stripped phases. This approach neglects that Case A type mass transfer can result in somewhat lower mass stripped stars (see e.g., Pols 1994). However, this evolutionary channel is responsible for less than a fourth of the stripped stars of intermediate or high mass. We therefore expect that the total effect is small for the ionizing emission. We interpolated the spectral models for stripped stars over the mass of the stripped star to find the best-fitting radiative representation for each stripped star. With the described approach, we assumed that the spectrum of a stripped star remains stationary throughout the stripped phase. This approximation is good for most of the time since stripped stars have roughly the same effective temperature and luminosity throughout the long lasting central helium burning phase (see Appendix C and Paper II). Our approach also assumes that stripped stars have the same properties as long as they are formed from a progenitor with the same initial mass and independent on the formation channel and initial orbital period and mass ratio. This assumption is also valid for most cases. The exception is for systems that interact on the main sequence, that have very low metallicity, or that have initially very long orbital periods (cf. Claeys et al. 2011; Yoon et al. 2017; Sravan et al. 2018).

We considered two different star formation histories. The first is an instantaneous starburst with initially 10⁶ M_⊙ of mass in stars. The second is a population in which stars form at a constant rate of 1 M_⊙ yr⁻¹. For the case of continuous star formation, we used the predictions for the co-eval stellar population and convolved the number of stripped stars, the spectral energy distribution and emission rates of ionizing photons over time. We performed the convolution every 1 Myr, which produces stochastic effects expected for a constant star formation rate of 1 M_⊙ yr⁻¹. In Appendix B we consider three additional, more realistic, star formation histories and we discuss how much the stripped stars in those populations affect the emission rates of ionizing photons.

2.3. Including the contribution from stripped stars to the model of a full stellar population

In the previous subsection, we described how we model the radiative contribution from stripped stars. To model the radiation from a realistic population, we also need to model the contribution of the remainder of the population, which includes primarily single stars and stars in binary systems that have not yet interacted.

STARBURST99 provides well-established models for the integrated spectra of stellar populations, including models of main-sequence stars, giant stars, and stars in more evolved stages of the stellar life (Leitherer et al. 1999, 2010). Most of the stars in a stellar population are main-sequence stars that have not yet interacted. Binary interaction primarily occurs at later evolutionary stages as the stars expand significantly after central hydrogen exhaustion and only mildly during the main-sequence. STARBURST99 thus constitutes a fair approximation for stars that have not interacted with a binary companion. However, including stripped stars implies that there should be fewer giant stars as a fraction of them have become stripped. Moreover, the companions to stripped stars are expected to have accreted material and thus become more massive and somewhat rejuvenated. Similar is expected for binary stars that have merged. This leads to a slight increase in the radiation in optical and UV wavelengths since the mass-gaining star, in most cases, is a main-sequence star (cf. Van Bever et al. 1999). We expect the total effect from mass-gainers, mergers, and the lack of giant stars on the emission in the optical and UV wavelengths to be small, compared to the total emission from the full stellar population.

We used the combination of models from STARBURST99 and our model for the contribution from stripped stars to represent the radiation of a full stellar population in which stripped stars are formed. We made our models publicly available on the STARBURST99 online interface², providing the addition from stripped stars to the spectral energy distribution, the emission rates of H I-, He I-, and He II-ionizing photons, and the high-resolution UV and optical spectra. In this study, we compare the contribution from stripped stars to the version of STARBURST99 that uses the initial mass function from Kroupa (2001) with mass limits of 0.1 M_⊙ and 100 M_⊙, that is, the same as what we employed for the stripped stars in the population. We chose to compare to the STARBURST99 models that use non-rotating stellar evolutionary models from the Geneva grids (Levesque et al. 2012) and spectral models from WM-BASIC for OB-stars (Pauldrach et al. 2001; Leitherer et al. 2010), CMFGEN for WR stars (Hillier & Miller 1998; Smith et al. 2002), and BASEL v3.1 for later type main sequence stars, cooler stars, and red supergiants (Lejeune et al. 1997). When comparing our models of various metallicity with STARBURST99, we used the STARBURST99 models of Z = 0.014, 0.008, 0.002, and 0.001 together with our models with metallicity of Z = 0.014, 0.006, 0.002, and 0.0002, respectively (see also Table A.1). We consider that the difference in metallicity between the two model sets is small and expect that the difference when using models with exactly the same metallicity is also small. We consider that the combination of the stripped star models and STARBURST99 is a good assumption for radiation with wavelengths shorter than ∼5000 Å. For the model to be accurate at longer wavelengths, we would need to decrease the radiation from giant stars to compensate for stars that we assume have become stripped. Giant stars emit their radiation primarily at wavelengths longer than ∼5000 Å. We expect the decrement of radiation at these long wavelengths to be at maximum about 30% (as also suggested by Eldridge et al. 2017, see e.g., their Figs. 5 and 15). This is approximately the fraction of massive stars that get stripped (Sana et al. 2012) and thus the fraction of giant stars that should be missing. This topic is beyond the scope of this paper, but we hope to address it more in detail at a later stage.

In several cases in this study, we compare our results to those from the code BPASS (version 2.2.1, or Tuatara, see Stanway & Eldridge 2018, but also Eldridge & Stanway 2009, 2012; Eldridge et al. 2017). We used the version of the code that accounts for binary interactions and employs a similar initial mass fraction as we do (α = −2.35 for M_init >  0.5 M_⊙ and α = −1.3 for M_init <  0.5 M_⊙, Kroupa et al. 1993) and the same mass limits as we do. In Tuatara, the binary fraction, period distribution, and mass ratio distribution have been adapted following the data presented in Moe & Di Stefano (2017). The BPASS code uses stellar evolutionary models computed with a modified version of the STARS code that accounts for binary interaction, and the spectral models that are used are from BASEL v3.1 and v2.2 for most part of the stellar evolution, WM-BASIC for O-stars, and POWR models for WR stars defined as when X_H, s <  0.4 and log₁₀T_⋆ ≥ 4.45. We therefore expect that the emission from stripped stars of masses above ∼0.5 M_⊙ is represented using models from POWR in BPASS. We present a comparison between our models and BPASS for individual stripped stars in Appendix D, where we also discuss differences in the modeling of stripped stars in our simulations and in BPASS. We compare our models with metallicity Z = 0.014, 0.006, 0.002, and 0.0002 with the models of BPASS with metallicity Z = 0.014, 0.006, 0.002, and 0.0001, respectively (see Table A.1).

3. The number of stripped stars in stellar populations

The model for a population of stripped stars gives estimates for the number of stripped stars that are present in stellar populations, what are their properties and how they were formed. The number of stripped stars in a population is expected to be much smaller than the number of main sequence stars since the helium-core burning phase is about ten times shorter than the main sequence evolution and only a fraction of stars loose their envelopes through binary interaction.

We show the number of stripped stars that are formed as a function of time for the two star formation histories and for solar metallicity in Fig. 1. For the case of the instantaneous starburst of 10⁶ M_⊙ (panel a), between ∼150 and ∼400 stripped stars are present in the population between 10 Myr and 2 Gyr after the starburst. Because we only consider progenitors with masses between 2 M_⊙ and 20 M_⊙, the formation of stripped stars start after about 10 Myr, when the 20 M_⊙ stars evolve off of the main sequence and can become stripped. The number of stripped stars in the population increases with time because stars of initially lower mass, that are also more common, evolve off of the main sequence and can get stripped.

Fig. 1. Number of stripped stars shown as function of time for a population of solar metallicity. Panel a: numbers for a co-eval stellar population with initially 106 M⊙ in stars. We highlight the three formation channels that we consider to lead to the formation of stripped stars: stable mass transfer initiated during the main sequence evolution of the donor star (Case A, RLOF in light green), stable mass transfer initiated when the donor star passes over the Hertzsprung gap (Case B, RLOF in medium green), and successful ejection of a common envelope that developed during interaction initiated during the donor’s Hertzsprung gap evolution (Case B, CEE). For comparison, we show the number of low-luminosity WR stars predicted by BPASS in light gray (Eldridge et al. 2017). We indicate the expected time of WR star formation through solely stellar winds and the time of subdwarf formation through low-mass stellar evolution above the panel. Panel b: case of continuous star formation at a rate of 1 M⊙ yr−1. We show the parts of the population of stripped stars that are created from progenitor stars with initial masses between 10 − 20 M⊙ in purple, 5 − 10 M⊙ in pink, 3 − 5 M⊙ in salmon, and 2 − 3 M⊙ in peach.

Panel a of Fig. 1 shows the fraction of stripped stars that are formed through the three considered formation channels: stable mass transfer initiated during the main sequence evolution of the donor star (labeled Case A, RLOF), stable mass transfer initiated during the Hertzsprung gap evolution of the donor star (labeled Case B, RLOF), and unstable mass transfer initiated during the Hertzsprung gap evolution of the donor star (labeled Case B, CEE). The figure shows that Case B type stable mass transfer is the most common evolutionary channel for stripped stars. This is expected since massive stars swell up significantly when the hydrogen shell is ignited. The expansion during the Hertzsprung gap decreases with lower stellar mass. This can be seen in the figure since the relative importance of Case A type mass transfer becomes a more important formation channel for stripped stars that form later and thus originate from lower mass stars. In our model, the contribution from common envelope evolution is relatively small and primarily occurring for stars more massive than ∼4 M_⊙. However, depending on uncertain details in the treatment of common envelope evolution (see Sect. 2.2) it is likely that the contribution from common envelope evolution is slightly different from our predictions. It is possible that our prediction for the fraction of stripped stars that form via Case A type mass transfer is slightly exaggerated. The reason is that we do not account for conservative mass transfer, which might lead the accreting star to expand sufficiently and form an overcontact binary that eventually leads to coalescence. However, this is expected to occur only for very short periods (see e.g., Almeida et al. 2015, for a massive example).

For reference, we also show the number of low-luminosity WR stars produced in BPASS in panel a of Fig. 1. These are defined as stars that have luminosities of log₁₀(L/L_⊙) < 4.9, temperatures of log₁₀(T_⋆/K)≥4.45 and surface abundance of hydrogen of X_H, surf ≤ 0.4 (see the entry for WNH in Table 3 of Eldridge et al. 2017, we note that the WNH entry in this table includes all low-luminosity WR stars). Since no main sequence stars should fit with these constraints and WR stars formed from wind mass loss likely have higher luminosity, it is plausible that these are stripped stars. We expect that the subdwarfs with lower mass than ∼0.5 M_⊙ are excluded since they have surface temperatures that are lower than 10^4.45 K. Such subdwarfs are primarily formed after ∼400 Myr, which could explain why the number of low-luminosity WR stars predicted by BPASS drops around that time. The predicted number of stripped stars is larger in BPASS than in our model for times before around 500 Myr. The reason could be that BPASS treats mass transfer and common envelope evolution different (cf. Eldridge et al. 2017). However, our predictions and the predictions from BPASS are of the same order with just a factor of about two difference at maximum.

For the case of continuous star formation (panel b of Fig. 1), the number of stripped stars steadily increases with time as lower mass stars start to interact with binary companions. The figure shows that a stellar population that has formed stars at a constant rate of 1 M_⊙ yr⁻¹ for 100 Myr contains roughly 10 000 stripped stars. About 3 000 of these originate from progenitor stars more massive than 10 M_⊙ and are therefore likely to be hotter than 70 kK and more luminous than 10⁴ L_⊙ (following Table 1 of Paper II). The figure shows that no or few subdwarfs are predicted to be present in such a young stellar population, but that the population of stripped stars is dominated by subdwarfs at later times. After about 1 Gyr of constant star formation, there is about half a million stripped stars present in the population, with the majority of these being subdwarfs (over 90%).

The effect of metallicity is small on the predicted number of stripped stars. We describe the results for lower metallicity environments in Appendix A.

4. The impact of stripped stars on the spectral energy distribution

In this section, we describe the effect of stripped stars on the total spectral energy distribution of a stellar population. Although the contribution from stripped stars to the total bolometric luminosity is small, the hard, ionizing spectra of stripped stars do significantly change the shape of the spectral energy distribution. The differences mainly occur in the extreme ultraviolet. We discuss the effects for co-eval stellar populations in Sect. 4.1 and for the case of continuous star formation in Sect. 4.2.

4.1. Predictions for co-eval stellar populations

Figure 2 shows how the spectral energy distribution is affected by the presence of stripped stars. The panels of the figure show snapshots taken 3, 11, 20, 50, 100, and 800 Myr after an instantaneous starburst of 10⁶ M_⊙, assuming solar metallicity.

Fig. 2. Spectral energy distribution of a co-eval stellar population. The contribution from stripped stars is highlighted with a blue line. The parts of the spectra that are H I-, He I-, and He II-ionizing are shaded in blue, while the UV is shaded in purple. For comparison, we show the spectral energy distribution of a population containing only single stars using STARBURST99 (green line), which can be interpreted as the contribution from the remaining stars in the stellar population. We also show the predictions from BPASS (gray line), where the effects of binary interaction are included. The panels correspond to different times after the instantaneous starburst of 106 M⊙. The model shown here assumes solar metallicity (see Appendix A for lower metallicity models).

Stripped stars dominate the ionizing output from the stellar population after about 10 Myr and up to at least 100 Myr. At 3 Myr after starburst, shown in Fig. 2a, the ionizing radiation originates from massive main-sequence stars. Stars stripped in binaries are not yet present, because it takes time for them to form. Their progenitors, the donor stars, need time to evolve and swell up to fill their Roche-lobe. Stripped stars are therefore formed with a delay corresponding roughly to the main-sequence lifetime of the progenitor star. For a 20 M_⊙ progenitor star, the main-sequence lifetime, and thus the delay with which stripped stars with such progenitors form, is about 10 Myr. Also, stars that initiate mass transfer during the main-sequence evolution result in stripped stars after a time delay corresponding to the main-sequence lifetime of the donor star. The reason for this is that the mass transfer rate is slow during the main-sequence evolution and it is not terminated until central hydrogen depletion is reached.

Stars stripped in binaries are created over an extended period of time. The main-sequence lifetime, and therefore the time delay with which stripped stars are created, varies with the mass of the progenitor star. The stars that form stripped stars 11, 20, 50, 100, and 800 Myr after starburst have initial masses of about 18, 12, 7, 5, and 2 M_⊙. The resulting stripped stars have masses of about 7, 4, 2, 1, and 0.5 M_⊙, respectively. The mass range of the stripped stars that are present at each point in time is small since the duration of the stripped phase is about 10% of the main-sequence lifetime of the progenitor star, which is the time delay with which the stripped stars are created. The temperature of stripped stars decreases with decreasing stellar mass as seen in Table 1 of Paper II, which shows that a 7 M_⊙ WR star has a temperature of 100 000 K, while a 1 M_⊙ subdwarf has a temperature of 40 000 K. This decrease in temperature causes their contribution to the integrated spectrum to become softer with time as the mass of the stripped stars that are present decreases.

As seen in Sect. 3 and Fig. 1, the number of stripped stars in a population increases as time proceeds. Despite the increase in their total numbers with time, we find that the total bolometric luminosity produced by stripped stars decreases. This is because the luminosity of individual stripped stars is a steep function of mass.

About 500 Myr after starburst, the stripped stars no longer significantly contribute with ionizing photons, according to our models. The reason for this is that the stripped stars that are still present at these late times are subdwarfs. These subdwarfs are affected by diffusion processes, which alter their surface composition and structure (for a discussion see Sect. 3 of Paper II). The result is an increase of the abundance of hydrogen at their surfaces, which creates a sharp cut-off of the spectral energy distribution at the Lyman limit. The integrated spectra are still significantly different from what is expected for a population of single stars, see panel f of Fig. 2. At these late times, we expect that white dwarfs contribute with ionizing radiation (Panagia & Terzian 1984). However, more detailed modeling is needed to further understand the relative contributions of ionizing photons in late starbursts.

For comparison, we also show the spectral energy distributions predicted by the BPASS models in Fig. 2. The BPASS predictions for the shape of the ionizing part of the spectral energy distribution match well with our predictions for populations younger than about 50 Myr. It is interesting to note that despite our simpler approach, similar results are reproduced for starbursts of these ages. After this time, the BPASS models predict that the ionizing radiation is harder than what we find in our simulations. The reason is likely the adopted atmosphere models for central stars in planetary nebulae in BPASS, which are represented by hot WR star models (see Gräfener et al. 2002, for the models used) and the implementation of emission from X-ray binaries (priv. comm. with JJ Eldridge). After about 300 Myr, the effects of diffusion on the stellar surface cause our stripped star models to appear slightly cooler, which could contribute to the difference between our models and the predictions of BPASS at late times. However, it is not likely that this is the full explanation, since the hard ionizing emission the BPASS predicts at 800 Myr shown in Fig. 2f requires very high temperatures, possibly similar to or higher than the massive WR stars.

The effect of metallicity on the shape of the spectral energy distribution is relatively small. At lower metallicity, stripped stars are slightly cooler (see Paper I for a discussion), which means that the ionizing part of the spectral energy distribution is slightly softer. The first notable differences occur at very low metallicities, when Z <  0.002, see Appendix A and Fig. A.2.

4.2. Predictions for continuous star formation

In Sect. 4.1 we showed that, for co-eval stellar populations, stripped stars make a distinct contribution to the ionizing spectra at late times, while massive main-sequence stars dominate at early times. In Fig. 3, we show the integrated spectrum for the idealized case of constant star formation. The spectrum is for a population in which stars have formed at a constant rate of 1 M_⊙ yr⁻¹, and for a prolonged period of time, here chosen to be 500 Myr. At this time, the ionizing spectrum has reached equilibrium. The figure shows that the spectrum is heavily dominated by single stars at almost the entire wavelength range. Also, the emission of H I-ionizing photons is dominated by massive stars. The contribution from stripped stars to the total bolometric luminosity is negligible. They only dominate the emission of the hardest He II-ionizing photons. We note that our predictions for this part of the spectrum are uncertain and depend on the treatment of the stellar winds.

Fig. 3. Spectral energy distribution of a stellar population in which star formation has taken place at a constant rate of 1 M⊙ yr−1 for 500 Myr. Otherwise similar to Fig. 2. We show models with solar metallicity.

Realistic stellar populations do not follow constant star formation, but are rather made up of a number of starbursts in star clusters that started to form stars at different times. The result is a star formation rate that varies with time. For such, more realistic, star formation histories, the relative contribution from stripped stars depends on the recent star formation activity. For populations that formed stars at a significant rate in the very recent past, ≲10 Myr, massive stars likely dominate the output of ionizing photons. For populations that did not form stars very recently, we expect stripped stars to play a significant role. Since stripped stars emit ionizing photons for more than ten times longer than massive stars, it is likely that they will significantly contribute with the emission of ionizing photons in-between starbursts in stellar populations, see Appendix B where we consider how the emission rate of ionizing photons is affected by the presence of stripped stars for three additional star formation histories.

5. Impact on the budget of ionizing photons

In this section, we discuss the emission rates of ionizing photons: Q_0, pop, Q_1, pop, and Q_2, pop for H I-, He I-, and He II-ionizing photons, respectively. We use the common definition of the emission rate of ionizing photons as the number of emitted photons with wavelengths shorter than the ionization threshold of the considered atom or ion, which can be calculated from the emitted luminosity, L_λ, in the following way:

$$\begin{array}{c}\hfill Q_i={\displaystyle \frac{1}{\mathit{hc}}}{{\displaystyle \int }}_0^{{\lambda }_i}\lambda L_{\lambda }\mathrm{d}\lambda ,\end{array}$$(1)

where h is the Planck’s constant, c is the speed of light, and the subscript i refers to the considered ion. This method provides a good approximation for the emission rates of ionizing photons if the surrounding medium is sufficiently dense. The probability that an ionizing photon will lead to ionization once it encounters an atom or ion is decreasing with increasing photon energy (Osterbrock & Ferland 2006). This decreasing ionization probability gives rise to a 100 times longer mean-free path for a photon with a wavelength of 228 Å compared to one of 912 Å in a medium containing only hydrogen (the mean-free path can be expressed as ⟨l⟩=1/(nσ), where n is the number density of the surrounding medium and σ is the ionization cross-section, which is wavelength dependent in the following way for hydrogen: σ = σ₀(λ₀/λ)³). In a typical density for H II regions of 10² cm⁻², we calculate that these mean-free paths are of order 0.0005 pc and 0.05 pc, which are negligible length-scales in terms of the size of star-forming regions. We, therefore, consider the definition of the emission rates of ionizing photons as a realistic assumption (see however McCandliss & O’Meara 2017 for cases when low densities are considered).

We compare the expected contribution from the stripped stars with that from the massive single stars in the same stellar population. For this comparison, we use STARBURST99 to represent the emission from the massive stars, as described in Sect. 2.3. We also show the predictions from the BPASS models.

5.1. Predictions for co-eval stellar populations

In Fig. 4a, we show the emission rate of H I-ionizing photons for a co-eval stellar population as a function of time after starburst. The massive stars in the stellar population primarily emit their H I-ionizing radiation during the first 5 Myr. Stripped stars play a significant role at later times.

Fig. 4. Emission rates of ionizing photons from stellar populations as a function of time. The contribution from stripped stars is shown in blue, while green represents the contribution from the massive main sequence and WR stars in the stellar population (models from STARBURST99, see Sect. 2.3). For reference, we show the predictions from BPASS models in which binary interactions are included using gray color. The solid lines correspond to the predictions from a population with solar metallicity, while the shaded regions of the same color represent the effect of lowering the metallicity, from solar down to Z = 0.0002 for the stripped stars (see Appendix A and Table A.1). Left column: emission rates from a co-eval stellar population with initially 106 M⊙ in stars, and right column: emission rates from a stellar population in which stars form at the constant rate of 1 M⊙ yr−1. Top, middle, and bottom rows: emission rates of H I-, He I-, and He II-ionizing photons, respectively (Q0, pop, Q1, pop, and Q2, pop).

At ∼10 Myr, when the first stripped stars have been created in our simulations, they emit ionizing photons with a rate of ∼10⁵¹ photons per second for a burst of 10⁶ M_⊙ formed stars. This is an emission rate of about a factor of ten higher than what massive single stars produce at that time. About 100 Myr after starburst, the emission rate of H I-ionizing photons from stripped stars has decreased to ∼10⁴⁹ s⁻¹, which corresponds to the typical emission rate of one WR-star or one massive O-star (Smith et al. 2002; Martins et al. 2005; Groh et al. 2014) and is several orders of magnitudes higher than expected from single stars with an age of 100 Myr. The emission of H I-ionizing photons keeps decreasing with time, and the trend suddenly steepens around 300 Myr. This is because the stripped stars that are present in the population are subdwarfs that have hydrogen-rich surfaces as a result of diffusion processes, as we discussed briefly in Sect. 4.1, see also Sect. 3 of Paper II. This results in a sharper Lyman cut-off in their spectra and thus a drop-off in the contribution of H I-ionizing photons. After ∼1 Gyr, the total ionizing emission rate is similar to that of one early B-type star (Q_0, pop ∼ 10⁴⁵ s⁻¹, cf. Smith et al. 2002).

Metallicity does not significantly affect the emission rate of H I-ionizing photons from stripped stars, as can be seen from the width of the shading in Fig. 4a. Because stripped stars have high temperatures (T_eff >  20 000 K) at all metallicities, the radiation peaks in the H I-ionizing wavelengths and Q_0, pop is, therefore, less dependent on metallicity. For stripped stars with progenitor masses that are lower than about 4 M_⊙, the emission rate of H I-ionizing photons increases with decreasing metallicity. However, because the time of stripping and duration of the stripped phase decreases with decreasing metallicity, the total emission rate of H I-ionizing photons from stripped stars in a population remains similar. The impact of metallicity on the H I-ionizing emission from main-sequence stars is larger than for stripped stars. The H I-ionizing emission rate from massive stars increases by up to a factor of five when metallicity is lowered. The boost of ionizing photons at late times due to the presence of stripped stars is independent of metallicity.

The presence of stripped stars appear to cause a sudden increase in H I-ionizing emission around 10 Myr. This is likely an artifact of our models not reaching higher progenitor masses than 20 M_⊙ and therefore do not account for stripped stars that are formed earlier than 10 Myr.

For comparison, we also show predictions from the BPASS models in Fig. 4a. They closely follow the predictions of STARBURST99 during the first 5 Myr but show a boost of ionizing photons at later times, similar to our models of stripped stars (see also Stanway et al. 2016 and Wofford et al. 2016). Our predictions for the contribution of stripped stars follow the trend of the BPASS models after 10 Myr. The BPASS models predict a shallower decrease after ∼300 Myr compared to our models. We expect that this is due to differences in how we treat the atmospheres of stars in binary systems.

Stripped stars emit He I-ionizing photons at a rate that is about five times lower than that of H I-ionizing photons. In Fig. 2 the emitted luminosities appear to be of similar order for H I- and He I-ionizing radiation. However, the difference in the emission rates of He I- and H I-ionizing photons is larger since photons with shorter wavelengths are more energetic. The difference in emission rates of He I- and H I-ionizing photons can be seen by comparing the panels b and a in Fig. 4. The diagrams show that the decline of Q_1, pop closely follows that of Q_0, pop. This is because the temperature of the stripped stars that are responsible for the ionizing emission does not change sufficiently to significantly affect the relative emission of He I- to H I-ionizing photons.

Once they are formed, stripped stars dominate the emission of He I-ionizing photons. The He I-ionizing emission rate predicted for single stars by STARBURST99 decreases much more steeply with time than it does for stripped stars. Comparing to the emission of He I-ionizing photons from the main-sequence and WR stars, we see that stripped stars boost the emission rate by four orders of magnitudes already 10 Myr after starburst.

Changing the metallicity does not significantly alter the emission rate of He I-ionizing photons from stripped stars. At late times, after around 300 Myr the evolution of Q_1, pop shows a feature, a temporary rise of the emission rate. The predictions for different metallicities deviate in this region as can be seen from the spread in the shading. This feature originates in the treatment of diffusion processes in surface layers of subdwarfs as discussed earlier.

The prediction of Q_1, pop from BPASS closely follows the predictions from stripped stars between 10 Myr and 100 Myr after starburst. At early times, ≲5 Myr, we see that the predictions from BPASS closely follow those from STARBURST99 at solar metallicity. However, for low metallicities, Z ≤ 0.004, the BPASS models show an increase at these early times. The increase is particularly prominent between 3 Myr and 20 Myr after starburst. This is a result of the treatment of stars that have gained mass or merged through binary interaction, causing them to rotate rapidly. Rapid rotation is thought to give rise to mixing processes in the stellar interior, causing the stars to evolve chemically homogeneously (Maeder 1987; Yoon & Langer 2005; Cantiello et al. 2007). Chemically homogeneous stars are thought to be very hot and bright and, if they exist, they can contribute with a significant amount of ionizing photons (Brott et al. 2011; Köhler et al. 2015; Szécsi et al. 2015; Kubátová et al. 2019). Other sources that could affect the emission of H I-ionizing photons from stellar populations are white dwarfs (Panagia & Terzian 1984), rejuvenated mass-gainers (Chen & Han 2009; de Mink et al. 2014), and central stars in planetary nebulae or post-AGB stars (e.g., Miller Bertolami 2016).

Our models predict that stripped stars are important contributors of the He II-ionizing photons emitted by stellar populations. Fig. 4c shows that stripped stars reach emission rates of ∼10^48.5 s⁻¹, which is similar to the emission rates from massive stars. Stripped stars maintain their high emission rate for longer because of their longer lifetimes and because they form from progenitors with a range of ages. The He II-ionizing emission is strongly temperature dependent. The decline of He II-ionizing emission with time is, therefore, steeper for Q_2, pop than for Q_0, pop or Q_1, pop. After about 50 Myr, the emission rate of He II-ionizing photons from stripped stars falls below to 10⁴⁵ s⁻¹.

Metallicity affects the emission rate of He II-ionizing photons from stripped stars, causing variations that span two orders of magnitude. The largest deviation occurs at the lowest metallicities, where we find that the stripping process fails to remove all the hydrogen, resulting in cooler stars (Paper I, Yoon et al. 2017). The emission rate of He II-ionizing photons is determined in the steep Wien-part of the stellar spectrum, meaning that a small shift in temperature results in a large shift in Q_2, pop. The emission rate of He II-ionizing photons is also dependent on parameters of the stellar wind and can increase by several orders of magnitudes if, for example, the mass-loss rate is just a few factors lower (Paper I).

A feature in the otherwise smooth decline of Q_2, pop for stripped stars is seen in Fig. 4c about 15 Myr after starburst. This feature occurs because of the change of the treatment of wind clumping in the atmospheres of the stripped stars (see Paper II for details). Future detailed spectral analysis of observed stripped stars is necessary to constrain the properties of the stellar winds.

The BPASS predictions for Q_2, pop roughly follow the predictions of STARBURST99 for the first ∼3 Myr. Between 10 Myr and 50 Myr, they are consistent with our predictions for stripped stars. After that, BPASS predicts an emission rate that is several orders of magnitude higher and that stays high for the remaining considered time, probably because of the adopted atmosphere models for central stars in planetary nebulae and the inclusion of the emission from X-ray binaries (priv. comm. with JJ Eldridge). The models of BPASS predict that metallicity variations give rise to a large range of values for Q_2, pop before about 10 Myr has passed. Similar to the case of Q_1, pop, we believe that these variations are due to rapidly rotating stars. Stars that likely play important roles as emitters of He II-ionizing photons are accreting white dwarfs and X-ray binaries (Chen et al. 2015; Fragos et al. 2013; Madau & Fragos 2017; Schaerer et al. 2019). These binaries reside in late evolutionary stages of interacting binaries, where one star has already died, possibly after interaction already occurred, and the second star now fills its Roche-lobe and transfers material to the compact object.

5.2. Predictions for continuous star formation

We show the predicted emission rates of ionizing photons for continuous star formation in Figs. 4d, e, and f. Once equilibrium is reached, stripped stars emit H I-ionizing photons at a rate of 10^51.8 s⁻¹. This is about 5% of the emission rate from massive stars (10^53.1 s⁻¹). As Fig. 4d shows, the massive stars dominate the emission of H I-ionizing photons in the case of continuous star formation.

Stripped stars emit He I-ionizing photons at a rate that is about five times lower than the emission rate from massive stars. The contribution from stripped stars corresponds to about 15% of the total stellar emission of He I-ionizing photons. Depending on the assumed metallicity, the contribution varies somewhat (between 10 and 20%). At lower metallicity, the massive main sequence stars are hotter and therefore contribute with a larger fraction compared to the stripped stars.

The emission of He II-ionizing photons shown in Fig. 4f is dominated by the contribution from stripped stars as they emit about five times more He II-ionizing photons than the massive stars. Stripped stars reach emission rates of He II-ionizing photons of 10⁴⁹ s⁻¹, while the massive stars only reach an emission rate of 10^48.3 s⁻¹. We note that the emission rate of He II-ionizing photons is sensitive to metallicity variations. At very low metallicity (Z ∼ 0.0002), the emission rate from stripped stars decreases by two orders of magnitude. Conversely, the emission rate from massive stars increases by a factor of five at very low metallicity.

6. Implications for observable quantities

In this section, we discuss the implications of accounting for stripped stars for various observable quantities commonly used to describe unresolved stellar populations. We summarize the values for the considered quantities in Table 1 at several snapshots after a starburst of 10⁶ M_⊙ and also for continuous star formation, taken after 500 Myr.

6.1. Diagnostics of the budget of ionizing photons

6.1.1. Production efficiency of ionizing photons, ξ_ion

The production efficiency of ionizing photons is a quantity that can be measured observationally. It relates the emission rates of ionizing photons to the UV luminosity and is, therefore, a parameter that describes the strength of the ionizing emission independent on the stellar mass or star formation rate of the population. The production efficiency of hydrogen-ionizing photons has already been measured for a large number of unresolved stellar populations (Robertson et al. 2013; Stark et al. 2015; Bouwens et al. 2016; Matthee et al. 2017; Shivaei et al. 2018).

The production efficiency of ionizing photons is defined as follows:

$$\begin{array}{c}\hfill {\xi }_{\mathrm{ion}}={\displaystyle \frac{Q_{\mathrm{pop}}}{L_{\nu }(1500\mathrm{\AA })}},\end{array}$$(2)

where Q_pop is the emission rate of ionizing photons and L_ν(1500 Å) is the luminosity at the wavelength 1500 Å in units of erg s⁻¹ Hz⁻¹. We use Q_0, pop, Q_1, pop, and Q_2, pop in (2) to calculate the production efficiencies of H I-, He I-, and He II-ionizing photons, which we refer to as ξ_ion, 0, ξ_ion, 1, and ξ_ion, 2, respectively. We average the UV luminosity between 1450 Å and 1550 Å to estimate the continuum luminosity and avoid fluctuations caused by spectral features.

For co-eval stellar populations, the ionizing radiation that stripped stars produce causes the production efficiencies of ionizing photons to remain at high values for much longer than what is predicted from single star populations. The integrated spectrum is harder if stripped stars are present, which can be recognized by comparing the production efficiency of H I-ionizing photons with either that of He I- or He II-ionizing photons. This is visualized in the hardness diagrams shown in Fig. 5.

Fig. 5. Diagrams showing the hardness of the ionizing part of the spectrum of stellar populations, constructed by the production efficiencies of H I- and He I- (panel a) or He II-ionizing photons (panel b) (ξion, 0, ξion, 1, ξion, 2). The colored regions represent the time spans indicated by the color bars for co-eval stellar populations and cover the predictions from all metallicities (from solar metallicity down to Z = 0.0002 for stripped stars, see Appendix A and Table A.1). Green shades show predictions for young, single star populations, while the purple and red shades represent the hardness of stellar populations in which stripped stars are included. The case of constant star formation is shown with solid lines and markers (pluses in black for single stars and circles in red for when stripped stars are included), taken after 500 Myr. Above the diagrams, we show the distribution of measured ξion, 0 for three samples of observed unresolved stellar populations from intermediate to high-redshift by Bouwens et al. (2016), Matthee et al. (2017), and Shivaei et al. (2018).

Figure 5a shows that the same range of values for the production efficiency of H I-ionizing photons can be produced by either massive stars in a stellar population younger than 10 Myr or stripped stars in a stellar population that is up to ten times older (cf. Wilkins et al. 2016). In the figure, a clear separation is visible between a population that contains stripped stars and one that does not. The difference is due to the harder ionizing spectra that stripped stars introduce, which shift the production efficiencies of helium-ionizing photons to higher values relative to what is expected from a single star population. In the case of ξ_ion, 2, the separation is several orders of magnitude. For continuous star formation, the role of stripped stars is small for the production efficiencies of ionizing photons but could be relevant for ξ_ion, 2, as can be seen from Fig. 5b and in Table 1.

For reference, we also show three distributions of measured ξ_ion, 0 in observational samples of galaxies on top of the diagrams in Fig. 5. These samples are of distant galaxies with various classifications and span a range of redshifts (Bouwens et al. 2016; Matthee et al. 2017; Shivaei et al. 2018, see also Tang et al. 2019 and Chevallard et al. 2018). The observed galaxies have a broad range of ξ_ion, 0. For consistency, we chose to show the samples with the dust-correction made in the same way when measuring the UV luminosity (assuming a Calzetti et al. 1994 slope for the dust extinction). We note that the method to account for the dust correction has an impact on the distribution of the estimated ξ_ion, 0 (Matthee et al. 2017, see also Hao et al. 2011; Murphy et al. 2011).

An interesting test of the underlying stellar population would be to measure the production efficiencies of helium-ionizing photons for observed galaxies, which would enable to place them individually in the hardness diagrams. We believe that this would provide a very valuable test for the models of stellar populations, and in particular for the impact of binary stellar evolution.

6.1.2. Ionization parameter, U

The ionization parameter is traditionally used to quantify the degree of ionization a stellar population causes on the surrounding nebula as it compares the flux of ionizing photons to the density of the surrounding gas (e.g., Osterbrock 1989). As a result, populations with the same ionization parameter typically show very similar nebular spectra even though they may have different stellar masses or star formation rates (e.g., Dopita et al. 2000; Nakajima & Ouchi 2014).

We follow the standard definition of the ionization parameter, U, (e.g., Shields 1990, see also Kewley et al. 2013), where the isotropically emitted ionizing radiation from a central source is compared to the density of the gas surrounding the source:

$$\begin{array}{c}\hfill U\equiv {\displaystyle \frac{Q_{0,\mathrm{pop}}}{4\pi R_{\mathrm{S}}^2{n}_{\mathrm{H}}c}}={\displaystyle \frac{{\alpha }_B^{2/3}}{3^{2/3}c}}{\left({\displaystyle \frac{Q_{0,\mathrm{pop}}{ϵ}^2{n}_{\mathrm{H}}}{4\pi }}\right)}^{1/3}·\end{array}$$(3)

In (3), n_H is the number density of hydrogen in the gas, ϵ is the volume filling factor of the gas, α_B is the recombination coefficient for hydrogen, and c is the speed of light. In the last equality of (3), we expanded the Strömgren radius, $R_{\text{S}}={[3Q_{\text{0},\text{pop}}/(4\pi n_{\text{H}}^2{\alpha }_Bϵ)]}^{1/3}$ (Strömgren 1939). When calculating the ionization parameter, we account for clumping in the nebula by assuming ϵ = 0.1, following Zastrow et al. (2013). We assume a typical gas temperature of 10 000 K, which leads α_B = 2.6 × 10⁻¹³ cm³ s⁻¹ (Case B type recombination, Osterbrock & Ferland 2006). We adopt the emission rates of ionizing photons presented in Sect. 5 and assume a range of gas densities from n_H = 10 to 10⁴ cm⁻³.

Our models show that stripped stars increase the ionization parameter for stellar populations in which star formation has stopped at least 10 Myr ago, as shown in Fig. 6. The stripped stars cause the ionization parameter to remain at high values for an extended time period. We note, however, that it is likely that the density of the surrounding gas changes with time, which will consequentially alter tthe effect of stripped stars on the ionization parameter. In the case of constant star formation, stripped stars increase U by about 2%.

Fig. 6. Top: ionization parameter computed for co-eval stellar populations as a function of time and for constant star formation, taken after 500 Myr. We show the predictions for stellar populations containing only single stars in gray shades and the model when stripped stars are included in purple shades. For constant star formation, the contribution from stripped stars is about 2% and therefore we do not show the markers for single star populations as they overlap. We show measurements of the ionization parameter for groups of observed galaxies to the right of the diagram (Moustakas et al. 2010; Nakajima & Ouchi 2014). Bottom: table explaining which gas density and stellar mass for co-eval stellar populations that correspond to which contour in the diagram using numbers as labels. In the case of constant star formation, the numbers are correlated with combinations of the gas density and star formation rate instead.

The ionization parameter depends on the gas density and the emission rate of ionizing photons. The latter is closely related to the stellar mass in case of co-eval stellar populations and the star formation rate in case of constant star formation. We, therefore, compute the evolution of the ionization parameter for a range of gas density, stellar mass, and star formation rate combinations, as seen in Fig. 6. For reference, we show the measured ranges of several samples of observed galaxies as vertical bars on the right side of the figure (Moustakas et al. 2010; Nakajima & Ouchi 2014). These galaxies are grouped roughly according to their properties and redshift. The local galaxy samples are the SDSS galaxies, the Green Pea galaxies (GPs), the Lyman Continuum (LyC) leakers, and the Lyman-Break Analogs (LBAs) summarized by Nakajima & Ouchi (2014) together with the local galaxy sample presented by Moustakas et al. (2010). The samples of distant galaxies are the z ∼ 1 galaxies, the Lyman Break Galaxies (LBGs) and the Lyman Alpha Emitters (LAEs) summarized in Nakajima & Ouchi (2014). The individual galaxies are likely to experience a bursty star formation rate that can resemble the behavior of single starbursts at times (see Appendix B, Stanway et al. 2014; Chisholm et al. 2019). It is, however, unlikely that they experience constant star formation as a very low star formation rate is inferred in most cases. As stripped stars prolong the time during which a stellar population can create ionizing photons, they could play an important role in shaping the distribution of ionization parameters found in the observed samples. Even though we likely have over-estimated the contribution from stripped stars to the ionization parameter because of the decreasing gas density, they remain interesting to consider since the observed galaxies show a large spread in measured ionization parameters. Assuming that the gas density has decreased by a factor of ten by the time stripped stars are formed, means that we have over-estimated the ionization parameter merely by about a factor of three, while the observed ionization parameters range over about two orders of magnitude.

Co-eval stellar populations, such as individual star clusters, have been observed to have ionization parameters that coincide with the distribution for the ionization parameters for the observed galaxies mentioned above. Xiao et al. (2019) show that ionization parameters for clusters are typically between −2.5 and −4. This suggests that the comparison with a co-eval stellar population for the observed galaxies is justified.

We do not expect that stripped stars contribute significantly to the ionizing emission from populations with very high measured values for the ionization parameter of log₁₀U  ≳   − 2 as observed by, for example, Erb et al. (2010), Berg et al. (2019), and Leitherer et al. (2018). When stripped stars dominate the ionizing emission, such high ionization parameters require that the galaxy is of high mass (≳10⁸ M_⊙) and that star formation has halted about 10 Myr ago (see Fig. 6). Galaxies with such high ionization parameters are observed to have ongoing star formation.

6.2. Impact on the UV luminosity and the UV continuum slope, β

The luminosity, L_ν, in the ultraviolet wavelengths has traditionally been used as a diagnostic for the star formation rate of stellar populations (Kennicutt 1998) as the wavelength range is dominated by the emission from young and massive stars. Despite stripped stars are very hot, they do not significantly impact the UV luminosity in stellar populations that form stars at a constant rate or in which star formation has halted less than 500 Myr ago. We display the results for L_ν at 1500 Å in Table 1.

The slope of the UV continuum, β, can be used to infer dust attenuation of stellar populations (e.g., Meurer et al. 1999). Similar to the UV luminosity, the slope of the UV continuum is not affected by the presence of stripped stars unless star formation has halted more than about 100 Myr ago. To quantify the effect, we take the common approach and define the UV continuum slope as the exponent in a power-law: L_λ ∝ λ^β, between 1250 Å and 2600 Å (Calzetti et al. 1994).

We find that stripped stars do not significantly affect the slope of the UV continuum if the slope is steep, β ≲ −0.5. For such cases, the UV is dominated by hot main-sequence stars. For shallower slopes of β ∼ −0.5, stripped stars change β by at least 0.1, and for even shallower slopes stripped stars can dominate the UV radiation. This resembles the observed phenomenon called the UV-upturn (Burstein et al. 1988), which has been considered to originate from subdwarfs that are formed late after star formation has ended, for example, through binary interaction in low-mass stars (Han et al. 2007). In our models, the UV slope becomes shallower with time for stellar populations in which stars are no longer forming. In the case of constant star formation, the UV is dominated by radiation from massive stars and the effect of stripped stars on the UV slope is negligible. The effect of metallicity is small on both the UV luminosity and the slope of the UV continuum. We show the results for lower metallicity in Tables A.2–A.4.

We conclude that both the UV luminosity and the slope of the UV continuum are un-affected by the presence of stripped stars in stellar populations in which star formation is ongoing or that are younger than about 100 Myr. Therefore, in such stellar populations, the method of inferring the star formation rate using the UV luminosity remains the same as well as inferring the dust attenuation using the UV continuum slope (cf. Reddy et al. 2018).

6.3. Impact on spectral features

We find that the stellar emission-line contribution from stripped stars is not distinguishable in the integrated spectrum of stellar populations. An emission feature that is observed to be very strong among stellar populations is He IIλ1640 (e.g., Kehrig et al. 2015) and the equivalent width of this line increases at most by 1 Å when stripped stars are included. This maximum increase occurs when the most massive stripped stars appear in co-eval and high-metallicity stellar populations. This also contrasts with the observed increase of the He IIλ1640 feature with decreasing metallicity (e.g., Shirazi & Brinchmann 2012; Senchyna et al. 2017, 2019). We note that higher wind mass-loss rates or slower stellar winds would increase the equivalent widths of the emission lines from stripped stars, as they are mainly formed by recombination in the stellar wind (see Sect. 2 and Paper II for a discussion). Eldridge & Stanway (2012) found that BPASS suggests that the He IIλ1640 feature could be significantly enhanced if rotation of massive stars leading to chemically homogeneous evolution is invoked. However, the observed He IIλ1640 emission is sometimes even stronger than what can be explained with rotating stars (Senchyna et al. 2017, 2019; Nanayakkara et al. 2019, see also Senchyna & Stark 2019).

The impact from stripped stars on the nebular spectrum is likely more interesting because their ionizing radiation affects the ionization state of the gas and thus also the spectral features emitted by the nebula. Nebular features are important diagnostics for the nature of observed galaxies and can, for example, be used to determine whether the ionizing source is stellar or quasar (Feltre et al. 2016; Gutkin et al. 2016, see also Stasińska et al. 2015). Here, we discuss likely effects that stripped stars have on the nebula based on simple considerations of the hardness of the ionizing spectrum. We assume that the gas density is sufficiently high for the ionizing emission from stripped stars to produce nebular emission lines. However, we note that the gas density likely has decreased when stripped stars are formed and the nebular features may therefore be relatively faint. Detailed modeling of the nebula is required to understand the strength of the nebular features caused by stripped stars. This is out of the scope of this study, but the topic of a forthcoming paper.

Figure 7 shows the shape of the ionizing part of the spectra of co-eval stellar populations compared to AGN (cf. Steidel et al. 2014; Stark et al. 2015; Feltre et al. 2016). We assume that the spectra of AGN are characterized by a power-law, L_ν ∝ ν^α, with a slope of −2 <  α <  −1.2 (e.g., Feltre et al. 2016). The figure shows that the spectra from single star populations always are softer than that of AGN as they are characterized by steep spectral slopes of α ≲ −2.5 (see also D’Aloisio et al. 2019). When stripped stars are present, the ionizing spectrum for photon energies lower than about 50 eV is harder. The slope is even close to flat for stellar populations younger than 50 Myr. The spectrum becomes softer than what young and massive stars can produce first after more than 100 Myr after a starburst.

Fig. 7. Shape of the ionizing part of the spectra of co-eval stellar populations compared to that of AGN, shown by normalizing the spectra at the H I-ionizing threshold at the photon energy 13.6 eV. We use gray for the spectra of stellar populations that only contain single stars and green for spectra of stellar populations in which stripped stars are included. Expected spectral slopes for AGN are indicated in transparent gray. We mark the ionization thresholds of different ionization stages for a variety of elements above the figure.

The ionizing part of the integrated spectrum is sufficiently hard when stripped stars are included for nebular oxygen to be ionized to O²⁺ and carbon to be ionized to C³⁺, while the ionizing radiation from massive stars favors lower ionization stages such as O⁺ and C²⁺. The high ionization of carbon has been detected at high-redshift (e.g., the C IVλ1548 Å feature, Stark et al. 2015, see also Berg et al. 2018 and Senchyna et al. 2019), making stripped stars interesting to consider as sources of such energetic photons. We also expect that stripped stars can give rise to high ratios of the nebular emission lines of O III to O II, which is a ratio used as an ionization measure and commonly labeled O32 (see e.g., Shapley et al. 2015; Strom et al. 2017, and Tang et al. 2019, for observed distributions of O32). Full ionization of nebular helium is favored for the hard ionizing radiation of AGN compared to that of stellar populations. However, when stripped stars are included more He II-ionizing photons are produced than when only single stars are considered (see Sect. 5), which potentially could give rise to nebular He II features, however likely not as strong as observed (cf. Kehrig et al. 2018). Promising sources responsible for the strong nebular emission of He II-ionizing photons are, for example, X-ray binaries (Schaerer et al. 2019).

We expect that stripped stars affect the location of stellar populations in the BPT diagram (Baldwin et al. 1981), which is a commonly used diagnostic diagram used to distinguish between star-forming galaxies and galaxies that host an AGN (Kewley et al. 2001; Kauffmann et al. 2003). For example, given the likely impact of stripped stars on the strength of O III lines, we expect that the ratio of [O III]/Hβ (the vertical axis of the BPT diagram) is affected. This is consistent with a recent study by Xiao et al. (2018, see also Stanway et al. 2014) who used BPASS models and found that populations that contain stripped stars may be located in the same region as small dwarf galaxies with strong ionizing emission are found (e.g., the Green Pea galaxies, Cardamone et al. 2009, see also the low-mass galaxies of Shapley et al. 2015). This questions the age estimate and the stellar content of such stellar populations and similar to H II regions.

For lower metallicity, the spectra of massive stars become harder and the spectra of stripped stars become softer. However, a distinct difference in the hardness between the spectra of massive stars and those that stripped stars introduce is visible for metallicities larger than Z = 0.0002. This is visible in Fig. A.4, which shows figures analogous to Fig. 7, but for lower metallicities. In the case of Z = 0.0002, stripped stars produce as hard ionizing emission as the very youngest massive stars, which is interesting as the impact from stripped stars is present for up to 100 Myr, while the massive stars die after about 5 Myr.

In the case of constant star formation, massive stars dominate and the result is a softer ionizing spectrum. Our models suggest that stripped stars primarily contribute with photons more energetic than 50 eV, but the increase is small and possibly difficult to distinguish.

Stripped stars produce photons capable of ionizing hydrogen and give rise to nebular Hα emission (Xiao et al. 2018). The contribution from stripped stars to nebular Hα is expected to be smaller than from massive stars for stellar populations in which stars form at a constant rate, since the line luminosity is related to the emission rate of H I-ionizing photons (see e.g., Leitherer & Heckman 1995; Schaerer 2003). However, in stellar populations in which star formation has halted, stripped stars can be responsible the formation of nebular Hα.

Nebular Hα is often used to infer the star formation rate of stellar populations as most H I-ionizing photons are produced by young and massive stars (Kennicutt 1998). Our models predict that this method is reliable unless star formation has stopped more than about 10 Myr ago.

7. Summary and conclusions

We analyzed the radiative contribution from stars stripped in binaries to the spectral energy distribution and ionizing emission of stellar populations. To do this, we simulated a stellar population including a realistic fraction of binaries to estimate at what time stripped stars are expected to be present. For the spectra, we used our previously computed grid of custom-made atmosphere models (Paper II).

Our main focus is the ionizing emission since stripped stars emit the majority of their radiation in the ionizing wavelengths. They thus constitute a stellar source of ionizing photons that is neglected in models that only account for single stars. We compare the emission rates of hydrogen and helium ionizing photons from stripped stars to the emission rates from massive single main-sequence stars and WR stars as predicted by STARBURST99. We also compare our results with predictions from the binary spectral synthesis code BPASS.

We quantify the effects of stripped stars on observable properties for co-eval stellar populations and continuous star formation. These include the production efficiencies of hydrogen and helium ionizing photons (ξ_ion, 0, ξ_ion, 1, and ξ_ion, 2), the ionization parameter (U), the UV luminosity, and the UV-continuum slope (β). We find that the presence of stripped stars increases the hardness of the ionizing part of the integrated spectrum. The harder spectrum may be distinguishable when both the production efficiency of hydrogen and helium ionizing photons are inferred. It is also likely that the hard ionizing spectrum gives rise to a different ionization structure of the surrounding nebula and thus a characteristic combination of nebular emission lines.

Our main findings are the following:

1. Stripped stars are important as an additional stellar source of ionizing radiation. For co-eval stellar populations, stripped stars dominate the ionizing output by several orders of magnitude once they are created, which in our simulations is at an age of ∼10 Myr. In the case of continuous star formation their emission rates of H I-, He I-, and He II-ionizing photons reach levels as high as 5, 15, and 500% compared to what is expected from massive stars alone.

2. The ionizing emission from stripped stars is primarily important in stellar populations in which star formation has recently halted. The reason is that stripped stars are created with a delay after the start of star formation. The emission of ionizing photons from stripped stars is significant after about 10 Myr and up to 100 Myr after star formation stopped.

3. Our models indicate that stripped stars only impact the integrated spectrum of a stellar population in the ionizing wavelengths.

4. Stripped stars introduce a characteristic hardness in the ionizing part of the spectral energy distribution. The effect on observable properties is not only an increase in the production efficiency of H I-ionizing photons (ξ_ion, 0), but a relatively larger increase in the production efficiency of He I- and He II-ionizing photons (ξ_ion, 1 and ξ_ion, 2). Current measurements of ξ_ion, 0 agree well with our predictions. We argue that future measurements of ξ_ion, 1 and ξ_ion, 2 will provide stringent tests for the theoretical models.

5. The presence of stripped stars in co-eval populations also affects the commonly used ionization parameter (U), causing it to remain between −3.5 ≲ log₁₀U ≲ −2.5 for a few 100 Myr. This is a range that often is observed for stellar populations. In the case of continuous star formation, stripped stars only affect the ionization parameter by up to 2%.

6. The ionizing radiation introduced by stripped stars is sufficiently hard to ionize the nebula to high ionization states, which potentially are visible via nebular emission lines if gas is still present at these times. Stripped stars likely have the largest relative impact on the nebular spectrum for stellar populations in which star formation recently halted. In these cases, it is likely that the nebula reaches high ionization states, such as O²⁺, C³⁺, and possibly also fully ionized helium. Because of the hard, integrated spectrum, we expect high ratios of nebular emission of O III to O II (O32). This means that the location of the stellar population in the BPT diagram is affected as well (cf. Xiao et al. 2018).

7. Stripped stars do not significantly affect the UV luminosity or the slope of the UV continuum (β) in populations in which stars form at a constant rate or in populations that are younger than ∼100 Myr. This suggests that in such populations, the UV luminosity can still be used as a reliable diagnostic for the star formation rate and the UV continuum slope can be used to make inferences about dust attenuation. In co-eval populations with ages between about 200 Myr and 1 Gyr, stripped stars significantly affect or even dominate the output of UV radiation, which significantly affects the UV continuum slope. However, at this point the slope is already close to flat, β ≳ −0.5. We note that we only consider the effect from stripped stars, while in reality other products of binary interaction may affect the ultraviolet part of the integrated spectral energy distribution.

8. For co-eval stellar populations, metallicity only modestly affects the emission of ionizing photons from stripped stars. However, in the case of continuous star formation, metallicity become more important. With lower metallicity, the stripped stars are cooler and their effect on the spectral hardness is therefore smaller. However, other products of binary interaction, such as rapidly spinning accretor stars, may still have large effects (cf. Eldridge et al. 2017).

Our results are mostly consistent with the results from the more complex binary population and spectral synthesis code BPASS for the parts that we expect are due to the presence of stripped stars (differences discussed in Appendix D). Our simulations have enabled an insightful analysis for the impact of stars stripped in binaries in particular. With our models, we show which differences that likely are due to stripped stars. We have discussed in detail (Sect. 2) and shown (Appendix C) that a simple model is sufficient for representing the emission from stripped stars. It has also allowed us to focus on refining detailed evolutionary models and spectral models especially made for stripped stars, to make sure that these objects are modeled according to our expectations.

Our models are publicly available on the CDS database, where we provide electronic tables with the contribution from stripped stars to the spectral energy distribution and to the emission rates of H I-, He I-, and He II-ionizing photons. Our models can also be obtained via the STARBURST99 online interface.

Acknowledgments

We thank the anonymous referee for a constructive report that helped improve the manuscript. The authors also acknowledge various people for helpful and inspiring discussion at various stages during the preparation of this manuscript, including Evan Bauer, Jared Brooks, Maria Drout, JJ Eldridge, Chris Evans, Rob Farmer, Miriam Garcia, Stephan Geier, Zhanwen Han, Edward van den Heuvel, Stephen Justam, Lex Kaper, Alex de Koter, Søren Larsen, Danny Lennon, Pablo Marchant, Colin Norman, Philipp Podsiadlowski, Onno Pols, Mathieu Renzo, Hugues Sana, Tomer Shenar, Nathan Smith, Elizabeth Stanway, Silvia Toonen, Jorick Vink, Manos Zapartas and the VFTS collaboration. This work was carried out on the Dutch national e-infrastructure with the support of SURF Cooperative. The authors acknowledges John Hillier for making his code, CMFGEN, publicly available. YG thank Martin Heemskerk for providing computing expertise and support throughout the project and Alessandro Patruno for allowing us to use the Taurus computer. SdM has received funding under the European Union’s Horizon 2020 research and innovation program from the European Research Council (ERC) (Grant agreement No. 715063). JHG acknowledges support from the Irish Research Council New Foundations Award 206086.14414 “Physics of Supernovae and Stars”. This work made use of v2.2.1 of the Binary Population and Spectral Synthesis (BPASS) models as last described in Eldridge et al. (2017) and Stanway & Eldridge (2018).

References

Appendix A: Effect of metallicity

A.1. Number of stripped stars

Envelope-stripping through binary interaction is not strongly dependent on metallicity (Paper I), meaning that stripped stars can be formed both at low and at high metallicity. This is not the case for wind-stripped WR stars and envelope-stripping has

therefore been considered a formation channel for WR stars in low-metallicity environments (e.g., Shenar et al. 2016).

Figure A.1 shows the number of stripped stars present in stellar populations at metallicities of Z = 0.006, 0.002, and 0.0002. The independence on metallicity for the formation of stars stripped in binaries is visible since both in the co-eval case and for continuous star formation, metallicity does not significantly change the number of stripped stars that are formed. In the case of solar metallicity (see Fig. 1) between around 150 and 400 stripped stars are present in the co-eval stellar population at the times when stripped stars are produced. In the case for very low metallicity (Z = 0.0002), the numbers are between 150 and 500. With lower metallicity, the contribution from Case B type mass transfer somewhat increases while the contribution from Case A type mass transfer somewhat decreases. The reason is that the expansion during the Hertzsprung gap is relatively larger for lower mass stars at lower metallicity. In agreement with the co-eval stellar population, the number of stripped stars in a population is very similar independent on metallicity for the case of continuous star formation.

Fig. A.1. Number of stripped stars present in a co-eval stellar population as function of time (left column) and for continuous star formation (right column). These panels are analogous to those in Fig. 1 but for populations with metallicity Z = 0.006 (top row), 0.002 (middle row), and 0.0002 (bottom row).

The number of low-luminosity WR stars computed in BPASS is also shown for the co-eval stellar populations in Fig. A.1. The number of these objects is larger for metallicities lower than solar compared to at solar metallicity. The reason could be that other evolutionary processes than envelope-stripping contribute to the formation of what is referred to as low-luminosity WR stars in BPASS.

A.2. Integrated spectrum and ionizing emission

Main-sequence stars are more luminous and hotter at low metallicity than at high metallicity. Their emission rates of ionizing photons therefore increase with decreasing metallicity. The emission rates of ionizing photons from stripped stars are also affected by their luminosity and effective temperature. The luminosity of stripped stars increases with lower metallicity, but the effective temperature decreases (Paper I; Paper II). The reason is primarily that envelope-stripping is less efficient at lower metallicity, which causes a larger amount of hydrogen to be left after mass transfer. The size of the stripped star is then dependent on how much leftover hydrogen there is, which leads them to be cooler at lower metallicity. Here, we combine the models for stripped stars at metallicities Z = 0.006, 0.002, and 0.0002, with the models from STARBURST99 for Z = 0.008, 0.002, and 0.001, respectively. We compare with the BPASS models with Z = 0.006, 0.002, and 0.0001. We summarize the metallicity combinations that we use for the different models in Table A.1.

We show the evolution of the contribution from stripped stars to the integrated spectrum of a co-eval stellar population in Fig. A.2. The figure shows that the emission from stripped stars becomes softer when metallicity decreases, while the emission from main-sequence stars becomes harder. As shown in Fig. 4, the total effect of metallicity on the emission from stripped stars is small for Q_0, pop and Q_1, pop, but large for Q_2, pop. The reason why the effect is large for Q_2, pop is that He II-ionizing photons are created in the steep Wien-part of the spectrum and the emission rate of He II-ionizing photons is therefore very sensitive to temperature variations. Since the stripped stars are cooler at lower metallicity, the emission rate of He II-ionizing photons drops.

Fig. A.2. Spectral energy distribution of a co-eval stellar population with initially 106 M⊙ stars, here shown for metallicities Z = 0.006, 0.002, and 0.0002 (horizontally) and for increasing time after starburst (vertically). The figure is analogous to Fig. 2. See Table A.1 for the metallicities of the STARBURST99 models.

The largest differences between our predictions and those from BPASS occur at low metallicity. Around 10 Myr after starburst and for Z ≤ 0.004, BPASS accounts for chemically homogeneous evolution for the accreting stars that were spun up during mass transfer (Eldridge et al. 2017). The result is that BPASS predicts higher emission rates of H I- and He I-ionizing photons than our models that only account for the ionizing emission from stripped stars at low metallicity.

For continuous star formation, we show the contribution from stripped stars to the spectral energy distribution in Fig. A.3 for the cases of lower metallicity. Comparing with Fig. 3, we find that the contribution from stripped stars is similar for Z ≳ 0.002, while the softer spectra of stripped stars are clearly visible at metallicity Z = 0.0002.

Fig. A.3. Spectral energy distribution in the case of constant star formation for metallicities Z = 0.006, 0.002, and 0.0002 (for the metallicities of the STARBURST99 models, see Table A.1). The models are for 1 M⊙ yr−1 and are taken 500 Myr after star formation started. The figures are analogous to Fig. 3.

The hardness of the ionizing part of the integrated spectrum affects the nebular ionization, as discussed in Sect. 6.3. We show the ionizing part of the spectra of co-eval stellar populations at low metallicity in Fig. A.4. The spectra of main-sequence stars are seen to become harder and those of stripped stars to become softer with lower metallicity. At Z = 0.0002, the hardness is similar for a population containing stripped stars and at an age of 20 Myr as for a population of only 2 Myr that contains massive main-sequence stars. However, we note that the duration for which massive stars give rise to such hard ionizing spectra is significantly shorter than the duration stripped stars emit ionizing radiation.

Fig. A.4. Ionizing part of the spectral energy distribution for co-eval stellar populations with metallicities Z = 0.006, 0.002, and 0.0002 from top to bottom. The spectra are normalized at the ionization threshold for hydrogen, 13.6 eV. The models from STARBURST99 have here Z = 0.008, 0.002, and 0.001 (see Table A.1). See also Fig. 7.

Fig. A.5. Ionizing part of the spectral energy distribution in the case of continuous star formation, shown after 500 Myr. We compare a population containing only single stars (gray) with one containing also stripped stars (green). These models have solar metallicity. See also Fig. 7.

As a complement to Fig. 7, we show the hardness of the ionizing part of the spectrum for a stellar population in which stars have formed at a constant rate for 500 Myr and with solar metallicity. The spectrum is only mildly affected by the presence of stripped stars. The largest differences from a population containing only single stars appear at high photon energies (≳40 eV). This could lead to stronger nebular emission of O III, C IV, and He II than what is expected from single star models.

We present our predictions for properties of stellar populations with metallicities of Z = 0.006, 0.002, and 0.0002 in Tables A.2–A.4. The general trends of metallicity are discussed in the sections that are mentioned in the tables.

Appendix B: Effect of star formation history

Stellar populations are in reality not co-eval or form stars at exactly a constant rate. More realistic stellar populations are either extended starbursts, composed of multiple stellar clusters that make up a bursty star formation rate or quenched of star formation. For simplicity, we considered the simplified cases of a starburst and a stellar population in which stars form at a constant rate in the main part of this study (see Sect. 5). If we want to better understand the effect of stripped stars in stellar populations, more realistic star formation histories are needed to compare observations with.

Here, we consider three additional star formation histories. First, an extended starburst, shaped with a Gaussian profile of width 40 Myr and in total 10⁶ M_⊙ of stars are created. This could be representative to a super star cluster similar to what are present in high-redshift galaxies and likely similar to 30 Dor in the Large Magellanic Cloud (e.g., Schneider et al. 2018). Second, a stellar population in which stars form at a constant rate for 200 Myr and then star formation is quenched and completely stops. This could be representative to aging stellar populations in which star formation has stopped (e.g., the retired galaxies of Stasińska et al. 2015). Third, several, randomly in time distributed, bursts of star formation with fluctuations in order of every 40 Myr and averaging a star formation rate of 1 M_⊙ yr⁻¹. This could be representative to a clumpy star-forming galaxy in which the emission comes from one or a few bursts at a time (see e.g., the models of Kimm & Cen 2014). This type of star formation is observed to be common for small galaxies that likely were in majority among galaxies at high redshift and could be the cause of cosmic reionization. The described star formation histories are visualized in the top row of panels in Fig. B.1.

Fig. B.1. Effect of different star formation histories on the emission rate of H I-ionizing photons and ionization parameter from a stellar population. We show the contribution from massive stars in green using STARBURST99 and the contribution from stripped stars in blue. From left to right, we show the effect from a Gaussian shaped burst of about 10 Myr, quenched star formation after 200 Myr, and a bursty star-forming environment. Top panels: adopted star formation history, middle panels: emission rate of H I-ionizing photons, bottom panels: ionization parameter calculated assuming a gas density of 102 cm−2.

We calculate the emission rates of ionizing photons for the described stellar populations by convolving the emission rate of ionizing photons from the co-eval stellar population (described in Sect. 5) over time and weighing with the specific star formation history. Since we consider stellar populations that have 10⁶ M_⊙ in stars or more, we do not expect the stochastic effects to be large. We also calculate the ionization parameter for the stellar populations following (3) and assuming a gas density of n_H = 100 cm⁻³. With time, the gas density likely decreases, causing slightly smaller ionization parameters than what we calculate. However, a decrease of a factor of ten in density leads to a smaller decrease in the ionization parameter of about a factor of three. We therefore consider this difference to be small. Here, we treat the models from STARBURST99 and the contribution from stripped stars separately.

We show the emission rate of H I-ionizing photons in the middle row of panels in Fig. B.1. The figures show that the massive stars dominate the ionizing output when star formation is actively ongoing, while stripped stars provide the majority of ionizing emission shortly after star formation has stopped. This is also visible in the ionization parameter, which is shown in the bottom panels.

For the burst, the emission rate of H I-ionizing photons is very similar to what is expected for the co-eval stellar population as seen when comparing the bottom-left panel of Fig. B.1 with panel a of Fig. 4. Since we consider the same mass produced in stars, the emission rates not only show the same shapes, but also the same numbers as described in Sect. 5. The result is that stripped stars appear to extend the duration of the starburst by providing ionizing emission for a longer time. However, we note that the massive stars in the stellar populations dominate the ionizing emission during the starburst itself.

For the stellar population in which star formation is quenched, it appears like a combination of the constant star formation population and the co-eval stellar population considered in the main manuscript. First, massive stars completely dominate the ionizing emission during the ongoing star formation. Once the star formation stops, the emission rate of ionizing photons from massive stars rapidly drops, leaving only stripped stars to emit ionizing photons at a significant rate. Their emission rate is lower, but lasts for longer, as seen in the co-eval stellar population.

In the case of the bursty star formation, the effect is similar; when star formation is ongoing the massive stars dominate the emission of ionizing photons, but once star formation stops their emission rapidly drops and the ionizing emission comes from stripped stars. For the bursty environment that we considered, the stripped stars prevent ionizing emission to ever completely die out, which contrasts with what is expected for massive stars. Stripped stars keep the emission rate of H I-ionizing photons to never drop below ∼10⁵¹ s⁻¹ for this particular case. However, we expect this to scale with the frequency of starbursts and the mean star formation rate.

The overall effect of the considered star formation history is similar for the ionization parameter compared to the emission rate of ionizing photons, since the two quantities are closely related (see Eq. (3)). After star formation has stopped, the stripped stars keep the ionization parameter at relatively high values for several hundred Myr. This can be seen for example in the case of the halted star formation, for which stripped stars cause the ionization parameter to remain larger than log₁₀U  ∼   − 3 for another ∼300 Myr after star formation has stopped. For the bursty environment that we considered here, stripped stars never allow the ionization parameter to decrease below log₁₀U  ∼   − 3, which is in contrast to what is predicted for single stars.

In summary, stripped stars appear to play an important role as contributors to the emission of ionizing photons in all the three considered star formation histories. Their effect is to extend emission of ionizing photons for a longer period in time. Their emission rate is smaller compared to what is expected for the massive stars formed at the same time.

Since stripped stars appear to affect the emission rate of ionizing photons significantly for a variety of different stellar populations, we conclude that they cannot be neglected when analyzing data from observed stellar populations. This is also implied from their effect on the ionization parameter.

Appendix C: Evolution during the stripped phase

Our approach in this study is simple, perhaps more simple than what is expected to reach sufficient accuracy in the results. A concern may be that one spectral model is not sufficient to represent the radiative emission from a stripped star during its entire stripped phase, which is a part of the method we employ and describe in Sect. 2.2. Here, we present additional spectral models computed at different times during the stripped phase to show that one spectral model per stripped star is sufficient to represent their emission and estimate their emission rate of ionizing photons. The models that we show here have solar metallicity.

C.1. Evolution in the Hertzsprung-Russell diagram during the stripped phase

Figure C.1 shows the evolution of five models in the Hertzsprung-Russell diagram computed using MESA (models presented in Paper II). The models are for donor stars in binary systems that initially have masses of 2.7, 5.45, 8.15, 12.17, and 18.17 M_⊙. The companion masses are decided such that the initial mass ratio is q = M_2, init/M_1, init = 0.8 and the initial period is set such that interaction is initiated early during the Hertzsprung gap evolution of the most massive star in the system. Figure C.1 shows that the stripped stars contract after Roche-lobe overflow is completed and become hotter than they were at the zero-age main sequence. The color shows the mass fraction of helium in the center of the stars. As the color changes, central helium burning progresses. The figure thus shows that the long-lasting central helium burning occurs in a very small part of the HR diagram, meaning that the luminosity and temperature of the stripped stars remain almost the same throughout the stripped phase.

Fig. C.1. Hertzsprung-Russell diagram showing the evolution of five stars that become stripped through mass transfer with a companion star. We show the evolution for models of donor stars with initial masses of 2.7, 5.45, 8.15, 12.17, and 18.17 M⊙, labeled as such. The tracks are gray for all but the central helium burning phase (defined as when the central helium mass fraction is between 0.9 and 0.1), which is colored with the central mass fraction of helium. The colored regions cover only a small part of the HR diagram and therefore indicates that the stripped stars have similar stellar parameters throughout the stripped phase. The models have solar metallicity. (We show here the surface temperature T⋆, which is referred to as the effective temperature in MESA).

The evolution after central helium burning is short in relation to the previous helium core burning. However, during this phase the stars expand and sometimes can reach helium giant phases. Even though the parameters change significantly both during this expansion phase and the earlier contraction after mass transfer, these two phases have short durations and only a small fraction of stripped stars would, therefore, reside in those phases.

Table C.1 displays the stellar parameters of the three most massive stripped stars shown in Fig. C.1. For nine times during the helium core burning of each star, we show the following parameters computed in MESA: mass (M_strip), surface radius (R_⋆), surface temperature (T_⋆), surface gravity (log₁₀ g_⋆), luminosity (L), wind mass loss rate (Ṁ_wind), and surface mass fraction of hydrogen and helium (X_H, s and X_He, s). These nine times are chosen when the central helium mass fraction is X_He, c = 0.9, 0.8, 0.7, 0.6, 0.5, 0.4, 0.3, 0.2, and 0.1, which spans almost the entire stripped phases. Most of the parameters are not varying significantly throughout the stripped phases. The stripped star with a progenitor of initially 8.15 M_⊙ shows least variation in the stellar parameters. The higher mass stripped stars loose the outermost hydrogen layers because of wind mass loss and appear somewhat hotter, smaller, and more compact.

From the small changes visible in the stellar parameters of the models in Table C.1, it is plausible that the changes in the emission rates of ionizing photons also are small. However, to estimate the change of the emission rates of ionizing photons during the stripped phases, modeled spectra give a more accurate result than the black-body spectra that can be created from the parameters computed in MESA.

C.2. Atmosphere models for different times during the stripped phase

To accurately estimate the emission rates of ionizing photons during the stripped phases, we create eight additional atmosphere models for each of the evolutionary models of the stripped stars with progenitors of initially 8.15, 12.17, and 18.17 M_⊙. These atmosphere models are created such that there is each a model when the central helium mass fraction is X_He, c = 0.9, 0.8, 0.7, 0.6, 0.5, 0.4, 0.3, 0.2, and 0.1, which corresponds to a close to constant temporal separation between the spectral models (see Table C.1). The technique for how to create the atmosphere models was described in detail in Paper I and Paper II.

C.2.1. Spectral energy distribution

In Fig. C.2, we show the spectral energy distribution for each of the stripped stars computed at the different times during the stripped phases. The figure shows that the shapes are very similar for most of the wavelength range throughout the stripped phases. The observable UV, optical and infrared spectra are close to identical for all but perhaps the model at X_He, c = 0.9, when the stars are still contracting and thus are slightly brighter. The H I- and He I-ionizing parts of the spectral energy distributions are also very similar, while the two most massive stripped stars show that the He II-ionizing part of the spectral energy distributions change with time. This is an effect of the outermost hydrogen layers being removed through wind mass loss, causing the star to be hotter. Since the wind mass loss rate from stripped stars is uncertain, the estimated He II-ionizing emission is also uncertain.

Fig. C.2. Spectral energy distributions modeled at nine different times during the evolution of stripped stars. The displayed models are for stars with initial masses of 8.15 (top), 12.17 (middle), and 18.17 M⊙ (bottom). The spectral models are created during the helium core burning at the time when the stripped stars have central helium mass fractions of XHe, c = 0.9, 0.8, ...0.1, which we highlight using purple, pink, and yellow colors. We use blue colors to mark the parts of the spectra that are H I-, He I-, and He II-ionizing. The spectral models are very similar during the stripped phases for the three stars, with the exception of the uncertain He II-ionizing part of the spectrum.

Figure C.2 shows that all but the uncertain He II-ionizing part of the spectral energy distributions seem mostly unaffected by the evolution during the stripped phase. This means that including multiple spectral models during the stripped phases of each star in the evolutionary grid would not significantly improve the estimate of how stripped stars affect the total spectral energy distribution of a stellar population (e.g., Fig. 2).

C.2.2. Emission rate of ionizing photons

We calculate the emission rates of H I-, He I-, and He II-ionizing photons at different times during the stripped phases using the additional spectral models. We present the resulting emission rates as a function of time in Fig. C.3 and give the values in Table C.2. The variation in the emission rate of H I-ionizing photons is a factor of ∼1.5 for all the three considered stripped stars. The same is true for the emission rate of He I-ionizing photons. The variation in the uncertain emission rate of He II-ionizing photons is from a factor of a few up to several ten times.

Fig. C.3. Emission rate of ionizing photons as function of time during the stripped phase of three models with initial masses of Minit = 18.17 (left column), 12.17 (middle column), and 8.15 M⊙ (right column). We show the emission rates of H I-, He I-, and He II-ionizing photons in the top, middle and bottom rows respectively. The constant rates that we assume when creating the population synthesis are shown as colored boxes (see Sect. 2.2), while the results from atmosphere models that are modeled at different times during the stripped phases are shown as colored dots connected with a line. The prediction from assuming blackbody radiation from the stars is shown in gray. We hatch the time when the star is not stripped.

We calculate the total number of emitted ionizing photons by the three stripped stars displayed in Fig. C.3 and in Table C.2 and compare whether multiple spectral models give rise to a significant improvement in the estimate for the total number of emitted ionizing photons. To do this, we first assume that the stripped phase lasts from when the central hydrogen mass fraction has reached X_H, c <  0.001 until the end of our computations. The stripped phase starts slightly later, but the difference in time is negligible. We then use either the emission rates as a function of time by using the multiple spectral models for different time bins or just the emission rates for the spectral model computed at X_He, c = 0.5.

The difference between the estimate from just one spectral model and from using a sequence of spectral models is remarkably small. Using just one single spectral model computed at X_He, c = 0.5 underestimates the total production of H I- and He I-ionizing photons with less than 20% for either of the stars. The emission rate of He II-ionizing photons is underestimated by at maximum 40% for the models we considered here.

In Fig. C.3 and Table C.2 we also show the evolution of the ionizing emission rates as predicted by blackbody radiation if we would only have known the stellar parameters from MESA. As we showed in Paper I, the difference between the predictions from blackbody and from CMFGEN is subtle for the emission rates of H I- and He I-ionizing photons, but very significant for the emission rates of He II-ionizing photons.

We conclude that modeling a sequence of spectra for each stripped star model will only marginally improve the estimated emission rate of ionizing photons. As Fig. C.3 and Table C.2 indicate, the stripped stars emit on average ionizing photons at a slightly higher rate than what the model computed half-way through central helium burning suggests. The difference is small and we consider therefore that the simple approach of using just one spectral model per stripped star is sufficient to estimate the emission rate of ionizing photons from stripped stars (as described in Sect. 2.2).

Appendix D: Comparison to models of BPASS

The code BPASS combines large grids of evolutionary models of interacting binaries and single stars with spectral libraries to simulate a stellar population (Eldridge et al. 2017). Stanway et al. (2016) show that BPASS predicts higher emission rate of ionizing photons and also a delayed emission, which is very similar to what we find in this study. In the main part of this manuscript (Sect. 5), we compare the emission rate of ionizing photons that we predict to come from stripped stars with the emission rate of ionizing photons that BPASS predicts. Most of our predictions are consistent with the predictions from BPASS, but there are parts where this is not the case. To address whether the reason is that stripped stars are treated in a way that is very different from ours or whether it actually comes from other types of stars, we go into detail and compare individual evolutionary and spectral models for stripped stars.

Here, we first compare the evolutionary models of stripped stars included in BPASS with our evolutionary models, using similar parameters for the initial binary systems. We then compare the spectral energy distribution that BPASS would choose for the individual stripped stars with the CMFGEN models that we have computed especially for stripped stars. We use the models from the BPASS version 2.2.1, Tuatara (Eldridge et al. 2017; Stanway & Eldridge 2018). These models are publicly available on the BPASS website³.

We use four models from our model grid with solar metallicity. These have initial masses of 18.17, 12.17, 8.15, and 5.45 M_⊙. We chose these masses both because they cover a large mass range and thus a diversity of stripped stars and also because these stripped stars are formed roughly at the times we display the integrated spectrum in Fig. 2 (11, 20, 50, and 100 Myr after a starburst).

D.1. Evolutionary models

Depending on which assumptions are made for stellar evolution models, the stars evolve slightly or very different. We show the evolutionary tracks of stars that undergo Case B mass transfer and get stripped in the left column of panels in Fig. D.1. Our models are shown in color while the models from BPASS with the most similar initial binary parameters are shown in pale green. As the figure shows, the stripped stars computed in the BPASS models are fainter than our models for stripped stars. The reason for the difference is likely the lower overshooting that is assumed in BPASS (Eldridge et al. 2017) and the resulting stripped stars therefore have lower masses in BPASS compared to in our models. This means that for the same population initial conditions, BPASS should predict fewer stripped stars than what our model does since they would reach the lower mass limit for central helium burning at a higher progenitor mass than what our models do.

To see whether BPASS can produce stripped stars that are similar to our models, we also show the evolution of a more massive star computed in BPASS in each the HR diagrams in Fig. D.1. The resulting stripped star in these cases appears in the same region of the HR diagram as our models and likely has very similar properties as our structure models for stripped stars.

Fig. D.1. Comparison between our models and the ones used in BPASS for the evolution (left) and spectra (right) of stripped stars. We use five example models of stars that lose their envelopes by Roche-lobe overflow initiated early on the Hertzsprung gap and that have initial masses of Minit, 1 = 5.5, 8.2, 12.2, and 18.2 M⊙ and initial mass ratios of q = Minit, 2/Minit, 1 = 0.8. We show the models we created in orange, yellow, blue, and green, and BPASS models in gray and pale green.

We also note that it is visible in Fig. D.1 that mass transfer is treated differently in BPASS than in MESA. The shape of the evolutionary tracks during mass transfer is different, and the systems computed in BPASS tend to tighten in contrast to our models. The reason for this could be that angular momentum is removed from the orbit of the binary system when it cannot be accreted by the accretor star in BPASS, while we assume that mass is lost with the specific angular momentum of the accretor star itself.

D.2. Spectral models

The code BPASS uses WR star models from POWR to represent the radiative emission from stripped stars with temperatures T_⋆ ≥ 4.45 and surface hydrogen mass fraction X_H, s <  0.4. The coolest subdwarfs are represented by B-star models from the BASEL grids. The POWR models are made in particular for WR stars and therefore are set to have high wind mass loss rates. The models also have a fixed terminal wind speed of v_∞ = 1000 km s⁻¹.

To better understand whether the choice of atmosphere model matters for the emission rate of ionizing photons from stripped stars, we compare the spectral model that we computed for a stripped star with what spectral model would have been assumed by BPASS for the same stripped star parameters (we follow the description for how this is done in Eldridge et al. 2017). We show the spectral energy distributions in the right column of Fig. D.1, where our models are shown in color and the model BPASS would have assumed for the same stellar parameters is shown in dark gray. In pale green, we also show the spectral energy distribution that BPASS would have assumed for a stripped star with similar initial binary parameters (marked in the corresponding HR diagram).

The figure shows that the spectral energy distributions of the various models are similar for most of the wavelength range. The exception is for the He II-ionizing part because the models assumed in BPASS show a steep drop in emission at the ionization limit of He⁺. The differences between our models (shown in colors) and the what is expected for stripped stars with similar progenitors (pale green) is also not very large, primarily a difference in brightness can be seen.

Since the stripped stars in BPASS are of lower luminosity than in our model, intuitively they should also emit less ionizing radiation. We measure the emission rates of ionizing photons in the spectral models of our models and compare to the emission rates for a stripped star with the same progenitor mass as modeled in BPASS. We find that the emission rate of H I-ionizing photons is 30% lower in the model assumed in BPASS compared to our model in the case of the stripped star with the progenitor of mass ∼18 M_⊙. With lower mass, the difference increases and for the stripped star with a progenitor of ∼5.5 M_⊙ the stripped star modeled in BPASS emits H I-ionizing photons at half the rate of our model. For the emission rate of He I-ionizing photons, the effect is similar, just slightly larger, spanning from 40% fainter for higher mass stripped stars in BPASS, to 75% fainter for the 5.5 M_⊙ progenitor. The emission rate of He II-ionizing photons is at highest 6% of the emission from our 18 M_⊙ model in the BPASS counterpart. In the other cases we consider, the He II-ionizing emission is about seven orders of magnitudes lower compared to the emission from our models. This questions the origin of the He II-ionizing emission from BPASS at later times (see Fig. 4).

In summary, the emission rate of ionizing photons is expected to be lower from the stripped stars in BPASS than the stripped stars we modeled with the same progenitor stars. Likely because of assumptions that affect how efficiently stripped stars are formed in BPASS, the predicted total emission rates of ionizing photons is higher in BPASS than in our model. The emission rate of He II-ionizing photons that BPASS assumes for stripped stars is close to negligible compared to what our models predict. This suggests that the hard ionizing emission originates from different objects in BPASS.

All Tables

All Figures

Fig. 1. Number of stripped stars shown as function of time for a population of solar metallicity. Panel a: numbers for a co-eval stellar population with initially 106 M⊙ in stars. We highlight the three formation channels that we consider to lead to the formation of stripped stars: stable mass transfer initiated during the main sequence evolution of the donor star (Case A, RLOF in light green), stable mass transfer initiated when the donor star passes over the Hertzsprung gap (Case B, RLOF in medium green), and successful ejection of a common envelope that developed during interaction initiated during the donor’s Hertzsprung gap evolution (Case B, CEE). For comparison, we show the number of low-luminosity WR stars predicted by BPASS in light gray (Eldridge et al. 2017). We indicate the expected time of WR star formation through solely stellar winds and the time of subdwarf formation through low-mass stellar evolution above the panel. Panel b: case of continuous star formation at a rate of 1 M⊙ yr−1. We show the parts of the population of stripped stars that are created from progenitor stars with initial masses between 10 − 20 M⊙ in purple, 5 − 10 M⊙ in pink, 3 − 5 M⊙ in salmon, and 2 − 3 M⊙ in peach.

In the text

Fig. 2. Spectral energy distribution of a co-eval stellar population. The contribution from stripped stars is highlighted with a blue line. The parts of the spectra that are H I-, He I-, and He II-ionizing are shaded in blue, while the UV is shaded in purple. For comparison, we show the spectral energy distribution of a population containing only single stars using STARBURST99 (green line), which can be interpreted as the contribution from the remaining stars in the stellar population. We also show the predictions from BPASS (gray line), where the effects of binary interaction are included. The panels correspond to different times after the instantaneous starburst of 106 M⊙. The model shown here assumes solar metallicity (see Appendix A for lower metallicity models).

In the text

	Fig. 3. Spectral energy distribution of a stellar population in which star formation has taken place at a constant rate of 1 M⊙ yr−1 for 500 Myr. Otherwise similar to Fig. 2. We show models with solar metallicity.
In the text

Fig. 4. Emission rates of ionizing photons from stellar populations as a function of time. The contribution from stripped stars is shown in blue, while green represents the contribution from the massive main sequence and WR stars in the stellar population (models from STARBURST99, see Sect. 2.3). For reference, we show the predictions from BPASS models in which binary interactions are included using gray color. The solid lines correspond to the predictions from a population with solar metallicity, while the shaded regions of the same color represent the effect of lowering the metallicity, from solar down to Z = 0.0002 for the stripped stars (see Appendix A and Table A.1). Left column: emission rates from a co-eval stellar population with initially 106 M⊙ in stars, and right column: emission rates from a stellar population in which stars form at the constant rate of 1 M⊙ yr−1. Top, middle, and bottom rows: emission rates of H I-, He I-, and He II-ionizing photons, respectively (Q0, pop, Q1, pop, and Q2, pop).

In the text

Fig. 5. Diagrams showing the hardness of the ionizing part of the spectrum of stellar populations, constructed by the production efficiencies of H I- and He I- (panel a) or He II-ionizing photons (panel b) (ξion, 0, ξion, 1, ξion, 2). The colored regions represent the time spans indicated by the color bars for co-eval stellar populations and cover the predictions from all metallicities (from solar metallicity down to Z = 0.0002 for stripped stars, see Appendix A and Table A.1). Green shades show predictions for young, single star populations, while the purple and red shades represent the hardness of stellar populations in which stripped stars are included. The case of constant star formation is shown with solid lines and markers (pluses in black for single stars and circles in red for when stripped stars are included), taken after 500 Myr. Above the diagrams, we show the distribution of measured ξion, 0 for three samples of observed unresolved stellar populations from intermediate to high-redshift by Bouwens et al. (2016), Matthee et al. (2017), and Shivaei et al. (2018).

In the text

Fig. 6. Top: ionization parameter computed for co-eval stellar populations as a function of time and for constant star formation, taken after 500 Myr. We show the predictions for stellar populations containing only single stars in gray shades and the model when stripped stars are included in purple shades. For constant star formation, the contribution from stripped stars is about 2% and therefore we do not show the markers for single star populations as they overlap. We show measurements of the ionization parameter for groups of observed galaxies to the right of the diagram (Moustakas et al. 2010; Nakajima & Ouchi 2014). Bottom: table explaining which gas density and stellar mass for co-eval stellar populations that correspond to which contour in the diagram using numbers as labels. In the case of constant star formation, the numbers are correlated with combinations of the gas density and star formation rate instead.

In the text

Fig. 7. Shape of the ionizing part of the spectra of co-eval stellar populations compared to that of AGN, shown by normalizing the spectra at the H I-ionizing threshold at the photon energy 13.6 eV. We use gray for the spectra of stellar populations that only contain single stars and green for spectra of stellar populations in which stripped stars are included. Expected spectral slopes for AGN are indicated in transparent gray. We mark the ionization thresholds of different ionization stages for a variety of elements above the figure.

In the text

	Fig. A.1. Number of stripped stars present in a co-eval stellar population as function of time (left column) and for continuous star formation (right column). These panels are analogous to those in Fig. 1 but for populations with metallicity Z = 0.006 (top row), 0.002 (middle row), and 0.0002 (bottom row).
In the text

	Fig. A.2. Spectral energy distribution of a co-eval stellar population with initially 106 M⊙ stars, here shown for metallicities Z = 0.006, 0.002, and 0.0002 (horizontally) and for increasing time after starburst (vertically). The figure is analogous to Fig. 2. See Table A.1 for the metallicities of the STARBURST99 models.
In the text

	Fig. A.3. Spectral energy distribution in the case of constant star formation for metallicities Z = 0.006, 0.002, and 0.0002 (for the metallicities of the STARBURST99 models, see Table A.1). The models are for 1 M⊙ yr−1 and are taken 500 Myr after star formation started. The figures are analogous to Fig. 3.
In the text

	Fig. A.4. Ionizing part of the spectral energy distribution for co-eval stellar populations with metallicities Z = 0.006, 0.002, and 0.0002 from top to bottom. The spectra are normalized at the ionization threshold for hydrogen, 13.6 eV. The models from STARBURST99 have here Z = 0.008, 0.002, and 0.001 (see Table A.1). See also Fig. 7.
In the text

	Fig. A.5. Ionizing part of the spectral energy distribution in the case of continuous star formation, shown after 500 Myr. We compare a population containing only single stars (gray) with one containing also stripped stars (green). These models have solar metallicity. See also Fig. 7.
In the text

Fig. B.1. Effect of different star formation histories on the emission rate of H I-ionizing photons and ionization parameter from a stellar population. We show the contribution from massive stars in green using STARBURST99 and the contribution from stripped stars in blue. From left to right, we show the effect from a Gaussian shaped burst of about 10 Myr, quenched star formation after 200 Myr, and a bursty star-forming environment. Top panels: adopted star formation history, middle panels: emission rate of H I-ionizing photons, bottom panels: ionization parameter calculated assuming a gas density of 102 cm−2.

In the text

Fig. C.1. Hertzsprung-Russell diagram showing the evolution of five stars that become stripped through mass transfer with a companion star. We show the evolution for models of donor stars with initial masses of 2.7, 5.45, 8.15, 12.17, and 18.17 M⊙, labeled as such. The tracks are gray for all but the central helium burning phase (defined as when the central helium mass fraction is between 0.9 and 0.1), which is colored with the central mass fraction of helium. The colored regions cover only a small part of the HR diagram and therefore indicates that the stripped stars have similar stellar parameters throughout the stripped phase. The models have solar metallicity. (We show here the surface temperature T⋆, which is referred to as the effective temperature in MESA).

In the text

Fig. C.2. Spectral energy distributions modeled at nine different times during the evolution of stripped stars. The displayed models are for stars with initial masses of 8.15 (top), 12.17 (middle), and 18.17 M⊙ (bottom). The spectral models are created during the helium core burning at the time when the stripped stars have central helium mass fractions of XHe, c = 0.9, 0.8, ...0.1, which we highlight using purple, pink, and yellow colors. We use blue colors to mark the parts of the spectra that are H I-, He I-, and He II-ionizing. The spectral models are very similar during the stripped phases for the three stars, with the exception of the uncertain He II-ionizing part of the spectrum.

In the text

Fig. C.3. Emission rate of ionizing photons as function of time during the stripped phase of three models with initial masses of Minit = 18.17 (left column), 12.17 (middle column), and 8.15 M⊙ (right column). We show the emission rates of H I-, He I-, and He II-ionizing photons in the top, middle and bottom rows respectively. The constant rates that we assume when creating the population synthesis are shown as colored boxes (see Sect. 2.2), while the results from atmosphere models that are modeled at different times during the stripped phases are shown as colored dots connected with a line. The prediction from assuming blackbody radiation from the stars is shown in gray. We hatch the time when the star is not stripped.

In the text

Fig. D.1. Comparison between our models and the ones used in BPASS for the evolution (left) and spectra (right) of stripped stars. We use five example models of stars that lose their envelopes by Roche-lobe overflow initiated early on the Hertzsprung gap and that have initial masses of Minit, 1 = 5.5, 8.2, 12.2, and 18.2 M⊙ and initial mass ratios of q = Minit, 2/Minit, 1 = 0.8. We show the models we created in orange, yellow, blue, and green, and BPASS models in gray and pale green.

In the text