# X-ray emission from star cluster winds in starburst galaxies

Annika Franeck  Astronomical Institute of the Czech Academy of Sciences, Boční
II 1401, CZ-14100 Prague, Czech Republic  Richard Wünsch Astronomical
Institute of the Czech Academy of Sciences, Boční II 1401, CZ-14100 Prague,
Czech Republic  Sergio Martínez-González CONACYT-Instituto Nacional de
Astrofísica, Óptica y Electrónica, AP 51, 72000 Puebla, Mexico  Ivana Orlitová
Astronomical Institute of the Czech Academy of Sciences, Boční II 1401,
CZ-14100 Prague, Czech Republic  Peter Boorman Astronomical Institute of the
Czech Academy of Sciences, Boční II 1401, CZ-14100 Prague, Czech Republic
Department of Physics & Astronomy, Faculty of Physical Sciences and
Engineering, University of Southampton, Southampton, SO17 1BJ, UK  Jiří
Svoboda Astronomical Institute of the Czech Academy of Sciences, Boční II
1401, CZ-14100 Prague, Czech Republic  Dorottya Szécsi I. Physikalisches
Institut, Universität zu Köln, Zülpicher Strasse 77, D-50937 Köln, Germany
Institute of Astronomy, Faculty of Physics, Astronomy and Informatics,
Nicolaus Copernicus University, Grudziądzka 5, 87-100 Toruń, Poland  Vanesa
Douna Instituto de Astronomía y Física del Espacio, CONICET-UBA, Ciudad
Autónoma de Buenos Aires, Argentina

###### Abstract

Inspired by the excess soft X-ray emission recently detected in Green Pea
galaxies, we model the soft X-ray emission (0.5 – 2.0 keV) of hot gas from
star cluster winds. By combining individual star clusters, we estimate the
soft X-ray emission expected from the typically unresolved diffuse hot gas in
starburst galaxies, devoid of competing emission from e.g., AGN or other
unresolved point sources. We use stellar models of sub-solar metallicities
($0.02$ $Z_{\odot}$ and $0.4$ $Z_{\odot}$), and take into account supernova
explosions for massive stars. For lower metallicities, we find that stellar
winds do not contribute significantly ($\lesssim 3$% of the mechanical energy)
to the observed soft X-ray emission of normal star forming galaxies. For
higher metallicities and possibly also for larger proportions of massive star
clusters in the simulated starburst galaxies, we reproduce well the observed
correlation between star formation rate and X-ray luminosity previously
reported in the literature. However, we find that no combination of model
assumptions is capable of reproducing the substantial soft X-ray emission
observed from Green Pea galaxies, indicating that other emission mechanisms
(i.e. unusually large quantities of High-/Low-Mass X-ray Binaries, Ultra-
Luminous X-ray sources, a modified initial mass function, Intermediate-Mass
Black Holes, or AGN) are more likely to be responsible for the X-ray excess.

galaxies, Green Pea galaxies — star clusters — hot gas — stellar winds — X-ray
emission

††journal: ApJ

##  1 Introduction

The hot gas within star forming regions is often invoked to contribute to the
X-ray energy budget in star-forming galaxies, yet only a limited number of
detailed studies of its contribution are available (see e.g. Silich et al.,
2005; Lopez et al., 2014; Rosen et al., 2014). More widely-studied sources of
X-ray emission in galaxies are believed to arise from X-ray coronae found in
accreting supermassive black holes (aka Active Galactic Nuclei; AGN) and X-Ray
Binaries. Here, X-ray Binaries refer to any compact object (typically a
stellar mass black hole or neutron star) accreting material from a High-Mass
or Low-Mass companion (aka High-Mass X-ray Binaries & Low-Mass X-ray Binaries,
respectively).

The relative proportions of each X-ray-emitting process to the total observed
X-ray luminosity depends on specific conditions within the galaxy. In the
nearby universe, the X-ray emission from individual starburst galaxies was
studied by e.g., Grimes et al. (2007); Otí-Floranes et al. (2012, 2014); Basu-
Zych et al. (2013, 2016), who attributed all the unresolved X-ray emission to
the hot gas. However, corrections were not applied for unresolved compact
sources. Mineo et al. (2012a, b) systematically studied the X-ray spectra of a
sample of nearby star-forming galaxies. To evaluate the X-ray contribution
from the hot gas, the resolved High-Mass X-ray Binary contributions to the
X-rays were removed, alongside an estimated unresolved High-Mass X-ray Binary
component. Contributions from other unresolved X-ray sources (including Low-
Mass X-ray Binaries, cataclysmic variables, AB stars, young stellar objects &
supernova remnants) were estimated and found to be relatively weak in the soft
0.5–2 keV band. Mineo et al. (2012b) then found the hot gas to contribute at a
level of $<30\%$ to the total observed X-ray emission, and the diffuse
emission in the 0.5 – 2 keV band to correlate with the star formation rate
(SFR) of the galaxy.

As with the X-ray emission from hot gas, the emission from High-Mass X-ray
Binaries also scales with SFR due to the rapid stellar evolution of higher-
mass stars. High-Mass X-ray Binaries are hence typically short-lived and thus
the contribution to the X-ray budget is broadly proportional to SFR with a
modification due to metallicity. In contrast, the X-ray contribution from
longer-lived Low-Mass X-ray Binaries scales with the mass of the galaxy.
Empirical laws can then be derived which connect SFR to the observed X-ray
luminosity of a given galaxy (Ranalli et al., 2003; Mineo et al., 2012a; Basu-
Zych et al., 2013; Douna et al., 2015; Brorby et al., 2016).

Unusually bright X-ray emission was recently reported from a class of
intermediate redshift ($z\sim 0.2-0.3$) compact starburst galaxies called
“Green Peas” (Svoboda et al., 2019). Green Pea galaxies have attracted
attention for their similarity to high-redshift galaxies near the Epoch of
Reionisation – namely low mass, low metallicity, vigorous star formation
(Cardamone et al., 2009), high ionization, and high fractions of escaping
ionizing ultraviolet radiation (e.g. Izotov et al., 2016, 2018; Kim et al.,
2020). Svoboda et al. analysed XMM-Newton observations of three Green Pea
galaxies, finding two to have X-ray luminosities a factor $\sim$4–6 times
brighter than predicted from empirical relations between the X-ray luminosity,
SFR and metallicity for star-forming galaxies. If the X-ray excess is
connected to starburst activity, possible explanations include: High-Mass
X-ray Binaries, Low-Mass X-ray Binaries, high numbers of Ultra-Luminous X-ray
sources, modified initial mass functions, intermediate-mass black holes, or a
substantial contribution from hot gas. The remaining intriguing possibility is
an X-ray flux contribution from a hidden AGN. Kawamuro et al. (2019) performed
relatively short exposure hard X-ray observations of the two XMM-Newton-
detected Green Pea galaxies with NuSTAR, finding both to be non-detections
$\geq$ 10 keV. Using an established mid-infrared vs. X-ray correlation for
local AGN-hosting Seyfert galaxies, Kawamuro et al. showed that if this
correlation holds for Green Pea galaxies, and the observed AGN-like mid-
infrared colours were powered by an AGN, the central engine would have to have
column densities in excess of $N_{\rm H}$ $\gtrsim$ 2$\times 10^{24}$ cm-2 to
explain the NuSTAR non-detections. Thus, none of the hypotheses have been
completely satisfactory to-date to explain the power source responsible for
Green Pea X-ray emission and none can currently be directly proven
observationally.

The contribution of X-ray emission from hot diffuse gas in star-forming
regions has been analysed and considered in a multitude of ways in previous
theoretical studies. Cerviño et al. (2002) computed the X-ray emission from
diffuse gas in a single stellar population heated by stellar winds and
supernovae. An efficiency parameter was used to represent the fraction of
output mechanical energy that heats the gas. However, the efficiency of this
process is very difficult to determine and many different values have been
used in the literature. Stevens & Hartwell (2003) modeled the properties of a
stellar cluster based on the analytical description of the wind by Chevalier &
Clegg (1985). The authors compared model radial temperature profiles and
diffuse emission of the gas to a sample of nearby clusters observed by
_Chandra_. Once again, disentangling the contribution of the diffuse gas from
that of other sources was not straightforward. Aiming to overcome this
problem, Oskinova (2005) followed a similar approach but additionally included
the star clusters’ evolution and emission of stars, before comparing the
predictions to an observed sample of clusters in the Large Magellanic Cloud
(LMC). The author attributed the differences between the observed and
predicted diffuse X-ray emission to the evolutionary stages of the clusters
and to the evolution of X-ray-emitting low-mass stars.

Marcolini et al. (2005) carried out a series of 2-D hydrodynamical simulations
of pc-scale conductive clouds immersed in a galactic-scale wind. The authors
found that as the clouds are shock-processed, the resulting soft X-ray
emission is dominated by cloud material even if the wind metallicity is
several times higher than that of the clouds.

To gain insight into the X-ray emission contribution from hot gas in galaxies
beyond that provided by empirical relations (as used in Svoboda et al.), here
we use a hydrodynamical approach that follows the thermalization of the
kinetic energy deposited by stellar winds and supernova explosions within
stellar clusters (Silich et al., 2004, 2011; Wünsch et al., 2007, 2011). For a
pre-defined population of stars in a given volume, our code – WINDCALC (Wünsch
et al., 2017) – computes the radial profiles of gas density and temperature,
while accounting for wind collisions from individual stars as well as the mass
and energy deposited by supernova explosions. The code has been tested against
other numerical codes and observations, yielding robust results (Silich et
al., 2005; Wünsch et al., 2007, 2011). Here we incorporate the X-ray emission
from such populations of stars. Our work is originally driven by the
unexplained X-ray luminosity of Green Pea galaxies, but is equally useful for
understanding gas emission in less extreme objects such as “ordinary” star-
forming galaxies as well as individual star clusters.

The paper is structured as follows: the numerical code that defines the set-up
of the star clusters and galaxies, as well as the computation of X-ray
luminosities is introduced in Section 2. The results of the modeling, the role
of various parameters, and the comparison with observational data are
presented in Section 3, and discussed in Section 4. Finally, the summary and
conclusions are provided in Section 5.

##  2 Numerical simulations of star clusters, X-ray emission, and galaxies

We aim to calculate the X-ray emission from shocked stellar winds in a
population of star clusters within starburst galaxies, such as those possible
in Green Pea galaxies. First we calculate the expected X-ray emission from
individual star clusters for a grid of cluster masses and cluster ages. We
then use those star clusters to construct a starburst galaxy according to the
SFR of the galaxy. Once the star cluster population has established mass and
age distributions, we integrate the X-ray spectra of the individual star
clusters to give the X-ray spectrum of the starburst galaxy. From the
resulting spectrum we calculate the X-ray luminosity $L_{\rm{X}}$ within the
soft 0.5 – 2.0 keV energy band. We then estimate $L_{\rm{X}}$ for starburst
galaxies as a function of SFR, and investigate how our results change with
different assumptions.

###  2.1 Simulation of a star cluster

####  2.1.1 Distribution of stars within a star cluster

