Open Access
Issue
A&A
Volume 682, February 2024
Article Number A69
Number of page(s) 19
Section Astrophysical processes
DOI https://doi.org/10.1051/0004-6361/202347880
Published online 02 February 2024

© The Authors 2024

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

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

1. Introduction

Ultra-luminous X-ray sources (ULXs) are extragalactic X-ray sources with luminosities exceeding 1039 erg s−1 and are off center from their host galaxies (Kaaret et al. 2017; Fabrika et al. 2021). Since their first detection by the Einstein Telescope (e.g., Fabbiano 1989), several hundred more ULX candidates have been observed (Swartz et al. 2004; Walton et al. 2011; Liu 2011; Earnshaw et al. 2019; Kovlakas et al. 2020; Soria et al. 2022; Bernadich et al. 2022; Salvaggio et al. 2023). It is the general consensus that these sources are bright X-ray binaries (XRBs) with a non-degenerate star transferring mass onto a neutron star (NS) or a black hole (BH; King et al. 2001; Zezas et al. 2002). The Eddington limit of a stellar-mass BH (with mass ∼10 M) is around 1039 erg s−1. Hence, ULXs are thought to be super-Eddington (assuming isotropic emission). The nature of the compact object (CO) involved in ULXs has been the subject of much speculation because of the apparent super-Eddington nature of ULXs. In order for the observed ULX luminosities to be explained by Eddington-limited accretion, one proposed explanation is the presence of accreting intermediate-mass BHs, that is, BHs with masses around 100–1000 M (Colbert & Mushotzky 1999; Miller 2006; Maccarone et al. 2007). However, it is now believed that the majority of ULXs might not belong to this class, as evidence for intermediate-mass BHs is scarce (e.g., Miller-Jones et al. 2012; Perera et al. 2017; Koliopanos et al. 2017; Tremou et al. 2018; Zocchi et al. 2019; Mann et al. 2019). Recently, it has been suggested that a ULX source (CXO J133815.6+043255) is an intermediate-mass BH based on study of its radio spectral energy distribution (Smith et al. 2023).

X-ray binaries are categorized into three types: low-mass X-ray binaries (LMXBs; donors with masses ≲2.0 M), intermediate-mass X-ray binaries (IMXBs; donors with masses in the range 2.0–8.0 M), and high-mass X-ray binaries (HMXBs; donors with masses ≳8.0 M). Generally, ULXs have been associated with active star-forming regions that are dominated by young and bright HMXBs (Roberts et al. 2002; Gao et al. 2003; Kaaret et al. 2004; Wolter & Trinchieri 2004; Zezas et al. 2007; Anastasopoulou et al. 2016; Wolter et al. 2018; Kovlakas et al. 2020). However, studies have found that intense mass-transfer phases in LMXBs and IMXBs (with accreting NSs) could produce ULX luminosities (King et al. 2001; Shao & Li 2015; Karino 2018; Wiktorowicz et al. 2019; Quast et al. 2019; Misra et al. 2020). Binaries with an expanding donor (after exhausting core hydrogen) with a radiative envelope can drive super-Eddington mass-transfer rates if the donor is slightly more massive than the accreting CO (within a certain limit so as to not initiate a dynamical instability). Also, some ULX observations have been linked to old stellar populations (that are dominated by LMXBs; Angelini et al. 2001; Colbert & Ptak 2002; Kim & Fabbiano 2004; Swartz et al. 2004; Fabbiano et al. 2006; Feng & Kaaret 2008).

The discovery of coherent X-ray pulsations (with a pulse period of 1.37 s) in M82 X-2 (Bachetti et al. 2014) confirmed that at least a fraction of ULXs have NS accretors, as X-ray pulses imply the presence of accreting NSs. Matter accreting onto the NS surface forms a hot spot that creates X-ray pulses with the spinning of the NS (Seward & Charles 1995). Since BHs do not have a well-defined surface, they cannot emit similar X-ray pulsations. Many more pulsating ULXs have been observed since M82 X-2 (Fürst et al. 2016; Israel et al. 2017a,b; Tsygankov et al. 2017; Carpano et al. 2018; Brightman et al. 2018; Doroshenko et al. 2018; Sathyaprakash et al. 2019; Zhang et al. 2019; Rodríguez Castillo et al. 2020; Quintin et al. 2021). The X-ray luminosities of pulsating ULXs are typically in the range of 1039–1041 erg s−1, exceeding the NS Eddington limit, which is around 1038 erg s−1, by a few orders of magnitude.

There are many explanations for the apparent breach of the Eddington limit. One suggested model explaining ULXs describes the formation of an accretion disk receiving matter at super-Eddington rates (Shakura & Sunyaev 1973; Lipunova 1999; Poutanen et al. 2007). Within the described accretion disk, at a given disk radius, strong outflows emerge perpendicular to the disk, taking away excess matter and angular momentum and thereby keeping the disk locally Eddington-limited while the total luminosity emitted from the disk exceeds the CO Eddington limit. The structure of the accretion disk changes from a thin disk to a geometrically thick disk. The outflowing wind carries angular momentum, which creates a hollow cone through which the outgoing X-ray emission escapes. This results in an additional increase in the perceived luminosity for an observer looking face-on at this source (King et al. 2001; King 2009). The degree of collimation of outgoing emission depends on the nature of the accreting object. Since the Eddington limit of NSs is an order of magnitude lower than that of BHs, their X-ray emission is expected to be more beamed to reach ULX-like luminosities. Additionally, various theoretical studies suggest that a high degree of collimation of emission from geometrically thick disks can explain ULX luminosities well enough (e.g., King & Lasota 2016, 2019, 2020; Middleton & King 2017; Lasota & King 2023). Alternatively, observations of certain sources suggest that the high luminosity from ULXs could be intrinsic to the accreting CO involved. For instance, studying the pulsed X-ray emission of the ULX source RX J0209.6-7427 (with an accreting pulsar) during its 2019 giant outburst, Hou et al. (2022) estimated that the dominant source of the pulsed X-ray emission is from the fan beam of the accretion columns and that the high luminosity of the source (∼ 1.11 × 1039 erg s−1) is intrinsic instead of beamed.

For accreting NSs, the presence of magnetic fields could add further complications. For instance, Chashkina et al. (2019) proposed an advection-dominated accretion disk when the matter is transferred at super-Eddington rates. The presence of strong magnetic fields would disrupt accretion disks within the magnetosphere, channeling all material onto the NS and modifying the critical mass-accretion rate. Magnetar-like fields (≳ 1014 G) have also been suggested to explain NS ULXs, as they reduce the electron scattering cross-section, leading to an increase in the Eddington limit (Basko & Sunyaev 1976; Mushtukov et al. 2015, 2017; Chashkina et al. 2017, 2019; Mushtukov et al. 2019). Various observational studies involving the first galactic pulsating ULX, Swift J0243.6+6124, indicate the presence of multi-pole magnetic field components (Doroshenko et al. 2020; Kong et al. 2022). The magnetic field strength estimated for the source is around 1.6 × 1013 G (based on study of the cyclotron resonance scattering features, or CRSF; Kong et al. 2022), which is an order of magnitude higher than that of typical accreting pulsars with a magnetic field assumed to have only dipole components. Hou et al. (2022) studied the pulsed emission from the ULX RX J0209.6-7427 and estimated the magnetic field strength to be about 4.8–8.6 × 1012 G or 1.7–2.2 × 1013 G, both estimates correspond to dipole and multi-pole magnetic fields of the NS, respectively. A similar result was reached by Liu et al. (2022) for this source and for SMC X-3 during its outburst. However, CRSF studies of some other ULXs have suggested pulsar-like magnetic fields in the range of 1011–1013 G, which is below the magnetic field strength required for magnetars (Walton et al. 2018; Brightman et al. 2018; Koliopanos et al. 2019). Recently, Lasota & King (2023) suggested that magnetar-like field strengths cannot explain ULX spin-up rates and inferred that beaming due to the accretion disk must be involved in ULXs.

While, observationally, non-pulsating ULXs dominate (so far) over pulsating ULXs in number, with only about ten of the more than 500 ULX candidates known to pulsate, the absence of observed X-ray pulsations does not negate the presence of an accreting NS. Hence, many more ULXs may contain NS accretors instead of BHs, and their X-ray pulses may simply not have been observed. This argument is supported by the fact that a majority of LMXBs (which are thought to contain accreting NSs) do not show measurable X-ray pulsations (Vaughan et al. 1994; Dib et al. 2005; Messenger & Patruno 2015; Patruno et al. 2018). There is also a source (namely, M51 ULX-8) confirmed to be a ULX by observed cyclotron resonance scattering features for which X-ray pulses have not been observed so far. The observed scattering features suggest the presence of an NS with a magnetic field of around 1012 G (Brightman et al. 2018; Middleton et al. 2019).

Israel et al. (2017a) discovered a pulsating ULX in NGC 5907, which is also one of the most luminous ULXs known at an X-ray luminosity of 1041 erg s−1, suggesting that other extremely bright ULXs might have accreting NSs. Pintore et al. (2017) studied the spectral properties of bright ULXs and concluded that most ULXs fit well with the X-ray spectra of accreting magnetic NSs in the Galaxy. In fact, many theoretical studies have suggested that NSs dominate over BHs as accretors in the total ULX population, and the low number of confirmed NS ULXs is an observational effect (e.g., Mushtukov et al. 2015; King & Lasota 2016, 2020; King et al. 2017; Middleton & King 2017). Walton et al. (2018) found spectral signatures in a large sample of bright ULXs (all so far non-pulsating) that are usually associated with pulsar-like emissions. Wiktorowicz et al. (2015) studied the evolutionary channels that lead to extreme mass-transfer phases in close XRBs using binary population synthesis. In their study, of the brightest ULXs that exceeded 1042 erg s−1 (these sources are also known as hyperluminous X-ray sources, or HLXs), half contained accreting NSs instead of BHs. While the authors focused on HLXs, the 50% contribution from NS XRBs to HLXs contrasts with the expectation that BH XRBs are always brighter than NS XRBs.

King & Lasota (2016, 2020) suggested that NS ULXs and BH ULXs are indistinguishable at some point during their individual evolution. Accretion would eventually align the NS spin and the accretion disk axes, obscuring the X-ray pulses. They also found that NS ULXs, when observed, would be brighter than BH ULXs, as they have stronger beaming due to their lower Eddington limit. This also suggests why HLXs are expected to have a significant contribution from NS ULXs. However, short durations of the intense mass-transfer phase and increased beaming would down-weight their observability. When studying ULX populations by taking into account their observability due to geometrically beamed emission, Wiktorowicz et al. (2019) also found that the majority of their simulated NS ULXs had beamed emission with the mass-transfer phase proceeding on a thermal timescale compared to BH ULXs (which were typically isotropically emitting and proceeded on a nuclear timescale). They found that NS ULXs dominated the total underlying population of simulated ULXs in populations with constant star formation and that the observed ULX population was populated by NS ULXs and BH ULXs equally. For starburst populations, they found that BH ULXs dominated both observed and total populations when the systems are young (≲ 1000 Myr) but that this dominance shifts to NS ULXs for older systems. As a consequence of geometrically beamed emission in NS ULXs, observations from flux-limited surveys would be dominated mostly by BH ULXs.

