High-Mass X-Ray Binaries with Be Donors as Ultraluminous X-Ray Sources

Since the detection of X-ray pulses from ultraluminous X-ray sources (ULXs) in 2014, neutron stars have been considered as their central objects. However, it remains unclear how neutron stars can be brighter than the Eddington luminosity, and no unified view exists on the magnetic field of neutron stars and the degree of beaming. Recent observations suggest that some X-ray pulsating ULXs have Be-type donors, and some of them occupy the same region as Be-type high-mass X-ray binaries (Be-HMXBs) on the Corbet diagram, which reveals the relation between spin and orbital periods. This suggests that at least some ULXs are special cases of Be-HMXBs. In this study, we use the framework of mass-accretion models for Be-HMXBs to investigate the conditions under which neutron stars achieve mass-accretion rates beyond the Eddington limit and become observable as ULXs. We show that a Be-HMXB may become a ULX if the magnetic field of the neutron star and the density of the Be disc meet certain conditions. We also show that, although a stronger magnetic field increases the brightness of a neutron star ULX with a Be donor, its brightness cannot exceed the Eddington limit by a more than a factor of ${\approx} 50$. Finally, we propose a scenario whereby some normal Be-HMXBs may evolve into ULXs as the donor evolves into a giant.


INTRODUCTION
Ultraluminous X-ray sources (ULXs) are X-ray sources located outside the galactic centre with luminosities above the Eddington luminosity of stellar-mass objects. A current controversy involving ULX sources debates whether they are stellar-mass black holes (BHs) with mass-accretion above the Eddington limit or intermediate-mass BHs with sub-Eddington mass-accretion rates (Kaaret et al. 2017). However, the detection of X-ray pulses from ULXs in M81 revealed that some ULX sources are in fact rotating neutron stars (NSs) (Bachetti et al. 2014). Since their first detection, several ULXs have been found to emit X-ray pulses, and it is now widely recognised that NSs constitute a certain fraction of ULXs (Carpano et al. 2018;Doroshenko et al. 2018;Fürst et al. 2016;Israel et al. 2017b,a;Tsygankov et al. 2017;Rodríguez Castillo et al. 2020). Despite the increasing number of NS ULXs being observed, it remains unclear how NSs emit X-rays far beyond the Eddington limit. In particular, whether the magnetic field of ULX NSs is strong or weak remains hotly debated (Eksi et al. 2015;Tong 2015;King & Lasota 2019;Mushtukov et al. 2015Mushtukov et al. , 2021Middleton et al. 2019;Brightman et al. 2018;Walton et al. 2018). In recent years, modelw that comprehensively treat spin, magnetic field and beaming of NS ULXs have been considered (Erkut et al. 2020;Mushtukov et al. 2021;Abarca et al. 2021).
In addition, as ever more NS ULXs are observed, a major question has become the fraction of ULXs that are of NS origin versus ★ E-mail: karino@ip.kyusan-u.ac.jp BH origin (Fragos et al. 2015;Shao & Li 2015;Misra et al. 2020).
To estimate the population of NS ULXs, we must understand the properties of binary systems in NS ULXs. Based on current observations, NS ULXs appear to be a type of high-mass X-ray binary (HMXB) objects with a relatively massive companion (Table 1 summarises the observed properties). Some NS ULXs have supergiants as companions, whereas other systems appear to Be-HMXBs with Be-type companions (see Table 1 in King & Lasota (2019)). To study of HMXBs, Corbet diagrams are often used to show how the NS orbital period is related to the rotation period (Corbet 1984(Corbet , 1986. In such diagrams, some NS ULXs occupy the same region as Be-HMXBs (see Fig. 1). If the donor star is a supergiant and the system has a short orbital period, mass-transfer from the donor to the NS will occur via Roche lobe overflow (RLOF), which leads to a large masstransfer rate (Bildsten et al. 1997;Frank et al. 2002). In this case, one can easily imagine that the mass-transfer rate exceeds the Eddington limit, allowing it to emit X-rays above the Eddington luminosity. However, it is another matter to determine whether the material flowing into the NS can be accreted onto the NS without being lost by outflow (Chashkina et al. 2019;Kosec et al. 2018;Lipunova 1999). Conversely, the X-ray luminosity of Be-HMXBs is generally less than that of HMXBs with supergiant donors (Bildsten et al. 1997). Be-HMXBs produce periodic X-ray outbursts with a period related to the orbital period of the Be-HMXB and sometimes produce giant outbursts that are abnormally bright (Rivinius et al. 2013) (the X-ray luminosity of a typical normal burst is about 10 36 erg s −1 , whereas a giant outburst it is over tenfold brighter). An interesting question is therefore whether such binary systems can produce X-ray luminosities above the Eddington luminosity (Okazaki et al. 2013;Karino & Miller 2016); in other words, are some NS ULXs special cases of Be-HMXBs? Given that Be-HMXBs represent a significant fraction of HMXBs (Bildsten et al. 1997;Coe & Kirk 2015), a significant number of NS ULXs may exist since they follow part of the regular evolution of Be-HMXBs (even if no pulses are detected) (Inoue et al. 2020;Kuranov et al. 2020;Mushtukov et al. 2021). The NS/BH fraction of the ULX population is important not only for calculating the initial mass function of massive stars but also for estimating the population of gravitational wave sources such as NS-NS binary stars (Mondal et al. 2020).
With these motivations, we examine in this study whether Be-HMXBs can indeed be ULXs. For an NS binary system with a Be donor to be a ULX, it must satisfy the following criteria: (1) The mass-transfer rate from the Be donor to the NS must exceed the Eddington limit.
(2) The material supplied from the donor must accrete onto the NS at a rate exceeding the Eddington limit.
In this study, we divide the mass-accretion process from the Be donor to the NS according to criteria (1) and (2) above and investigate under what conditions Be-HMXBs can become ULXs.
In the next section, we discuss under which conditions a Be-HMXB can become a ULX. Based on the conditions obtained in Section 2, we consider in Section 3 what type of binary system can become a ULX, discuss its evolutionary path, and compare the results with observations. Finally, Section 4 is devoted to the conclusion.

BINARY MODELS AND MASS-TRANSFER PROCESSES IN THE SYSTEMS
This section starts by giving the assumptions used for describing Be-NS binary systems, and then discusses the models.

Binary Model
In this study, the NS mass NS and radius NS are fixed at 1.4 ⊙ and 10 km, respectively. In addition, we assume a Be-type donor star with mass d = 10 ⊙ , which is almost the median value of the typical mass range of Be-type stars (3 ⊙ -18 ⊙ ) (Townsend et al. 2004). We further assume that the Be star ejects its mass through a deccretion disc. The mass loss rate of a single 10 ⊙ B-type star in its main-sequence phase is estimated to be <10 −8 ⊙ yr −1 (Vink et al. 2001), and even a Be star in a binary system is unlikely to exceed this rate by many orders of magnitude. Therefore, the donor mass decreases only negligibly on the timescale of 10 million years, which is the typical main-sequence lifetime of a 10 ⊙ star. Conversely, the donor radius d varies rather strongly during the main-sequence phase. Figure 2 shows the approximate variation of the star radius over time, calculated as per Hurley et al. (2000). In this figure, the end of the main-sequence corresponds to the point where the stellar radius expands rapidly over a short time. These results show that the stellar radius generally increases from ≈3 ⊙ to 10 ⊙ during the main-sequence phase. The observed orbital period orb of Be-HMXB ranges from a few days to several hundred days, and the orbital eccentricity ranges from near-circular to highly distorted (Townsend et al. 2011;Klus et al. 2014;Coe & Kirk 2015). In this study, we consider the orbital period and orbital eccentricity as the main parameters of a binary system and examine various combinations thereof over a wide range of parameter space.
The NS magnetic field is typically about 10 12 G in most HMXB systems, as determined from their spin rates and cyclotron resonant scattering features (Christodoulou et al. 2017). In contrast, the NS magnetic field of ULXs remains debated (Eksi et al. 2015;Tong 2015;Mushtukov et al. 2021). The NS magnetic field initially decays on the timescales of Ohmic decay, Hall dissipation, and mass-accretion (Aguilera et al. 2008;Zhang & Kojima 2006). Unfortunately, the age of NSs in NS ULXs, the initial magnetic field strength, and the cumulative accreted mass all remain unknown. Recent studies argue that the magnetic field strength of NSs in ULXs are almost on a par with that of normal HMXBs: NS = 10 11 -10 12 G Middleton et al. 2019;Walton et al. 2018). Thus, we focus herein on these typical values and vary the magnetic field in several ways.

Mass-Transfer Rate
In Be-HMXBs, the compact object receives mass from the Be stellar disc. Okazaki et al. (2013) shows that, in the case of a large outburst (type II outburst), the disc material is accreted in a Bondi-Hoyle-Littleton (BHL) process. According to the BHL process, when an NS interacts with the disc material, it captures a part of the Be disc inside the accretion radius: (1) In this case, the mass-transfer rate from the Be disc to the NS is (Bondi & Hoyle 1944;Hoyle & Lyttleton 1939;Edgar 2004) where rel is the relative velocity between the orbital velocity of the NS and the rotational velocity of the Be disc material, which varies depending on the inclination of the Be disc axis and the NS orbital  ( Figure 2. Donor radius. Time evolution of radii of Be stars in typical mass range (6 ⊙ -12 ⊙ ). Stars in this mass range have a radius of 4-10 solar radii during the main-sequence, although the radius increases rapidly above 10 ⊙ at the end of the main-sequence.
axis. is the density of the Be disc and is modelled as where 0 is the density of the inner edge of the Be disc. The density of the Be disc and the rotational velocity of the disc material both vary with distance from the Be star. However, we do not take into account the gradient of the physical quantities of the Be disc but represent them all at the periastron point = (1 − ), where is the semi-major axis of the orbit and can be obtained from Kepler's law by specifying d , NS , and orb . In order to consider BHL accretion from a medium with a gradient of physical quantities, the integration of physical quantities within the accretion radius must be properly considered (Wang 1981;Karino et al. 2019). The accretion of interstellar matter to a single NS is given by the mass-accretion rate (2); however, in a binary system, the gravity of the donor must also be taken into account. As a result, the NS gravity exceeds donor gravity only inside the Roche radius: where is the mass ratio (Eggleton 1983). Therefore, the mass-transfer rate from the Be disc to the NS neighbourhood in Be-HMXBs becomes (Okazaki et al. 2013). However, if the Roche radius of the Be star is less than the donor radius at the periastron, instead of accretion from the disc, RLOF occurs from the donor itself. In this case, a massive mass-transfer may occur, in which case the orbit would rapidly approach a circular orbit and/or, in some cases, a common envelope may form (Eggleton 2006). In such a case, the binary system would be unstable and may change rapidly, so the system would not be a steady ULX.

Model of Accretion Disc around Neutron Star
The material transferred from the Be disc to the NS forms an accretion disc and accretes onto the NS (such accretion discs have been studied in detail by Okazaki et al. (2013)). If the accretion disc is a standard one (Shakura & Sunyaev 1973), the accretion proceeds on a viscous timescale, which is exceedingly long, so the material stored in the disc falls slowly into the NS. In this case, even if a large amount of mass is transferred from the Be disc, accretion onto the NS is gradual, and the resulting X-ray curve will not be as bright and sharp as the observed giant outbursts seen from Be-HMXBs. Therefore, the disc is likely to be a slim disc, which permits thermally stable, advectiondominated, supercritical accretion flow (Abramowicz et al. 1988). To achieve supercritical accretion above the Eddington limit, the magnetic radius mag of the NS must be less than the photon-trapping radius trap (Okazaki et al. 2013), which means that the mass-transfer rate exceeds a certain value. Assuming the conditions given by Okazaki et al. (2013), massaccretion above the Eddington limit can be achieved only when the mass-transfer rate is supercritical: T > sc = 1.7 × 10 19 7/9 4/9 30 .
Here we use the normalised magnetic moment of the NS, where is a parameter related to the accretion geometry ( = 1 for spherically symmetric accretion and = 0.5 for disc accretion). Although the geometry of the slim disc may be intermediate between spherical and disc accretion, we use = 0.5 in this work. Any system with a binary parameter set that does not satisfy the above condition for the mass-transfer rate [eq. (6)] is rejected as a candidate for a NS ULX.

Mass-Accretion Rate
Material transferred from the Be disc to the neighbourhood of the NS will shift further through the slim disc to the vicinity of the NS if the condition (6) is satisfied. Eventually, when disc material is trapped in the magnetic field of the NS, it falls through the magnetic field lines and forms an accretion column near the magnetic poles, where it accretes and emits extremely bright X-rays (Pringle & Rees 1972). However, if the mass-transfer rate exceeds the Eddington limit, radiation pressure becomes dominant somewhere within the accretion disc, causing a significant thickening of the disc and outflows due to radiation (Shakura & Sunyaev 1973;Lipunov 1982;Lipunova 1999;Poutanen et al. 2007). The radius sph at which radiation pressure becomes dominant can be written as (King et al. 2017), where Edd is the Eddington critical massaccretion rate. As matter accretes beyond this radius, its accretion rate decreases as once it becomes trapped in the NS magnetic field. Here, the magnetic radius mag is Eqs. (9) and (10) give acc , which is the mass-accretion rate of NS ULX candidates with Be-type donors. From this we obtain the luminosity which we can compare with the observed ULX luminosities. Since sph ∝ , the mass-accretion rate onto the NS, given in eq. (9), does not depend on and only depends on mag . After some calculations, the mass-accretion rate can be written as which implies that the luminosity of the NS ULX is determined only by the magnetic field strength of the NS (Lipunov 1982).

NS ULX MODEL
The previous section considered the mass-transfer process from the Be disc to the vicinity of the NS and calculated how much of the transferred material is accreted onto the NS. The result indicates that the mass-transfer rate must exceed a certain level for super-Eddington mass-accretion to proceed in the accretion disc. In addition, the RLOF must be absent from Be stars to produce a stable ULX. To summarise, to become a ULX, a NS that accretes mass from a Be-type donor must satisfy the following three conditions: The luminosity of a NS when these conditions are satisfied is expressed as Figure 3 shows the NS luminosity with the orbital period and eccentricity serving as variables and shows the areas in parameter space that satisfy the conditions (13)-(15). The colour contour shows the luminosity normalised by the Eddington luminosity. The mass and radius of the donor are 10 ⊙ and 6 ⊙ , respectively. These results reveal a slightly inclined system between the plane of the Be disc and the orbital plane ( = /8, see Section 4). Shown are the results given with nine combinations of NS magnetic field and base density of the Be disc. The former varies from 10 11 and 10 13 G and the latter varies from 10 −11 to 10 −10 g cm 3 .

Possible Neutron Star ULX Systems
In figure 3, the white area is the region where Be-HMXBs cannot be a ULX because it does not satisfy the conditions mentioned above. The lower-right region is mainly where the NS orbit does not pass through the high-density region of the Be disc, which means that the Be disc does not transfer enough matter in this region of parameter space. Conversely, the upper-left region is where the NS orbit is too close to the donor, causing RLOF (Eggleton 2006). In such a system, if the outer layer of the donor is convective, a common envelope forms and the system will not emit X-rays. Moreover, the mass ratio becomes large for Be-NS binary systems, in which case a common envelope inevitably forms (van den Heuvel et al. 2017). Although conditions for the formation of the common envelope remain under debate, the binary system causing RLOF may not be a stable system because the orbital parameter changes on a short time scale (Ivanova et al. 2013).
In addition, if the magnetic field is weak, the mass-accretion rate of the NS will be insufficient to transform it into a ULX because the outflow from the disc will dominate.
The base density of the Be disc has a significant effect on the formation of ULXs. Specifically, a Be disc with a density less than 10 −11 g cm 3 cannot produce a sufficiently high-mass-transfer rate to produce a ULX. However, Be stars with densities greater than 5 × 10 −11 g cm −3 are observed not infrequently (Rivinius et al. 2013;Touhami et al. 2011;Draper et al. 2014).
As discussed in the previous section, the luminosity of ULXs increases with NS magnetic field strength. For a strong magnetic field (≈10 13 G), the proposed model shows that the source can be up to 50 times brighter than the Eddington luminosity. However, the brightest systems are limited to narrow regions in parameter space corresponding to short orbital perimeters and small eccentricities. Even in the parameter region where ULX is possible, the ULX is only moderately bright in many case (i.e. less than 20 times the Eddington limit), which is consistent with the observational suggestion that NS ULXs with Be donors have relatively small luminosity (see Table 1). Very bright NS ULXs, exceeding 10 40 erg s −1 , are probably not Be-HMXBs, but rather RLOF systems with supergiant donors, which means that at least two subgroups of NS ULX possibly exist: brighter systems with supergiant donors, and dimmer systems with Be donors.

Warp of Be disc
The misalignment of the orbital plane of the binary system with the Be disc is suggested to be more favourable for the accretion onto the NS in Be-HMXB (Okazaki et al. 2013) because, when these planes  Figure 3. Region in orb -plane where a binary system with a Be donor can be a ULX. The colour contours show the maximum possible X-ray luminosity normarized by Eddington luminosity, and the white areas correspond to the areas where the Eddington limit is not reached. All Be donors in the binary are fixed at 10 ⊙ and 6 ⊙ , and an inclusion angle is also fixed as = /8. From left to right in the nine figures, the base density of the Be disc is 0 = 2 × 10 −11 , 5 × 10 −11 , and 1 × 10 −10 g cm −3 . The magnetic field of the NS is NS = 1 × 10 13 , 10 12 , and 10 11 G, starting from the top line. Given that the determination of 2 is rather worse, we have drawn three lines in the diagram to represent the conditions (21) and which correspond to 2 = 2.7, 1, and 0.5, from the top to bottom, respectively. The region to the lower left of these lines is where the warped Be disc forms. are aligned, the tidal truncation reduces the size of the Be disc. As a result, the NS is less likely to gain mass from the Be disc. If the plane of the Be disc and the orbital plane are misaligned, the outer edge of the Be disc is warped by tidal action. In addition, the interaction of an NS with this warped region is suggested to increase the mass-transfer rate and lead to a giant outburst in Be-HMXB (Okazaki et al. 2013;Reig et al. 2007;Negueruela & Okazaki 2001;Moritani et al. 2011).
For such mass-transfer from the warped disc to occur, the Be disc radius must be larger than the typical radius at which warping is induced by tidal action. The radius tr of the truncated disc is taken to be 0.4 to 0.5 times the periastron distance, (Paczynski 1977), so we assume tr = 0.5 (1 − ).
Conversely, the position at which the disc warps due to tidal action is determined by the ratio of the shear viscosity in the azimuthal direction of the Be disc to that in the vertical direction. Following the work of Martin et al. (2011), which gives the warping conditions for a Be disc in a binary system, the radius at which the disc warps can be written as warp = 4.91 × 10 11 1 − 2 3/4 1/2 2 ℎ orb 1d where ℎ is the ratio of the scale height of the disc thickness to the Be stellar radius. Here, we set ℎ = 0.04 (Rivinius et al. 2013), and 2 is the coefficient of the vertical shear viscosity and is given by where 1 is the dimensionless viscosity parameter (Shakura & Sunyaev 1973;Ogilvie 1999). Numerical results from Lee et al. (1991) suggest that 1 > 0.1 because an excessively small value of 1 does not provide sufficient mass to the Be disc. For this reason, 1 = 0.3 is adopted by Martin et al. (2011), for which eq. (19) gives 2 = 2.7. Conversely, Cheng et al. (2014), who attempt to explain the existence of giant outbursts in Be-HMXBs based on the warp condition of the Be disc, argue that 2 = 0.5-1 explains the observations. In this study, we consider 2 = 0.5, 1, and 2.7.
The warping condition of the Be disc in Be-HMXBs is written as Using the Keplerian law in Eq. (17) allows the condition (20) to be rewritten as orb < 6.91 × 10 2 In this equation, we set the NS mass to NS = 1.4 ⊙ . Figure 3 shows this condition for three different viscosities represented by three different curves that correspond to 2 = 0.5, 1, and 2.7, from top to bottom. For each viscosity, the lower part of the curve is where the condition (21) is satisfied (i.e. where the warped disc forms). These figures show that, when the vertical viscosity is small ( 2 = 0.5), this condition does not significantly affect the conditions for outbursts in the ULX class, and only systems with exceptionally large orbital periods and eccentricities are excluded. Conversely, for large 2 (=2.7), as in the case of Martin et al. (2011), this condition significantly limits the occurrence of giant outbursts, which means that systems with orbital periods of more than a few tens of days and orbital eccentricities of more than 0.3-0.5 are excluded. The nearby NS ULX systems Swift J0243.6+6124 and SMC X-3, which likely have Be donors, have orbital periods of 28 and 45 d, respectively. Thus, a large 2 is unlikely to work, and a viscosity of about 2 = 1, as adopted by Cheng et al. (2014), may be reasonable (see Section 4).

DISCUSSION
In this study, we examine several conditions for Be-HMXBs to occur by mass-accretion beyond the Eddington limit and find the parameter region in which these binary systems become ULXs. Figure 3 shows the general results. Roughly speaking, only Be-HMXBs with small orbital eccentricity and short orbital period can become ULXs. In the orb -plane, a system with a small orbital period suffers RLOF, and such systems are not known to be stable and may form a common envelope (Eggleton 2006). Conversely, a system with a long orbital period cannot pass through the high-density part of the Be disc if the eccentricity is small, so the mass-transfer rate does not become sufficiently large. In addition, a system with a long orbital period and a large eccentricity cannot satisfy the warp condition of the Be disc (note that the warp condition depends on the viscosity parameters, which are highly indefinite). However, in the region where both the orbital period and the eccentricity are large, the probability that a binary system becomes a ULX is relatively low because the allowed ULX region narrows.
When the orbital eccentricity is small, the period of the binary system depends on the base density of the circumstellar disc of the Be donor (Touhami et al. 2011;Draper et al. 2014;Rivinius et al. 2013). In other words, the larger the base density of the Be disc, the more material the NS can acquire and the more likely it is to become a super-Eddington source. Observationally, the density of the Be disc is about 10 −10 -10 −12 g cm −3 , and if the disc density is near the upper limit, the mass supply from the Be disc suffices for the NS to become a ULX, and even a binary system with a period of a few tens of days can become a ULX.
The size of the ULX parameter region also depends on the magnetic field of the NS. If the magnetic field is strong, the condition becomes more stringent for the slim disc to be captured by the magnetic field around the NS (6), and the ULX parameter region shrinks. Conversely, if the magnetic radius is large and the accretion disc is captured by the magnetic field of the NS before significant mass loss occurs due to radiation-pressure-induced outflow, both the massaccretion rate and the X-ray luminosity increase. Therefore, the X-ray luminosity increases when the NS has strong magnetic field. As a result, a system with sph ≈ mag is likely to be a ULX (King et al. 2017;King & Lasota 2019. However, there is an upper limit to the X-ray luminosity because the ULX parameter region narrows as the magnetic field becomes stronger. Under the conditions used in our models, the ULX parameter region almost disappears when the neutron star has a magnetar-level magnetic field ( NS ≈ 10 14 G). Therefore, NS ≈ 10 13 G, shown in Fig. 3, is almost the upper limit of the magnetic field allowed for a Be-HMXB to become a ULX, and its brightness is limited to about 50 times the Eddington luminosity. Table 1 summarises the donor information and the X-ray luminosities of observed NS ULXs. The X-ray luminosities of the ULXs with Be-type donors are relatively small, and none of them exceeds 25 times the Eddington limit. Thus, the proposed model is consistent with observations on this point.

Inclination Angle
A giant outburst is facilitated when the Be disc plane is inclined (i.e. nonparallel) with respect to the orbital plane because the disc is truncated in this case and becomes smaller (Okazaki et al. 2013). However, the inclination angle between these two planes affects the relative velocity rel of the NS relative to the matter in the Be disc (where rel is the orbital velocity of the NS at the periastron). Here we assume that the matter in the Be disc rotates around the Be star with Keplerian velocity Kep . Thus, depending on the inclination angle, the relative velocity between the NS and the Be disc matter can vary from orb − Kep to orb + Kep . We introduce the parameter to represent the inclination angle between the Be disc plane and the orbital plane, and express the relative velocity as 2 rel = 2 orb + 2 Kep − 2 orb Kep cos ( ) .
Thus, is the angle between the orbital velocity vector and the Keplerian velocity vector of the disc matter, so = 0 corresponds to the coplanar case and = /2 corresponds to the disc being perpendicular to the orbital plane. We now compare the results as a function of inclination angle . Figure 4 shows the results for NS = 10 12 G and 0 = 5 × 10 −11 g, cm −3 . Here, as in Fig. 3, all Be-type donors are fixed to 10 ⊙ and 6 ⊙ . From left to right, figure 4 shows the results for = 0, /8, and /4. As increases, the relative velocity increases, and the BHL accretion rate decreases, so the ULX parameter region narrows rapidly and disappears from our parameter set. Conversely, as the relative velocity decreases, the mass-transfer rate due to BHL accretion increases, so the potential ULX parameter region increases. Okazaki et al. (2013) argue that their simulations show that the brightest burst occurs in the coplanar case, despite the largest burst theoretically occurring when the two planes are moderately offset.

Conditions Concerning Mass-Transfer Rate
For a Be-HMXB to become a ULX, this study imposes two conditions concerning the mass-transfer rate: mag < sph and T > sc  Figure 4. The same diagram as in Fig. 3 but with different inclusion angles: all Be donors are fixed at 10 ⊙ and 6 ⊙ , the base density of the Be disc is 0 = 5 × 10 −11 g cm −3 , and the NS magnetic field is NS = 1 × 10 12 G. From left to right, we set = 0, /8, and /4. As increases, the relative velocity increases, and the BHL accretion rate decreases, so the ULX parameter region narrows rapidly and disappears from the parameter set.  Figure 5. The ULX parameter region when mag < sph and T > sc are relaxed. The parameters are the same as those in the centre of Fig. 3. When the condition is relaxed, the ULX parameter region expands, but the expanded area corresponds to very dim ULX systems. (King et al. 2017;Okazaki et al. 2013). The former condition is stronger than the latter, so when the former is imposed, the latter can be relaxed, and even if we set T > Edd , this relaxation will have a negligible effect on the results. For most of NS ULXs where pulsations have been observed, mag < sph (King et al. 2017;. Therefore, given the observational constraints on spin-up rates and the existence of outflows, we consider it reasonable to impose the condition mag < sph here (Kosec et al. 2018;Vasilopoulos et al. 2021).
The results remains essentially unchanged when the magnetic field is moderate (10 12 G) and T > sc , even if mag > sph . However, if the magnetic field grows above 10 12 G, the ULX parameter region extends significantly. Figure 5 shows the results for the standard parameters (i.e. with these two conditions relaxed). Comparing this with the corresponding results in in the centre of figure 3, we see that the ULX parameter region is extended. However, the majority of the extended parameter region corresponds to a luminosity less than twice the Eddington luminosity of 1.4 ⊙ objects and is therefore only very faint as a ULX. In fact, the number of extra-galactic ULXs exceeds 10 Edd and, even were the number of fainter objects (<10 Edd ) to increase, they would not be observed. Therefore, for comparison with observations, the results are independent of the conditions relaxed above.

Evolution of Binary Systems
The parameters of a binary system evolve over time. In particular, due to its relatively high-mass, the radius of a Be donor changes rapidly (<1 Myr). The ULX parameter region shown in Fig. 3 depends on the radius of the Be donor. Figure 2 shows the evolution of the ULX parameter region as the radius of the Be donor varies from d = 6 to 10 ⊙ . As the Be donor evolves and its radius increases, the potential ULX region clearly broadens and shifts toward the longorb side. Toward the end of the main-sequence, when the radius reaches about 10 ⊙ , even a system with an orbital period of about 100 days can become a ULX. Even a 10 ⊙ donor loses only about 2% of its mass during the 10 Myrs of the main-sequence phase, and any change in mass can be neglected during the HMXB phase (Vink et al. 2001;Hurley et al. 2000). Although the maximum accretion rate onto the NS is large, the intermittent nature of accretion in a Be-HMXB means that the cumulative mass-accretion is not large, so changes in NS mass and orbit can be ignored.
The fact that the ULX parameter region expands as the donor evolves means that even a normal Be-type HMXB with a luminosity below the Eddington luminosity can become ULX as the Be donor evolves. As seen from the Corbet diagram in figure 1, many Be-type HMXBs exist with orbital periods of around 100 d (Klus et al. 2014;Townsend et al. 2011). Therefore, a certain fraction of Be-HMXBs may become ULXs during their evolution. However, numerous unknown features remain, such as whether late-main-sequence Be-type stars retain their discs and how the NS magnetic field evolves.
Additionally, note that HMXBs containing NSs could be progenitors of NS-NS binary systems. Passing through HMXBs with moderate orbital periods is especially important for the formation of NS-NS binaries that merge within the age of the Universe (Belczynski et al. 2002;Tauris et al. 2015). In this case, Be-HMXBs with orbital periods ranging from a few days to several hundred days could play an important role as candidate progenitors. If Be-HMXBs become ULXs that can be observed even outside the galaxy, their population will provide valuable information on the abundance and merger rate of NS-NS binary systems.

Comparison with Observed Systems
The results obtained are compared with the observed number of NS ULXs. Although the magnetic field of the NS and the base density and viscosity of the Be disc may differ in each binary system, we assume the following parameters for the binary systems: 0 = 1 × 10 −10 g cm −3 , d = 10 ⊙ , and 2 = 1.0. The parameter region in which a Be-HMXB can exceed the Eddington luminosity and become  Figure 6. Same diagram as Fig. 3, but with different donor radii: the Be-type donor is 10 ⊙ , the base density of the Be disc is 0 = 5 × 10 −11 g cm −3 , the NS magnetic field is NS = 1 × 10 12 G, and = /8. The donor radii are 6 ⊙ , 8 ⊙ , and 10 ⊙ , from left to right. As the donor evolves and the radius increases, the ULX parameter region widens, and even systems with large orbital periods could become ULXs. a ULX is shown by the area enclosed by the curves in figure 7. The two solid lines indicate the edges of the ULX region in the orb plane (as in Fig. 3) and come from the RLOF condition and the warped disc condition. The three dashed lines correspond to the different X-ray luminosities: 15 Edd , 10 Edd , and 5 Edd , from the left to the right. Figure 7 also plots the observed NS ULX systems with known orbital periods and which are suggested to have Be donors (Tsygankov et al. 2017;Doroshenko et al. 2018;King & Lasota 2019).
The luminosity of the NS ULX depends on NS and, for HMXBs, we assume a typical magnetic field strength of NS = 10 12 G. Although NS ULXs were initially thought to have a strong magnetic fields (Bachetti et al. 2014;Eksi et al. 2015;Mushtukov et al. 2015), recent evidence suggests that this is not necessarily the case (King & Lasota 2019;Middleton et al. 2019;Mushtukov et al. 2019). Observational evidence for SMC X-3, shown here as an example, suggests that the magnetic field strength is NS = (1 to 5) × 10 12 G (Tsygankov et al. 2017).
Among the compared systems, SMC X-3 and Swift J0243.6+6124 are in the ULX region. In the case of these two systems, the BHL process from the warped Be disc transports the material to the vicinity of the NS at a mass-transfer rate above the Eddington rate, and most of this material accretes into the NS magnetosphere. However, in the case of NGC2404 ULX, when 2 = 1, the unknown eccentricity must be less than about 0.4 to satisfy the warp condition. NGC2404 ULX is also likely to have an orbital period close to the lower limit. Conversely, if the NS magnetic field is small, the ULX parameter region is extended, making it possible for NGX2404 ULX to enter the ULX parameter region even if the orbital period is not so small.
We now consider the X-ray luminosity of these systems when they produce a giant outburst. Assuming that the magnetic field NS = 10 12 G, the luminosity of Be-HMXB as a ULX is expected to be at most 10-20 times the Eddington luminosity (cf. Fig. 3 and Fig. 7). Even if the magnetic field were a little stronger, by analogy with Fig. 3, the X-ray luminosity would be unlikely to exceed 20 times the Eddington luminosity since these sources are located near the right edge of the ULX parameter region. However, the luminosity here does not consider how beaming affects the observed luminosity. The observed luminosity should thus be the intrinsic luminosity divided by the beaming factor < 1 (King et al. 2017). Some recent studies have argued that beaming is not very strong, so we assume that the effect of beaming is a factor of order unity (Takahashi & Ohsuga 2017;Erkut et al. 2020;Mushtukov et al. 2021) 1 . In particular, analyses based on various models suggest that beaming is weak in Swift 1 Recent numerical simulation indicates that beaming could be stronger when mag is much smaller than sph (Abarca et al. 2021).
J2043.6+6124 (Erkut et al. 2020). The X-ray luminosity of observed Be donor systems is 10-20 times the Eddington luminosity (see Table 1), which is consistent with this prediction when we consider weak beaming. Figure 3 also shows that the X-ray luminosity of NS ULXs with Be donors does not exceed ≈50 Edd without assuming unreasonable conditions. This is also consistent with the fact that the X-ray luminosity of ULXs with Be donors is less than that of ULXs with supergiant donors, as shown in Table 1.
Although the number of observed systems is not large, Table 1 shows that most bright NS ULXs have orbital periods of a few days. Conversely, systems with Be donors have X-ray luminosities less than 20 Edd and orbital periods of more than 10 days. Note that, for another reason, it is unlikely that Be donors exist among the luminous systems with orbital periods of less than 10 days,. That is, in Be-HMXBs, the NS captures the disc matter as it passes through the Be disc. Therefore, the donor star must supply the captured mass to the disc. If the orbital period of the system is too short, it may be impossible to supply and store sufficient mass in the circumstellar disc. Although an uncertainty exists, the mass supply rate from the Be donor to the circumstellar disc is (Rivinius et al. 2013) disc ≈ 10 18 g s −1 . Considering that most of the mass is ejected from the outer edge of the disc into interstellar space, it is not surprising that it takes more than 10 days to store sufficient mass in the disc to cause a giant outburst. A future challenge is to quantify how much mass is transferred from the donor to the disc and from the disc to the NS and interstellar space.

CONCLUSION
Recently, several ULXs have been found to have NSs as emitters. Although the number of observations remains small, some of these ULXs are known to have Be-type donors. In this paper, we consider the following conditions for a Be-HMXB to become a ULX: i. The donor does not suffer from RLOF near the periastron of the binary system.
ii. The mass capture rate from the Be disc undergoing the BHL process exceeds the Eddington limit.
iii. Mass-accretion above the Eddington rate proceeds to the NS magnetosphere before massive outflow from the accretion disc.
In addition to these criteria, the condition of the Be disc, although highly indefinite, is also considered: iv. The Be disc warps inside the truncation radius. We consider that a Be-HMXB satisfying all these conditions should become a NS ULX, and find the requisite parameter region.
The results indicate that a NS should receive mass-accretion above the Eddington rate for a band-like parameter region extending from the lower left to the upper right in the orb -plane. However, the warp conditions of the Be disc do not allow for excessively large eccentricity and orbital period. The size of the ULX parameter region increases with the base density of the Be disc. Although the maximum luminosity of the ULX increases as the magnetic field of the NS increases, the ULX parameter region itself decreases as the magnetic field increases, placing an upper limit on the luminosity of NS ULXs with Be-type donors. Therefore, NS ULXs much brighter than the Eddington luminosity likely form a separate subclass that is subject to RLOF from supergiant donors.
The ULX parameter region also depends on the radius of the Be donor: the larger the radius, the wider the ULX parameter region. Therefore, as the Be donor evolves and expands, even Be-type HMXBs, whose luminosity is currently below the Eddington luminosity, may evolve into ULXs. Since Be-HMXBs are quite common, it would be useful to know how many of them become ULXs and for how long during their evolution, so that we can estimate the fraction of ULXs in the observed HMXB population. If the ULX population can be deduced from the HMXB population, it will also be useful to estimate how many non-pulsating ULXs are of NS origin. Finally, a general evolutionary model should be constructed for Be-type HMXBs and NS ULXs with Be donors in combination with binary-evolution calculations.

DATA AVAILABILITY
The data underlying this article will be shared on reasonable request to the corresponding author.