We follow the semi-analytic approach of Wünsch et al. (2017) to model the
hydrodynamics of the shocked stellar winds in individual star clusters, with
masses ranging from $10^{3}$ to $10^{8}$ $M_{\odot}$. For each cluster, we
draw a distribution of stars according to the standard initial mass function
with an index of $\alpha=-2.3$ at the highest mass end (Kroupa, 2001), and
with star masses $m_{\rm{star}}$ between 0.01 and 120 $M_{\odot}$. The stars
are created instantaneously, and are distributed according to a Schuster
profile of $\rho_{\rm{star}}\propto[1+(r/R_{\rm{C}})^{2}]^{-b}$ where
$R_{\rm{C}}$ is the core radius and $b=3/2$ (Palouš et al., 2013; Tenorio-
Tagle et al., 2015). The distance from the centre is denoted by $r$.
Furthermore, we use the cutoff radius $R_{\rm{SC}}$ to define the radius of
the star cluster. If not indicated otherwise, all clusters have a core radius
of $R_{\rm{C}}=1$ pc and a cutoff radius of $R_{\rm{SC}}=3$ pc.

####  2.1.2 Stellar evolution models

We follow the evolution of individual stars with stellar models taken from the
BoOST project (Szécsi et al., 2020). These models were computed with the
‘Bonn’ Code (e.g. Szécsi et al., 2015b, and references therein) assuming two
different metallicities $Z=0.02\,Z_{\odot}$ as for the dwarf galaxy IZw18, and
$Z=0.4\,Z_{\odot}$ as for the Large Magellanic Cloud (LMC). We consider low
metallicities for consistency with the measured sub-solar values for Green Pea
galaxies (e.g., Izotov et al. 2011). Additionally, in Section 3.5 another set
of stellar models from the BoOST project are tested, which incorporate
chemically-homogeneous evolution. We interpolate linearly between the stellar
models using the code synStars (e.g. sect. 4 of Szécsi et al., 2020) to get
synthetic populations that consist of smoothly interpolated stellar tracks.
Additionally, our treatment of massive stars exploding as supernovae is
similar to the procedure used in starburst99 (Leitherer et al., 1999). In low
metallicity environments, (e.g., for the $Z=0.02\,Z_{\odot}$ case), stars with
masses between 8 and 40 $M_{\odot}$ explode as regular supernovae (Heger et
al. 2003, also cf. the detailed discussion in Section 3.5 of Szécsi & Wünsch
2019 about this question), while for the higher LMC-like metallicity we assume
that all stars with masses $m_{\rm{star}}\geq 8\,M_{\odot}$ explode. This
results in different starting times of the supernovae in the framework of the
cluster’s lifetime. Since massive stars live for shorter periods of time, the
supernovae in the higher-metallicity environment ($Z=0.4\,Z_{\odot}$) arise at
$t_{\rm{clus}}\sim 2$ Myr. In the lower metallicity environment
($Z=0.02\,Z_{\odot}$) supernovae start from $t_{\rm{clus}}\sim 4.5$ Myr. Each
supernova contributes an energy of $10^{51}$ erg to the gas and returns a mass
of $m_{\rm{SN}}=m_{\rm{star}}-m_{\rm{remnant}}$ to the free gas available in
the cluster, where the remnant mass is fixed to $m_{\rm{remnant}}=1.4$
$M_{\odot}$.111We note that this assumption is unlikely to hold for the
highest-mass supernovae in our simulations, but the predicted difference in
observed X-ray emission from the resulting clusters is expected to be small.
See Szécsi et al., in prep. for more details.

We distribute the supernovae in a way that the supernova rate is constant
until the end of each stellar track, but changing to the next stellar track
results in a different supernova rate. Thus the supernova rate is in general
variable throughout a full simulation.

At the state of the art one cannot predict in which type of supernova a star
will explode based on its stellar evolution track, and such considerations are
outside the scope of this work. Therefore, predictions about chemical changes
due to supernovae cannot be well constrained. Thus we restrict ourselves to
use the final surface chemical composition of each star for calculating X-ray
spectra.

####  2.1.3 Radial profiles of the star clusters

Figure 1: Radial profiles of the gas temperature (upper left), and gas density
(upper right) for clusters of two different cluster masses (solid and dashed
lines, respectively) and metallicities ($0.02$ $Z_{\odot}$ and $0.4$
$Z_{\odot}$ in black and green, respectively). The profiles are shown for a
cluster evolutionary time of $t=1.5$ Myr. In the lower left panel we present
the corresponding X-ray spectra by integrating up to an aperture radius
$R_{\mathrm{ap}}=100$ pc. Higher densities yield brighter X-ray fluxes whilst
line emission is more pronounced at higher metallicities. The orange shaded
area marks the soft X-ray energy band (0.5 – 2.0 keV), used in this study. The
lower right panel shows the cumulative luminosity of this energy band from the
xspec components as a function of the aperture $R_{\mathrm{ap}}$ of the star
cluster. For all cases, the dense, hot gas in the centre of the cluster
dominates the X-ray emission of the whole cluster. The grey shading in the
upper left panel indicates temperatures that do not contribute to the X-ray
emission in the apec model. Hence, the different curve lengths in the lower
right panel – see the text for details.

We follow the evolution of each star cluster up to an age of
$t_{\rm{clus,max}}=20$ Myr. Our model provides the energy loss, mass loss and
chemical composition of the gas as a function of time.

Assuming spherical symmetry, we calculate the radial profiles of the gas
temperature and gas density in the star cluster up to a distance of
$R_{\rm{max}}=100$ pc from the cluster centre. Figure 1 shows examples of the
radial profiles for two different cluster masses and the two metallicities at
a cluster’s evolutionary time of $t_{\rm{clus}}=1.5$ Myr, before the onset of
supernovae. We find the temperatures in the cluster’s centre to depend mainly
on the gas metallicity. At larger radii, the temperature in the higher mass
star cluster drops faster compared to the more massive cluster, due to the
high mass cluster being in a strongly radiative regime (see Wünsch et al.,
2007; Silich et al., 2004, for further details).

Figure 1 also shows the gas density profiles of the clusters. The mass of the
cluster mainly determines the density of the gas due to the fixed cutoff
radius $R_{\rm{SC}}$– i.e. while the volume of the clusters remains fixed, the
amount of stars grows proportionally with $m_{\rm{clus}}$. Thus the
contribution to the stellar wind increases and the gas density scales linearly
with $m_{\rm{clus}}$. Further, we find the gas density of two clusters with
the same mass to be larger for higher metallicities due to the higher mass
loss of stars at higher metallicity.

###  2.2 X-ray spectrum of a star cluster

####  2.2.1 Calculation of the X-ray emission with xspec

Based on the temperature and density profiles of the star clusters (see top
panels in Figure 1), we calculate the X-ray spectra by integrating up to
$R_{\mathrm{ap}}=100$ pc. We use the model apec (Astrophysical Plasma Emission
Code, v.12.10.1 Smith et al., 2001) within the xspec package (Arnaud, 1996) to
calculate the continuum and line emission from optically thin collisionally-
ionized gas. Assuming spherical symmetry, for each cluster, we take into
account all gas within the aperture, i.e. $r<R_{\mathrm{ap}}$ (this neglects
the contribution from the gas at large radii projected into the aperture,
however, in Section 2.2.2 it is shown that this contribution is very small).
We divide the sphere into 300 shells distributed equidistantly in the
logarithm of temperature. Each of the 300 shells is considered as a single
apec component, whose X-ray contribution is calculated based on four input
parameters:

  * (i)

the gas plasma temperature $k_{B}T$, where we use the mass weighted
temperature of each shell, and the Boltzmann constant $k_{B}$,

  * (ii)

the chemical abundances derived from the stellar models (see Section 2.1.2),

  * (iii)

the redshift, which we set to $z=0$ to obtain the spectrum in the intrinsic
rest frame.

  * (iv)

a normalization parameter, which is the emission measure scaled by the
distance, as $10^{-14}/\\{4\pi[D_{\rm{A}}\,(1+z)^{2}]\\}\times\int
n_{\rm{e}}n_{\rm{H}}\rm{dV}$. Here, $D_{\rm{A}}\,(1+z)^{2}$ corresponds to the
luminosity distance $D_{\rm{L}}$ which we set to the value for SDSS
J074936.77+333716.3 (aka ‘GP1’); $D_{\rm{L}}=1437.9$ Mpc (Svoboda et al.,
2019)222The distance was calculated from redshift $z$ = 0.2733 with
cosmological parameters: $H_{0}$ = 67.8 km s-1 Mpc-1, $\Omega_{\rm M}$ =
0.308, $\Lambda_{0}$ = 0.69.. Despite step (iii), we applied a representative
Green Pea distance to the normalization to obtain appropriate flux levels that
could be directly compared to real measurements from Green Peas (Svoboda et
al., 2019). The number densities of electrons and hydrogen ($n_{\rm{e}}$,
$n_{\rm{H}}$, respectively) result from the stellar models (see Section
2.1.2).

Finally, we sum the X-ray emission of the individual shells to obtain the
X-ray emission of the whole star cluster.

####  2.2.2 X-ray spectra and luminosity of the star clusters

For the modelled star clusters introduced in Section 2.1, we show in the lower
left panel of Figure 1 the corresponding X-ray spectra. The spectra are shown
in the rest frame. Only, the flux level corresponds to a source that would be
located at $z\approx 0.3$, typical for Green Pea galaxies. In general, the
X-ray flux in high mass clusters is larger compared to that of the lower mass
clusters due to the higher densities that ensure more mass is available for
radiating X-ray emission. The temperature of the gas – determined by the gas
metallicity – spreads the continuum to higher energies. Therefore, the spectra
with $Z=0.4\,Z_{\odot}$ drop at higher energies compared to the lower
metallicity cases (green vs. black lines in Figure 1).

To quantify the X-ray emission from the hot gas, we calculate the X-ray
luminosity in the soft X-ray band of 0.5 – 2.0 keV (shaded area in the lower
left panel of Figure 1) by

$L_{\rm{X}}=4\pi d^{2}\int_{0.5\,\rm{keV}}^{2.0\,\rm{keV}}F(E)\rm{d}E,$ (1)

with the distance $d$ set to the same value for GP1 and the photon energy, $E$
= $h\nu$ where $h$ is the Planck constant and $\nu$ is the photon frequency.
In Figure 1 we show the cumulative luminosity $L_{\rm{X}}$ as a function of
the aperture. We find the inner part of the cluster to dominate the X-ray
emission for all cluster masses and metallicities. The apec model only
includes contributions with plasma temperatures $k_{\rm{B}}T\gtrsim 8.0\times
10^{-3}$ keV, (corresponding to $T\sim 9.2\times 10^{4}$ K). Therefore, the
lines in the Figure end once the temperature is beyond this value, and the
temperatures not considered for the X-ray emission calculation are shaded in
grey in the temperature profiles in Figure 1. Further, we use an aperture
radius of $R_{\mathrm{ap}}=100$ pc to provide us with a reasonable safety
factor that ensures we accounted for the majority of the X-ray emission.

####  2.2.3 Evolution of $L_{\rm{X}}$ and the energy insertion rate of the
star clusters