The relative importance of NSs versus BHs in the intrinsic ULX populations depends not only on their accretion physics but also on the properties of the host galaxy, including age and metallicity. Shao & Li (2015) used a combination of parametric population synthesis calculations (modified BPS code; Hurley et al. 2002; Shao & Li 2014) with detailed binary evolution models (TWIN version of the Eggleton code; Eggleton 1971, 1972) and found that accreting NSs dominate ULX populations over BHs in galaxies like M82 and the Milky Way when comparing relatively young and old populations, respectively, with constant star-formation rates (SFRs). Wiktorowicz et al. (2017) found that BHs dominate only at early times in starbursts, while NSs dominate at later times in starbursts as well as in constant star-formation scenarios (for solar metallicity regions). However, for galaxies with a constant SFR in sub-solar metallicity environments, BH ULXs were still more abundant, as a lower metallicity leads to more compact, massive stellar cores that collapse to more massive COs. Wiktorowicz et al. (2019) found that young populations (with stellar ages ≲ 10 Myr) that were dominated by BH ULXs had more beamed BH ULXs than in older populations, while the majority of NS ULXs were beamed at all epochs. Consequently, the question of the dominant type of accretor involves the considerations of many assumptions that affect the underlying stellar populations.

In this work, we study populations of ULXs formed at starbursts of different ages, investigating the effects of age on ULXs. We also compare two sets of model assumptions that differ in certain physical properties, such as the natal kick velocities. Starburst populations give us valuable insights into the dependence of ULXs on the CO accretion and on the age of the stellar population and how these dependencies affect the sub-populations of NS ULXs and BH ULXs. Alternatively to synthetic starburst populations, simulations with continuous star-formation scenarios give us information about the general properties of ULX populations in large samples of galaxies while obscuring the effects of age on the populations. To carry out this study, we used the newly developed binary population synthesis code POSYDON (Fragos et al. 2023), an open-source framework that allows for population studies where the entire life of a binary is modeled using detailed and self-consistent stellar structure and binary evolution calculations. Our ULX population models are based on the earlier study by Misra et al. (2023), who modeled the X-ray luminosity function of extragalactic populations of HMXBs. Even though more ULXs per SFR have been observed in lower than solar metallicity environments (for instance, Prestwich et al. 2013; Basu-Zych et al. 2016; Kovlakas et al. 2020), the majority of observed ULXs are in solar metallicity environments (Kovlakas et al. 2020). Hence, we limited this study to solar metallicity (Z = 0.0142; Asplund et al. 2009).

In Sect. 2 we discuss the population synthesis code we employed and present the model assumptions used to carry out the population synthesis study and the calculations of X-ray luminosities of accreting COs. In Sect. 3, we present the various burst populations corresponding to a range of ages spanning from 5 to 1000 Myr. We discuss the demographics of the ULX populations at various ages, including the nature of the accretors and the donors, and explore the effect of accretion on observable X-ray pulses. In Sect. 4, we compare our results to those in the literature and present our concluding remarks.

2. Methods

We used POSYDON (Fragos et al. 2023, code version v1), a newly developed binary population synthesis code that incorporates detailed stellar structure and binary evolution tracks computed using the stellar evolution code Modules for Experiments in Stellar Astrophysics (code version 11701; Paxton et al. 2011, 2013, 2015, 2018, 2019; Jermyn et al. 2023). Accurate information about the structure of the stars via the detailed evolutionary sequences is expected to lead to a physically accurate estimation of the mass-transfer rate stability (particularly for thermal timescale mass-transfer), stellar rotation, and transport of angular momentum within and between the binary stars. Since all the previous studies of synthetic populations of NS and BH ULXs have been carried out, at least in part, using parametric population codes, it is imperative to study these populations by including information from detailed stellar structure calculations. Parametric codes, or rapid population codes, approximate the stellar structure using fitting formulas (e.g., BSE; Hurley et al. 2000, 2002) or look-up tables (e.g., ComBine; Kruckow et al. 2018) constructed from single stellar evolutionary tracks that can introduce systematic biases when binary interaction is involved.

2.1. Parameters of models used

In this study, we focus on the sub-population dominating ULX populations that are young, bright HMXBs (with high mass-transfer rates), and found in active star-forming galaxies. (Roberts et al. 2002; Gao et al. 2003; Kaaret et al. 2004; Wolter & Trinchieri 2004; Zezas et al. 2007; Anastasopoulou et al. 2016; Wolter et al. 2018; Kovlakas et al. 2020). Misra et al. (2023) studied the effects of different combinations of physics parameters on the synthetic X-ray luminosity function (XLF) of HMXBs and found models that best matched the observed XLFs of HMXBs from star-forming galaxies (taken from Lehmer et al. 2019). We took two models from that study, which we indicate as A and B (64 and 44, respectively, in Misra et al. 2023), with the model parameters described by Table 1. These models correspond to the ones that best fit the observations from two different approaches: first, matching the slope of the XLFs (here model A) and, secondly, matching the normalization of the XLFs (here model B).

Table 1.

Physical parameters corresponding to the best-fitting models A and B from Misra et al. (2023).

The primary difference between the two models is in the normalization of the supernova (SN) kicks during the formation of a BH. In model A, BH kicks are normalized using the mass of the newly formed BH, while in model B, the BH kicks are not normalized, and they receive the same kick velocities as NSs (drawn from a Maxwellian distribution with a velocity dispersion of 265 km s−1). The strengths of the kicks have an observable effect on XRB populations. Stronger kicks, similar to kicks received by NS binaries, result in the disruption of wide binaries (with orbital periods ≳ 100 days) during the SN event (Misra et al. 2023). Another difference between the two models is in the efficiency of the common-envelope (CE) phase, or αCE, which is 1.0 for model A (corresponding to full orbital energy being available for envelope ejection) and 0.3 (corresponding to 30% of the orbital energy being available for envelope ejection) for model B. The αCE parameter affects XRBs with He-rich donors, with lower efficiency leading to a higher merger rate of binaries in the CE phase. The formation of most XRBs, however, especially those that have H-rich donors and BH accretors, do not go through a CE, and hence the αCE parameter does not affect their evolution. As we will see later, XRBs with He-rich donors that consist of accreting NSs are affected by αCE, and this effect is more apparent when NSs start to dominate the ULX populations. Finally, there is an additional difference between the models, namely, in the observability of wind-fed accreting binaries. However, this parameter affects the XRBs at luminosities lower than 1038 erg s−1 and is hence below our limit for ULXs and has a negligible effect in our study.

2.2. Properties of the initial binary population

In order to study the effect of the age of the host population on the properties of ULXs, we simulated burst populations of 107 binaries, initially at zero-age main-sequence (ZAMS), with the range of burst ages being 5 Myr, 10 Myr, 40 Myr, 100 Myr, 300 Myr, and 1000 Myr. The reason we chose this range of burst ages is because theoretical studies of burst populations using rapid population codes have shown that variations in ULX demographics cease beyond ∼ 1 Gyr, which is when the populations are composed of only NS ULXs. From this point onward, these NS ULXs only decrease in number as their donors reach the end of their evolutionary lifetimes (Wiktorowicz et al. 2017, 2019). Each of our burst populations averages at around 2.7 × 108 M, and for the same initial total population mass, synthetic populations with younger burst ages sample binaries for a Milky Way-like galaxy better (which would have an SFR of about 54 M yr−1 for a starburst age of 5 Myr) as compared to older burst ages (which would have an SFR of about 0.27 M yr−1 for a starburst age of 1000 Myr). Over recent history, the SFR of the Milky Way is estimated to be less than 2 M yr−1 (Chomiuk & Povich 2011; Licquia & Newman 2015). Hence, the generated synthetic populations sufficiently sample the young stellar populations in star-forming regions where ULXs are prominent.

The initial binaries were generated by drawing the primary stellar mass from the Kroupa initial mass function (Kroupa 2001) in the mass range [7.0, 120.0] M and the secondary mass from a uniform mass-ratio distribution (Sana et al. 2013). The initial orbits were taken to be circular, with periods drawn from Sana et al. (2013) and with the same range as the POSYDON grids from 0.35 to 103.5 days (extrapolated down to 0.35, as the distribution from Sana et al. 2013 is limited to 1 day; Fragos et al. 2023). We normalized the populations using a correction factor of ≈5.89 (including a binary mass function of 0.7; Sana et al. 2012) to account for the unsampled region of the initial mass function, 0.08–7.0 M (Bavera et al. 2020). As mentioned previously, all the simulations were carried out at solar metallicity (0.0142; Asplund et al. 2009).

2.3. X-ray luminosity calculation

The X-rays from an XRB are generated from stellar material being captured by the CO. There are three types of mass-transfer phases in XRBs: mass transfer from the inner Lagrangian point in semi-detached Roche lobe overflowing binaries, mass transfer from stellar winds leaving the donor surface that are captured by the gravitational pull of the accretor, and mass transfer from a decretion disk of the donor (which is a highly-spinning B star, i.e., a Be star) when the CO interacts with the disk in an eccentric orbit. Roche-lobe overflow (RLO) occurs when the donor star starts to fill out the volume of the gravitational equipotential surface passing through the inner Lagrangian point, and material leaves the donor surface accreting onto the CO companion. Massive stars tend to have very high rates of stellar-wind loss, a fraction of which gets captured by the accretor. To account for wind mass transfer, we used the mechanism described by Bondi & Hoyle (1944).

The third type of mass transfer for Be XRBs is treated separately. Since the modeling of the decretion disk was not carried out in POSYDON, we identified Be XRBs and assigned them X-ray luminosities (using the same treatment as Misra et al. 2023). Abdusalam et al. (2020) carried out population synthesis studies for a binary population with a constant SFR of 3 M yr−1, including the effect of geometrical beaming of the X-ray emission. They found that pulsating ULXs tend to have Be- or intermediate-mass donors. The Be XRB luminosities in our populations reach up to ∼ 1038 erg s−1, and none appear as ULXs (see Misra et al. 2023). While there are a few Be XRBs that have been observed to have luminosities reaching ∼ 1039 erg s−1 during outbursts (Skinner et al. 1982; Wilson-Hodge et al. 2018; Tao et al. 2019; Vasilopoulos et al. 2020; Chandra et al. 2020), so far no Be XRBs have been observed to be persistently emitting ULX-like luminosities.

When the mass-transfer rate (tr) is below the Eddington limit, mass accretion onto the CO is assumed to be fully conservative: All the mass lost by the donor via RLO is accreted onto the CO. We used the following equation to estimate the X-ray luminosity for sub-Eddington mass-transfer rates, including RLO and wind-fed accretion,

L X RLO / wind = η M ˙ acc c 2 , if m ˙ 1.0 , $$ \begin{aligned} L^\mathrm{RLO/wind}_{X} = \eta \dot{M}_{\rm acc} c^2, \ \ \ \text{ if}\ \dot{m} \le 1.0, \end{aligned} $$(1)

where acc is the mass-accretion rate, η is the radiative efficiency of accretion, c is the speed of light, and tr/Edd is the Eddington ratio. The radiative efficiency of accretion is estimated by using the properties of the accretor,

η = G M acc R acc c 2 , $$ \begin{aligned} \eta = \frac{GM_{\rm acc}}{R_{\rm acc} c^2}, \end{aligned} $$(2)