Figure 2: Mechanical energy insertion rate $L_{\rm{SC}}$ (left), X-ray
luminosity $L_{\rm{X}}$ in the energy band of 0.5 – 2.0 keV (middle), and
fraction of the inserted energy transformed into X-ray luminosity from the
X-ray-emitting hot gas (right). In line with the scaling used in previous
works, we normalize $L_{\rm{SC}}$ by the cluster mass $m_{\rm{clus}}$, and
$L_{\rm{X}}$ by $m_{\rm{clus}}^{2}$ in the left and middle panels,
respectively (see Sections 2.2.3 and 3.2.1). The scaling $L_{\rm{X}}\propto
m_{\rm{clus}}^{2}$ is only an approximation, but holds for $m_{\rm{clus}}$
between $10^{3}$ and $1.23\times 10^{7}$ $M_{\odot}$ (solid lines). The dashed
lines in the middle and right panels show the relations for a higher cluster
mass of $m_{\rm{clus}}=8.11\times 10^{7}$ $M_{\odot}$, with the lower
metallicity (black dashed line) showing a disagreement t $\leq$ 4.5 yr before
the onset of supernovae – see Section 2.2.3 for details. The shaded regions in
the right panel show the evolution of $L_{\rm{X}}$ / $L_{\rm{SC}}$ with time
for clusters with masses between 103 (bottom edge) and 1.23 $\times$ 107
$M_{\odot}$ (top edge). Note the zig-zag patterns in the left panel are an
artefact caused by interpolation at the end of stellar tracks when the energy
from supernovae drops abruptly to zero. In all panels, metallicities of Z =
0.02 Z⊙ and Z = 0.4 Z⊙ are represented by black and green colors,
respectively.

Due to stellar winds and supernovae, mechanical energy is transferred to the
gas. This can be measured with the energy insertion rate $L_{\rm{SC}}$ (also
referred as the mechanical luminosity of the cluster), and can be compared to
the X-ray luminosity $L_{\rm{X}}$.

Figure 2 shows the evolution of $L_{\rm{SC}}$ and $L_{\rm{X}}$ as a function
of the cluster age (left and middle panel, respectively). The gas density
$\rho$ scales linearly with the energy insertion rate $L_{\rm{SC}}$ (left
panel), and since $\rho\propto m_{\rm{clus}}$ (see top right panel of Figure
1), we have $m_{\rm{clus}}\propto L_{\rm{SC}}$ (Silich et al., 2004, 2011;
Palouš et al., 2013). Before the onset of supernovae ($t_{\rm{clus}}\leq 2$
Myr for $Z=0.4$ $Z_{\odot}$ and $t_{\rm{clus}}\leq 4.5$ Myr for $Z=0.02$
$Z_{\odot}$), $L_{\rm{SC}}$ is determined by stellar winds, which depend on
the metallicity. Stars with higher metallicities have stronger stellar winds,
and thus, their energy insertion rate $L_{\rm{SC}}$ is larger. Once supernovae
have set in, the energy and mass release are dominated by supernovae, and the
metallicity of the gas is only of minor importance to $L_{\rm{SC}}$.

The X-ray luminosity $L_{\rm{X}}$ (Figure 2 middle panel) scales with
$\rho^{2}$, as $L_{\rm{X}}\propto n_{\rm{e}}n_{\rm{H}}\Lambda\propto\rho^{2}$
(see Silich et al., 2005; Tenorio-Tagle et al., 2015, Equation 3 therein).
Since $\rho\propto m_{\rm{clus}}$, we find approximately $L_{\rm{X}}\propto
m^{2}_{\rm{clus}}$. This scaling holds for $m_{\rm{clus}}$ between $10^{3}$
and $1.23\times 10^{7}$ $M_{\odot}$ in our calculations. For larger cluster
masses ($m_{\rm{clus}}=8.11\times 10^{7}$ $M_{\odot}$, dashed lines), it
breaks down for low metallicity environments at cluster ages between 2 and 4.5
Myr, since the gas inside the clusters become thermally unstable and part of
it cools to low temperatures and forms dense warm/cold clumps that do not
radiate in X-rays. With the onset of supernovae, all clusters in the
considered mass range stabilize, and $L_{\rm{X}}$ becomes comparable for both
metallicities. The maximum of $L_{\rm{X}}$ is reached shortly after the onset
of supernovae. This behaviour is consistent with previous results reported in
the literature (e.g. Cerviño et al., 2002; Oskinova, 2005).

The right panel of Figure 2 shows which fraction of the inserted energy
$L_{\rm{SC}}$ is transformed into X-ray luminosity $L_{\rm{X}}$. The largest
efficiency of transforming mechanical energy into X-rays is achieved at $t\sim
3$ Myr for $Z=0.4$ $Z_{\odot}$, at the same time when $L_{\rm{X}}$ has its
maximum. For a star cluster at this time with a mass of
$m_{\rm{clus}}=10^{7}\,M_{\odot}$, we obtain $L_{\rm{X}}$ / $L_{\rm{SC}}$
$\sim 3\times 10^{-2}$. Thus, 3% of the mechanical energy is transformed into
soft X-ray emission for this cluster mass. For the low metallicity
environment, the peak values of $L_{\rm{X}}$ are smaller. A cluster with
$m_{\rm{clus}}=10^{7}\,M_{\odot}$ at $Z=0.02$ $Z_{\odot}$ would transform a
maximum 0.7% of its mechanical energy into X-rays.

####  2.2.4 X-ray emission of star clusters for a grid of masses and ages

Figure 3: X-ray luminosity $L_{\rm{X}}$ of the soft X-ray band (0.5 – 2.0 keV;
color coded) for star clusters in a grid of different masses and ages with
metallicities 0.02 $Z_{\odot}$ (left) and 0.4 $Z_{\odot}$ (right). We use the
individual X-ray spectra of clusters to construct the X-ray spectrum of a
starburst galaxy (see Section 2.3). In general, the hot gas X-ray emission is
dominated by the younger, more massive star clusters.

We calculate the X-ray spectra for a grid of cluster masses ($10^{3}$ –
$10^{8}$ $M_{\odot}$, $\Delta\log_{10}(m_{\rm{clus}}/M_{\odot})=0.133$) and
ages (0 – 20 Myr, $\Delta t_{\rm{clus}}=1$ Myr). This is the basis for
constructing the X-ray spectra of starburst galaxies. Figure 3 shows
$L_{\rm{X}}$ for the spectra of all clusters for both metallicities (left and
right panels, respectively). Young and massive clusters ($t_{\rm{clus}}\sim 3$
– 6 Myr, $m_{\rm{clus}}\geq 10^{7}\,M_{\odot}$) show the largest X-ray
luminosities $L_{\rm{X}}$ up to $10^{41}$ erg s-1. In line with the results
presented in Figure 2, $L_{\rm{X}}$ becomes comparable at later times
($t_{\rm{clus}}\gtrsim 4.5$ Myr) for both metallicities – i.e. once the
supernova explosions have set in. Before this, the X-ray emission of star
clusters with higher metallicity is larger, since the stellar winds of the
high-metallicity stars contain more mass.

We use this grid of X-ray spectra to estimate the X-ray spectrum of a whole
starburst galaxy, as described in the following section.

###  2.3 Simulation of starburst galaxies composed of star clusters

Next, we expand our study to starburst galaxies by assuming a given starburst
galaxy is composed of star clusters, and estimate the mass & age distributions
of star clusters within such a complex. The mass distribution of star clusters
in galaxies can be described by a cluster mass function, and was recently
studied in observations of the LEGUS dwarf galaxy sample by Cook et al.
(2019). The authors found a mass function of
${\rm{dN}}/{\rm{d}}m_{\rm{clus}}\propto m_{\rm{clus}}^{\beta}$, with
$\beta=-2$333Note that the index of the cluster mass function $\beta$ is not
the same as the index of the initial mass function for stars in the clusters,
denoted $\alpha$ in this work., in agreement with earlier studies by e.g.,
Zhang & Fall (1999); de Grijs et al. (2003); Gieles & Bastian (2008).

We simulate a starburst galaxy cluster population with a fixed SFR as a
function of time up to $t_{\rm{gal,max}}=300$ Myr. For each time step, we draw
star clusters according to the cluster mass function. With proceeding time,
the star clusters increase in age, until $t_{\rm{clus,max}}=20$ Myr, at which
point the clusters are discarded. The minimum time step within the galaxy’s
evolution is set to d$t_{\rm{gal,min}}=10^{3}$ yr. In combination with the
SFR, this determines the minimum mass of clusters that are drawn in each time
step, as

$m_{\rm{min}}={\rm{d}}t_{\rm{gal,min}}\times\rm{SFR}.$ (2)

Thus, we can draw several clusters within one time step, until the sum of
their masses $\sum_{i}m_{{\rm{clus}},i}\geq m_{\rm{min}}$. To fulfill the SFR,
we draw new clusters according to the mass drawn in the last time step
$\sum_{i}m_{{\rm{clus}},i}$ and pause it for subsequent time steps if
necessary. In this case, the clusters currently present in the starburst
galaxy only increase in age.

Note, that the time evolution depicts the age of the starburst galaxy for
fixed SFR, but does not correspond to a real evolutionary time step of the
galaxy. Instead, each time step of the galaxy represents a possible
realisation of the starburst galaxy, and is therefore of statistical
relevance.

Figure 4: Time evolution of (left) the starburst mass (solid lines) and their
cumulative masses (dashed lines), (middle) the number of clusters for
starburst galaxies with different SFRs. We assume here star clusters within a
mass range of $10^{3}$ – $10^{7}\,M_{\odot}$ and a power law index of the
cluster mass function of $\beta=-2$. (Right) the X-ray luminosity $L_{\rm{X}}$
for $Z=0.02$ $Z_{\odot}$. The dashed lines and the shaded areas indicate the
median values and the range between the 25th and 75th percentile for each
starburst galaxy (see Sections 2.3 an 2.4 for details).

Figure 4 shows in the left and middle panels how the mass of the starburst
galaxy and the number of included clusters changes with time, assuming a
constant SFR (color coded). We draw clusters from a cluster mass function with
a power law index of $\beta=-2$, while the cluster masses range from $10^{3}$
– $10^{7}\,M_{\odot}$. In the first 20 Myr the full composition of the
starburst galaxy is still building up, and its mass increases. After this
initial time, we define the starburst galaxy mass as the final mass
accumulated within 20 Myr, as set by

$m_{\rm{gal,max}}=t_{\rm{clus,max}}\times\rm{SFR}.$ (3)

Note, that this mass only contains the contribution of young star clusters,
with ages below 20 Myr. The total cumulative stellar mass formed by the
starburst is shown by dashed lines in the left panel of Figure 4. Comparing it
to stellar masses of galaxies in the Mineo et al. (2012b) sample ($3\times
10^{8}$ – $6\times 10^{10}$ $\mathrm{M}_{\odot}$), we see that for instance
the lowest considered SFR = 0.1 $\mathrm{M}_{\odot}$ yr-1would have to last
$3$ Gyr to explain the whole stellar mass of the lowest mass galaxy from the
sample. On the other hand, the highest considered SFR = 90
$\mathrm{M}_{\odot}$ yr-1needs 3 Myr and 670 Myr to explain the whole stellar
mass of the lowest and the highest mass galaxy, respectively.

We further see in Figure 4 fluctuations in the mass of the starburst galaxy
and the number of clusters, for example for the galaxy with $\rm{SFR}=1$
$M_{\odot}\,\rm{yr}^{-1}$. These fluctuations appear once a massive star
cluster was drawn that exceeds the minimum mass $m_{\rm{min}}$ (see Equation
(2)). Thus, to fulfill the SFR, no further star clusters are added to the
starburst galaxy in the following time steps, while the age of the clusters
continues to increase. The starburst mass as well as the number of star
clusters then decreases, causing the fluctuations. As can be seen in the left
and middle panel of Figure 4, the dips occur after the peaks in the mass of
the starburst galaxy.

###  2.4 X-ray spectrum of starburst galaxies

Figure 5: Distribution of star cluster mass and age in one starburst galaxy
realization with ${\rm{SFR}}=30$ $\mathrm{M}_{\odot}$ yr-1 (left). The color
bar indicates the number of present clusters. The right panel shows the
corresponding X-ray spectra for starburst galaxies with each metallicity
considered ($Z=0.02\,Z_{\odot}$ and $0.4\,Z_{\odot}$ in black and green,
respectively). The X-ray luminosity in the (shaded) soft X-ray band 0.5 – 2.0
keV for each spectrum is indicated in the legend.

We calculate the X-ray spectrum of one starburst galaxy by summing the X-ray
spectra of the individual star clusters present. Figure 5 shows the mass and
age distribution of star clusters for one starburst galaxy realization with
$\rm{SFR}=30$ $\mathrm{M}_{\odot}$ yr-1 at a time step of $t_{\rm{gal}}=100$
Myr (see Figure 4 for comparison). According to the cluster mass function,
most star clusters have masses $m_{\rm{clus}}\lesssim 10^{5}$ $M_{\odot}$, and
only a few are of higher masses. As we have not applied further constraints
with respect to the age distribution of the clusters, the mass distribution is
fulfilled for all cluster ages. Further, there are age bins with less clusters
present. These appear subsequently after we have drawn a high mass cluster in
the previous time step (see Section 2.3).

In the right panel of Figure 5 we present the X-ray spectra of the
corresponding starburst galaxies for both metallicities ($0.02\,Z_{\odot}$ and
$0.4\,Z_{\odot}$ in black and green, respectively). In accordance with the
spectra of the individual clusters, we find more pronounced emission lines for
higher metallicity (see Figure 1). The luminosities of the soft X-ray band
(shaded area in Figure 5) for these two cases are $L_{\rm{X}}\sim 3\times
10^{38}$ erg s-1 for $0.02\,Z_{\odot}$ and $L_{\rm{X}}\sim 1\times 10^{40}$
erg s-1 for $0.4\,Z_{\odot}$. Thus, $L_{\rm{X}}$ of the starburst galaxy is of
a similar order of magnitude as $L_{\rm{X}}$ for individual young, massive
star clusters (see Figure 3). We therefore conclude that the X-ray emission of
the hot gas from the whole starburst galaxy can be dominated by the emission
from the young, massive star clusters.

We calculate the X-ray spectra of the starburst galaxies for each time step
until $t_{\rm{gal,max}}=300$ Myr. For the examples of starburst galaxies
presented in Figure 4, we show the corresponding values of $L_{\rm{X}}$ for a
metallicity of $0.02\,Z_{\odot}$. As a result of the fluctuations in number of
clusters and the changing X-ray emission of the clusters at different ages,
$L_{\rm{X}}$ of the starburst galaxy likewise fluctuates with time. To
characterize $L_{\rm{X}}$ of one starburst galaxy statistically, we estimate
the median value of $L_{\rm{X}}$ in a range of $t_{\rm{gal}}$ between 20 and
300 Myr. We further use the 25th and 75th percentiles of $L_{\rm{X}}$ to
determine the range of fluctuation. These values are shown as shaded
horizontal bands in Figure 4. In the following, we use these measurements for
describing $L_{\rm{X}}$ in our study.

##  3 X-ray emission from hot gas with various assumptions

Table 1: Fit of the X-ray luminosity $L_{\rm{X}}$ as a function of the SFR for all models. As a fit function we use $L_{X}=C\times{\rm{SFR}}^{\zeta}$. We fit the data of our models separately for a lower (SFR $<30$ $\mathrm{M}_{\odot}$ yr-1) and upper (SFR $\geq 30$ $\mathrm{M}_{\odot}$ yr-1) regime, giving ($C_{1}$, $\zeta_{1}$) and ($C_{2}$, $\zeta_{2}$), respectively. We set the limits between the regimes for each model with the break point in the SFR.  Met. | $\beta$ | Mass range | Break. point | $C_{1}$ | $\zeta_{1}$ | $C_{2}$ | $\zeta_{2}$
---|---|---|---|---|---|---|---
$Z_{\odot}$ |  | $M_{\odot}$ |  $M_{\odot}$ yr-1 | $({\rm erg\,s}^{-1})/(M_{\odot}\,{\rm yr}^{-1})^{\zeta}$ |  | $({\rm erg\,s}^{-1})/(M_{\odot}\,{\rm yr}^{-1})^{\zeta}$ |
0.02 | $-2.0$ | 103 \- 106 | 3 | $1.50\times 10^{36}$ | $1.30$ | $(2.20\pm 0.03)\times 10^{36}$ | $1.01\pm 0.00$
0.02 | $-2.0$ | 103 \- 107 | 30 | $(2.00\pm 0.20)\times 10^{36}$ | $1.70\pm 0.03$ | $(2.20\pm 0.19)\times 10^{37}$ | $0.99\pm 0.02$
0.02 | $-2.0$ | 103 \- 108 | 150 | $(4.50\pm 1.40)\times 10^{36}$ | $1.44\pm 0.08$ | $2.50\times 10^{35}$ | $1.96$
0.02 | $-2.0$ | 104 \- 107 | 30 | $(4.50\pm 0.16)\times 10^{36}$ | $1.54\pm 0.01$ | $(1.90\pm 0.11)\times 10^{37}$ | $1.07\pm 0.01$
0.02 | $-2.0$ | 105 \- 107 | 30 | $(7.10\pm 0.37)\times 10^{36}$ | $1.52\pm 0.02$ | $(3.50\pm 0.26)\times 10^{37}$ | $1.02\pm 0.01$
0.02 | $-1.5$ | 103 \- 107 | 10 | $(8.10\pm 0.37)\times 10^{36}$ | $1.88\pm 0.02$ | $(6.50\pm 0.22)\times 10^{37}$ | $1.00\pm 0.01$
0.02 | $-1.2$ | 103 \- 107 | 10 | $(9.00\pm 0.33)\times 10^{36}$ | $1.93\pm 0.02$ | $(8.40\pm 0.26)\times 10^{37}$ | $0.99\pm 0.01$
0.4 | $-2.0$ | 103 \- 106 | 3 | $1.00\times 10^{37}$ | $1.27$ | $(1.70\pm 0.04)\times 10^{37}$ | $1.00\pm 0.00$
0.4 | $-2.0$ | 103 \- 107 | 30 | $(1.20\pm 0.11)\times 10^{37}$ | $1.69\pm 0.04$ | $(1.40\pm 0.18)\times 10^{38}$ | $1.02\pm 0.02$
0.4 | $-2.0$ | 103 \- 108 | 150 | $(2.30\pm 0.47)\times 10^{37}$ | $1.67\pm 0.05$ | $8.40\times 10^{36}$ | $1.80$
0.4 | $-2.0$ | 104 \- 107 | 30 | $(2.30\pm 0.30)\times 10^{37}$ | $1.60\pm 0.05$ | $(1.70\pm 0.14)\times 10^{38}$ | $1.04\pm 0.02$
0.4 | $-2.0$ | 105 \- 107 | 30 | $(3.90\pm 0.56)\times 10^{37}$ | $1.60\pm 0.05$ | $(2.60\pm 0.22)\times 10^{38}$ | $1.02\pm 0.02$
0.4 | $-1.5$ | 103 \- 107 | 10 | $(4.80\pm 0.56)\times 10^{37}$ | $1.80\pm 0.07$ | $(4.80\pm 0.07)\times 10^{38}$ | $1.00\pm 0.00$
0.4 | $-1.2$ | 103 \- 107 | 10 | $(5.20\pm 0.21)\times 10^{37}$ | $1.99\pm 0.02$ | $(6.40\pm 0.28)\times 10^{38}$ | $0.99\pm 0.01$

###  3.1 $L_{\rm{X}}$ as a function of the SFR

Figure 6: X-ray luminosity ($L_{\rm{X}}$) of the starburst galaxies as a
function of their SFR, assuming different metallicities (color coded). The
starburst galaxies consist of star clusters with masses in the range of
$10^{3}$ – $10^{7}$ $M_{\odot}$, distributed according to a cluster mass
function with power law index $\beta=-2$. In general, $L_{\rm{X}}$ is larger
for higher metallicities. The dashed-dotted lines are the fits to our data
(see text). For comparison, the dashed line shows the observationally obtained
relation between $L_{\rm{X}}$ and the SFR for star-forming galaxies from Mineo
et al. (2012b), and the green symbols indicate $L_{\rm{X}}$ for the Green Pea
galaxies from Svoboda et al. (2019).

As a first approach, we analyse the soft X-ray emission of starburst galaxies
composed of star clusters with masses between $10^{3}$ to $10^{7}$
$M_{\odot}$, and a cluster mass function with power law index $\beta=-2$.
Figure 6 shows $L_{\rm{X}}$ of these starburst galaxies as a function of the
SFR for different metallicities. The points indicate the median values of
$L_{\rm{X}}$, and the shaded areas the range between the 25th and 75th
percentiles (see Section 2.4). In general, $L_{\rm{X}}$ increases with
increasing SFR and in line with the previous results, we find $L_{\rm{X}}$ to
be larger for starburst galaxies at higher metallicities. For comparison, we
plot the observed relation between hot gas X-ray luminosity and SFR from Mineo
et al. (2012b); $L_{0.5-2{\rm keV}}/\textrm{SFR}\sim(5.2\pm 0.2)\times
10^{38}$ $({\rm erg\,s}^{-1})/(M_{\odot}\,{\rm yr}^{-1})^{\zeta}$ which was
obtained by fitting X-ray spectra with the mekal model inside xspec. In
contrast to the results of Mineo et al., our results cannot be described by
one single power law. Instead, our relation between $L_{\rm{X}}$ and the SFR
consists of two ranges, in which the range of smaller SFRs ($<30$
$\mathrm{M}_{\odot}$ yr-1) has a steeper index, and the range of larger SFRs
($\geq 30$ $\mathrm{M}_{\odot}$ yr-1) a flatter one. We therefore fit our data
separately for these two ranges, each by a function of

$L_{X}=C\times{\rm{SFR}}^{\zeta}$ (4)