where Macc and Racc are the mass and radius of the accretor (for BHs, Racc is the spin-dependent innermost stable circular orbit; Podsiadlowski et al. 2003), and G is the gravitational constant.

As the mass-transfer rate approaches and exceeds the Eddington limit, mass accretion onto the CO is limited by the Eddington limit (LEdd), and we followed the accretion disk model by Shakura & Sunyaev (1973) in order to model the accretion disk that receives material at a super-Eddington rate. This was combined with the prescription from King et al. (2001) and King (2009) for the additional effect of geometrical beaming. The isotropic-equivalent observed X-ray luminosity is defined as follows,

L X , isotropic RLO / wind = L Edd b ( 1 + ln m ˙ ) , if m ˙ > 1.0 , $$ \begin{aligned} L^\mathrm{RLO/wind}_{{X,\mathrm {isotropic}}} = \frac{L_{\mathrm{Edd} }}{b}(1 + \ln {\dot{m}}), \ \ \ \ \text{ if} \dot{m} > 1.0, \end{aligned} $$(3)

where is the ratio of the mass-transfer rate and the Eddington rate, and b is the beaming factor reflecting the geometrical collimation of the emission from the thick disk. The approximate value of b is given by King (2009),

b = { 73 m ˙ 2 , if m ˙ > 8.5 , 1 , otherwise . $$ \begin{aligned} b= \left\{ \begin{array}{ll} \frac{73}{\dot{m}^2},&\text{ if}\ \dot{m}> 8.5,\\ 1,&\text{ otherwise}. \end{array}\right. \end{aligned} $$(4)

For very high mass-transfer rates, the prescription above might lead to extremely strong beaming (b ≪ 10−3). We followed Wiktorowicz et al. (2017) in setting a lower limit for b at 3.2 × 10−3, which approximately corresponds to an opening angle of 9°. The beaming factor was used to down-weight all the synthetic XLFs presented in this work. While geometrical beaming reduces the number of binaries observed, it also leads to an increase in the luminosity of the brightest ULXs, as the binaries are more strongly beamed.

2.4. Suppression of pulses caused by accretion onto a neutron star

There is at least one confirmed NS ULX that does not (yet) show X-ray pulsations, and that is M51 ULX-8 (Brightman et al. 2018; Middleton et al. 2019). Mushtukov et al. (2017) suggested the lack of observed X-ray pulses is caused by an optically thick envelope that smears and obscures them. The reprocessing of the emitted hard radiation turns it into a blackbody-like emission, which depends on the optical thickness of the envelope and makes the cyclotron scattering lines disappear (which is characteristic of magnetic NSs). The observation of pulses, therefore, depends on the viewing angle with respect to the NS spin and the envelope.

King & Lasota (2020) postulated that the observations of pulses in accreting NSs are due to the misalignment of the NS spin axis and the beaming axis (refer to Fig. 2 in their paper for more clarity). They showed that sinusoidal pulsed light curves are created when either of two conditions is achieved: (1) the spin and accretion disk axes are strongly misaligned and the magnetic axis is aligned close to the spin, or (2) the spin and disk axes are strongly aligned and the magnetic axis is misaligned with the spin. For both of these configurations, the escaping pulses are different. If the spin and disk axes are strongly misaligned, pulses escape the beaming tunnel (as the magnetic axis, and hence the hot spot, follow the spin axis) and can be observed. However, accretion tends to rapidly align the spin axis of an NS with the orbital axis. When aligned, as pulsed light curves are created when the magnetic axis is misaligned with the NS spin, most pulses are scattered by the beaming tunnel and pulsations are negligible.

The magnetic field of an NS can affect its accretion phase if it is strong enough (about ≳ 1012 G). A strong magnetic field has a large Alfvén radius, allowing the accreted angular momentum to efficiently spin up the NS. However, a weak magnetic field allows the accretion disk to extend all the way to the NS surface. To investigate how accretion-induced alignment of the NS axes would result in observable pulses, we performed an order of magnitude calculation of the angular momentum of the incoming matter (Jacc) at the magnetosphere radius (RM) and at the NS surface (RNS) and equated it to the angular momentum of a fiducial NS (JNS) in order to find the amount of mass that should be accreted for the pulses to be buried,

J acc = J NS I acc ω acc = I NS ω NS , $$ \begin{aligned} J_{\rm acc} = J_{\rm NS} \Leftrightarrow I_{\rm acc} \omega _{\rm acc} = I_{\rm NS} \omega _{\rm NS}, \end{aligned} $$(5)

where INS is the accretor moment of inertia and ωNS is the accretor angular velocity. When the condition in Eq. (5) is fulfilled, the NS is assumed to have gained enough angular momentum to align its spin and orbit axes, and there will be no observable pulses. Using the Keplerian velocity for the velocity of matter at radius r, we got ω acc = G M NS / r 3 $ \omega_{\mathrm{acc}} = \sqrt{GM_{\mathrm{NS}}/r^3} $ and the moment of inertia for accreted matter Δm as Iacc = Δmr2. Thus, we arrived at an equation describing the accreted mass,

Δ m = I NS ω NS G M NS r . $$ \begin{aligned} \Delta m = \frac{I_{\rm NS}\omega _{\rm NS}}{\sqrt{GM_{\rm NS}r}}. \end{aligned} $$(6)

We assumed INS = 1045 g cm2 (approximate value for a 1.3 M NS) and took ωNS = 2π/Pspin, with Pspin = 1s (similar to the spin period of the most well-studied pulsating ULX M82 X-2, which is 1.37 s; Bachetti et al. 2014). For the matter being accreted at the surface of the NS (hence, r is taken to be the NS radius, which is RNS = 12.5 km; Most et al. 2018; Miller et al. 2019; Riley et al. 2019; Landry et al. 2020; Abbott et al. 2020; Kim et al. 2021; Biswas 2021; Raaijmakers et al. 2021), we could get an estimate for the amount of matter that needs to be accreted for pulses to be unobservable,

Δ m ( R NS ) = 2.07 × 10 4 ( M NS 1.4 M ) 1 / 2 ( R NS 12.5 km ) 1 / 2 ( P spin 1 s ) 1 M . $$ \begin{aligned} \Delta m (R_{\rm NS}) = 2.07\times 10^{-4} \bigg (\frac{M_{\rm NS}}{1.4\,M_{\odot }}\bigg )^{-1/2} \bigg (\frac{R_{\rm NS}}{12.5\mathrm {km}}\bigg )^{-1/2} \bigg (\frac{P_{\rm spin}}{1\,\mathrm s}\bigg )^{-1} \,M_{\odot }. \end{aligned} $$(7)

The case described above would correspond to the case when there is a weak magnetic field, and the accretion disk can be approximated to extend down to the NS surface. In an NS with a strong magnetic field, the magnetosphere radius is defined as the radius where the magnetic stresses dominate the accretion flow (Frank et al. 2002), and it is defined as

R M = ( 2 G M NS ) 1 / 7 M ˙ NS 2 / 7 μ 4 / 7 , $$ \begin{aligned} R_{\rm M} = (2GM_{\rm NS})^{-1/7} \dot{M}_{\rm NS}^{-2/7} \mu ^{4/7}, \end{aligned} $$(8)

where μ is the NS magnetic dipole moment (which has a typical value of 1030 G cm3 that corresponds to an NS surface magnetic field strength of B NS μ / R NS 3 10 12 G $ B_{\mathrm{NS}} \sim \mu / R_{\mathrm{NS}}^{3} \approx 10^{12}\,\rm G $). By simplifying Eq. (8) and rewriting it in convenient units, we obtained

R M = 1.83 × 10 8 ( M NS 1.4 M ) 1 / 7 ( M ˙ NS 10 8 M yr 1 ) 2 / 7 ( μ 10 30 G c m 2 ) 4 / 7 cm . $$ R_{\rm M} = 1.83\times 10^{8} \bigg (\frac{M_{\rm NS}}{1.4\, M_{\odot }} \bigg )^{-1/7} \bigg (\frac{\dot{M}_{\rm NS}}{10^{-8}\, M_{\odot } \,\mathrm{yr}^{-1}} \bigg )^{-2/7} \bigg (\frac{\mu }{10^{30}\,\mathrm G\, cm^{2}}\bigg )^{4/7} \mathrm {cm}. $$

Using Eqs. (9) and (6), we could then estimate the mass required to spin up the NS to Pspin if all the matter at RM is accreted due to the magnetic field of the NS,

Δ m ( R M ) = 1.71 × 10 5 ( M NS 1.4 M ) 3 / 7 ( M ˙ NS 10 8 M yr 1 ) 1 / 7 ( μ 10 30 G c m 2 ) 2 / 7 ( P spin 1 s ) 1 M . $$ \begin{aligned}&\Delta m (R_{\rm M}) = 1.71\times 10^{-5} \bigg (\frac{M_{\rm NS}}{1.4\, M_{\odot }}\bigg )^{-3/7} \bigg (\frac{\dot{M}_{\rm NS}}{10^{-8}\,M_{\odot }\, \mathrm{yr}^{-1}} \bigg )^{1/7} \nonumber \\&\qquad \qquad \bigg (\frac{\mu }{10^{30}\,\mathrm G\, cm^{2}}\bigg )^{-2/7} \bigg (\frac{P_{\rm spin}}{1\,\mathrm s}\bigg )^{-1} \, M_{\odot }. \end{aligned} $$(9)

After the NS has accreted enough matter so that the axis of the NS will be aligned with the disk axis, pulsations will not be observed. Hence, we used Eqs. (7) and (9) to identify ULXs with NSs that have not yet accreted enough matter in order to suppress the emitted X-ray pulses. The default populations presented in our work do not automatically take into account the suppression of beamed emission in the detectability of NS ULXs. We present the results including this effect separately in Sect. 3.4.

We note that the presented model assumes that the X-ray radiation comes from the hot spot created on the NS surface that follows the magnetic fields lines. If the emission is from the walls of the accretion column with a fan-beamed pattern, as postulated for the ULX source RX J0209.6-7427 (Hou et al. 2022), the assumed model might not be accurate. In this case, the inferred accretion rate is also super-Eddington. Since the deciding factor in the aligning of the spin and beaming axes is the amount accreted by the NS (defined by Eqs. (6) and (9)), for accretion rates higher than the standard scenario, the alignment would be much quicker and even fewer NS ULXs would be seen.

3. Results

Stars of various masses evolve differently with time, leading to different demographics of the ULX populations from progenitors formed at different ages. Using population synthesis, we investigated the effect of age on ULXs. Figure 1 shows the distributions of isotropic-equivalent X-ray luminosities (calculated using the description in Sect. 2.3) of ULXs in the burst populations (for both the models A and B) without down-weighting the geometrically beamed ULXs. At high luminosities in an observed population, the sources would be beamed, and only a fraction of them would be observed. Overall, as the burst age increases, the number of ULXs decreases. This is because older populations are dominated by LMXBs that are comparatively less likely to reach ULX-like luminosities than HMXBs. For both models, the oldest population (with the burst age of 1000 Myr) has no ULXs and, therefore, does not appear in the figures presented. The highest luminosity reached for XRBs with age 1000 Myr for model A is 6 × 1038 erg s−1, while for model B, it is 5 × 1037 erg s−1 (well below our X-ray threshold for ULXs).

thumbnail Fig. 1.

Distribution of X-ray luminosities of ULXs for all the burst ages (see legend) for the populations following models A and B (left and right panel, respectively). Geometrical beaming is accounted for in the X-ray luminosity, but the observability of a ULX is not. In the observed population, the brighter luminosities would be slightly suppressed, as they would be beamed.

Going from model A to B, the number of ULXs reduces overall, while the peak luminosities remain similar (see Fig. 1). The effect is stronger for younger populations (Nmodel A/Nmodel B ≈ 3, 4, 7, corresponding to burst ages of 5 Myr, 10 Myr, and 40 Myr, respectively), while older populations do not show the same trend (Nmodel A/Nmodel B ≈ 0.7). Younger populations are dominated by BH accretors, which are affected more by the kick normalization assumed with respect to older populations. For model A, BH kicks are normalized using the mass of the newly formed BH, causing heavy BHs to receive negligible kicks, while the kicks for model B are not normalized and receive NS-like kicks. As seen in the parameter study by Misra et al. (2023), different kick normalizations result in a reduction in the number of wind-fed XRBs with increasing kick strengths. The corresponding effect on RLO XRBs is not as straightforward. With increasing kick strength, though more binaries are disrupted, the number of SN surviving binaries that undergo RLO in close orbits increases since more post-SN binaries have high eccentricities. Therefore, even though we do see a slight reduction in ULXs with stronger kicks, there is not much change in the peak X-ray luminosity.

Figure 2 shows the different burst populations sub-divided into the types of XRBs, specifically for model A (correspondingly, the burst populations for model B are shown in Fig. A.1), extending down to luminosities around 1038 erg s−1, which is below our ULX threshold. As expected, younger populations are fully dominated by HMXBs (below 10 Myr), and with age, the prevalence shifts to IMXBs (at ages 40–100 Myr) and then to LMXBs (older than 100 Myr). Populations older at 1000 Myr do not contain ULXs. The dominance of IMXBs in ULXs lasts for a shorter duration (being prominent only at 40 Myr) than LMXBs since IMXBs have shorter RLO phases. Their intermediate-mass donors cannot drive winds as high as those of high-mass donors until the RLO provides them with enough of a mass-transfer rate to emit X-rays. The behavior of the populations is similar between models A and B with respect to the burst ages. Therefore, IMXBs are viable candidates to explain ULX luminosities at these ages, reaching ∼ 1040 erg s−1, as was also suggested by Tauris et al. (2017) and Misra et al. (2020).

thumbnail Fig. 2.

Synthetic XLFs of burst populations for model A. The panels show the different types of XRBs, namely, HMXBs, IMXBs, and LMXBs, as described in the legend. The combined XLF is shown by the dotted-gray line.

3.1. Nature of the accretors

To investigate the nature of accretors present in ULXs, we looked at the population of accreting BHs with respect to the total XRB populations. Figure 3 shows the percentage of BH XRBs (and BH ULXs) present in the total XRBs (and ULXs) across all the burst ages. The populations with the burst age of 1000 Myr (for both models A and B) have no ULXs. Therefore, they are not included in this figure. At 5 Myr, the populations have no NSs since the population is too young for NS progenitors to evolve out of the main sequence (MS), and they are 100% dominated by accreting BHs. The stars that would form NSs have initial masses in the approximate range of 8–25 M, with the IMF favoring the lower masses. The MS lifetime for the stars within this range is ∼1 to 10 Myr. Since lower masses dominate the populations, many of the NS progenitor stars are still in their MS. The number of BH XRBs (and BH ULXs) decreases with time as NSs increase in number because the donors present in BH XRBs are mostly in the intermediate- to high-mass range. With time, an increasing number of these stars would undergo core collapse. Due to stronger BH kicks in model B, the number of BH binaries surviving the SN is less than in model A.

thumbnail Fig. 3.

Percentage of BH binaries present in the simulated ULXs (solid lines) and in the simulated XRBs (dashed lines) following the parameters of models A (shown by orange) and B (shown by blue).

3.2. Properties of typical ultra-luminous X-ray sources in each population

To get a better idea of the properties of ULXs in the synthetic burst populations, we looked at the synthetic XLFs of the populations. Figure 4 shows the sub-populations of each of the simulated populations, describing the type of the XRB (RLO or wind), type of accretor (BH or NS), and type of donor (H-rich or He-rich) for model A. The same for model B is shown in Fig. A.2. Again, the figures extend to luminosities of ∼ 1038 erg s−1, which is below our ULX threshold. In general, all ULXs are dominated by XRBs undergoing RLO, with some contribution from wind-fed XRBs at wide orbits for model A. There are significantly fewer wind-fed ULXs in the ULX populations when following model B due to a disruption of wide orbits from strong SN kicks (see Fig. A.2). Additionally, the distributions of the orbital properties for all burst ages and both models are presented in Figs. 5 (for accretor masses), 6 (for donor masses), and 7 (for orbital periods). The orbital properties for the simulated ULXs are not significantly different between models B and A; hence, we primarily discuss the ULXs for model A and point out the differences with model B when significant.

thumbnail Fig. 4.

Synthetic XLFs of burst populations for model A. The populations are further split into the type of the XRBs (RLO or wind), type of accretors (BH or NS) and type of donors (H-rich or He-rich). The combined XLF is shown by the dotted-gray line.

thumbnail Fig. 5.

Distribution of accretor masses for ULXs in the synthetic populations across ages 5–300 Myr.

thumbnail Fig. 6.

Distribution of donor masses for ULXs in the synthetic populations across ages 5–300 Myr.

3.2.1. Burst age of 5 Myr

We utilized the aforementioned figures to derive the orbital properties of typical ULXs for all burst ages, described by the peaks of the respective distributions. For the population with a burst age of 5 Myr, a typical ULX has a BH in the mass range of about 10–15 M (the minimum BH mass being 10.25 M), with a massive 28 M H-rich MS donor undergoing RLO in an orbit of about 10 days. The binary, after the primary SN event, evolves further until the massive donor fills its Roche lobe, at around 10 days, and it starts an intense mass-transfer phase while it is still on the MS, leading to its appearance as a ULX. As the stronger kicks in model B are more prone to disrupting wider orbits, the number of ULXs drops off very quickly with an increasing orbital period compared to model A (see Fig. 7). The distribution of the donor masses differs between the two models (see Fig. 6). Model B has a flatter distribution compared to model A because the differing binaries between the two models are the ones that were disrupted by the strong kicks in model B. These binaries correspond to a wide range of orbital periods (many hundreds of binaries with ≳3 days) and massive donor masses (≳25 M).

thumbnail Fig. 7.

Distribution of orbital periods for ULXs in the synthetic populations across ages 5–300 Myr.

3.2.2. Burst age of 10 Myr

A population at 10 Myr has a typical ULX, with a 9 M BH, a 10 M H-rich MS donor, and an orbital period of 4 days undergoing RLO. The distribution of the orbital parameters is similar to the ULXs at 5 Myr (see Sect. 3.2.1), although with the difference of a shift in the donor mass to less massive (see Fig. 6) and more compact stars, and hence the difference also corresponds to narrower orbits, as can be seen in Fig. 7. The radii of donors at 5 Myr peaks at around 36 R, while for 10 Myr the peak is at 12 R. Additionally, there is an emergence of NS accretors at this age.

3.2.3. Burst age of 40 Myr

At 40 Myr, a typical ULX has a similar CO mass as that at 10 Myr, with a 9 M accreting BH and an intermediate-mass (around 4.5 M) H-rich MS donor in an RLO orbit with a period of 10 days. The peak of the orbital period distribution is larger than that at 10 Myr because at this age (see Fig. 7), even though the typical donor star is less massive than before, the mass is being transferred from a lower-mass star to a more massive CO, leading to orbital expansion, which contrasts with the process at the younger age when the mass ratios were more even. There is another peak in the ULX population at orbital periods ≳1000 days (only for model A) corresponding to H-rich post-MS donors that have expanded to large radii and are transferring matter effectively via wind-fed accretion. This population of wind-fed BH XRBs can already be seen in model A at 10 Myr. However, as can be seen when comparing Figs. 2 and A.1, the wind-fed HMXBs disappear completely in the population following model B due to the strong BH kicks disrupting wide orbits that would have led to wind-fed accretion. Additionally, NS accretors begin to gain noticeable numbers at this age (see Fig. 5), with about 5% and 40% of NSs present in the XRBs (for models A and B, respectively), as can be seen in Fig. 3.

3.2.4. Burst age of 100 Myr

At 100 Myr, a typical ULX has an NS accretor with a mass of 1.2 M and a low-mass stripped-helium (stripped-He) donor with 1.2 M in an orbit of 2 days. After the primary SN event, the donors are H-rich MS stars and are in wide orbits (centered around 1000 days). Consequently, when the donors fill their Roche lobes, they have depleted their core-He. Since they are rapidly expanding He-giants, the binary interaction is more prone to instability, and a CE phase occurs that strips the outer H-rich envelope of the donor. The surviving NS+He stripped binary then evolves until the donor fills its Roche lobe again, starting its XRB phase. The peak of the orbital period distribution for model B is at a much smaller orbit, at 0.5 days, compared to model A (which is at 2 days), as shown in Fig. 7. Since the effect of different BH natal kicks is not applicable to a majority of the population at this burst age, the secondary difference between the models becomes more apparent, such as CE efficiency or αCE. The value of the αCE is approximately proportional to the final orbital separation at the end of the CE phase and determines whether this interaction is stable or not. Going from model A to model B, the αCE decreases from 1.0 to 0.3, which would intuitively lead to more binaries with narrower post-CE orbits and higher instability (refer to Misra et al. 2023, for a detailed study on the effect of various assumptions of the CE phase on XRBs). However, in the context of ULXs, the model B ULXs experience an increase in the fraction of binaries with narrower orbits after the CE in the surviving binaries. Hence, the difference in the distributions of orbital periods in Fig. 7.

In the distribution of the ULX donor masses shown in Fig. 6, apart from the He-rich donors, there is another population of intermediate-mass stars (around 5 M). The same binaries result in the second peak in the orbital periods seen for model B in Fig. 7. These ULXs correspond to post-MS donors that previously underwent stable RLO and did not lose their H envelopes. For these binaries (and for both sets of model assumptions), RLO begins after the first SN when the donor is an H-giant star (in the H shell-burning phase). The relative expansion of the donor radius with respect to the Roche-lobe radius is not as high as when the donors are He-giants. The stability of the mass transfer is linked to the evolutionary state of the donor (Misra et al. 2020), as stars at certain later evolutionary stages have large radii with deep convective regions that expand rapidly on mass loss. If the rapid mass loss exposes an underneath radiative layer, the radius contracts, stabilizing the binary. Otherwise, rapidly expanding donors would lead to a CE.

3.2.5. Burst age of 300 Myr

A typical ULX at 300 Myr has an NS accretor with 1.2 M and a low-mass post-MS H-rich donor with 0.6 M in an orbit of 8 days. For model B, there is an additional XRB population with stripped-He donors that is not present in model A. However, their luminosities do not cross our ULX threshold of 1039 erg s−1. These XRBs dominate the population at luminosities ≲ 1038 erg s−1 at this age while appearing as ULXs at 100 Myr.

3.3. Beamed X-ray emission

An aspect that would affect the observed versus underlying population of ULXs is the effect of geometrical beaming. The stronger the geometrical beaming of the X-ray emission, the smaller is the chance of observing the ULX. In our synthetic populations, we defined beamed emission using the factor b, which is defined in Eq. (4), with beamed emission as b < 1 and unbeamed emission as b = 1. Hence, the observed luminosities from super-Eddington accretion disks are enhanced by the factor b.

Figure 8 shows the average beaming factor for each of the ULX populations, with the populations categorized by the type of accretor (NS or BH). We note that the populations with no ULXs (at 1000 Myr), no NS ULXs (at 5 Myr), or no BH ULXs (at 100 and 300 Myr for model B) are not shown. The figure also shows the fraction of beamed ULXs in each population. As expected, NSs are more strongly beamed than BHs, with their beaming factors less than approximately 0.25 across all burst ages and both models. The BHs have relatively less beamed emission and are visible more than 50% of the time whenever they are present in the population. With age, the average beaming factor for BH ULXs increases, that is, their emission becomes more isotropic on average, in agreement with Wiktorowicz et al. (2019). This effect is stronger for model A, unlike model B where BHs disappear from the population by 100 Myr. For the populations with no NS ULXs at 5 Myr, the beamed-to-unbeamed ratio is almost equal, with a slight dominance for unbeamed ULXs in model B, as strong kicks in model B reduce the number of close BH binaries that would undergo RLO and have high mass-transfer rates.

thumbnail Fig. 8.

Average beaming factor ( b ¯ ) $ (\bar{b}) $ for the ULX populations spanning all burst ages, distinguishing between the type of accretor (NS, shown by circles, or BH, shown by squares). Populations that did not have any ULXs are not shown. The figure also shows the fraction of beamed (with b < 1.0) to total ULXs in each population, shown as dashed lines. The blue colored symbols and lines correspond to the results for model A, and the orange-colored symbols and lines are for model B.

At 10 Myr, the beamed sources decrease to about 20% for both models while still being dominated by BH ULXs. The reason for the decrease in the number of beamed ULXs is the shift to less massive donor masses as the population age increases (the peak of the donor mass distribution goes from ∼25 M to ∼10 M; see Fig. 6). These less massive donors have longer MS lifetimes and smaller radii compared to the massive donors at 5 Myr. Therefore, a smaller fraction of donors has super-Eddington mass-loss rates that require considerable overflowing of the Roche lobe during RLO, resulting in decreased beaming (for mass transfer rates below 8.5 times the Eddington limit, the beaming factor is one, following the model from King 2009). Less massive MS donors typically do not drive mass-transfer rates to values as high as more massive MS donors. This is also the reason for the decrease in geometrical beaming for BH ULXs at 10 Myr. The small population of NS ULXs that appears at this burst age is fully beamed.

For model A at the burst age of 40 Myr, the ratio of beamed versus unbeamed ULXs remains similar to the one at 10 Myr, but the ratio is reversed for model B, as most ULXs become beamed (∼71%). The reason for this discrepancy is the change in the demographics of the beamed ULXs between the two models. Due to the difference in BH kicks, the ULX population for model B is dominated by NSs, which are always beamed, owing to their Eddington limits that are an order of magnitude lower than BH XRBs. For ages greater than 40 Myr, the fraction of beamed ULXs is around ≳94% in model A and 100% in model B since both populations are dominated by NS ULXs. As can be seen in Fig. 8, NS ULXs, whenever they are present, are strongly beamed, while BH ULXs are not as beamed across all the burst ages.

3.4. Suppression of pulses caused by accretion onto a neutron star

It has been suggested that the accretion phase would affect the observations of pulses seen from an NS accretor (Mushtukov et al. 2017; King & Lasota 2020). Using Eqs. (7) and (9) to define the accreted matter (at various radii) required to align the NS spin axis to the orbital axis, we modeled the suppression of pulses and investigated the relative contribution of pulsed versus non-pulsed systems in ULX populations as well as the distribution of the types of COs (BHs and NSs). In our investigation, the term “non-pulsed” simply refers to NSs that have accreted matter greater than the amounts described by Eqs. (7) and (9), as their X-ray pulsations are suppressed. Figure 9 shows the number of NS ULXs per stellar mass for models A and B, depicting the total population and the NS ULXs that would have observable pulses for angular momentum accreted at two radii. The two radii in question are the NS radius (RNS) and the magnetospheric radius (RM). The populations that are 10 Myr old and have a very small number of NS ULXs present (∼10−9/M) all have observable pulses. However, since the populations at 10 Myr are dominated by BH ULXs (almost 100%), the chance to observe these NS ULXs would be insignificant.

thumbnail Fig. 9.

The number of NS ULXs for models A and B (left and right panels, respectively) per stellar mass in the burst populations. With green filled circles, we show the number of NS ULXs with observable pulses using the NS radius (RNS) as the accretion radius. With pink circles we show the number of NS ULXs with observable pulses using the magnetospheric radius (RM), and the mustard circles show the total number of NS ULXs.

For other ages, we observed that a small fraction of the total NS ULXs emit observable pulses, especially for populations with ages from 40 Myr to 100 Myr (when NSs start to dominate the ULX population) with pulses from NS ULXs observable for 30% to 50% of the binaries. An important feature is present at 100 Myr, where the pulses from about 30% of NS ULXs may be observed when the accretion occurs at the RNS and ≲2% when accretion is at RM. The matter accreted at RM carries more angular momentum compared to RNS due to its larger distance. Therefore, less material needs to be accreted for the pulses to be suppressed, and the number of NS ULXs showing pulses is significantly lower for accretion at RM. For a typical NS (with surface magnetic field ∼ 1012 G), the magnetospheric radius is around 2 × 108 cm. However, for an NS with higher field strength, such as ∼ 1014 G (i.e., a magnetar), the magnetospheric radius is larger by a factor of around 14. This increased accretion radius would lead to a faster suppression of observable pulses, as the angular momentum transferred from the accreted material would be larger than the accretion from a smaller radius. Consequently, magnetar ULXs with pulses are more difficult to observe. At 300 Myr, the fraction of NS ULXs with observed pulses drops down to 5% for accretion at RNS and to 0% for accretion at RM since the accreted matter increases with time, eventually crossing the threshold set by Eqs. (7) and (9).

Using Eqs. (7) and (9), we estimated the timescales of the accretion phases in the ULXs. Varying the NS spin period between 0.1 and 10 s, the mass between 1.2 and 2.0 M, and the magnetic field strength between 108 and 1014 G, the resulting timescales are between 1 and 105 yr (over both estimates of the accretion radius). The medians for the timescale distributions are around 500 yr (for accretion at RNS) and around 100 yrs (for accretion at RM). While these periods of time imply that most NS ULXs should not be seen, Figure 9 clearly shows about 30% of the observed NS ULXs at 100 Myr (for model A). This unexpectedly high number can be explained by the fact that most of the NS ULXs (between the ages 40 and 100 Myr) that have any observable pulses have H-giant donors, whereas the rest of the hidden population is dominated by either stripped-He or H-rich MS donors. Stable mass transfer for giants donor stars occurs over thermal timescales, and for MS donors stars, it occurs over nuclear timescales (a few orders of magnitude longer than the thermal timescales). Hence, even though the mass-transfer rates of the mass transfer from giant stars (with deep convective envelopes) go to high values (often exceeding the Eddington limit of the accretor), the phase itself is short lived. At 10 Myr, the few NS ULXs present have pulsed emission, as the systems are young and not enough time has passed to fully align the spin and beaming axes. At 300 Myr, both sets of NS ULXs (pulsed and non-pulsed) have H-giant donors, with higher masses (∼3 M) tending toward the pulsed NS ULXs and lower masses (≲2 M) going toward non-pulsed.

Comparing the NS ULX populations between the two models in Fig. 9, we found an increase in the total number of NS ULXs when going from model A to model B, most notably at the burst age 100 Myr, where the total NS ULX number increases by about 1.5 times. The reason is the difference in the values of the αCE between the models (see Table 1). Lower values of αCE lead to narrower orbits of binaries that survive the CE phase. These narrow binaries lead to higher mass-transfer rates during RLO and therefore a higher number of ULXs. However, this increase in the number of NS ULXs with decreasing αCE was not as distinctly reflected in the number of ULXs where escaped pulses could be seen. This occurs because even though more NS binaries are in the ULX phase, due to a high degree of super-Eddington mass transfer, more of them are able to accrete enough material to fulfill the conditions described by Eqs. (7) and (9). Hence, the number of NS ULXs with observable pulses increases only marginally.

3.5. Hyperluminous X-ray sources

In our synthetic populations, HLXs are identified by their extremely high X-ray luminosities of ≳ 1041 erg s−1, which is much brighter than typical ULXs. Even though most observed HLXs are believed to be accreting intermediate-mass BHs (Farrell et al. 2011; Servillat et al. 2011; Lasota et al. 2011; Davis et al. 2011; Yan et al. 2015), it is worth exploring extreme luminosities in the context of accreting NSs and stellar-mass BHs. Wiktorowicz et al. (2015) found a large contribution from NS ULXs to HLXs, as accreting NSs would tend to exhibit stronger geometrical beaming (their Eddington luminosity is much lower compared to BHs). In our ULX populations, all the ULXs reaching these extreme luminosities are geometrically beamed (having b < 1). However, the amount of beaming differs at different burst ages.

We refer to Figs. 4 and A.2 to study the kinds of binaries that appear as HLXs in our POSYDON populations. In these figures, we also show the fraction of the different accretor types that received the strongest beaming (with b = 0.0032; our lower-threshold for b) for our populations of HLXs in Fig. 10. For young populations at burst ages 5 Myr and 10 Myr, BH XRBs with massive H-rich donors undergoing RLO easily reach extreme luminosities via geometrical beamed emission. Since the population at 5 Myr only has BH accretors, all the strongest beamed HLXs are BH accreting binaries. At 10 Myr, even though the HLXs are still dominated by BH XRBs, HLXs with the strongest geometrical beaming begin to include some NS ULX systems (≲5%). At 40 Myr, different populations increase to comprise the HLX population, including wind-fed BH XRBs and RLO NS XRBs, all of which have H-rich donors (model B does not have the population of wind-fed BH ULXs present in model A). At this age in model B, the strongly beamed HLXs are fully dominated by the NS ULXs (at 100% of the HLXs), as strong kicks have reduced the number of surviving BH binaries. All the strongly beamed HLXs in the older populations are fully dominated by NS accretors with different types of donors (H-rich or He-rich).

thumbnail Fig. 10.

Percentage of HLXs (having X-ray luminosities ≳ 1041 erg s−1) with maximally beamed emission (with b = 3.2 × 10−3) for all the burst ages, distinguishing between the type of accretor (NS, shown by circles, or BH, shown by squares). The blue symbols show the results for model A, and the orange symbols are for model B.

4. Discussion and conclusions

To summarize, we investigated the effect of stellar age on populations of ULXs by generating various populations at fixed burst ages. We used POSYDON, a new binary population synthesis code that utilizes detailed calculations of binary interactions, to generate populations of ULXs at different burst ages, following two of the best-fitting models from Misra et al. (2023). The main differences between the two models are listed in Table 1. We ran simulations of 107 ZAMS binaries per population, each population corresponding to burst ages of 5 Myr, 10 Myr, 40 Myr, 100 Myr, 300 Myr, and 1000 Myr. Using the comparison of the resulting ULX populations, we reached the following key conclusions:

  • The treatment of the BH kicks greatly affects the total amount of XRBs observed. The higher kicks in model B (where BH kicks are not normalized) led to a higher disruption of binaries, reducing the overall binaries with ULX-like luminosities (decrease of up to seven times at the starburst age of 40 Myr). As this affects the demographic of the populations, it also results in a relative dominance of NS ULXs in populations with stronger kicks (model B).

  • We find that geometrically beamed emission is affected by the burst age as well as the BH kick prescription, and on average, less than 50% of BH ULXs and all of NS ULXs are beamed. Populations with stronger kicks (where NSs dominate) have a higher fraction of ULXs with beamed emission.

  • There is an inverse correlation between the number of ULXs and the age of the burst, with the largest number of ULXs at 5 Myr and none at 1000 Myr. We find that younger ULX populations (≲ 10 Myr) are dominated by HMXBs, populations with ages ∼ 40 Myr are dominated by IMXBs, and older populations (100–300 Myr) are dominated by LMXBs.

  • As the populations age, the ratio of BH to NS accretors decreases, which affects the other orbital properties of ULXs. The type of ULX donor changes, with younger populations (< 100 Myr) having massive H-rich MS donors, while older populations (> 100 Myr) are dominated by low-mass H-rich post-MS donors. At 100 Myr, stripped He-rich donors are prominent in ULXs.

  • Accretion plays a significant role in the observation of pulsating ULXs, affecting the fraction of NS ULXs inferred in a population. Depending on the radius at which the accretion disk begins (at the NS radius or the magnetospheric radius), the number of NS ULXs where pulses are observed varies (from 30% to none of the total NS ULX population depending on the stellar age), implying many more NS ULXs might be present in the ULX populations compared to those that show coherent pulsations. Additionally, NS ULXs with strong magnetic fields (≳ 1012 G) are increasingly difficult to observe.

  • Hyperluminous X-ray sources (ULXs with luminosities ≳ 1041 erg s−1) are present at all the burst ages explored (except at 1000 Myr), with an inverse correlation with age. All the HLXs, at all ages, have geometrically beamed emissions.

In addition to the effects of age on ULX populations, we also observed a reflection of the physical assumptions made about binary evolution. We find that BHs dominate ULXs in younger starburst populations, and NSs start to occupy more space in the demographic after about 100 Myr, confirming the results of Wiktorowicz et al. (2017). Wiktorowicz et al. (2017) also found that the production of BH ULXs is limited to up to 200 Myr. Similarly, in our POSYDON burst populations, BH ULXs appear before 100 Myr (dominating the populations up to 40 Myr), after which NS ULXs are the only main contributors, with negligible contributions from BH ULXs. The specific contribution from BH ULXs also depends on the natal kicks because model A (with weaker BH kicks) has more BH accretors than model B (see Fig. 3 and Sect. 3.1). Since we generated ULX populations from star-forming bursts of specific ages, the results presented in this study are not directly comparable to observations, but the populations can be convolved with star-formation histories to do so.

Our results regarding the geometrically beamed emission in ULXs are similar to those presented by Wiktorowicz et al. (2019). They found that while most NS ULXs are beamed whenever they are present, the fraction of BH ULXs that have beamed emission decreases with age. This also corroborates the work by King & Lasota (2016), where they suggested that NS ULXs are more likely to be beamed compared to BH ULXs due to the lower NS Eddington limit. However, the calculations involved are heavily dependent on the simplistic model used to define X-ray luminosity from super-Eddington accretion disks (described in Sect. 2.3) and might show different results with a different model. King & Lasota (2020) introduced the idea that many pulsating ULXs have a strong misalignment between the NS spin axis and the accretion disk axis, which enables the X-ray pulses to escape and be observed. As the NS accretes matter and hence angular momentum, the axes align and X-ray pulses are further suppressed. Generally, we found in our study that X-ray pulses are suppressed for 70%–100% of the NS ULXs, depending on the population age and assuming a certain accretion model.

While the ULX populations were generated at solar metallicity, ULXs have been observed in excess in low-metallicity regions compared to solar metallicity, when normalized by the SFR (Soria et al. 2005; Soria 2007; Mapelli et al. 2009, 2010; Linden et al. 2010; Kaaret et al. 2011; Prestwich et al. 2013; Basu-Zych et al. 2016; Kovlakas et al. 2020). This behavior is attributed to the fact that at lower metallicities, stars have more compact and massive stellar cores (Mapelli et al. 2009; Zampieri & Roberts 2009) that can initiate RLO in close orbits, leading to a higher production of HMXBs per SFR (Majid et al. 2004; Dray 2006; Soria 2007; Mapelli et al. 2009; Linden et al. 2010; Fragos et al. 2013a,b; Basu-Zych et al. 2016), which has been reproduced by theoretical studies (Zuo et al. 2014; Wiktorowicz et al. 2017). In the future, we will carry out a study of ULX populations with detailed population synthesis techniques, as used in this work, at low-metallicities, as it will further our understanding of ULXs in galaxies where the aforementioned excess is observed, as well as in the high-redshift Universe.

This study was carried out to observe the effects of age on ULX populations and the general trends in the populations while also investigating the effects of certain model assumptions. We find that while the populations overall followed similar trends between different models of physical assumptions, comparative studies between observations and simulations could shed light on certain physical processes. Additionally, extending this study to include a wide range of metallicities and star-forming histories would be the next step in answering many questions related to ULXs and the field of binary evolution as a whole.

Acknowledgments

The authors thank the anonymous referee for their con- structive comments that helped improve the manuscript. The POSYDON project is supported primarily by two sources: the Swiss National Science Foundation (PI Fragos, project numbers PP00P2_211006 and CRSII5_213497) and the Gordon and Betty Moore Foundation (PI Kalogera, grant award GBMF8477). DM acknowledges support from the Swiss National Science Foundation (project number PP00P2_176868) and the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (grant agreement No. 101002352). KK acknowledges support from the Federal Commission for Scholarships for Foreign Students for the Swiss Government Excellence Scholarship (ESKAS No. 2021.0277), and the Spanish State Research Agency, through the María de Maeztu Program for Centers and Units of Excellence in R&D, No. CEX2020-001058-M. S.S.B., T.F., and Z.X. were supported by the project number PP00P2_211006. S.S.B. was also supported by the project number CRSII5_213497. Z.X. acknowledges support from the Chinese Scholarship Council (CSC). A.D., K.A.R., P.M.S. and M.S. were supported by the project number GBMF8477. EZ acknowledges funding support from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (Grant agreement No. 772086) as well as from the Hellenic Foundation for Research and Innovation (H.F.R.I.) under the “3rd Call for H.F.R.I. Research Projects to support Post-Doctoral Researchers” (Project No: 7933). The computations were performed at Northwestern University on the Trident computer cluster (funded by the GBMF8477 award) and at the University of Geneva on the Yggdrasil computer cluster.

References

  1. Abbott, B. P., Abbott, R., Abbott, T. D., et al. 2020, ApJ, 892, L3 [NASA ADS] [CrossRef] [Google Scholar]
  2. Abdusalam, K., Ablimit, I., Hashim, P., et al. 2020, ApJ, 902, 125 [NASA ADS] [CrossRef] [Google Scholar]
  3. Anastasopoulou, K., Zezas, A., Ballo, L., & Della Ceca, R. 2016, MNRAS, 460, 3570 [NASA ADS] [CrossRef] [Google Scholar]
  4. Angelini, L., Loewenstein, M., & Mushotzky, R. F. 2001, ApJ, 557, L35 [NASA ADS] [CrossRef] [Google Scholar]
  5. Asplund, M., Grevesse, N., Sauval, A. J., & Scott, P. 2009, ARA&A, 47, 481 [NASA ADS] [CrossRef] [Google Scholar]
  6. Bachetti, M., Harrison, F. A., Walton, D. J., et al. 2014, Nature, 514, 202 [NASA ADS] [CrossRef] [Google Scholar]
  7. Basko, M. M., & Sunyaev, R. A. 1976, MNRAS, 175, 395 [Google Scholar]
  8. Basu-Zych, A. R., Lehmer, B., Fragos, T., et al. 2016, ApJ, 818, 140 [NASA ADS] [CrossRef] [Google Scholar]
  9. Bavera, S. S., Fragos, T., Qin, Y., et al. 2020, A&A, 635, A97 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  10. Bernadich, M. C., Schwope, A. D., Kovlakas, K., Zezas, A., & Traulsen, I. 2022, A&A, 659, A188 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  11. Biswas, B. 2021, ApJ, 921, 63 [NASA ADS] [CrossRef] [Google Scholar]
  12. Bondi, H., & Hoyle, F. 1944, MNRAS, 104, 273 [Google Scholar]
  13. Brightman, M., Harrison, F. A., Fürst, F., et al. 2018, Nat. Astron., 2, 312 [CrossRef] [Google Scholar]
  14. Carpano, S., Haberl, F., Maitra, C., & Vasilopoulos, G. 2018, MNRAS, 476, L45 [NASA ADS] [CrossRef] [Google Scholar]
  15. Chandra, A. D., Roy, J., Agrawal, P. C., & Choudhury, M. 2020, MNRAS, 495, 2664 [NASA ADS] [CrossRef] [Google Scholar]
  16. Chashkina, A., Abolmasov, P., & Poutanen, J. 2017, MNRAS, 470, 2799 [Google Scholar]
  17. Chashkina, A., Lipunova, G., Abolmasov, P., & Poutanen, J. 2019, A&A, 626, A18 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  18. Chomiuk, L., & Povich, M. S. 2011, AJ, 142, 197 [Google Scholar]
  19. Colbert, E. J. M., & Mushotzky, R. F. 1999, ApJ, 519, 89 [NASA ADS] [CrossRef] [Google Scholar]
  20. Colbert, E. J. M., & Ptak, A. F. 2002, ApJS, 143, 25 [NASA ADS] [CrossRef] [Google Scholar]
  21. Davis, S. W., Narayan, R., Zhu, Y., et al. 2011, ApJ, 734, 111 [NASA ADS] [CrossRef] [Google Scholar]
  22. Dib, R., Ransom, S. M., Ray, P. S., Kaspi, V. M., & Archibald, A. M. 2005, ApJ, 626, 333 [NASA ADS] [CrossRef] [Google Scholar]
  23. Doroshenko, V., Tsygankov, S., & Santangelo, A. 2018, A&A, 613, A19 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  24. Doroshenko, V., Zhang, S. N., Santangelo, A., et al. 2020, MNRAS, 491, 1857 [Google Scholar]
  25. Dray, L. M. 2006, MNRAS, 370, 2079 [NASA ADS] [CrossRef] [Google Scholar]
  26. Earnshaw, H. P., Roberts, T. P., Middleton, M. J., Walton, D. J., & Mateos, S. 2019, MNRAS, 483, 5554 [NASA ADS] [CrossRef] [Google Scholar]
  27. Eggleton, P. P. 1971, MNRAS, 151, 351 [CrossRef] [Google Scholar]
  28. Eggleton, P. P. 1972, MNRAS, 156, 361 [NASA ADS] [Google Scholar]
  29. Fabbiano, G. 1989, ARA&A, 27, 87 [NASA ADS] [CrossRef] [Google Scholar]
  30. Fabbiano, G., Kim, D. W., Fragos, T., et al. 2006, ApJ, 650, 879 [NASA ADS] [CrossRef] [Google Scholar]
  31. Fabrika, S. N., Atapin, K. E., Vinokurov, A. S., & Sholukhova, O. N. 2021, Astrophys. Bull., 76, 6 [NASA ADS] [CrossRef] [Google Scholar]
  32. Farrell, S. A., Servillat, M., Wiersema, K., et al. 2011, Astron. Nachr., 332, 392 [NASA ADS] [CrossRef] [Google Scholar]
  33. Feng, H., & Kaaret, P. 2008, ApJ, 675, 1067 [NASA ADS] [CrossRef] [Google Scholar]
  34. Fragos, T., Lehmer, B., Tremmel, M., et al. 2013a, ApJ, 764, 41 [NASA ADS] [CrossRef] [Google Scholar]
  35. Fragos, T., Lehmer, B. D., Naoz, S., Zezas, A., & Basu-Zych, A. 2013b, ApJ, 776, L31 [CrossRef] [Google Scholar]
  36. Fragos, T., Andrews, J. J., Bavera, S. S., et al. 2023, ApJS, 264, 45 [NASA ADS] [CrossRef] [Google Scholar]
  37. Frank, J., King, A., & Raine, D. 2002, Accretion Power in Astrophysics (Cambridge University Press) [Google Scholar]
  38. Fürst, F., Walton, D. J., Harrison, F. A., et al. 2016, ApJ, 831, L14 [Google Scholar]
  39. Gao, Y., Wang, Q. D., Appleton, P. N., & Lucas, R. A. 2003, ApJ, 596, L171 [NASA ADS] [CrossRef] [Google Scholar]
  40. Hirai, R., & Mandel, I. 2021, PASA, 38, e056 [NASA ADS] [CrossRef] [Google Scholar]
  41. Hou, X., Ge, M. Y., Ji, L., et al. 2022, ApJ, 938, 149 [NASA ADS] [CrossRef] [Google Scholar]
  42. Hurley, J. R., Pols, O. R., & Tout, C. A. 2000, MNRAS, 315, 543 [Google Scholar]
  43. Hurley, J. R., Tout, C. A., & Pols, O. R. 2002, MNRAS, 329, 897 [Google Scholar]
  44. Israel, G. L., Belfiore, A., Stella, L., et al. 2017a, Science, 355, 817 [NASA ADS] [CrossRef] [Google Scholar]
  45. Israel, G. L., Papitto, A., Esposito, P., et al. 2017b, MNRAS, 466, L48 [NASA ADS] [CrossRef] [Google Scholar]
  46. Jermyn, A. S., Bauer, E. B., Schwab, J., et al. 2023, ApJS, 265, 15 [NASA ADS] [CrossRef] [Google Scholar]
  47. Kaaret, P., Alonso-Herrero, A., Gallagher, J. S., et al. 2004, MNRAS, 348, L28 [NASA ADS] [CrossRef] [Google Scholar]
  48. Kaaret, P., Schmitt, J., & Gorski, M. 2011, ApJ, 741, 10 [NASA ADS] [CrossRef] [Google Scholar]
  49. Kaaret, P., Feng, H., & Roberts, T. P. 2017, ARA&A, 55, 303 [Google Scholar]
  50. Karino, S. 2018, MNRAS, 474, 4564 [NASA ADS] [CrossRef] [Google Scholar]
  51. Kim, D.-W., & Fabbiano, G. 2004, ApJ, 611, 846 [NASA ADS] [CrossRef] [Google Scholar]
  52. Kim, M., Kim, Y.-M., Sung, K. H., Lee, C.-H., & Kwak, K. 2021, A&A, 650, A139 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  53. King, A. R. 2009, MNRAS, 393, L41 [NASA ADS] [Google Scholar]
  54. King, A., & Lasota, J.-P. 2016, MNRAS, 458, L10 [Google Scholar]
  55. King, A., & Lasota, J.-P. 2019, MNRAS, 485, 3588 [Google Scholar]
  56. King, A., & Lasota, J.-P. 2020, MNRAS, 494, 3611 [NASA ADS] [CrossRef] [Google Scholar]
  57. King, A. R., Davies, M. B., Ward, M. J., Fabbiano, G., & Elvis, M. 2001, ApJ, 552, L109 [NASA ADS] [CrossRef] [Google Scholar]
  58. King, A., Lasota, J.-P., & Kluźniak, W. 2017, MNRAS, 468, L59 [Google Scholar]
  59. Koliopanos, F., Vasilopoulos, G., Godet, O., et al. 2017, A&A, 608, A47 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  60. Koliopanos, F., Vasilopoulos, G., Buchner, J., Maitra, C., & Haberl, F. 2019, A&A, 621, A118 [EDP Sciences] [Google Scholar]
  61. Kong, L.-D., Zhang, S., Zhang, S.-N., et al. 2022, ApJ, 933, L3 [NASA ADS] [CrossRef] [Google Scholar]
  62. Kovlakas, K., Zezas, A., Andrews, J. J., et al. 2020, MNRAS, 498, 4790 [Google Scholar]
  63. Kroupa, P. 2001, MNRAS, 322, 231 [NASA ADS] [CrossRef] [Google Scholar]
  64. Kruckow, M. U., Tauris, T. M., Langer, N., Kramer, M., & Izzard, R. G. 2018, MNRAS, 481, 1908 [CrossRef] [Google Scholar]
  65. Landry, P., Essick, R., & Chatziioannou, K. 2020, Phys. Rev. D, 101 [CrossRef] [Google Scholar]
  66. Lasota, J.-P., & King, A. 2023, MNRAS, 526, 2506 [NASA ADS] [CrossRef] [Google Scholar]
  67. Lasota, J. P., Alexander, T., Dubus, G., et al. 2011, ApJ, 735, 89 [NASA ADS] [CrossRef] [Google Scholar]
  68. Lehmer, B. D., Eufrasio, R. T., Tzanavaris, P., et al. 2019, ApJS, 243, 3 [NASA ADS] [CrossRef] [Google Scholar]
  69. Licquia, T. C., & Newman, J. A. 2015, ApJ, 806, 96 [NASA ADS] [CrossRef] [Google Scholar]
  70. Linden, T., Kalogera, V., Sepinsky, J. F., et al. 2010, ApJ, 725, 1984 [Google Scholar]
  71. Lipunova, G. V. 1999, Astron. Lett., 25, 508 [NASA ADS] [Google Scholar]
  72. Liu, J. 2011, ApJS, 192, 10 [NASA ADS] [CrossRef] [Google Scholar]
  73. Liu, J., Vasilopoulos, G., Ge, M., et al. 2022, MNRAS, 517, 3354 [CrossRef] [Google Scholar]
  74. Maccarone, T. J., Kundu, A., Zepf, S. E., & Rhode, K. L. 2007, Nature, 445, 183 [NASA ADS] [CrossRef] [Google Scholar]
  75. Majid, W. A., Lamb, R. C., & Macomb, D. J. 2004, ApJ, 609, 133 [NASA ADS] [CrossRef] [Google Scholar]
  76. Mann, C. R., Richer, H., Heyl, J., et al. 2019, ApJ, 875, 1 [Google Scholar]
  77. Mapelli, M., Colpi, M., & Zampieri, L. 2009, MNRAS, 395, L71 [NASA ADS] [CrossRef] [Google Scholar]
  78. Mapelli, M., Ripamonti, E., Zampieri, L., Colpi, M., & Bressan, A. 2010, MNRAS, 408, 234 [NASA ADS] [CrossRef] [Google Scholar]
  79. Messenger, C., & Patruno, A. 2015, ApJ, 806, 261 [NASA ADS] [CrossRef] [Google Scholar]
  80. Middleton, M. J., & King, A. 2017, MNRAS, 470, L69 [Google Scholar]
  81. Middleton, M. J., Brightman, M., Pintore, F., et al. 2019, MNRAS, 486, 2 [Google Scholar]
  82. Miller, J. M. 2006, in New Horizons in Astronomy: Frank N. Bash Symposium, eds. S. J. Kannappan, S. Redfield, J. E. Kessler-Silacci, M. Landriau, & N. Drory, ASP Conf. Ser., 352, 121 [NASA ADS] [Google Scholar]
  83. Miller, M. C., Lamb, F. K., Dittmann, A. J., et al. 2019, ApJ, 887, L24 [Google Scholar]
  84. Miller-Jones, J. C. A., Wrobel, J. M., Sivakoff, G. R., et al. 2012, ApJ, 755, L1 [Google Scholar]
  85. Misra, D., Fragos, T., Tauris, T. M., Zapartas, E., & Aguilera-Dena, D. R. 2020, A&A, 642, A174 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  86. Misra, D., Kovlakas, K., Fragos, T., et al. 2023, A&A, 672, A99 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  87. Most, E. R., Weih, L. R., Rezzolla, L., & Schaffner-Bielich, J. 2018, Phys. Rev. Lett., 120 [CrossRef] [Google Scholar]
  88. Mushtukov, A. A., Suleimanov, V. F., Tsygankov, S. S., & Poutanen, J. 2015, MNRAS, 454, 2539 [Google Scholar]
  89. Mushtukov, A. A., Suleimanov, V. F., Tsygankov, S. S., & Ingram, A. 2017, MNRAS, 467, 1202 [NASA ADS] [Google Scholar]
  90. Mushtukov, A. A., Ingram, A., Middleton, M., Nagirner, D. I., & van der Klis, M. 2019, MNRAS, 484, 687 [NASA ADS] [CrossRef] [Google Scholar]
  91. Patruno, A., Wette, K., & Messenger, C. 2018, ApJ, 859, 112 [Google Scholar]
  92. Patton, R. A., & Sukhbold, T. 2020, MNRAS, 499, 2803 [NASA ADS] [CrossRef] [Google Scholar]
  93. Paxton, B., Bildsten, L., Dotter, A., et al. 2011, ApJS, 192, 3 [Google Scholar]
  94. Paxton, B., Cantiello, M., Arras, P., et al. 2013, ApJS, 208, 4 [Google Scholar]
  95. Paxton, B., Marchant, P., Schwab, J., et al. 2015, ApJS, 220, 15 [Google Scholar]
  96. Paxton, B., Schwab, J., Bauer, E. B., et al. 2018, ApJS, 234, 34 [NASA ADS] [CrossRef] [Google Scholar]
  97. Paxton, B., Smolec, R., Schwab, J., et al. 2019, ApJS, 243, 10 [Google Scholar]
  98. Perera, B. B. P., Stappers, B. W., Lyne, A. G., et al. 2017, MNRAS, 468, 2114 [NASA ADS] [CrossRef] [Google Scholar]
  99. Pintore, F., Zampieri, L., Stella, L., et al. 2017, ApJ, 836, 113 [Google Scholar]
  100. Podsiadlowski, P., Rappaport, S., & Han, Z. 2003, MNRAS, 341, 385 [NASA ADS] [CrossRef] [Google Scholar]
  101. Poutanen, J., Lipunova, G., Fabrika, S., Butkevich, A. G., & Abolmasov, P. 2007, MNRAS, 377, 1187 [NASA ADS] [CrossRef] [Google Scholar]
  102. Prestwich, A. H., Tsantaki, M., Zezas, A., et al. 2013, ApJ, 769, 92 [NASA ADS] [CrossRef] [Google Scholar]
  103. Quast, M., Langer, N., & Tauris, T. M. 2019, A&A, 628, A19 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  104. Quintin, E., Webb, N. A., Gúrpide, A., Bachetti, M., & Fürst, F. 2021, MNRAS, 503, 5485 [NASA ADS] [CrossRef] [Google Scholar]
  105. Raaijmakers, G., Greif, S. K., Hebeler, K., et al. 2021, ApJ, 918, L29 [NASA ADS] [CrossRef] [Google Scholar]
  106. Riley, T. E., Watts, A. L., Bogdanov, S., et al. 2019, ApJ, 887, L21 [Google Scholar]
  107. Roberts, T. P., Warwick, R. S., Ward, M. J., & Murray, S. S. 2002, MNRAS, 337, 677 [NASA ADS] [CrossRef] [Google Scholar]
  108. Rodríguez Castillo, G. A., Israel, G. L., Belfiore, A., et al. 2020, ApJ, 895, 60 [Google Scholar]
  109. Salvaggio, C., Wolter, A., Belfiore, A., & Colpi, M. 2023, MNRAS, 522, 1377 [NASA ADS] [CrossRef] [Google Scholar]
  110. Sana, H., de Mink, S. E., de Koter, A., et al. 2012, Science, 337, 444 [Google Scholar]
  111. Sana, H., de Koter, A., de Mink, S. E., et al. 2013, A&A, 550, A107 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  112. Sathyaprakash, R., Roberts, T. P., Walton, D. J., et al. 2019, MNRAS, L104 [Google Scholar]
  113. Servillat, M., Farrell, S. A., Lin, D., et al. 2011, ApJ, 743, 6 [Google Scholar]
  114. Seward, F. D., & Charles, P. A. 1995, Exploring the X-Ray Universe (Cambridge, UK: Cambridge University Press) [Google Scholar]
  115. Shakura, N. I., & Sunyaev, R. A. 1973, A&A, 24, 337 [NASA ADS] [Google Scholar]
  116. Shao, Y., & Li, X.-D. 2014, ApJ, 796, 37 [Google Scholar]
  117. Shao, Y., & Li, X.-D. 2015, ApJ, 802, 131 [Google Scholar]
  118. Skinner, G. K., Bedford, D. K., Elsner, R. F., et al. 1982, Nature, 297, 568 [NASA ADS] [CrossRef] [Google Scholar]
  119. Smith, K. L., Magno, M., & Tripathi, A. 2023, ApJ, 956, 3 [NASA ADS] [CrossRef] [Google Scholar]
  120. Soria, R. 2007, Ap&SS, 311, 213 [NASA ADS] [CrossRef] [Google Scholar]
  121. Soria, R., Cropper, M., Pakull, M., Mushotzky, R., & Wu, K. 2005, MNRAS, 356, 12 [NASA ADS] [CrossRef] [Google Scholar]
  122. Soria, R., Kolehmainen, M., Graham, A. W., et al. 2022, MNRAS, 512, 3284 [NASA ADS] [CrossRef] [Google Scholar]
  123. Swartz, D. A., Ghosh, K. K., Tennant, A. F., & Wu, K. 2004, ApJS, 154, 519 [Google Scholar]
  124. Tao, L., Feng, H., Zhang, S., et al. 2019, ApJ, 873, 19 [Google Scholar]
  125. Tauris, T. M., Kramer, M., Freire, P. C. C., et al. 2017, ApJ, 846, 170 [Google Scholar]
  126. Tremou, E., Strader, J., Chomiuk, L., et al. 2018, ApJ, 862, 16 [Google Scholar]
  127. Tsygankov, S. S., Doroshenko, V., Lutovinov, A. A., Mushtukov, A. A., & Poutanen, J. 2017, A&A, 605, A39 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  128. Vasilopoulos, G., Ray, P. S., Gendreau, K. C., et al. 2020, MNRAS, 494, 5350 [NASA ADS] [CrossRef] [Google Scholar]
  129. Vaughan, B. A., van der Klis, M., Wood, K. S., et al. 1994, ApJ, 435, 362 [NASA ADS] [CrossRef] [Google Scholar]
  130. Walton, D. J., Roberts, T. P., Mateos, S., & Heard, V. 2011, MNRAS, 416, 1844 [Google Scholar]
  131. Walton, D. J., Fürst, F., Heida, M., et al. 2018, ApJ, 856, 128 [NASA ADS] [CrossRef] [Google Scholar]
  132. Wiktorowicz, G., Sobolewska, M., Sądowski, A., & Belczynski, K. 2015, ApJ, 810, 20 [NASA ADS] [CrossRef] [Google Scholar]
  133. Wiktorowicz, G., Sobolewska, M., Lasota, J.-P., & Belczynski, K. 2017, ApJ, 846, 17 [NASA ADS] [CrossRef] [Google Scholar]
  134. Wiktorowicz, G., Lasota, J.-P., Middleton, M., & Belczynski, K. 2019, ApJ, 875, 53 [NASA ADS] [CrossRef] [Google Scholar]
  135. Wilson-Hodge, C. A., Malacaria, C., Jenke, P. A., et al. 2018, ApJ, 863, 9 [Google Scholar]
  136. Wolter, A., & Trinchieri, G. 2004, A&A, 426, 787 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  137. Wolter, A., Fruscione, A., & Mapelli, M. 2018, ApJ, 863, 43 [NASA ADS] [CrossRef] [Google Scholar]
  138. Yan, Z., Zhang, W., Soria, R., Altamirano, D., & Yu, W. 2015, ApJ, 811, 23 [NASA ADS] [CrossRef] [Google Scholar]
  139. Zampieri, L., & Roberts, T. P. 2009, MNRAS, 400, 677 [NASA ADS] [CrossRef] [Google Scholar]
  140. Zezas, A., Fabbiano, G., Rots, A. H., & Murray, S. S. 2002, ApJ, 577, 710 [NASA ADS] [CrossRef] [Google Scholar]
  141. Zezas, A., Fabbiano, G., Baldi, A., et al. 2007, ApJ, 661, 135 [CrossRef] [Google Scholar]
  142. Zhang, Y., Ge, M., Song, L., et al. 2019, ApJ, 879, 61 [NASA ADS] [CrossRef] [Google Scholar]
  143. Zocchi, A., Gieles, M., & Hénault-Brunet, V. 2019, MNRAS, 482, 4713 [Google Scholar]
  144. Zuo, Z.-Y., Li, X.-D., & Gu, Q.-S. 2014, MNRAS, 437, 1187 [NASA ADS] [CrossRef] [Google Scholar]

Appendix A: Synthetic XLFs for model B

We present the figures describing the populations for model B. In Figure A.1, we show the synthetic XLFs of burst populations at different burst ages divided into sub-populations of different types of XRBs (namely, HMXBs, IMXBs, and LMXBs; similar to Figure 2). Figure A.2 further splits the populations by type of accretion (RLO or wind), accretor (BH or NS), and donor (H-rich or He-rich).

thumbnail Fig. A.1.

Synthetic XLFs of burst populations for model B. The panels show the different types of XRBs, namely, HMXBs, IMXBs, and LMXBs, as described in the legend. The combined XLF is shown by the dotted-gray line.

thumbnail Fig. A.2.

Synthetic XLFs of burst populations for model B. The populations are further split into the type of the XRBs (RLO or wind), type of accretors (BH or NS), and type of donors (H-rich or He-rich). The combined XLF is shown by the dotted-gray line.

All Tables

Table 1.

Physical parameters corresponding to the best-fitting models A and B from Misra et al. (2023).

All Figures

thumbnail Fig. 1.

Distribution of X-ray luminosities of ULXs for all the burst ages (see legend) for the populations following models A and B (left and right panel, respectively). Geometrical beaming is accounted for in the X-ray luminosity, but the observability of a ULX is not. In the observed population, the brighter luminosities would be slightly suppressed, as they would be beamed.

In the text
thumbnail Fig. 2.

Synthetic XLFs of burst populations for model A. The panels show the different types of XRBs, namely, HMXBs, IMXBs, and LMXBs, as described in the legend. The combined XLF is shown by the dotted-gray line.

In the text
thumbnail Fig. 3.

Percentage of BH binaries present in the simulated ULXs (solid lines) and in the simulated XRBs (dashed lines) following the parameters of models A (shown by orange) and B (shown by blue).

In the text
thumbnail Fig. 4.

Synthetic XLFs of burst populations for model A. The populations are further split into the type of the XRBs (RLO or wind), type of accretors (BH or NS) and type of donors (H-rich or He-rich). The combined XLF is shown by the dotted-gray line.

In the text
thumbnail Fig. 5.

Distribution of accretor masses for ULXs in the synthetic populations across ages 5–300 Myr.

In the text
thumbnail Fig. 6.

Distribution of donor masses for ULXs in the synthetic populations across ages 5–300 Myr.

In the text
thumbnail Fig. 7.

Distribution of orbital periods for ULXs in the synthetic populations across ages 5–300 Myr.

In the text
thumbnail Fig. 8.

Average beaming factor ( b ¯ ) $ (\bar{b}) $ for the ULX populations spanning all burst ages, distinguishing between the type of accretor (NS, shown by circles, or BH, shown by squares). Populations that did not have any ULXs are not shown. The figure also shows the fraction of beamed (with b < 1.0) to total ULXs in each population, shown as dashed lines. The blue colored symbols and lines correspond to the results for model A, and the orange-colored symbols and lines are for model B.

In the text
thumbnail Fig. 9.

The number of NS ULXs for models A and B (left and right panels, respectively) per stellar mass in the burst populations. With green filled circles, we show the number of NS ULXs with observable pulses using the NS radius (RNS) as the accretion radius. With pink circles we show the number of NS ULXs with observable pulses using the magnetospheric radius (RM), and the mustard circles show the total number of NS ULXs.

In the text
thumbnail Fig. 10.

Percentage of HLXs (having X-ray luminosities ≳ 1041 erg s−1) with maximally beamed emission (with b = 3.2 × 10−3) for all the burst ages, distinguishing between the type of accretor (NS, shown by circles, or BH, shown by squares). The blue symbols show the results for model A, and the orange symbols are for model B.

In the text
thumbnail Fig. A.1.

Synthetic XLFs of burst populations for model B. The panels show the different types of XRBs, namely, HMXBs, IMXBs, and LMXBs, as described in the legend. The combined XLF is shown by the dotted-gray line.

In the text
thumbnail Fig. A.2.

Synthetic XLFs of burst populations for model B. The populations are further split into the type of the XRBs (RLO or wind), type of accretors (BH or NS), and type of donors (H-rich or He-rich). The combined XLF is shown by the dotted-gray line.

In the text

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

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

Initial download of the metrics may take a while.