The results of all fits are summarized in Table 1. We indicate by the breaking
point the limit between the two SFR ranges, which is 30 $\mathrm{M}_{\odot}$
yr-1in this case. We find $\zeta_{1}\sim 1.7$ for the range of lower SFRs for
both metallicities, and $\zeta_{2}\sim 1$ for the range of larger SFRs. Thus,
our normalization $C$ can only be compared to the results of Mineo et al.
(2012b) for the range of larger SFRs. Here our data results in $C_{2}\sim
2.20\times 10^{37}$ and $\sim 1.40\times 10^{38}$ $({\rm
erg\,s}^{-1})/(M_{\odot}\,{\rm yr}^{-1})^{\zeta}$ for the metallicities of
$Z=0.02\,Z_{\odot}$ and $0.4\,Z_{\odot}$, respectively. These indices are
about a factor of $24$ and $4$ times smaller than the index found by Mineo et
al. (2012b). We additionally note that for SFRs above 17 $\mathrm{M}_{\odot}$
yr-1, the relation of Mineo et al. that we compare to is an extrapolation.

In Figure 6, the green points indicate $L_{\rm{X}}$ as measured for the three
Green Pea galaxies in Svoboda et al. (2019). As the observed X-ray emission is
even higher for the Green Pea galaxies than the relation by Mineo et al.
(2012b), our model predicts even less of the X-ray emission to originate from
the hot gas. However, our results are obtained on a basis of several
assumptions/restrictions, such as the star cluster’s mass range and the power
law index of the cluster mass function $\beta$, as well as the fixed core and
cutoff radii of the star clusters. In addition, the inclusion of chemically
homogeneous evolving stars could have an impact on our results. Thus in the
following sections, we test how these factors influence the value and scaling
of $L_{\rm{X}}$.

###  3.2 Dependence on star cluster masses

Figure 7: Similar to Figure 6, we show how different assumptions in
constructing the starburst galaxies influence the relation between
$L_{\rm{X}}$ and SFR. We vary the upper (left) and lower mass limit (middle)
of the cluster’s mass range with the power law index of the cluster mass
function $\beta=-2$. In contrast we vary $\beta$ in the (right) panels whilst
$m_{\rm{clus}}\sim 10^{3}$ – $10^{7}$ $M_{\odot}$. The rows represent the
metallicities of $Z=0.02\,Z_{\odot}$ (upper) and $Z=0.4\,Z_{\odot}$ (lower).
$L_{\rm{X}}$ increases if the number of massive star clusters increases and/or
star clusters with larger masses are present. This is the case for a mass
range shifted towards higher cluster masses, and for flatter indices $\beta$
(see Section 3.2).

The star cluster’s mass range in galaxies is an ongoing matter of debate.
Young massive clusters and super star clusters are observationally found to
have masses between $10^{5}$ and $10^{8}$ $M_{\odot}$ (e.g. Whitmore, 2003;
O’Connell, 2004; Bastian et al., 2006; Reines et al., 2008). Recently, Norris
et al. (2019) provided a catalogue of compact stellar systems and argue in
favour of a fundamental maximum star cluster mass of $\sim 5\times
10^{7}\,M_{\odot}$, requiring an initial cluster mass of $10^{8}\,M_{\odot}$.
The low mass end of the star cluster’s mass range is hard to constrain due to
the lower luminosities (see e.g. Figure 5). However, such clusters enclose the
mass in the galaxy and thus have an influence on the amount of mass left for
more massive star clusters. Various observational studies investigated the
distribution of cluster masses; Adamo et al. (2010a) find in Haro 11 more than
30% of the clusters to have masses $>10^{5}\,M_{\odot}$, while the fraction of
low-mass and very young clusters ($<10^{4}\,M_{\odot}$, 1 – 3 Myr) is missing.
Likewise, in SBS 0335 only clusters with high masses ($10^{5}$ to
$10^{6}\,M_{\odot}$) are present (Adamo et al., 2010b). Thus next we study the
dependence of cluster mass range and power law index on the distribution
functions we derive (see Figure 7). Whilst the individual clusters in Green
Pea galaxies are unresolved, we assume the structure could be similar to Haro
11. Furthermore some of the parameters found by Adamo et al. appear extreme
compared to regular starburst galaxies but should be explored in the case of
Green Pea galaxies which also appear extreme.

####  3.2.1 Varying the upper limit of the star cluster mass range

In Figure 7 we show in the left panels the relation between $L_{\rm{X}}$ and
the SFR, assuming different upper limits on the cluster mass range, while
keeping the lower limit constant at $10^{3}$ $M_{\odot}$. For comparison with
observations, we plot the relation found from Mineo et al. (2012b), and the
X-ray luminosity measurements of the Green Pea galaxies from Svoboda et al.
(2019). If we decrease the upper mass of the cluster’s mass range to
$10^{6}\,M_{\odot}$ (light green), the X-ray luminosity $L_{\rm{X}}$ becomes
smaller compared to the cases with larger upper masses for all SFRs. This is a
result of the absence of high mass star clusters in these starburst galaxies.
Although the total mass of the starburst galaxy is fixed by the SFR (see
Equation (3)), the distribution of cluster masses changes with the upper limit
of the cluster mass range. Therefore, as for a single star cluster we find
approximately $L_{\rm{X}}\propto m_{\rm{clus}}^{2}$ (see Section 2.2.3), as
the absence of high mass star clusters in the starburst galaxy leads to a
decrease of $L_{\rm{X}}$. Similarly, $L_{\rm{X}}$ increases when we set the
star cluster mass range to $10^{3}$ – $10^{8}$ $M_{\odot}$.

Further, $L_{\rm{X}}$ fluctuations decrease with increasing SFR (denoted by
the shaded area around the points). The higher the SFR, the more star clusters
at the high mass end of the mass range that can coexist. Since the
fluctuations arise from the appearance of massive clusters, once the number of
high mass clusters stabilizes, the fluctuations decrease. This happens faster
for a lower upper cluster mass of $10^{6}$ $M_{\odot}$.

For a cluster mass range of $10^{3}$ – $10^{8}$ $M_{\odot}$ we only calculate
$L_{\rm{X}}$ for $\rm{SFR}\geq 10$ $\mathrm{M}_{\odot}$ yr-1, as otherwise the
mass of the most massive cluster is larger than the mass of the total
starburst galaxy $m_{\rm{gal,max}}$ (see Equation (3)). In this particular
case, a starburst galaxy at $\rm{SFR}=1$ $\mathrm{M}_{\odot}$ yr-1, has
$m_{\rm{gal,max}}=2\times 10^{7}\,M_{\odot}$ according to Equation (3), and
thus, cannot be represented when including cluster masses up to $10^{8}$
$M_{\odot}$. This is why the points of low SFRs are missing for this cluster
mass range in Figure 7.

We compare the two metallicities in the upper and lower rows of Figure 7.
Qualitatively, $L_{\rm{X}}$ scales similarly with the SFR for both
metallicities. Quantitatively, $L_{\rm{X}}$ is systematically below the
observational scaling by Mineo et al. (2012b), and the X-ray emission of the
Green Pea galaxies for the panels with $Z\sim 0.02$ $Z_{\odot}$. For $Z=0.4$
$Z_{\odot}$ and a cluster mass range of $10^{3}$ – $10^{8}$ $M_{\odot}$
however $L_{\rm{X}}$ is consistent with the observational relation by Mineo et
al. (2012b) within an order of magnitude, which appears to be driven by
increasing the number of higher-mass clusters in a given simulation.

As above, we fit the relation between $L_{\rm{X}}$ and SFR with a broken power
law for two ranges of SFR (Equation (4)). The breaking points separating these
two regimes differ with the cluster mass ranges, and are 3
$\mathrm{M}_{\odot}$ yr-1 for a cluster mass range of $10^{3}$ –
$10^{6}\,M_{\odot}$, and 150 $\mathrm{M}_{\odot}$ yr-1 for a cluster mass
range of $10^{3}$ – $10^{8}\,M_{\odot}$. In the range of lower SFR, the
relation between SFR and $L_{\rm{X}}$ is determined by the number of clusters,
and by the mass of the most massive clusters. In the range of larger SFR, the
relation depends only on the number of clusters present, as several of the
most massive clusters are at the upper mass limit. Thus the breaking point
shifts with the upper limit of the cluster mass range. For all plots shown in
the left panels, we find the range of low SFR $\zeta_{1}>1$, while in the
range of high SFR $\zeta_{2}\sim 1$ (see Table 1). Thus, only the upper range
of SFR is comparable with the relation by Mineo et al. (2012b).

For the mass range of $10^{3}$ – $10^{6}\,M_{\odot}$ we find normalizations
$C_{2}$ of $2.2\times 10^{36}$ and $1.7\times 10^{37}$ $({\rm
erg\,s}^{-1})/(M_{\odot}\,{\rm yr}^{-1})^{\zeta}$ for the metallicities of
$0.02$ and $0.4\,Z_{\odot}$, respectively – i.e. a factor of 236 and 30 times
smaller than the index reported by Mineo et al. (2012b).

####  3.2.2 Varying the lower limit of the star cluster mass range

In Figure 7 we show the relation between $L_{\rm{X}}$ and SFR by varying the
lower limit of the cluster mass range, and keeping the upper limit at
$10^{7}\,M_{\odot}$ (middle panels). By increasing the lower mass limit from
$10^{3}$ to $10^{5}\,M_{\odot}$, the X-ray luminosity $L_{\rm{X}}$ increases,
as the available mass in the starburst galaxy is distributed to clusters of
higher mass. The increase of $L_{\rm{X}}$ is smaller compared to the case
where we varied the upper mass limit, as in the individual star clusters the
X-ray emission scales with $m_{\rm{clus}}^{2}$ (see Sections 2.2.3 and 3.2.1).
Although for the higher metallicity we find larger values for $L_{\rm{X}}$
(bottom middle panel), we do not reach the values expected from observations
by Mineo et al. (2012b) nor by Svoboda et al. (2019). For the fits with
Equation (4) the breaking points between the two ranges in SFR do not change
in this case, because the upper mass of the mass range is constant, and the
mass distributed on lower mass scales remains similar. For the upper SFR
range, $\zeta_{2}\sim 1$ for all cases, and thus, comparable with the relation
found by Mineo et al. (2012b). The indices of the fits in the upper SFR
regimes are between $\sim$2.0 and $\sim$3.5$\times 10^{37}$ $({\rm
erg\,s}^{-1})/(M_{\odot}\,{\rm yr}^{-1})^{\zeta}$ for a metallicity of 0.02
$Z_{\odot}$, and between $\sim$1.4 and $\sim$2.6$\times 10^{38}$ $({\rm
erg\,s}^{-1})/(M_{\odot}\,{\rm yr}^{-1})^{\zeta}$ for a metallicity of 0.4
$Z_{\odot}$ (see Table 1). Comparing with Mineo et al. (2012b), these values
correspond to $\leq 1/15$ and $\leq 1/2$ of the index found in their work for
the two metallicities.

###  3.3 Varying the power law index of the cluster mass function

As introduced in Section 2.3, the distribution of cluster masses within the
given range is described by the cluster mass function with a power law index
$\beta$. In accordance with e.g. Zhang & Fall (1999); de Grijs et al. (2003);
Gieles & Bastian (2008) we set $\beta=-2$ in the analyses described above.
However, it is debated whether $\beta$ varies with the type of galaxy or SFR.
For example in the LEGUS dwarf galaxy sample with SFRs between 0.005 to $0.5$
$\mathrm{M}_{\odot}$ yr-1 Cook et al. (2019) find $\beta=-2$, whereas in NGC
4449 with ${\rm{SFR}}\sim 1.5$ $\mathrm{M}_{\odot}$ yr-1 Annibali et al.
(2011) report a flatter power law index of $\beta\sim-1.5$. Likewise, Adamo et
al. (2010a) find the index is found to be flatter in Haro 11 with a current
SFR of $22\pm 3$ $\mathrm{M}_{\odot}$ yr-1, as observationally indicated by
the luminosity function from the star clusters (Adamo et al., 2010a). With the
luminosity function, the cluster mass function is expected to change, as they
are directly connected to the star formation processes (Adamo et al., 2010a).

Since we construct starburst galaxies at higher SFRs in our simulations, we
analyse how a flatter power law index of the cluster mass function influences
the scaling between $L_{\rm{X}}$ and SFR. The results are shown in the right
panels of Figure 7. We assume a cluster mass range of $10^{3}$ –
$10^{7}\,M_{\odot}$, and vary the cluster mass function $\beta$ between
$-2.0$, $-1.5$, $-1.2$. For flatter indices we find an increase in
$L_{\rm{X}}$ as well as smaller fluctuations of $L_{\rm{X}}$ around the median
value. This is a result of the increase in number of high mass clusters. The
flatter the power law index, the more high mass clusters are present, which
dominate the X-ray luminosity of the whole starburst galaxy. Comparing results
for the two metallicities (upper and lower right panels), we find $L_{\rm{X}}$
to be in agreement with the observational relation of Mineo et al. (2012b) for
$Z=0.4$ $Z_{\odot}$ and $\beta\geq-1.5$. The fit in the upper SFR range for
$\beta=-1.2$ has a normalization of $C_{2}\sim 8.4\times 10^{37}$ and $\sim
6.4\times 10^{38}$ $({\rm erg\,s}^{-1})/(M_{\odot}\,{\rm yr}^{-1})^{\zeta}$
for the two metallicities, respectively, which corresponds to 1/6 and 1/23 of
the index found by Mineo et al. (2012b). However, none of our simulations
reproduces the luminosities observed for the Green Pea galaxies as reported in
Svoboda et al. (2019).

###  3.4 Varying the core and cutoff radii

Figure 8: Similar to Figure 6, we show for a metallicity of $0.02$ $Z_{\odot}$
the relation between $L_{\rm{X}}$ and SFR for starburst galaxies with smaller
cutoff and core radii $R_{\rm{SC}}$ and $R_{\rm{C}}$, as found by Espinoza et
al. (2009) in the Arches cluster (blue symbols). The black symbols represent
the same values as shown in Figure 6. Smaller radii lead to an increased X-ray
luminosity, but $L_{\rm{X}}$ is still beyond the observed values.

So far we have assumed the core radius $R_{\rm{C}}$ and cutoff radius
$R_{\rm{SC}}$ of the star cluster to be $R_{\rm{C}}=1$ pc, and $R_{\rm{SC}}=3$
pc, respectively (see Section 2.1.1). However, some star clusters have been
found to be more compact than this. Espinoza et al. (2009) studied the Arches
cluster with a mass of $\sim 2\times 10^{4}$ $M_{\odot}$, finding
$R_{\rm{C}}=0.14$ pc, and $R_{\rm{SC}}=0.7$ pc (see their Table 8). With these
values, the Arches cluster is one of the densest known young clusters in the
Milky Way. The smaller radii lead to higher concentrations of mass in the
cluster, which leads to increased gas density, temperature and hence X-ray
emission.

We study the influence of $R_{\rm{C}}$ and $R_{\rm{SC}}$ on our results by
repeating the calculation of $L_{\rm{X}}$ as a function of SFR for starburst
galaxies with cluster masses $10^{3}$ – $10^{7}$ $M_{\odot}$, a cluster mass
function power law index of $\beta=-2$, and $R_{\rm{C}}$ and $R_{\rm{SC}}$ set
to the values found by Espinoza et al. (2009). Figure 8 shows the resulting
relation derived after decreasing the cluster radii, which shifts the X-ray
luminosity on average to higher values. However even with these assumptions
$L_{\rm{X}}$ is still smaller than the observational results from Mineo et al.
(2012b) and Svoboda et al. (2019).

###  3.5 Including chemically homogeneous stars

Figure 9: Similar to Figure 6, we show for a metallicity of $0.02$ $Z_{\odot}$
the relation between $L_{\rm{X}}$ and the SFR for starburst galaxies assuming
that 20% of the stars are chemically homogeneous evolving stars (CHE, blue
symbols). The black symbols represent the same values as shown in Figure 6. By
including chemically-homogeneously evolving stars, $L_{\rm{X}}$ increases
slightly. However, it does not reach the X-ray emission expected from
observations.

Star clusters with gas at low-metallicity can host fast rotating stars, which
may evolve quasi chemically-homogeneously (Yoon et al., 2006; Brott et al.,
2011; Szécsi et al., 2015b, a; Szecsi, 2017; Garcia et al., 2019). The
temperatures of such stars are higher than of regular stars, and some may have
strong winds (e.g., Kubátová et al. 2019). Here we test whether the hot
stellar winds of these chemically-homogeneously evolving stars would
contribute significantly to the X-ray emission. We apply the hot chemically-
homogeneously evolving stellar models from the BoOST project (Szécsi et al.,
2020) with low metallicity (Z = 0.02 $Z_{\odot}$), to obtain the mass loss in
each star cluster. For these models, a correction for the wind optical depth
is taken into account, following the method presented in Chapter 4.5.1 of
Szécsi (2016). Based on observational comparisons to the star-forming dwarf
galaxy IZw18, Szécsi et al. (2015b) find that $\sim$20% of the stellar
population is likely evolving homogeneously. Therefore, for our low-
metallicity cluster (Z = 0.02 $Z_{\odot}$) we combine mass loss from 80%
regular and 20% chemically-homogeneously evolving stars when constructing our
star clusters, and calculate their X-ray emission for a grid of masses
($10^{3}$ – $10^{7}$ $M_{\odot}$) and ages ($\leq 20$ Myr) in the same way as
in Section 2.3. From Figure 9, we find the chemically-homogeneously evolving
stars to increase $L_{\rm{X}}$ slightly, but still significantly below the
X-ray luminosities measured for the Green Pea galaxies.

##  4 Discussion

###  4.1 Comparison with previous work

The X-ray emission from hot gas in the interstellar medium from starburst
galaxies was previously simulated in a population synthesis model by Cerviño
et al. (2002). The authors assumed a single starburst galaxy and derived the
energy output from every star using stellar evolution models. Our starburst
galaxies are composed of several star clusters. For modelling the stars,
Cerviño et al. assumed a Salpeter initial mass function, and considered stars
with masses between 2 and 120 $M_{\odot}$, parabolic interpolation between the
Geneva stellar models (Schaller et al., 1992; Meynet et al., 1994) and
supernova explosions. The X-ray emission is then a result of the contribution
from supernova remnants and the hot gas. However, their approach of
calculating the hot gas properties of a starburst galaxy differs significantly
from ours. Instead of computing the wind dynamics, Cerviño et al. define an
efficiency parameter $\epsilon_{\rm{eff}}^{\rm{x}}$ describing the variable
fraction of total mechanical energy that contributes to the heating of the
diffuse interstellar medium.

Cerviño et al. (2002) and Otí-Floranes et al. (2012, 2014) use this model to
find a good set of parameters to explain observations by first constraining
the age of the observed object, then the corresponding kinetic energy released
to the gas at this time, and adopt the mass transformed into stars and
efficiency parameter $\epsilon_{\rm{eff}}^{\rm{x}}$ required to fit the
observation. The authors typically find efficiencies from 0.4% up to 10% for
star-forming objects (see e.g. Summers et al., 2001, 2004; Mas-Hesse et al.,
2008; Otí-Floranes et al., 2012). This is in line with our models; as we have
shown in Section 2.2.3, $\lesssim 3$% of the energy insertion rate in a star
cluster is transformed into X-rays. The actual value depends on the
metallicity as well as the mass and evolutionary state of the star cluster.

###  4.2 Comparison with observations

The X-ray emission of a sample of nearby star-forming galaxies was studied in
a series of papers by Mineo et al. (2012a, b, 2014). The galaxies considered
had SFRs of $\sim 0.1$ – $17$ $\mathrm{M}_{\odot}$ yr-1. To infer the gas
emission, Mineo et al. (2012b) measured the total X-ray flux of each galaxy
and subtracted the flux of resolved and unresolved faint High-Mass X-ray
Binaries. In addition, the authors attempt to correct for unresolved sources
of X-ray emission from Low-Mass X-ray Binaries, young stellar objects and
cataclysmic variables. However, the contributions of both are found to be
negligible and the models estimating their contribution to be too uncertain to
perform an accurate subtraction. Therefore, the X-ray luminosity of the
diffuse interstellar medium (0.5 – 2.0 keV) presented by Mineo et al. (2012b)
needs to be considered as an upper limit.

Mineo et al. (2012b) distinguish two components in the unresolved emission:
the diffuse X-ray emission and the X-ray emission of the thermal component
(modelled with the mekal model in xspec). Here, the diffuse emission refers to
the contribution of sources with hard spectra (e.g., unresolved faint compact
sources and/or hotter gas) that could not be reproduced well from the lower-
temperature mekal components. The authors find linear relations with indices
of $L_{0.5-2\,\rm{keV}}^{\rm{diff}}/{\rm SFR}\sim(8.3\pm 0.1)\times
10^{38}\,$$({\rm erg\,s}^{-1})/(M_{\odot}\,{\rm yr}^{-1})^{\zeta}$ and
$L_{0.5-2\,\rm{keV}}^{\rm{mekal}}/{\rm SFR}\sim(5.2\pm 0.1)\times
10^{38}\,$$({\rm erg\,s}^{-1})/(M_{\odot}\,{\rm yr}^{-1})^{\zeta}$. We
compared the relation between the X-ray luminosity and the SFR of our work to
the mekal component of Mineo et al. (2012b). In general, the X-ray luminosity
of our models varies when assuming different ranges of cluster masses within
the starburst galaxy, and assuming different power law indices $\beta$ of the
cluster mass function. For the widely reported value of $\beta=-2$ and
$m_{\rm{clus}}$ in the range of $10^{3}$ – $10^{7}$ $M_{\odot}$ we find our
X-ray emission to be smaller compared to the emission expected from e.g.,
Mineo et al.. However, if we allow more massive clusters to be present by
taking a larger mass range ($m_{\rm{clus}}\sim 10^{3}$ – $10^{8}$ $M_{\odot}$)
or by adopting a flatter power law index of $\beta\geq-1.5$, as suggested by
some studies of dwarf galaxies (e.g. Adamo et al., 2010a; Annibali et al.,
2011), the X-ray emission of our models for $Z=0.4$ $Z_{\odot}$ (see Figure 7)
is in line with the expectations by Mineo et al. (2012b).

###  4.3 Limitations of our model

The model we designed for constructing starburst galaxies assumes stars form
in individual star clusters within the starburst galaxy. We demonstrated how
the star cluster mass range and cluster mass function power law index $\beta$
influences the X-ray luminosity of the hot gas from the colliding stellar
winds. Here we discuss additional aspects influencing the X-ray emission from
the hot gas.

When calculating the X-ray emission of the hot gas in the interstellar medium,
we have accounted for the collision of nearby stellar winds, but the
interaction with the interstellar medium surrounding the star clusters has
been neglected. We have also not modeled the presence of clouds exposed to the
galactic-scale wind, as done by Marcolini et al. (2005). Such effects are
mainly important for massive clusters and high density environments, and could
be taken into account in future studies.

For each star cluster, we assume the same core radius $R_{\rm{C}}$ and cutoff
radius $R_{\rm{SC}}$, independent of the cluster mass (see Section 2.1).
Therefore, the volume of each cluster is the same, and the gas density in the
clusters increases with increasing cluster mass. We have chosen moderate
values of $R_{\rm{C}}=1$ pc and $R_{\rm{SC}}=3$ pc, but setting the radii to
smaller values increases the gas density and thus increases $L_{\rm{X}}$ (see
Section 3.4). However, from observations we expect the star clusters to
populate a range in mass and size (e.g. Baumgardt & Hilker, 2018), such that
density would not scale directly with cluster mass. Since our model does not
account for such a spread in cluster sizes, we highlight that $L_{\rm{X}}$ may
be underestimated for low mass star clusters, whilst overestimated for high
mass clusters if the cluster sizes are likewise larger. However, there is only
a weak systematic relation between star cluster mass and radius (see e.g.
Larsen, 2004; Krumholz et al., 2019).

We construct starburst galaxies out of star clusters according to their SFR.
To do so, we add star clusters to the starburst galaxies so that the SFR is
fulfilled with ${\rm{SFR}}=m_{\rm{min}}/{\rm{d}}t$ (see Section 2.3 and
Equation (2)). Since we only keep the star clusters up to an age of 20 Myr, we
thus limit the total mass of the starburst galaxy and likewise the maximum
mass of the most massive star cluster that can be present (see Equation (3)).
For example, a starburst galaxy with ${\rm{SFR}}=1$ $\mathrm{M}_{\odot}$ yr-1
(${\rm{SFR}}=0.1$ $\mathrm{M}_{\odot}$ yr-1) can only contain star clusters up
to $m_{\rm{clus}}\sim 10^{7}$ $M_{\odot}$ ($m_{\rm{clus}}\sim 10^{6}$
$M_{\odot}$). Thus the mass range of star clusters limits the range of SFR we
can study with our model.

We construct starburst galaxies in a way that the SFR is constant, averaged
over time scales of $\leq 20$ Myr. In observations, the formation of stars can
take place in several episodes, and thus, the SFR varies over time. When
estimated observationally, the SFR is then an average over time spans between
a few (from H$\alpha$ measurements) and several hundred (from FUV
measurements) Myr. For comparing our models with observations, we ideally need
to compare to recent SFRs, since the young massive star clusters are those
that dominate the X-ray emission of the starburst galaxy.

###  4.4 Interpretation of Green Pea galaxies

The motivation of this study was to test whether the excess X-ray emission in
Green Pea galaxies studied by Svoboda et al. (2019) can be explained by the
emission from hot gas in star clusters. According to our models, the stellar
winds do not significantly contribute to the X-ray emission with Green Pea-
like metallicities, and are hence unlikely to explain the X-ray excess
observed in these sources. We found the X-ray emission to increase when
considering larger cluster masses, smaller core and cutoff radii of the star
clusters, and a flatter power law index $\beta$ of the cluster mass function.
In general, young massive star clusters dominate the X-ray spectrum of the
starburst galaxy in our model. To obtain X-ray luminosities of $L_{\rm{X}}\sim
10^{41}$ – $10^{42}$ erg s-1 as found for the Green Pea galaxies, the
starburst galaxy could host massive star clusters. The most massive star
clusters our model calculates have $m_{\rm{clus}}\sim 10^{8}$ $M_{\odot}$,
which result in an X-ray luminosity of $L_{\rm{X}}\sim 10^{40}$ – $10^{41}$
erg s-1 depending on the metallicity (see Sections 2.2.3 and 2.2.4). The Green
Pea galaxies have SFRs in the range of $\sim 40$ – $60$ $\mathrm{M}_{\odot}$
yr-1, which in our model would correspond to $m_{\rm{gal,max}}\sim 10^{9}$
$M_{\odot}$ (see Equation (3)). Thus if the complete mass $m_{\rm{gal,max}}$
of a Green Pea galaxy is enclosed in massive star clusters of $10^{8}$
$M_{\odot}$, we could reproduce the measured luminosity. However, we find it
more likely that there are other sources of X-ray emission in the X-ray bright
Green Pea galaxies discovered by Svoboda et al..

In low metallicity environments such as Green Pea galaxies, it is expected
that massive stars lose less mass during their lives and thus form heavier
black holes than in high-metallicity environments. Therefore, High-Mass X-ray
Binaries could be more common in low metallicity environments (Schaerer et
al., 2019), and are thus another possible source for X-ray emission that can
be tested in future work.

##  5 Summary

We calculated the soft X-ray emission in the energy band 0.5 – 2.0 keV from
colliding stellar winds in star clusters for a grid of cluster masses and
ages, assuming sub-solar metallicities ($0.02$ and $0.4$ $Z_{\odot}$) and
including supernova explosions of massive stars. We use these star clusters to
construct starburst galaxies with a given SFR defined via the mass function of
the star clusters, and to estimate the X-ray emission in the starburst
galaxies from the hot gas in star clusters. How we construct the starburst
galaxies is based on various assumptions such as the range of cluster masses
and the cluster mass function power law index $\beta$. For $Z=0.02$
$Z_{\odot}$, we find the X-ray emission of our model to underpredict the X-ray
luminosities observed in star-forming galaxies (see e.g. Mineo et al., 2012b).
For $Z=0.4$ $Z_{\odot}$, we are in good agreement with Mineo et al. (2012b)
once we include star clusters up to $10^{8}$ $M_{\odot}$, or choose a flatter
cluster mass function index of $\beta\gtrsim-1.5$. We then varied model
parameters to study the general behaviour of the X-ray emission under
different model assumptions.

  * (i)

We find the X-ray emission to be higher in star clusters and starburst
galaxies with larger metallicities. Stars at higher metallicity have stronger
mass loss, and therefore the gas density and ultimately the X-ray emission in
these clusters is larger.

  * (ii)

The X-ray emission of the starburst galaxies is dominated by the young massive
star clusters, and varies with the assumed cluster mass range. Increasing the
range of considered cluster masses towards higher values yields increased
X-ray luminosity, as well as the X-ray luminosity fluctuations.

  * (iii)

When excluding low mass clusters in the considered mass range, the mass of the
starburst galaxies is distributed among clusters with higher masses, forcing
the X-ray emission of the starburst galaxies higher, and reducing the
fluctuations.

  * (iv)

The power law index $\beta$ of the cluster mass function influences the
distribution of the star clusters, and therefore the estimated X-ray emission.
Varying the cluster mass function $\beta$ from the fiducial value of $-2$ to a
flatter index ($\gtrsim-1.5$) allows a larger number of high mass clusters to
be present in the constructed starburst galaxies and thus the X-ray luminosity
increases.

  * (v)

Assuming smaller cluster radii would increase the gas density and cause more
X-ray emission.

  * (vi)

Including chemically-homogeneous evolving stars also leads to an increased
X-ray luminosity, though it is not sufficient to close the gap between our
estimated luminosities and those observed for e.g., the Green Pea galaxies.

Interpreting these results in the context of the enhanced soft X-ray emission
measured in the Green Pea galaxies by Svoboda et al. (2019), we conclude that
the wind material does not significantly contribute to the X-ray emission.
Thus, other sources of X-ray emission must be tested in combination with
deeper multi-wavelength observations to disentangle the elusive emission
mechanism(s) responsible for such bright X-ray emission in Green Pea galaxies.
Possible remaining mechanisms for the Green Pea X-ray excess include:
unusually high quantities of High- and/or Low-Mass X-ray Binaries, large
quantities of Ultra-Luminous X-ray sources, modified initial mass functions,
Intermediate-Mass Black Holes or AGN.

We thank the referee for constructive comments which helped improve the paper.
A.F. thanks Angela Adamo and Nathalie Webb for fruitful discussions about star
clusters and Green Pea galaxies. This work has been supported by the
institutional project RVO:67985815. AF & RW acknowledge financial support from
the Czech Science Foundation project No. 19-15008S. S.M.G acknowledges support
from CONACYT-México research grant A1-S-28458 and Cátedra n.482. P.B. & J.S.
acknowledge financial support from the Czech Science Foundation project No.
22-22643S. D.Sz. has been supported by the Alexander von Humboldt Foundation
and was funded in part by the National Science Center (NCN), Poland under
grant number OPUS 2021/41/B/ST9/00757. For the purpose of Open Access, the
author has applied a CC-BY public copyright license to any Author Accepted
Manuscript (AAM) version arising from this submission. The authors thankfully
acknowledge the computer resources, technical expertise and support provided
by the Laboratorio Nacional de Supercómputo del Sureste de México, CONACYT
member of the network of national laboratories.

## References

  * Adamo et al. (2010a) Adamo, A., Östlin, G., Zackrisson, E., et al. 2010a, MNRAS, 407, 870, doi: 10.1111/j.1365-2966.2010.16983.x
  * Adamo et al. (2010b) Adamo, A., Zackrisson, E., Östlin, G., & Hayes, M. 2010b, ApJ, 725, 1620, doi: 10.1088/0004-637X/725/2/1620
  * Annibali et al. (2011) Annibali, F., Tosi, M., Aloisi, A., & van der Marel, R. P. 2011, AJ, 142, 129, doi: 10.1088/0004-6256/142/4/129
  * Arnaud (1996) Arnaud, K. A. 1996, in Astronomical Society of the Pacific Conference Series, Vol. 101, Astronomical Data Analysis Software and Systems V, ed. G. H. Jacoby & J. Barnes, 17
  * Bastian et al. (2006) Bastian, N., Saglia, R. P., Goudfrooij, P., et al. 2006, A&A, 448, 881, doi: 10.1051/0004-6361:20054177
  * Basu-Zych et al. (2016) Basu-Zych, A. R., Lehmer, B., Fragos, T., et al. 2016, ApJ, 818, 140, doi: 10.3847/0004-637X/818/2/140
  * Basu-Zych et al. (2013) Basu-Zych, A. R., Lehmer, B. D., Hornschemeier, A. E., et al. 2013, ApJ, 774, 152, doi: 10.1088/0004-637X/774/2/152
  * Baumgardt & Hilker (2018) Baumgardt, H., & Hilker, M. 2018, MNRAS, 478, 1520, doi: 10.1093/mnras/sty1057
  * Brorby et al. (2016) Brorby, M., Kaaret, P., Prestwich, A., & Mirabel, I. F. 2016, MNRAS, 457, 4081, doi: 10.1093/mnras/stw284
  * Brott et al. (2011) Brott, I., de Mink, S. E., Cantiello, M., et al. 2011, A&A, 530, A115, doi: 10.1051/0004-6361/201016113
  * Cardamone et al. (2009) Cardamone, C., Schawinski, K., Sarzi, M., et al. 2009, MNRAS, 399, 1191, doi: 10.1111/j.1365-2966.2009.15383.x
  * Cerviño et al. (2002) Cerviño, M., Mas-Hesse, J. M., & Kunth, D. 2002, A&A, 392, 19, doi: 10.1051/0004-6361:20020785
  * Chevalier & Clegg (1985) Chevalier, R. A., & Clegg, A. W. 1985, Nature, 317, 44, doi: 10.1038/317044a0
  * Cook et al. (2019) Cook, D. O., Lee, J. C., Adamo, A., et al. 2019, MNRAS, 484, 4897, doi: 10.1093/mnras/stz331
  * de Grijs et al. (2003) de Grijs, R., Anders, P., Bastian, N., et al. 2003, MNRAS, 343, 1285, doi: 10.1046/j.1365-8711.2003.06777.x
  * Douna et al. (2015) Douna, V. M., Pellizza, L. J., Mirabel, I. F., & Pedrosa, S. E. 2015, A&A, 579, A44, doi: 10.1051/0004-6361/201525617
  * Espinoza et al. (2009) Espinoza, P., Selman, F. J., & Melnick, J. 2009, A&A, 501, 563, doi: 10.1051/0004-6361/20078597
  * Garcia et al. (2019) Garcia, M., Evans, C. J., Bestenlehner, J. M., et al. 2019, arXiv e-prints, arXiv:1908.04687.  https://arxiv.org/abs/1908.04687
  * Gieles & Bastian (2008) Gieles, M., & Bastian, N. 2008, A&A, 482, 165, doi: 10.1051/0004-6361:20078909
  * Grimes et al. (2007) Grimes, J. P., Heckman, T., Strickland, D., et al. 2007, ApJ, 668, 891, doi: 10.1086/521353
  * Heger et al. (2003) Heger, A., Fryer, C. L., Woosley, S. E., Langer, N., & Hartmann, D. H. 2003, ApJ, 591, 288, doi: 10.1086/375341
  * Izotov et al. (2011) Izotov, Y. I., Guseva, N. G., & Thuan, T. X. 2011, ApJ, 728, 161, doi: 10.1088/0004-637X/728/2/161
  * Izotov et al. (2016) Izotov, Y. I., Orlitová, I., Schaerer, D., et al. 2016, Nature, 529, 178, doi: 10.1038/nature16456
  * Izotov et al. (2018) Izotov, Y. I., Worseck, G., Schaerer, D., et al. 2018, MNRAS, 478, 4851, doi: 10.1093/mnras/sty1378
  * Kawamuro et al. (2019) Kawamuro, T., Ueda, Y., Ichikawa, K., et al. 2019, ApJ, 881, 48, doi: 10.3847/1538-4357/ab2bf6
  * Kim et al. (2020) Kim, K., Malhotra, S., Rhoads, J. E., Butler, N. R., & Yang, H. 2020, ApJ, 893, 134, doi: 10.3847/1538-4357/ab7895
  * Kroupa (2001) Kroupa, P. 2001, MNRAS, 322, 231, doi: 10.1046/j.1365-8711.2001.04022.x
  * Krumholz et al. (2019) Krumholz, M. R., McKee, C. F., & Bland -Hawthorn, J. 2019, ARA&A, 57, 227, doi: 10.1146/annurev-astro-091918-104430
  * Kubátová et al. (2019) Kubátová, B., Szécsi, D., Sander, A. A. C., et al. 2019, A&A, 623, A8, doi: 10.1051/0004-6361/201834360
  * Larsen (2004) Larsen, S. S. 2004, A&A, 416, 537, doi: 10.1051/0004-6361:20034533
  * Leitherer et al. (1999) Leitherer, C., Schaerer, D., Goldader, J. D., et al. 1999, ApJS, 123, 3, doi: 10.1086/313233
  * Lopez et al. (2014) Lopez, L. A., Krumholz, M. R., Bolatto, A. D., et al. 2014, ApJ, 795, 121, doi: 10.1088/0004-637X/795/2/121
  * Marcolini et al. (2005) Marcolini, A., Strickland, D. K., D’Ercole, A., Heckman, T. M., & Hoopes, C. G. 2005, MNRAS, 362, 626, doi: 10.1111/j.1365-2966.2005.09343.x
  * Mas-Hesse et al. (2008) Mas-Hesse, J. M., Otí-Floranes, H., & Cerviño, M. 2008, A&A, 483, 71, doi: 10.1051/0004-6361:20078398
  * Meynet et al. (1994) Meynet, G., Maeder, A., Schaller, G., Schaerer, D., & Charbonnel, C. 1994, A&AS, 103, 97
  * Mineo et al. (2014) Mineo, S., Gilfanov, M., Lehmer, B. D., Morrison, G. E., & Sunyaev, R. 2014, MNRAS, 437, 1698, doi: 10.1093/mnras/stt1999
  * Mineo et al. (2012a) Mineo, S., Gilfanov, M., & Sunyaev, R. 2012a, MNRAS, 419, 2095, doi: 10.1111/j.1365-2966.2011.19862.x
  * Mineo et al. (2012b) —. 2012b, MNRAS, 426, 1870, doi: 10.1111/j.1365-2966.2012.21831.x
  * Norris et al. (2019) Norris, M. A., van de Ven, G., Kannappan, S. J., Schinnerer, E., & Leaman, R. 2019, MNRAS, 488, 5400, doi: 10.1093/mnras/stz2096
  * O’Connell (2004) O’Connell, R. W. 2004, in Astronomical Society of the Pacific Conference Series, Vol. 322, The Formation and Evolution of Massive Young Star Clusters, ed. H. J. G. L. M. Lamers, L. J. Smith, & A. Nota, 551.  https://arxiv.org/abs/astro-ph/0405038
  * Oskinova (2005) Oskinova, L. M. 2005, MNRAS, 361, 679, doi: 10.1111/j.1365-2966.2005.09229.x
  * Otí-Floranes et al. (2012) Otí-Floranes, H., Mas-Hesse, J. M., Jiménez-Bailón, E., et al. 2012, A&A, 546, A65, doi: 10.1051/0004-6361/201219318
  * Otí-Floranes et al. (2014) —. 2014, A&A, 566, A38, doi: 10.1051/0004-6361/201323069
  * Palouš et al. (2013) Palouš, J., Wünsch, R., Martínez-González, S., et al. 2013, ApJ, 772, 128, doi: 10.1088/0004-637X/772/2/128
  * Ranalli et al. (2003) Ranalli, P., Comastri, A., & Setti, G. 2003, A&A, 399, 39, doi: 10.1051/0004-6361:20021600
  * Reines et al. (2008) Reines, A. E., Johnson, K. E., & Goss, W. M. 2008, AJ, 135, 2222, doi: 10.1088/0004-6256/135/6/2222
  * Rosen et al. (2014) Rosen, A. L., Lopez, L. A., Krumholz, M. R., & Ramirez-Ruiz, E. 2014, MNRAS, 442, 2701, doi: 10.1093/mnras/stu1037
  * Schaerer et al. (2019) Schaerer, D., Fragos, T., & Izotov, Y. I. 2019, A&A, 622, L10, doi: 10.1051/0004-6361/201935005
  * Schaller et al. (1992) Schaller, G., Schaerer, D., Meynet, G., & Maeder, A. 1992, A&AS, 96, 269
  * Silich et al. (2011) Silich, S., Bisnovatyi-Kogan, G., Tenorio-Tagle, G., & Martínez-González, S. 2011, ApJ, 743, 120, doi: 10.1088/0004-637X/743/2/120
  * Silich et al. (2005) Silich, S., Tenorio-Tagle, G., & Añorve-Zeferino, G. A. 2005, ApJ, 635, 1116, doi: 10.1086/497532
  * Silich et al. (2004) Silich, S., Tenorio-Tagle, G., & Rodríguez-González, A. 2004, ApJ, 610, 226, doi: 10.1086/421702
  * Smith et al. (2001) Smith, R. K., Brickhouse, N. S., Liedahl, D. A., & Raymond, J. C. 2001, ApJ, 556, L91, doi: 10.1086/322992
  * Stevens & Hartwell (2003) Stevens, I. R., & Hartwell, J. M. 2003, MNRAS, 339, 280, doi: 10.1046/j.1365-8711.2003.06184.x
  * Summers et al. (2001) Summers, L. K., Stevens, I. R., & Strickland, D. K. 2001, MNRAS, 327, 385, doi: 10.1046/j.1365-8711.2001.04722.x
  * Summers et al. (2004) Summers, L. K., Stevens, I. R., Strickland, D. K., & Heckman, T. M. 2004, MNRAS, 351, 1, doi: 10.1111/j.1365-2966.2004.07749.x
  * Svoboda et al. (2019) Svoboda, J., Douna, V., Orlitová, I., & Ehle, M. 2019, ApJ, 880, 144, doi: 10.3847/1538-4357/ab2b39
  * Szécsi (2016) Szécsi, D. 2016, PhD thesis, Mathematisch-Naturwissenschaftlichen Fakultät der Universität Bonn
  * Szecsi (2017) Szecsi, D. 2017, in XII Multifrequency Behaviour of High Energy Cosmic Sources Workshop (MULTIF2017), 65.  https://arxiv.org/abs/1710.05655
  * Szécsi et al. (2015a) Szécsi, D., Langer, N., Sanyal, D., et al. 2015a, in Wolf-Rayet Stars, ed. W.-R. Hamann, A. Sander, & H. Todt, 189–192
  * Szécsi et al. (2015b) Szécsi, D., Langer, N., Yoon, S.-C., et al. 2015b, A&A, 581, A15, doi: 10.1051/0004-6361/201526617
  * Szécsi & Wünsch (2019) Szécsi, D., & Wünsch, R. 2019, ApJ, 871, 20, doi: 10.3847/1538-4357/aaf4be
  * Szécsi et al. (2020) Szécsi, D., Wünsch, R., Agrawal, P., & Langer, N. 2020, arXiv e-prints, arXiv:2004.08203.  https://arxiv.org/abs/2004.08203
  * Tenorio-Tagle et al. (2015) Tenorio-Tagle, G., Silich, S., Martínez-González, S., Terlevich, R., & Terlevich, E. 2015, ApJ, 800, 131, doi: 10.1088/0004-637X/800/2/131
  * Whitmore (2003) Whitmore, B. C. 2003, in A Decade of Hubble Space Telescope Science, ed. M. Livio, K. Noll, & M. Stiavelli, Vol. 14, 153–178
  * Wünsch et al. (2017) Wünsch, R., Palouš, J., Tenorio-Tagle, G., & Ehlerová, S. 2017, ApJ, 835, 60, doi: 10.3847/1538-4357/835/1/60
  * Wünsch et al. (2007) Wünsch, R., Silich, S., Palouš, J., & Tenorio-Tagle, G. 2007, A&A, 471, 579, doi: 10.1051/0004-6361:20077282
  * Wünsch et al. (2011) Wünsch, R., Silich, S., Palouš, J., Tenorio-Tagle, G., & Muñoz-Tuñón, C. 2011, ApJ, 740, 75, doi: 10.1088/0004-637X/740/2/75
  * Yoon et al. (2006) Yoon, S. C., Langer, N., & Norman, C. 2006, A&A, 460, 199, doi: 10.1051/0004-6361:20065912
  * Zhang & Fall (1999) Zhang, Q., & Fall, S. M. 1999, ApJ, 527, L81, doi: 10.1086/312412