Free Access
Issue
A&A
Volume 626, June 2019
Article Number A85
Number of page(s) 28
Section Astrophysical processes
DOI https://doi.org/10.1051/0004-6361/201834350
Published online 18 June 2019

© ESO 2019

1. Introduction

Post-starburst galaxies are systems in which star-formation appears to have been rapidly quenched after a violent starburst episode (Dressler & Gunn 1983; Couch & Sharples 1987). Such systems are seen both at low redshift (Dressler & Gunn 1983; Couch & Sharples 1987; French et al. 2015, 2018; Rowlands et al. 2015; Alatalo et al. 2016) and high redshift (Watson et al. 2015; Hashimoto et al. 2018; Owen et al. 2019) and the abundant reservoirs of interstellar molecular gas often found in these objects would ordinarily enable vigorous star-formation to continue – however, it instead appears to be quenched (French et al. 2015; Rowlands et al. 2015; Alatalo et al. 2016). Researchers have therefore started to explore alternative mechanisms that would explain the observed quiescence. Proposed mechanisms include: starvation (Gunn & Gott 1972; Cowie & Songaila 1977; Nulsen 1982; Moore et al. 1999; Boselli & Gavazzi 2006), strangulation (Larson et al. 1980; Boselli & Gavazzi 2006; Peng et al. 2015), feedback from the starburst and associated stellar end-products (Goto 2006; Kaviraj et al. 2007; French et al. 2015), and the removal of dense molecular gas in outflows (Narayanan et al. 2008; potential evidence of which would be low-ionisation emission lines in the central regions of post-starburst outflows, for example Yang et al. 2006; Yan et al. 2006; Tremonti et al. 2007).

The presence of the large molecular reservoir, dust and a substantial interstellar medium (ISM) in many of these systems disfavours the action of outflows (Roseboom et al. 2009; Rowlands et al. 2015; Alatalo et al. 2016; French et al. 2018; Smercina et al. 2018). However, a deficiency of dense molecular gas pockets recently observed in low-redshift post-starburst systems (French et al. 2018) implies the operation of processes able to provide additional pressure support to the molecular gas reservoir against fragmentation and collapse. There is debate over the exact process which may be at work here but, if it is feedback from the starburst episode itself, a promising candidate would be the action of high energy (HE) cosmic rays (CRs): these are able to directly provide additional pressure support (e.g. Salem & Bryan 2014; Zweibel 2016), as well as drive interstellar heating effects (e.g. Begelman 1995; Zweibel 2013, 2016; Owen et al. 2018) so as to also increase thermal pressure support. Indeed, such CR heating would be particularly focussed in denser regions of the host galaxy’s ISM. This is where increased thermal pressure support would be particularly required to hamper gravitational collapse and halt the formation of the dense pockets required for star-formation.

Not only do CRs provide the required feedback channels, they are also known to be abundant in starburst galaxies – evidence of which includes the bright γ-ray glow of, for example, Arp 220 (e.g. Griffin et al. 2016; Peng et al. 2016; Yoast-Hull et al. 2017), NGC 253 and M82 (Acero et al. 2009; VERITAS Collaboration 2009; Abdo et al. 2010a; Rephaeli & Persic 2014) and M31 (Abdo et al. 2010b) among others (e.g. Abdo et al. 2010c,d; Hayashida et al. 2013; Tang et al. 2014; Rojas-Bravo & Araya 2016) which results from CR interactions. The source of these CRs is thought to be stellar end products, for example supernova (SN) remnants (Hillas 1984; Ackermann et al. 2013; Kotera & Olinto 2011; Morlino 2017), which can accelerate CRs to relativistic energies by, for example, Fermi processes (Fermi 1949). The high star-formation rates of starburst protogalaxies leads to frequent core-collapse SN explosions from their rapidly ageing population of massive, low-metallicity Type O and B stars. As such, starburst environments are effective CR factories (Karlsson 2008; Lacki et al. 2011; Lacki & Thompson 2012; Wang & Fields 2014; Farber et al. 2018), and it would therefore not be surprising for CRs to have an important role in feedback and quenching, both in high-redshift protogalaxies as well as their more local counterparts.

This study particularly considers the role of CRs in high-redshift protogalactic environments, when these systems were forming their first populations of stars. We organise this paper as follows: in Sect. 2 we consider the indications of bursty star-formation histories and the presence of quenching in starburst and post-starburst systems observed at high redshift. We then introduce a parametric model of a characteristic protogalaxy environment in Sect. 3, in which we consider the injection and interactions of CRs. In Sect. 4, we model CR propagation in protogalaxies, accounting for their interactions, secondaries, cooling and the ambient magnetic field. In Sect. 5, we then consider how CR heating arises and can be driven by two processes – either by direct Coulomb thermalisation, or via inverse-Compton radiative emission – and we specify a computational scheme to calculate the CR heating power in each case. We apply our heating model to a sample of observed high-redshift starburst and post-starburst systems in Sect. 6, and discuss the implications of the results for galactic evolution and star-formation histories. We summarise the work in Sect. 7.

2. Quenching of star-formation in high-redshift protogalaxies

Recent studies have shown substantial star-forming activity in galaxies at high redshifts (e.g. Watson et al. 2015; Hashimoto et al. 2018). These galaxies can be separated into two groups: those currently undergoing a starburst, and those where the starburst appears to have been quenched. An example of such quiescence is MACS1149-JD1, in which there is evidence of a second, older population of stars alongside those which would have been born more recently during the observed ongoing star-formation phase. Thus, the onset of star-formation had to arise in this system at a time long before the epoch at which it was observed (z = 9.11), possibly as early as z ≈ 15, only 250 Myr after the Big Bang (Hashimoto et al. 2018; Owen et al. 2019). Consideration of the observed star-formation rate, stellar mass and the spectroscopic best-fit ages of the young and old stellar populations in this system leads to an intriguing conclusion: to attain its estimated stellar mass by z = 9.11, MACS1149-JD1 must have formed stars at a much higher rate on average than the measured star-formation rate. Hence, the observed star-formation must be relatively suppressed compared to the earlier episode during which the bulk of the stellar population would have formed. Moreover, the star-formation history (SFH) must have been burst-like but subsequently quenched – not unlike the SFHs inferred from closer post-starburst systems (French et al. 2015, 2018; Rowlands et al. 2015).

A number of similar high-redshift galaxies were found to exhibit burst-like SFHs, and we have identified 11 other examples (see Table 1). These systems tend to show a relatively high stellar mass combined with a relatively low inferred star-formation rate, being high-redshift galaxies in a post-starburst, quenched stage of their evolution. If enforcing such an assumption on their SFH, sufficient information is available in these 12 cases to estimate the timescale over which their rapid star-formation episode progressed, if also adopting an upper estimate for the redshift of formation of the first stars in these galaxies, zf. Recent studies (e.g. Hashimoto et al. 2018) estimate this redshift to be around (see also Thomas et al. 2017 which suggests a similar but slightly lower redshift of formation of zf = 14.8), which would correspond to a time when the age of the Universe was around 250–300 Myr (Wright 2006). It is likely that this would be an extreme case for the formation of the earliest galaxies, so could reasonably be regarded as a lower age limit, with many galaxies forming their first stars at later times – for example, evidence suggests a later zf ≈ 9 for A1689-zD1 is more appropriate (Mainali et al. 2018).

Table 1.

High-redshift systems with a comparatively suppressed rate of star-formation, but show evidence of a developed stellar population indicative of an earlier starburst episode.

It is also informative to consider certain ongoing starbursts at high-redshift where observations allow for the star-formation rate to be found, and the characteristic stellar population age to be estimated to give indications as to the duration of the active starburst. We have been able to identify four examples where this is already possible from the literature (see Table 2). This provides us with information about pre-quenched systems as well as their quenched counterparts.

Table 2.

Sample of high-redshift starburst galaxies identified in the literature which appear to be undergoing a starburst event during the epoch at which they are observed.

2.1. Post-starburst protogalaxies and indications of quenching

The 12 examples of post-starburst high-redshift candidates are listed in Table 1 (from lowest to highest inferred burst-phase star-formation rates), with literature references indicated. In all cases, redshift, star-formation rate, stellar mass, stellar population age and dynamical timescale were measured or could be inferred from the literature. Estimates for the upper limit of the starburst period length, the lower limit for the star-formation rate during this starburst period, and upper limit of the dynamical timescale, τdyn, for each system are also shown. τdyn is estimated either from the velocity dispersion of emission lines in the spectrum (Binney & Tremaine 2008; Förster Schreiber et al. 2009; Epinat et al. 2009; Gnerucci et al. 2011) if available1, or by taking the stellar mass as a lower bound on the galaxy’s dynamical mass (presumably there would be components in the galaxy other than stars, for example interstellar gas, which means the true dynamical mass would be larger) to allow a lower bound on τdyn to be found (e.g. Binney & Tremaine 2008). We note that, physically, τdyn represents the minimum time required for the effects of a feedback process to manifest themselves (whether by strangulation, heating or other mechanisms) and impact on the progression of subsequent star-formation.

It is possible to crudely estimate the timescale over which the starburst phase would have been active using

(1)

which is the difference between the approximate time since the starburst ended (tend ≈ τ*, the characteristic stellar population age), and the estimated lower limit of the formation time of the galaxy, at redshift zf. This gives a rough upper limit to the plausible timescale over which the starburst proceeded, and may be justified when considering that tSB is considerably smaller than t(zobs)−t(zf) such that the starburst is bound to be triggered near the t(zf) end of the timeline. Indeed, we later find that this approximation is not the main source of uncertainty in our model, and so is sufficient for our purposes. This allows a lower bound for the starburst phase star-formation rate to be estimated as

(2)

which holds if the vast majority of the galaxy’s stellar mass forms within the starburst period and does not significantly evolve after the end of the starburst episode (while some fraction of the more massive stars in these systems will have gone through their evolutionary lifecycle by the time these galaxies are observed, the number of stars with masses sufficiently large to yield a SN within a few hundred Myr would only comprise a small fraction of their total stellar population and would not have a large impact on the overall galaxy stellar mass). The quantities and tSB are calculated in Table 1 for each of the 12 post-starburst galaxies, which begins to give some insight into the physical processes operating within them: For instance, it can be seen that quenching does not occur immediately after the onset of star-formation. Instead, it takes several dynamical timescales to reduce star-formation to a relatively quiescent state, suggesting feedback is a gradual process rather than instantaneous.

One outlier in this sample would appear to be A1689-zD1. This has a particularly long estimated starburst period, with a correspondingly low starburst star-formation rate. Moreover, it has also been found to harbour an abundance of interstellar dust – at a level of around 1% of its stellar mass (Watson et al. 2015; Knudsen et al. 2017). Dust in galaxies is thought to be built up over time as stellar populations evolve, being injected by, for example, AGB stars (Ferrara et al. 2016), but a large fraction would be destroyed by the frequent SN events in a starburst (Bianchi & Schneider 2007; Nozawa et al. 2007; Nath et al. 2008; Silvia et al. 2010; Yamasawa et al. 2011). In this case, it would seem that the SN rate is insufficient to destroy forming dust on a competitive timescale to its formation, given that the dust destruction time would be of the order 0.1–1 Gyr for a star-formation rate of ℛSF ≈ 5.0 M yr−1 (Temim et al. 2015; Aoyama et al. 2017).

2.2. Protogalaxies with ongoing starburst activity

The four examples of high-redshift systems with an ongoing starburst episode are listed in Table 2 in order of their star-formation rates, with literature references indicated. Stellar ages, star-formation rates, the timescales over which each starburst has been ongoing as well as associated dynamical timescales are shown. These are found by spectral fitting methods (as outlined in e.g. Robertson et al. 2010; Stark et al. 2013; Oesch et al. 2014; Mawatari et al. 2016). The length of the starburst phase of these systems can be estimated as:

(3)

where M* is the stellar mass of the system, and ℛSF is the star-formation rate. M is the mass of the Sun.

From Table 2, the estimated stellar population ages in the two most actively star-bursting systems (GN-108036 and SXDF-NB1006-2) are substantially lower than their dynamical timescales. This means that there has not been sufficient time for the impacts of any feedback or quenching process to have taken hold, and so is consistent with their ongoing high star-formation rates. By contrast, the minimum dynamical timescale for the other two star-forming systems (GN-z11 and EGS-zs8-1) is shorter than the estimated timescale over which their observed starburst appears to have been ongoing, τdyn <  tSB. Thus, in these cases it is possible that feedback effects may be starting to influence their ongoing star-forming activity, perhaps accounting for their comparatively lower star-formation rates. We note that ISM cooling is relatively fast, arising on timescales of around

(4)

for ne as the electron (or ionised gas) number density and Te as the electron (gas) temperature. This means that, for star-formation to remain quenched, the underlying heating process must continue to operate for hundreds of Myr to reconcile the observations. It is clear that stellar radiation would be substantially reduced after the end of the starburst episode, and so is not a plausible candidate to maintain ongoing heating. While heating due to CRs would be more powerful (see Owen et al. 2018), the CR escape timescale (via diffusion) is too short to allow direct CR heating to maintain prolonged quenching,

(5)

where ℓ is the characteristic length-scale of the protogalaxy (around 1 kpc), and D is the characteristic diffusion coefficient (around 3 × 1028 cm2 s−1, see Sect. 4 and, for example Berezinskii et al. 1990; Aharonian et al. 2012; Gaggero 2012), which gives a CR escape timescale of just a few Myr. This would be even shorter in the presence of a fast galactic outflow (Owen et al. 2019). Given that other usual candidates for driving ongoing quenching are not present in these starburst systems (e.g. AGN), a second feedback mechanism must be operating to account for the long-lasting effect. Owen et al. (2019) proposed a strangulation mechanism (Larson et al. 1980; Boselli & Gavazzi 2006; Peng et al. 2015), where cold filamentary inflows of gas able to reignite star-formation are stifled and prevented from returning to the galaxy for hundreds of Myr, which would be consistent with the quenching timescales seen in some of the post-starburst systems considered here. In such a picture, we may reconcile the gas cooling with the molecular gas reservoirs by noting that the reservoir density would be substantially lower than the inner ISM where star-formation would be concentrated. This would allow for much longer cooling timescales and prolonged pressure support of the reservoir, but requires feedback activity operating on two scales: locally, with direct heating able to quench star-formation in the ISM, and more broadly with some larger scale heating or ionisation process required to strangulate star-formation for a substantial period of time.

3. Cosmic rays in protogalaxies

As actively star-forming environments, high redshift starburst protogalaxies are expected to be rich in CRs. The high star formation rates of such galaxies (see e.g. Di Fazio et al. 1979; Solomon & Vanden Bout 2005; Pudritz & Fich 2012; Ouchi et al. 2013; Knudsen et al. 2016) yield frequent SN events which produce the environments such as SN remnants known to be CR accelerators (see e.g. Dar & de Rújula 2008; Kotera & Olinto 2011), and these can efficiently accelerate CRs up to PeV energies (Bell 1978; Hillas 1984; Kotera & Olinto 2011; Schure & Bell 2013; Bell et al. 2013) via, for example diffusive shock acceleration processes such as Fermi acceleration (Fermi 1949; Berezhko & Ellison 1999; Allard et al. 2007). We specify an appropriate parametric model for a protogalaxy to assess the impacts of these CRs. The model is comprised of a radiation and density field with which the CRs interact, together with a magnetic field through which the CRs propagate (see also Owen et al. 2018, which adopts a similar parametrisation).

3.1. Modelling protogalactic environments

3.1.1. Radiation field

A protogalactic radiation field is comprised of photons from stellar radiation as well as the cosmological microwave background (CMB). The CMB is modelled as a spatially uniform black body of a characteristic temperature determined by the redshift, for which the photon number density is given by

(6)

where λC is the Compton wavelength of an electron (h/mec), Γ(...) is the gamma function and ζ(...) is the Riemann zeta function. Further, Θ(z)=kBT(z)/mec2, and T(z)=T0(1 + z) is the CMB temperature with T0 = 2.73 K as its present value (Planck Collaboration XIII 2016). At a typical redshift of z = 7, the energy density in CMB photons is around UCMB ≈ 1100 eV cm−3. The stellar radiation field may be modelled as the distributed emission from an extended ensemble of sources, weighted according to the protogalactic density profile. In a high-redshift rapidly star-forming galaxy, we may assume a stellar population dominated by low-metallicity, massive Type O and B stars with short lifetimes. A system of N = 106 such stars can be shown to give a total protogalactic luminosity NL* ≈ 2.8 × 1044 erg s−1, and this would correspond to a galaxy with ℛSN = 0.1 yr−1, or ℛSF ≈ 16 M yr−1 (e.g. Owen et al. 2018). The stellar lifetime may be considered as independent of ℛSN, so it would follow that this stellar luminosity should scale with ℛSN as stars would build-up in proportion to their birth rate. The energy density of the stellar radiation field then follows as U ≈ 150 (ℛSF/16 M yr−1) (rgal/1 kpc)−2 eV cm−3 where rgal is the characteristic radius of the galaxy, meaning that the stellar radiation is dominated by the CMB until the star-formation rate is well in excess of ℛSF ≈ 120 M yr−1 for a 1 kpc sized galaxy.

3.1.2. Density and magnetic field

We model the density field as an over-density on a uniform background intra-cluster medium (ICM) using the profile

(7)

(Tremaine et al. 1994; Dehnen 1993), with ξc = (r/rc)2 and ξh = (r/rh)2 as the dimensionless form of the parameters rc = 1 kpc and rh = 2 kpc as the protogalaxy core and halo radius respectively. We adopt values of nb, 0 = 10 cm−3 for the central value of the interstellar density profile and np = 10−3 cm−3 for the background density. We find the exact choice of these parameters, if reasonable, bears little influence on the conclusions drawn from the results presented in this study.

Observational studies favour the rapid development of magnetic fields in protogalaxies, reaching strengths comparable to the Milky Way within a few Myr of their formation (Bernet et al. 2008; Beck et al. 2012; Hammond et al. 2012; Rieder & Teyssier 2016; Sur et al. 2018). The development of a magnetic field in young galaxies is often attributed to a turbulent dynamo mechanism driven by the frequent SN explosions during a starburst phase (see e.g. Rees 1987; Balsara et al. 2004; Beck et al. 2012; Latif et al. 2013; Schober et al. 2013), which allows for the emergence of an ordered μG magnetic field in the required timescales. The saturation level of such a mechanism, for example as that introduced in the model of Schober et al. (2013), may be approximated by invoking equipartition with the kinetic energy of the turbulent gas, BL, sat = 𝒮(r)f where

(8)

with ρ as the local density and νf as the fluctuation velocity (νf ≈ Rgal(2πρG/3)1/2 for the protogalaxy, if adopting a pressure with gravity equilibrium approximation; see Schober et al. 2013). f represents a deviation from exact equipartition to account for the efficiency of energy transfer from the turbulent kinetic energy to magnetic energy, and simulation work estimates this to be around 10% (see, e.g. Federrath et al. 2011; Schober et al. 2013).

3.2. Cosmic ray injection

The power of injected primary CR protons is directly related to the total power of the SNe injecting them. Hence, following Owen et al. (2018), a function describing the injection of hadronic CRs at a position r may be described as a product of two separable components,

(9)

where the volumetric CR injection rate is SN(r), and the power-law energy spectrum component, typical of that arising from stochastic acceleration processes, has the normalisation

(10)

where we adopt a power-law index for the protons of Γ = 2.1, similar to regions where CRs are expected to be freshly injected (see, e.g. Allard et al. 2007; Kotera et al. 2010; Kotera & Olinto 2011; Owen et al. 2019). We take a reference energy E0 of 1 GeV and a maximum energy of interest of 1 PeV. LCR, eff is the total power in the CR protons, which can be expressed in terms of SN event rate ℛSN:

(11)

We use the parameters ESN as the total energy generated per SN event (around 1053 erg for core-collapse Type II P SNe), ε as the fraction of SN power converted into to CR power (0.1 is chosen as a conservative value2) and ξ as the fraction of SN energy available after accounting for losses to neutrinos (we use a value of 1% for the energy retained3), α ≈ 0.05 as the fraction of stars able to yield a core-collapse SN event, which is set by the nature of the initial mass function (IMF) of the stars developing in a galaxy and the upper cut-off mass MSN = 50 M (see, e.g. Fryer 1999; Heger et al. 2003) of a star to ultimately produce a SN event (use of the upper mass limit here ensures a conservative estimate for the CR luminosity).

3.3. Cosmic ray interactions in protogalaxies

Above a threshold energy of around a GeV (e.g. Kafexhiu et al. 2014), HE CR particles can undergo hadronic interactions with the interstellar gases and also with the radiation fields permeating the galaxy (see Dermer & Menon 2009; Owen et al. 2018). Hadronic interactions allow high energy CRs to deposit energy much more rapidly than their lower-energy counterparts and this enables them to drive a heating effect, although the eventual thermalisation of this energy is governed by the secondary particles of these interactions and the processes they undergo. These secondaries are comprised of further hadrons, charged leptons, photons and neutrinos (Pollack & Fazio 1963; Gould & Burbidge 1965; Stecker et al. 1968; Berezinsky & Grigor’eva 1988; Mücke et al. 1999). Neutrinos and high energy photons interact only minimally with their environment and so these effectively free-stream away from the vicinity of the interaction, taking their energy with them. The leptons injected from the pion decays typically each inherit around 3–5% of the energy of the original CR primary and can thermalise into their medium to drive a heating effect, or cool through other mechanisms (e.g. Blumenthal 1970; Rybicki & Lightman 1979; Dermer & Menon 2009; Schleicher & Beck 2013) – notably, radiative emission: at high-redshift, the increased number density of CMB photons allow CR secondary electrons to undergo inverse-Compton scattering (Schober et al. 2015). This may also arise in intensely star-forming environments where the energy density in the stellar radiation field is sufficiently high. This process up-scatters relatively low-energy CMB photons to keV energies, creating an X-ray glow (Lacki & Beck 2013; Schober et al. 2015). This X-ray emission can then heat the environment by scattering off ambient free electrons.

The injection of secondary CR electrons is discussed in more detail in Appendix A (specifically Appendices A.1 and A.3) where these interactions have been presented in general terms. Of interest to this study is the role these processes have on the protogalactic environment specified in Sects. 3.1.1 and 3.1.2. In the following, we consider the propagation path lengths in the centre of the protogalaxy model where CR losses would be most severe (as the density and magnetic fields are the strongest) to assess the relative importance of each process in cooling and absorbing the CR primary protons and secondary electrons. We define the energy loss path length of freely-streaming particles undergoing interactions (either cooling or absorption) as

(12)

where β is the velocity of the particle normalised to the speed of light (at high energies β ≈ 1) and τint ≈ [int]−1 is the timescale of an interaction (the inverse of the rate at which interaction events occur). γe is the electron Lorentz factor (proportional to its energy) and b(γe) is the electron cooling rate at the energy Ee = γemec2. This can be used to determine the relative importance of different processes in a given environment. The interactions which are experienced most strongly by a particle occur more rapidly and therefore have a shorter associated path length. While a beam of particles will experience all relevant absorption or cooling mechanisms, their relative path lengths give the weighted contribution of each process – in many cases, one will dominate.

The effective average path lengths experienced by CR protons and electrons due to each of the cooling and absorption processes discussed in Appendices A.1 and A.3 are shown in Fig. 1, with proton losses plotted in the left panel and electron losses in the right panel. These are calculated for an infinite uniform environment of conditions equivalent to those found at the centre of the protogalaxy model, at a redshift of z = 7 and with an ambient magnetic field strength of 5 μG which itself is not assumed to perturb or influence the propagation of the CR particles (an idealistic scenario such that the synchrotron cooling rate can be reasonably compared to other processes). This shows that, for the CR primary protons, losses are mainly dominated by the pp interaction above its threshold energy. Only at energies above 1019 eV, CMB photo-pion (dashed blue line) losses begin to become important, however the CR flux associated with such high energies would typically be negligible. At the lowest energies, below the pp interaction threshold, the free-streaming CRs are essentially only cooled by the adiabatic expansion of the Universe (if freely streaming). This indicates that the vast majority of secondary electrons are injected into the protogalactic medium by the pp pion-production interaction, and so the absorption parameter α* used in Eq. (A.14) is dominated by the term. For the secondary CR electrons, losses are dominated by Coulomb losses and inverse-Compton losses with the CMB (additionally stellar photons may contribute to inverse-Compton losses in extremely active star-forming systems), with other mechanisms having a comparatively negligible effect.

thumbnail Fig. 1.

Left: CR proton losses in terms of effective interaction path lengths. This results from a protogalaxy model with target density of np = 10 cm−3, a CMB radiation field specified at a redshift of z = 7. The thick black line indicates the total effect that would be experienced by a CR proton at the indicated energy accounting for all energy loss contributions. The lines, as labelled, are (1) pp Pion-production (pp, pπ±π0); (2) CMB Photo-pion (pγ, π±π0); (3) CMB Photo-pair (pγ, e±); (4) adiabatic losses due to cosmological expansion; (5) total. Right: secondary CR electron loss lengths with the same model. The lines, as labelled, are (1) adiabatic losses due to cosmological expansion; (2) inverse-Compton (radiative) cooling due to the CMB with (2a) as the estimated path length in the Thomson limit, and (2b) being the estimated path length if accounting for the transition to Klein-Nishina behaviour at higher energies (we find this does not significantly impact our later results, so the Thomson limit is used in our main calculations); (3) Coulomb losses (responsible for directly heating the ISM plasma); (4) synchrotron (radiative) cooling in a 5 μG magnetic field, being an appropriate strength for this model; (5) Triplet Pair-Production (eγ, e±); (6) free–free (Bremsstrahlung) losses; (7) total.

Open with DEXTER

The cooling and absorption mechanisms are also redshift dependent. Figure 2 shows the total CR losses in the protogalactic environment at redshifts of z = 0, 2, 4, 6, 8 and 10, indicated by the solid grey, dot-dashed black, dashed black, dotted black, dashed grey and solid black lines respectively. The protons (left panel) are largely unaffected by the redshift evolution. This is because their losses are dominated by the pp-pion production channel which is governed by the local density field of the galaxy (this does not evolve with redshift) and is therefore effectively decoupled from cosmological expansion. The electrons (right panel) are more affected by the expanding Universe in that their losses become much more influenced by inverse-Compton cooling at high redshift. This dominates over Coulomb cooling at lower energies as redshift increases, with the turn-over (the peak in the cooling curve) occurring above 10 GeV at z = 0, falling to 0.1 GeV by z = 10. This is consistent with the increasing photon number density of the CMB at higher redshifts.

thumbnail Fig. 2.

Interaction loss lengths accounting for all relevant processes (i.e. the solid line of Fig. 1), but calculated at redshifts of z = 0 − 10 as shown by the legend. Left: path lengths for protons; right: path lengths for electrons.

Open with DEXTER

4. Propagation of cosmic rays

The path lengths calculated according to Eq. (12) assume that the CR particles freely stream from their point of emission, propagating in a straight line. However, in general the propagation of charged CRs in ambient galactic (and protogalactic) magnetic fields is more complicated. A high energy (relativistic) particle of charge q in a uniform magnetic field describes a curved path, of size characterised by the gyro-radius (or Larmor radius):

(13)

This can be used to set the scale of a phenomenological description of CR propagation in terms of their diffusion through a turbulent interstellar magnetic field. For this, we may invoke the transport equation:

(14)

where processes such as advection of CRs in bulk flows, winds and galactic outflows have been neglected4. Here, Q(E, r) is the particle source term, S(E, r) is the absorption term, b(E, r) is the particle cooling rate (being the sum of all relevant cooling processes) and D(E, r) is introduced as the diffusion coefficient which takes a parametric form of

(15)

where ⟨|B|⟩|r = |B(r)| is the characteristic mean strength of the magnetic field at some position r, and D is normalised to D0 = 3.0 × 1028 cm2 s−1, a value comparable to empirical measurements of the Milky Way ISM diffusion coefficient (see, e.g. Berezinskii et al. 1990; Aharonian et al. 2012; Gaggero 2012)5 for a 1 GeV CR proton in a 5 μG magnetic field (with corresponding Larmor radius rL, 0). δ is introduced to encode the interstellar magnetic turbulence. We adopt a value of 1/2 here (e.g. Berezinskii et al. 1990; Strong et al. 2007). This is appropriate for a Kraichnan-type turbulence spectrum following a power law of the form P(k) dk ≈ k−2 + δ, which is thought to be a reasonable description for the turbulence in an ISM (see Yan & Lazarian 2004; Strong et al. 2007). For incompressible Kolmogorov turbulence δ = 1/3, but previous studies suggest that this would not allow for the required scattering and diffusion of CRs to be consistent with Milky Way observations (see Goldreich & Sridhar 1995; Lazarian & Pogosyan 2000; Stanimirović & Lazarian 2001; Brandenburg & Lazarian 2013). We solve the transport equation to model steady-state distributions of CR protons and electrons in a protogalaxy.

4.1. Proton transport

Figure 1 shows that CR protons are predominantly absorbed rather than cooled. We may therefore simplify the transport Eq. (14) to

(16)

where we have used the contraction np = np(E, r). This is the (differential) number density of protons per energy interval between Ep and Ep + dEp,

(17)

for Np as the number of protons and dV as a differential volume element. The source term Qp is the local injection rate of protons per unit volume by SN events: Eq. (9) gives the energy injection rate of CR protons at a specific energy Ep; since we require the (differential) injected rate of CRs per energy interval between Ep and Ep + dEp, we use

(18)

The sink term is dominated by pp pion production over the energies of interest, with the CR absorption rate given as ; see Eq. (A.8). Given that our model protogalaxy is axisymmetric, Eq. (16) may be written in spherical coordinates as:

(19)

It is useful to rewrite Eq. (19) as an initial value problem for a single proton injection event from a source region of size 𝒱S:

(20)

where r′ is introduced as the radial distance between the injection point ri and some general location r, given by |rri|. This inhomogeneous partial differential equation (PDE) can be converted to a homogenous form by use of the substitution,

(21)

where we invoke the approximation that the ISM density, nb is locally uniform. While this is not strictly true, it is a close enough approximation for our purposes: we will model the CR propagation only in the inner regions of the protogalaxy ISM where radial density variations are small. Beyond this, matters such as how the protogalactic magnetic field connects to that of the circumgalactic and intergalactic medium have an important effect on CR propagation, and must be considered. These are complicated matters worthy of a dedicated investigation, and robust modelling of this interfacing region lies outside of the scope of the present work. Here, the element c dt = ds′ represents a displacement along the path actually traversed by a CR proton. This is the quantity that determines the amount of attenuating material through which a proton traverses. The propagation of CRs is macroscopically modelled diffusively, but the true speed of these relativistic particles remains close to c. Their frequent scatterings in the magnetic field causes them to exhibit a random walk such that their movement away from their injection point is substantially less than c when considered in terms of their large-scale displacement. Thus, we have introduced the coordinate s′ as the proton path coordinate (i.e. along the path followed by a CR as it scatters). The relation between s′ and the spatial dimension r′ is thus ds′=ψ(r′) dr′, where ψ(r′) quantifies the magnetic scattering (calculated locally at r′). Following Owen et al. (2018), ψ(r′) can be determined from the ratio of the free-streaming (i.e. the free attenuation length due to CR interactions ℓ) and diffusive path lengths of the CRs,

(22)

where we have used the short-hand notation D = D(Ep, r′), and . It then follows that , where we have substituted for the absorption coefficient (in the case of pp pion production dominating the CR absorption), as defined in Eq. (A.14). From this, we may define the CR attenuation factor in terms of Eq. (21) as:

(23)

which quantifies the level of attenuation experienced by a beam of CR protons between a source at location ri, and some general location r.

The substitution Eq. (21) now allows Eq. (20) to be written as the homogeneous PDE:

(24)

where an injection episode at r = ri is taken as an initial condition and provides a normalisation constant ξ0. The solution is well known and can be found by method of Green’s functions to give:

(25)

(Owen et al. 2018) where t′ is the time elapsed since the injection took place. m is introduced as a geometrical parameter which, in the case of a three-dimensional spherical system takes the value m = 3. The normalisation ξ0 depends on the strength of the injection (being the product of the volumetric injection rate Qp(Ep) and a characteristic assigned source size, 𝒱S) and the time Δt over which it was active, ξ0 = Qp(Ep,ri) 𝒱S Δt. Transposing this back to a solution of np yields:

(26)

which practically separates the diffusive evolution from the attenuation term, thus allowing the attenuation to be calculated separately from the diffusion of CRs, and combined afterwards. The solution for CRs injected in multiple bursts from a single source follows as the convolution of individual, instantaneous injections. However, for analytical tractability, we may consider the discrete injection episodes from each source as a continuous process such that the convolution becomes a time-integral. The solution follows as:

(27)

for each source at location ri, where

(28)

The result evidently tends rapidly towards the steady-state solution which, as t → ∞, becomes:

(29)

where the upper incomplete Gamma function has been evaluated as

(30)

and x → 0 as t → ∞, such that . CR protons rapidly settle into this steady state distribution, dnp/dt = 0 (compared to galactic timescales), and so remain in effective equilibrium while CR injection is ongoing. The full galactic steady-state solution can be found by convolving this continuous single-point solution over a weighted distribution of injecting sources. This lends itself naturally the Monte-Carlo (MC) numerical scheme outlined in Sect. 5.3.

4.2. Electron transport

Unlike CR protons, HE electrons predominantly interact with their environment by continuous energy transfer, which may be regarded as cooling. As a result, Eq. (14) is instead reduced to:

(31)

where the contraction ne = ne(E, r) has been used for the CR electron number density per energy interval between Ee and Ee + dEe,

(32)

and b = b(Ee, r)=|dEe/dt| is the magnitude of the sum of the relevant cooling processes experienced by the CR electrons. From Fig. 1, it can be seen this is dominated by inverse-Compton and Coulomb processes. The electron source term Qe(Ee, r) is the secondary CR electron injection rate per unit volume per energy interval due to CR primary interactions, as determined by Eq. (A.24). Like the CR protons, invoking an axisymmetric system allows Eq. (31) to be written in spherical coordinates as:

(33)

We require the CR electron distribution in the protogalaxy in its steady-state, when dne/dt = 0, because this represents the persisting condition into which the system settles. After the proton distribution has stabilised to a steady-state, the injection of electrons will do so too as it is balanced by cooling and diffusion. This enables Eq. (33) to be written as:

(34)

As with Eq. (19), we may restate this as a boundary value problem for a single electron injection event from a source region of size 𝒱s at some location ri, which yields the homogenous PDE:

(35)

Here, r′ is again introduced as the distance between an injection point ri and an observation location r with the injection term giving the boundary condition as the steady-state number of electrons at the source, given by the product of their total injection rate and cooling rate in an energy interval between Ee and dEe, that is Ne, 0(Ee, ri)=Qe(Ee, ri) 𝒱Sτcool(Ee). This may be thought of as the number of electrons injected by the source point over a cooling timescale. Such a boundary condition invokes the assumption that the contribution from higher-energy electrons cooling into a given energy interval is negligible compared to fresh injections, which we argue is a reasonable approximation due to the underlying power-law nature of the electron injection spectrum.

Figure 1 (right panel) shows that the electron cooling term, b, is dominated by inverse-Compton processes at higher energies, with Coulomb losses becoming important at lower energies. According to Eq. (A.28), the inverse-Compton component depends on the energy density of the ambient radiation field. This is governed by the spatial distribution of photons, including both those from the CMB and stellar radiation fields. Section 3.1.1 indicates that, at the redshifts of interest, low to moderate star-formation rates result in the radiation energy density being easily dominated by the contribution from CMB photons, which thus govern the inverse-Compton cooling of electrons. The CMB is isotropic and homogeneous and, as such, the inverse-Compton cooling function is likewise. However, in systems with a relatively high star-formation rate, the stellar radiation photons can begin to have an important effect on the electron cooling: for example, at a redshift of z = 7, if the star-formation rate of a galaxy with a characteristic size of 1 kpc were to exceed 120 M yr−1, the stellar radiation field would instead provide a dominant contribution to the photon energy density and thus principally govern the inverse-Compton cooling function. Such a stellar radiation field would not be spatially independent, and it is therefore important to assess whether this contribution needs to be taken into account for this work: When considering the sample of galaxies in Tables 1 and 2, some systems have a remarkably high star-formation rate. However, their compactness is also relevant in calculating their stellar radiation energy density: for example, while a 1 kpc radius galaxy with a star-formation rate of 120 M yr−1 would offer a comparable energy density in stellar photons to the CMB at z = 7, another galaxy with the same star-formation rate and at the same redshift, but with a size of 10 kpc would have a negligible energy density compared to that of the CMB. When accounting for their size, we find that none of the star-forming galaxies in our sample exhibit an energy density in their stellar radiation fields which would exceed (or would even be comparable to) that of the CMB (see Table 3).

Table 3.

Peak energy densities of the radiation fields for all 16 of our sample of star-forming galaxies during their starburst phase, compared to the energy density of the CMB radiation field at their observed redshift.

However, these estimates assume that our sample of galaxies exhibit a relatively even distribution of star-formation throughout their volume. There are suggestions that, instead, high-redshift star-formation may be clumpy (e.g. Elmegreen et al. 2009; Förster Schreiber et al. 2011; Murata et al. 2014; Tadaki et al. 2014; Kobayashi et al. 2016). It is therefore useful to consider an extreme case in which star-formation is clumped in the inner 1 kpc of a galaxy, in order to assess the impact this would have on the peak stellar radiation energy density compared to the energy density of the CMB. In the most active case of our sample, SXDF-NB1006-2 with ℛSF ≈ 350 M yr−1, if all the star-formation were concentrated in the central 1 kpc of this system, its stellar radiation energy density would reach around U ≈ 3300 eV cm−3. This is substantially greater than the contribution from the CMB in this case (UCMB ≈ 1200 eV cm−3) and so it cannot immediately be argued that the stellar radiation is unimportant. We find that the overall impact is still not particularly substantial: in Fig. 3, we show various (normalised) electron spectra when injected from a single source and subjected to inverse-Compton (and Coulomb) cooling process6. Of particular note is the lower panel: The spectra shown in black correspond to the case where there is no significant contribution from stars to the ambient radiation field (i.e. it is dominated by the CMB energy density). In blue, the spectra are for the case when accounting for the presence of stellar photons, with their energy density consistent with a 1 kpc radius galaxy violently forming stars at a rate of ℛSF = 350 M yr−1. The different lines in each case relate to different redshifts (and thus CMB energy densities), as noted in the caption. It is clear that, while there is a notable difference at higher energies between the two cases, the underlying power-law nature of the electron spectrum means that the energy density of the electrons in this part of the distribution is relatively small compared to lower energies nearer to the spectral peak. The contribution from the part of the spectrum above 109 eV would account for less than 0.1% of the total electron energy density. The middle plot confirms this assertion: here the black lines correspond to the spectra at 0.2 kpc from the discrete electron source, while those in blue show the result further away, at 0.5 kpc. The different lines in these cases relate to different star-formation rates associated with the stellar radiation field (none for the solid line, then respectively 3.5, 35 and 350 M yr−1 for the dashed, dot-dashed and dotted lines). This shows that CR electrons must traverse through a substantial fraction of their host galaxy before their spectra differ considerably between the case where a strong stellar radiation field is present and the case where inverse-Compton cooling is entirely dominated by the CMB. In our later calculations (see Sect. 5.3), we will consider a distribution of injecting sources – by the time electrons have propagated far enough for any spectral difference between the two cases to be apparent, their contribution to the overall galactic electron spectrum would have been “drowned out” by another source emitting electrons nearby, which would be well within a 0.2 kpc radius. We therefore argue that it is reasonable for us to ignore the contributions from stellar photons in determining the radiative cooling rates of our CR secondary electrons, and duly resort to a CMB-only model which is spatially independent.

thumbnail Fig. 3.

Electron spectra following their injection from a single discrete source (normalised to the injection level at the minimum energy indicated by the dashed red line, ne, 0), accounting for variations in the energy densities of the radiation fields governing inverse-Compton cooling as well as electron diffusion. In all panels, the minimum energy at which CR secondary electrons may be injected is indicated by the dashed red vertical line – this would correspond to the secondaries produced when a 0.28 GeV CR proton primary undergoes a pp interaction. The injection spectral index would be approximately Γ = −2.1, as indicated in blue in the upper panel. Indices of subsequent steady-state spectra are steeper than this due to the energy dependence of the electron cooling rate. Upper panel: spectral variation with distance: dashed light grey line shows the steady-state spectrum expected at the injection point (we note that this is not the injection spectrum); remaining lines (solid black through to solid grey) denote the spectrum at 0.1–0.5 kpc in 0.1 kpc increments when adopting a CMB energy density at z = 7 and no stellar radiation field. The steepening of the high-end of the spectrum with distance is due to the energy-dependence of the radiative cooling function and of electron diffusion (with higher energy particles diffusing and undergoing radiative cooling more quickly). Middle panel: as above, but spectra calculated at 0.2 kpc (black lines) and 0.5 kpc (blue lines) from their injection location, with CMB at z = 7 and the spectral lines indicating evolution according to: no stellar radiation field (solid lines), stellar radiation field consistent with ℛSF = 3.5 M yr−1 (dashed lines), ℛSF = 35 M yr−1 (dot-dashed lines) and ℛSF = 350 M yr−1 (dotted lines) where star-formation is concentrated in a 1 kpc volume. Lower panel: as above, but with all lines calculated at 0.2 kpc from their injection source location, lines in black corresponding to a CMB-only radiation field and lines in blue corresponding to a strong radiation field with ℛSF = 350 M yr−1, shown at different redshifts over a similar range to our high-redshift galaxy sample: z = 6 (solid line), 7 (dashed line), 8 (dot-dashed line) and 9 (dotted line).

Open with DEXTER

Coulomb CR electron cooling depends on the local number density and ionisation state of the interstellar gas. Although the ionisation state is assumed to be fixed, the density profile does exhibit some spatial dependence in our model (see Eq. (7)). The gas density varies significantly on kpc scales, but can be assumed to be comparatively uniform more locally. The length-scale over which electron Coulomb cooling operates is:

(36)

(Owen et al. 2018) where τC is the electron Coulomb cooling timescale. For a 40 MeV secondary CR electron diffusing through an interstellar medium of density nb = 10 cm−3, this gives a length-scale of ℓC  ≈  0.2 kpc, being substantially less than the size of the protogalaxy. As a first approximation, we therefore model the Coulomb cooling to be spatially independent locally. As this is also the case with the radiative cooling function discussed earlier, we take b to be spatially independent in general for each individual injection location in the following treatment (we note that this is only a local approximation, and our model does account for variation in the cooling function at different source injection locations). Thus, using the substitution ζ = bne, Eq. (35) may be expressed as:

(37)

where we have invoked the spatial independence of b. As this PDE has a similar structure to Eq. (24), it may be solved by the same method. This gives:

(38)

Here, the normalisation constant ζ0 = bNe, 0 (from the boundary condition), and we note that τcool(Ee)=Eeb−1. Back-substituting for ζ = bne yields the steady-state result for electrons continuously injected at a single point source:

(39)

5. Cosmic ray heating

Previous studies have considered a range of different mechanisms by which CR heating may operate. These include damping of magnetohydrodynamical (MHD) waves (e.g. Wentzel 1971; Wiener et al. 2013; Pfrommer 2013), direct collisions and ionisations (see Schlickeiser 2002; Schleicher & Beck 2013) as well as Coulomb interactions between charged CRs and a plasma (e.g. Guo & Oh 2008; Ruszkowski et al. 2017). Some of these mechanisms can yield heating rates7 of up to 10−27 erg cm−3 s−1 in appropriate conditions (Wiener et al. 2013; Walker 2016), however the scales over which heating power may be imparted depends on the mechanism at work. Some operate chiefly within (Field et al. 1969; Walker 2016) their host galaxy (including within molecular clouds via ionisation, see Papadopoulos 2010; Papadopoulos et al. 2011; Juvela & Ysard 2011; Galli & Padovani 2015), but others are more effective in heating the medium around it (e.g. Loewenstein et al. 1991; Guo & Oh 2008; Ruszkowski et al. 2017, 2018) or far beyond (e.g. Nath & Biermann 1993; Samui et al. 2005; Sazonov & Sunyaev 2015; Leite et al. 2017). Hadronic interactions are one such proposed mechanism (e.g. Enßlin et al. 2007, 2011; Guo & Oh 2008; Ruszkowski et al. 2017), and could operate within the ISM of a protogalactic host. This mechanism is particularly important in star-forming galaxies due to the relative abundance of CR protons above the hadronic interaction threshold energy, perhaps accounting for as much at 99% of their total energy density (e.g. Benhabiles-Mezhoud et al. 2013).

In general, hadronic heating of a medium by CRs arises after the (pp) interaction of the CR primary injects secondary electrons (see Appendix A). These secondary electrons can then cool and thermalise into the ISM. As discussed previously in Sect. 4.2, the cooling of these electrons can proceed through two channels: via Coulomb interactions with the ambient ISM plasma or, at higher energies, via inverse-Compton interactions with (predominantly) CMB photons. The Coulomb channel leads to direct thermalisation of the electron energy into the ISM. We refer to this hereafter as direct Coulomb heating, or DC heating. Conversely, the inverse-Compton channel does not itself lead to direct energy transfer to the ambient gases: instead (and particularly in high-redshift environments where the CMB is of higher energy density), CR electrons up-scatter CMB photons into (principally) the X-ray band, causing the host galaxy to glow (Schleicher & Beck 2013; Schober et al. 2015). It is thought that the X-ray luminosity in starburst galaxies at high-redshift could become very intense: calculations in Schober et al. (2015) suggest luminosities of above 1039 erg s−1 in the 0.5–10 keV band8 could be achieved by z = 7 for systems with ℛSN ≈ 1.0 yr−1, and scaling in proportion to the SN event rate (or star-formation rate), however Appendix B shows this value could be even higher. These keV X-rays can then propagate and scatter in the ionised or semi-ionised interstellar and circumgalactic medium to drive a heating effect. We refer to this process as Indirect X-ray heating, or IX heating hereafter.

5.1. Direct coulomb heating

The volumetric heating rate due to the DC mechanism at a location r is given by:

(40)

(Owen et al. 2018) where Emin = 0.28 GeV is the minimum energy required for a hadronic interaction (e.g. Kafexhiu et al. 2014), and Emax ≈ 1 PeV is the approximate maximum energy a CR could be accelerated to in environments expected inside the host system, for example SN remnants (see Kotera & Olinto 2011; Schure & Bell 2013; Bell et al. 2013). np(Ep, r) is the local number density of primary CR protons (as calculated in Sect. 4.1), α*(Ep, r) is the absorption coefficient (given by Eq. (A.14)) and fC(Ep) is the Coulomb CR heating efficiency factor, which encodes the fraction of energy transferred from the interacting CR primary proton into this particular heating channel. It is defined as:

(41)

which is the rate of energy transfer from the secondary electrons to the ISM by Coulomb exchanges as a fraction of the total electron cooling rate by all processes, specified at an electron energy Ee. This is related to the CR primary energy by , which suggests that the energy of a typical CR secondary electron is around 40 MeV for a 1 GeV primary proton. The term is introduced to account for electron production multiplicity (denoted ℳ) and energy transfer efficiency to the secondary electrons () – see Appendix A.2 for details. This approach implies that local and instantaneous thermalisation arises. This is justified given that we have previously shown the thermalisation length-scale for 40 MeV electrons (from a 1 GeV primary) is ℓC ≈ 0.2 kpc, being substantially smaller than the 1 kpc characteristic size of the host galaxy. Instantaneous thermalisation can also be argued over galactic timescales: a 0.4 Myr thermalisation timescale is negligible compared to the dynamical timescale of the host system – typically 10s of Myr.

5.2. Indirect X-ray heating

The total inverse-Compton cooling rate of a single electron with energy Ee = γemec2 in a radiation field of energy density Uph is given by

(42)

(Rybicki & Lightman 1979) where σT = 6.65 × 10−25 cm2 is introduced as the Thomson cross-section. This is the rate at which energy is transferred from a CR electron to photons in the radiation field, and so dictates the inverse-Compton power due to the single electron,

(43)

The characteristic energy Eph to which low-energy photons in the radiation field can be up-scattered by the high energy electrons (of characteristic energy ) is

(44)

where may be estimated from the peak of the local electron spectrum (e.g. those shown in Fig. 3). Θ(z) retains its earlier definition of kBT(z)/mec2. This would suggest a characteristic emitted up-scattered photon energy of around 10–100 keV, in the X-ray band. For this energy range, we note that the inverse-Compton Klein-Nishina cross section does not exhibit a strong energy dependence, and is well approximated by the Thomson scattering cross-section to within 10%, |σKN − σT|/σKN <  0.1. This allows us to make a substantial simplification in our calculations hereafter, in that the specific frequency of the up-scattered photons is not of particular importance when considering their subsequent scatterings and interactions. As such, in the following, we will only need to consider the total X-ray inverse-Compton emission energy – not the full spectrum of this radiation.

To calculate the local intensity of inverse-Compton X-rays produced by the high energy CR electrons, we may adopt a similar approach to that which was used in solving the transport equation. In general, the intensity of the X-ray radiation field due to some discrete source at a position ri and spectral luminosity LIX(E) may be written as:

(45)

where r′ retains its earlier definition as the distance between the source position ri and some general location r, and xi is the local gas ionisation fraction such that the number of charged particles in the ambient medium is given by xinb. The exponential term encodes the attenuation effect experienced by the X-ray radiation as it propagates from its source9. Invoking our approximation whereby the spectral distribution of the X-rays is not important, this reduces to:

(46)

where ℒIX is now the total inverse-Compton X-ray luminosity of the source, related to spectral luminosity by

(47)

It is calculated as:

(48)

where 𝒱S retains its earlier definition as a characteristic source size, which contains Ne = ne 𝒱S electrons. The resulting total X-ray intensity ℐIX(r) through a distribution of such sources at a set of positions ri can be found by convolution of the distribution of points with the result given by Eq. (46). From this, the X-ray heating rate would then simply follow as:

(49)

5.3. Computational scheme

Equations (29) and (39) give the steady-state diffusion solution for CR protons and electrons emitted by a single source of some characteristic volume 𝒱S. Further, Eq. (46) gives the X-ray intensity due to a single point-like source. In our calculations, we must extend these single point solutions to model a spatially extended distribution of CR and radiative emission. For CR protons, this distribution would presumably follow the underlying density profile of the protogalaxy, since CRs originate in stellar end products: stars form by gravitational collapse of their underlying matter distribution so it follows that their end products should correlate with this too. We thus sample the NS points throughout a spherical volume by a Monte-Carlo (MC) method, weighted by the density profile. In this case, we specify a maximum radius for the distribution of 2 kpc, beyond which no further points are placed. Each point is assigned a characteristic volume, , where rgal is the characteristic size of the protogalactic system (typically 1 kpc) and NS is the number of individual sources sampled. The resulting distribution of points is then convolved with the steady-state proton distribution for a single source to yield the full galactic CR proton steady-state profile.

The CR electron profile is calculated in a similar manner – however, this must be a two-step process to properly account for the varying density profile of the galaxy in determining the level of electron injection at each position (in accordance with Eqs. (A.25) and (A.14)). The algorithm proceeds as follows: first, sample a points distribution according to the underlying density profile (as was done for the CR protons). Next, calculate the steady-state proton distribution at each point independently and sample a further NS, e points by the same spherical MC sampling method, this time weighted by the single-source steady-state proton solution. The resulting distribution is then convolved with the steady-state electron solution for a single source to give the resulting CR secondary electron profile for a single primary CR injection point. This result is then convolved with the original CR proton source point distribution to give the full galactic profile for electrons. We show the spectrally-integrated result as the CR electron energy density (blue line) together with that for the CR protons (black line) in Fig. 4, when adopting a SN event rate in the protogalaxy of ℛSN = 0.1 yr−1. Here (and hereafter) we have used NS = 10, 000 and Ne = 1000 (the exact electron distribution is less instrumental to our final results), which we find gives an acceptable signal to noise ratio.

thumbnail Fig. 4.

Steady-state protogalactic CR proton (primary) and electron (secondary) profile integrated over all energies for a SN rate ℛSN = 0.1 yr−1. The blue line shows the profile for electrons, while the black line corresponds to protons. The solid part of the black line is the proton distribution within the interstellar medium of the protogalaxy, where the magnetic field model is more robust. Outside the protogalaxy, the magnetic field is less well understood and so our result (the dashed part of the black line) should be regarded as an upper limit.10

Open with DEXTER

The scheme used to calculate the X-ray intensity profile proceeds as follows: a further NS, X points are sampled by the MC method, but now are weighted according to the full galactic CR electron profile. The X-ray intensity profile for each point is determined using Eq. (46), and this is then convolved with the weighted points distribution. The total galactic inverse-Compton X-ray intensity profile then results. We again adopt a value of NIX = 10 000. This, together with the CR electron profile, is used to determine the CR heating power via the two mechanisms according to Eqs. (40) and (49).

5.4. Heating power

In Fig. 5, we show the total power associated with the two CR heating processes: the direct Coulomb (DC) mechanism in black, and the indirect inverse-Compton X-ray (IX) mechanism in red. This uses a model protogalaxy at redshift z = 7, characteristic size 1 kpc, central density nb,0 = 10 cm−3, and considers SN-rates of ℛSN = 0.01, 0.1 and 1.0 yr−1 (broadly corresponding to star-formation rates of ℛSF ≈ 1.6, 16 and 160 M yr−1). We show that a heating rate as high as 10−24 erg cm−3 s−1 could be sustained by the DC channel when ℛSN = 1.0 yr−1 as CRs become trapped by the magnetic field11, with a z = 7 heating IX power of around 10−25 erg cm−3 s−1. This would correspond to a total X-ray luminosity of around 1041 erg s−1 (see also Appendix B). With regard to the IX mechanism, we also indicate shaded regions in Fig. 5 to illustrate the level that would be reached if only stellar photons were up-scattered. Thus, if the IX heating power due to up-scattering of CMB photons were to dominate (the plotted red lines), it must remain above this shaded area.

thumbnail Fig. 5.

CR heating rates of the protogalaxy ISM by (1) direct Coulomb heating, in black, up to 0.6 kpc, and (2) indirect X-ray heating, in red. The solid lines are the results calculated for ℛSN = 1 yr−1, dotted lines for ℛSN = 0.1 yr−1 and dot-dashed for ℛSN = 0.01 yr−1. ℛSN determines the energy budget of the system, and so the CR heating power in both channels scales directly. The shaded regions show the stellar radiation “floor” for the IX process: any evolution that brings the IX line into this area will instead be dominated by the up-scattering of stellar-radiation by CR electrons, preventing the X-ray heating power from falling into this region. The lower region accounts for a stellar radiation field when ℛSN = 0.01 yr−1, the intermediate case above, for ℛSN = 0.1 yr−1, and the highest region arises for ℛSN = 1.0 yr−1. We note that the DC heating rate falls to negligible levels outside the protogalaxy (not shown) due to the magnetic containment of the CRs.

Open with DEXTER

It is clear that up-scattering of CMB photons would be more important than that for stellar photons for any reasonable starburst model at z = 7. At low redshifts the picture is different. Consider Fig. 6, which shows the redshift evolution of the two CR heating effects between z = 0 and z = 12. The starlight provides an effective “floor” to the inverse-Compton emission (and resulting IX heating) at lower redshifts. For the particular model that we have adopted in Figs. 5 and 6, a cut-off redshift of approximately z = 4 arises (we note that this is model-dependent). Below this, the inverse-Compton scattering of stellar-photons dominates over CMB up-scattering, with the final IX heating channel attaining a power roughly 103 times smaller than the corresponding DC channel, making the IX process relatively unimportant by the time such a cut-off redshift has been reached. In general, the IX heating would scale with ℛSN when in the starlight-dominated regime. Thus, we argue the IX mechanism could only operate at a competitive level with the DC mechanism at high-redshift (which itself is largely redshift independent – as indicated by the coinciding black DC heating lines in Fig. 6). Given the scales over which these processes deposit their energy, this is intriguing: the DC mechanism predominantly operates within the ISM of the host galaxy – it is governed by the structure and extent of the magnetic field of the host, and so would drop to negligible levels as low as 10−40 erg cm−3 s−1 (Owen et al. 2018) outside the protogalaxy12. By contrast, the IX mechanism would operate on larger scales – the emitted X-rays can propagate into the circumgalactic medium to drive heating effects around the host, rather than just within it. The redshift evolution of the IX effect thus means that CR heating would be split comparably between the interior and exterior of a galaxy at high-redshift, but would become increasingly focussed onto the ISM over cosmic time. This means that the resulting CR feedback processes at work in and around starbursts could be fundamentally different in the early Universe compared to the present epoch (see also Sect. 6).

thumbnail Fig. 6.

CR heating rates via the direct Coulomb (solid black lines) and indirect X-ray mechanisms (red dashed lines) for a protogalaxy with nb, 0 = 10 cm−3 and ℛSN = 0.1 yr−1 for redshifts z = 0, 2, 4, 6, 8, 10 and 12. The black DC lines (only one line is visible as they overlap) do not show any discernible evolution over this redshift range, but those for the IX effect show significant variation due to the evolving CMB. The central heating rate at z = 12 by both mechanisms is comparable at a level of around 10−25 erg cm−3 s−1. The shaded orange region represents the stellar radiation “floor” for the IX process: any redshift evolution that brings the IX line into this area (i.e. below z ≈ 4) will instead be dominated by the up-scattering of stellar-radiation, preventing the IX power from falling into this region. The DC heating power becomes negligible outside the protogalaxy (not shown).

Open with DEXTER

5.5. Comment on cosmic ray heating in interstellar and circumgalactic structures

At high-redshifts, previous studies have suggested that the ISM of galaxies is clumpy: a single galaxy could be split into just 5–10 large kpc scale clumps (van den Bergh et al. 1996; Cowie et al. 1996; Elmegreen et al. 2004a,b, 2007; Elmegreen & Elmegreen 2005; Bournaud & Elmegreen 2009; Krumholz & Dekel 2010; Wardlow et al. 2017), resulting from the high instability of the cool molecular-gas rich environment (Bournaud & Elmegreen 2009; Krumholz & Dekel 2010) fuelled for example by inflows (Dekel et al. 2009a), with recent observations with ALMA pointing towards a range of different ISM conditions, from relatively smooth to highly clumpy (see Gullberg et al. 2018). Moreover, the circumgalactic medium (CGM) may have a complex structure, comprised of low-density H II regions together with cold dense inflows (e.g. Ribaudo et al. 2011; Sánchez Almeida et al. 2014; Peng et al. 2015; Owen et al. 2019) which themselves could be ionised and beaded with dense, neutral clumps (Fumagalli et al. 2011). As such, more careful consideration of how CRs heat and propagate through different densities of medium is important to properly account for their interactions with their environment and the way in which CR feedback may operate in and around protogalaxies. Consider Fig. 7, which shows how the efficiency of the two CR heating channels depends on ambient density (assuming full ionisation of the medium) for a 40 MeV CR secondary electron (resulting from the hadronic interaction of a 1 GeV proton). This is estimated from the relative timescales of the key loss processes of the CR electrons responsible for driving the CR heating effects, as are indicated by the lines in grey (inverse-Compton, Coulomb losses and synchrotron cooling). Figure 7 points towards a two-phase heating mechanism: CRs effectively heat dense pockets of gas (if exhibiting a reasonably high ionisation fraction) by Coulomb scattering, but are strongly impacted by the galactic magnetic fields and their substructures, effectively containing them to the ISM. Conversely, the IX channel is very efficient in hot diffuse gas – although it can still operate in higher density systems on large scales when the ionisation fraction is sufficiently high. This implies that the DC heating mechanism may be important for quenching star-formation, whereas the IX mechanism could play a role in strangulating the host system by heating inflowing filaments13.

thumbnail Fig. 7.

Fraction of 40 MeV CR secondary electron energy passed to the inverse-Compton X-ray heating channel (red solid line) compared to the direct Coulomb heating channel (black solid line) as a function of density. The left hand y-axis corresponds to the relative fraction, and this indicates a strong preference for the IX heating channel in low density media, or the DC channel in higher density media (if ionised). Timescales of respective processes are shown in grey to illustrate the underlying physics: the solid line gives the total energy loss timescale for CR electrons. This is calculated by combining the contributions from inverse-Compton (dashed horizontal line), Coulomb (dotted line) and free–free losses (dot-dashed line). The respective timescale y-axis is shown to the right of the plot.

Open with DEXTER

5.6. Parametrisation

Our calculated CR heating values can be reasonably scaled for both channels, if conditions are comparable. This allows us to parameterise our simulation results, to make them more easily applicable to our sample of observed high-redshift galaxies (see Sect. 6). The direct Coulomb CR heating rate is principally governed by the input CR power (i.e. via the star-formation rate ℛSF of the system) and the interstellar density, which provides the target for the underlying hadronic interactions. As such, it may reasonably be scaled as

(50)

where reference values are ℋC, 1 = 1.7 × 10−25 erg cm−3 s−1, ℛSF,1 = 16 M yr−1 and n1 = 10 cm−3. Subscript “2” denotes the scaled values. This requires that ℋC, 2 and ℋC, 1 are computed for systems where the ISM is of comparable density such that the CR heating efficiency (see Fig. 7) is largely unchanged – this would usually be the case for our purposes, given that this CR heating mechanism operates within the host galaxy and is unaffected by the variations that might arise in the external circumgalactic environment.

Similarly, the internal IX heating effect can be estimated using a parametrised scaling – however, this requires an additional dependence on redshift. The parametrisation in this case would be:

(51)

where ℋX, 1 = 6.4 × 10−27 erg cm−3 s−1, z1 = 7 and other quantities retain their earlier definitions. We note that, if the interstellar medium of a system were clumpy rather than relatively smooth, this scaling may not be strictly valid. For instance, if the ISM gas was comprised of clumps of around 102 cm−3 interspersed with low-density H II regions of densities around 10−1 cm−3, it would be quite feasible for the inverse-Compton X-ray heating to dominate over the Coulomb CR heating at high-redshift in the majority of the volume of the galaxy due to the relative efficiencies of the two processes in different density media (see Fig. 7). While it is beyond the scope of the current paper, this effect could have important implications for star-formation and is being investigated in a dedicated study.

For external CR heating, it is only necessary to consider IX heating, as the DC process is confined to the host. At distances where the circumgalactic medium density is not greatly influenced by the density profile of the host protogalaxy, the distance from the host is sufficiently large that we may consider it to be a point source. At such distances, the above parametrisation (51) additionally has a distance dependence, however the circumgalactic medium has been assumed to be uniform and so is normalised to a much lower level (consistent with the profile seen for the external heating power in Fig. 5), as given by

(52)

where ℋX, 1 = 7.9 × 10−37 erg cm−3 s−1. We introduce the distance dependence here, with r1 = 100 kpc and r2 as the distance of the required point from the centre of the protogalaxy. Furthermore, we use ne, 1 = xin1 = 10−2 cm−3, as the number density of electrons, with xi being the ionisation fraction. A check on these parameterisations is presented in Appendix C.

6. Thermodynamics in high-redshift protogalaxies

The two CR heating channels impact the dynamics and evolution of their host systems and environments in different ways: the DC channel can raise the temperature of the ISM gases to a level where star-formation cannot continue because of the increased pressure support provided by the heated gas against gravitational collapse (e.g. French et al. 2015) – this may be regarded as a quenching mechanism; the IX channel is more closely attributed to strangulation – cutting off the infall supply of cold gas from its source. If star-formation is driven by the supply of cold gas to the host galaxy through inflowing filaments (see Ribaudo et al. 2011; Sánchez Almeida et al. 2014), external IX heating could be sufficient to heat them as they approach the galaxy. These two processes are illustrated in the context of a galaxy hosting an inflow and active, ongoing star-formation in the schematic in Fig. 8.

thumbnail Fig. 8.

Schematic illustration of the CR heating processes operating in and around a starburst protogalaxy. At A, internal interstellar heating is dominated by the DC channel as CRs and their secondaries are not able to easily escape from the galactic magnetic field of their host system. This increases thermal gas pressure and leads to quenching (French et al. 2015). At B, the IX channel heating mechanism is more important, as this does not rely on CR escape. The X-rays emitted from the galaxy begin to heat and evaporate any cold inflowing filaments, and can do so on a timescale shorter than the inflow time. Such filamentary inflows are required to drive star-formation by providing cold gas which can more easily undergo fragmentation and collapse than ambient interstellar gas (Ribaudo et al. 2011; Sánchez Almeida et al. 2014), so their destruction could strangulate the galaxy (Ribaudo et al. 2011; Sánchez Almeida et al. 2014; Peng et al. 2015).

Open with DEXTER

6.1. Feedback processes

6.1.1. Internal heating and quenching

The condition for quenching can be defined in terms of the virial temperature, being that required to prevent gravitational collapse. It is defined as:

(53)

(Binney & Tremaine 2008), where M is the mass of the system and rh is the characteristic size of the system. For a high-redshift protogalaxy, this would strictly be the mass and size of proto-stellar clouds. However, since observations suggest galactic-wide reservoirs of cool gas being supported against gravitational collapse (French et al. 2015; Rowlands et al. 2015; Alatalo et al. 2016), we may apply the virial theorem to the galaxy as a whole to find the virial temperature for the full system. Thus, any system where the average ISM temperature exceeds this virial temperature would certainly be unable to collapse, and quenching would ensue. With this cut-off temperature Tvir, we can then assign a quenching timescale τQ as the time required for this temperature to be attained by a gas subjected to a given heating rate. However, the full CR heating power only develops after the magnetic field of the host galaxy has saturated. To estimate τQ, we must take this magnetic field evolution into account such that the quenching time would be the magnetic saturation time τmag plus the time taken for entrenched CR heating to raise the ISM temperature above Tvir. We may estimate the magnetic saturation time as:

(54)

(Schober et al. 2013; Owen et al. 2018). The time taken for the virial temperature to be reached due to the DC heating channel then follows as

(55)

where the second term gives the timescale for this temperature to be reached given a CR Coulomb heating power ℋDC, and nb is the local ISM gas density. We apply our model to the observed high-redshift galaxies introduced in Sect. 2 to find the DC heating rate for each system (column ℋDC in Table 4). From this, we use the above analysis to estimate the associated quenching timescale in Myr, shown by column τQ in Table 4.

Table 4.

Derived information for high-redshift quenched post-starburst and starburst systems, as introduced in Tables 1 and 2.

6.1.2. External heating and strangulation

The IX channel operates on larger scales with its effects being principally felt in the circumgalactic medium and the structures therein, where the DC process cannot reach – in particular, inflowing filaments Dekel et al. (2009b), Stewart et al. (2013), Goerdt & Ceverino (2015), Goerdt et al. (2015) which supply the cold gas required to feed star-formation in the host galaxy (Dijkstra et al. 2006; Sancisi et al. 2008; Ribaudo et al. 2011; Sánchez Almeida et al. 2014; Peng et al. 2015) via, for example, cold mode accretion (Birnboim & Dekel 2003; Kereš et al. 2005). If sufficiently strong, heating and evaporation of these filaments can lead to to strangulation, if they are prevented from delivering a supply of cold gas to the ISM of the host (see Fig. 8). The power with which the IX heating channel operates depends both on the local density of the filament being irradiated by X-rays, and the ionisation fraction of that filament. The ionisation fraction is very uncertain, however we adopt a fiducial value of xi = 0.9. This is in line with suggestions from simulation work (and may even be conservative) which suggest a largely ionised structure beaded with denser, cooler and less ionised clumps (see Fumagalli et al. 2011; Faucher-Giguère & Kereš 2011). We may estimate the characteristic density of a filament using mass continuity, if assuming that inflowing filaments are the principal driver of the star-formation in a protogalaxy at high-redshift (e.g. Kereš et al. 2005; Fumagalli et al. 2011; Sánchez Almeida et al. 2014). The star-forming efficiency of the supplied gas from these inflows is around 30% (i.e. ϵ = 0.3) – a value roughly consistent with simulations of high-redshift starbursts14. We set the characteristic inflow velocity to be νLyα ≈ 400 km s−1 (based on the Lyman-α velocity offset observed in Hashimoto et al. 2018) and adopt a total covering fraction of the inflow(s) into the host system of 10% (this is towards the conservative end of the range plausible values, with literature estimates varying between 5 and 40%; see, e.g. Ribaudo et al. 2011; Fumagalli et al. 2011; Faucher-Giguère & Kereš 2011). This yields an effective inflow radius of around 30% of the radius of the host galaxy. From this prescription, we may estimate the number density in the inflow filament as:

(56)

where is its approximate cross sectional area. This gives reasonable values within the expected range, that is nb ≈ 20 cm−3 for a system with ℛSF = 16 M yr−1 (values of the order nb ≈ 10−1 − 101 cm−3 would arguably be sensible; see, e.g. Kereš et al. 2009; Ceverino et al. 2010; Goerdt et al. 2012; Goerdt & Ceverino 2015; Falgarone et al. 2017). To analyse the impact of IX heating on the dynamics, we now require a suitable approximation for the virial temperature of a cold filament: consider the temperature required to prevent a segment of filamentary cloud from collapsing. This will give the minimum temperature above which the cloud begins to disperse (giving a lower limit for the strangulation timescale for a given heating rate). Thus, we consider a spherical region of cloud with radius (this follows from the 10% covering fraction of the filament), and use this to estimate the mass of a filamentary region, . For the example considered above, with ℛSF = 16 M yr−1, this would indicate a filamentary region mass of around 6 × 107 M. By Eq. (53), the corresponding virial temperature is Tvir ≈ 1.2 × 104 K. This is roughly consistent with the expectation for the temperature range of these structures (i.e between 104 − 105 K; see Dekel & Birnboim 2006; Dekel et al. 2009b). The strangulation timescale is calculated in a similar manner to the quenching timescale in the previous section, i.e the sum of the magnetic containment time (required for the X-ray glow to develop) and the heating timescale: τS = τmag + τIX. The IX heating timescale , where nb is the filament number density such that and xinb gives the number density of electrons, as required for X-ray scattering driven heating. We have introduced as the effective IX heating rate throughout the filament. The true IX heating power along the length of the filament would be subject to the inverse square law along its length, meaning its value would vary substantially. However, adopting an effective heating power for the inflow allows us to approximate the timescales associated with the heating process. This is estimated as follows: varies as r−2 away from the galaxy. Thus, consider a point a at which the total (integrated) heating rate up to r = a is equal to that above r = a, being an effective midpoint:

(57)

where rin is the innermost radius of the filament, and R is its outermost extent. Taking R = 50 kpc beyond which inflows are not clearly evident in modelling with simulations (e.g. Dekel et al. 2009b; Stewart et al. 2013; Goerdt & Ceverino 2015), and rin = 1 kpc gives a = 1.96 kpc. For an inflowing filament, the gas would experience equal amounts of heating on either side of a. Here, the proximity of a to rin suggests a rapid evaporation of the inflow filament as it nears the galaxy, with much more gentle heating arising further out. Using the earlier parameterisations (see Sect. 5.6), the characteristic IX heating power at a is estimated by applying an inverse-square law scaling to Eq. (51), yielding15 for a system with ℛSF = 16 M yr−1 and nb = 20 cm−3 at z = 7. This gives a strangulation timescale of around 200 Myr (when also accounting for the magnetic containment time) which, for an inflow traversing a distance of 50 kpc at a velocity of 400 km s−1, compares to an infall timescale of over 100 Myr. The filament’s dynamical timescale is around τdyn ≈ 20 Myr, so this paints a picture of the cold inflowing filament being destroyed by the X-ray glow of the host galaxy, providing a reasonable expectation that the system could be strangulated by the IX heating process. We apply these calculations to determine the associated timescales for the observed high-redshift systems in Table 4.

6.2. Quiescence, cooling and inflow

For star-formation to resume in a galaxy after quenching and strangulation has occurred, either the heated interstellar gases must cool sufficiently, or inflowing cold filaments must be allowed to resume their supply of gas into the system. The timescale for the cooling of the ambient gases can be estimated using Eq. (4). Quenching by CR heating can be sustained for some time after the cessation of burst activity in the host. Typically, CRs will be absorbed or will escape from the galaxy by diffusion. Those absorbed will contribute to the heating effect, otherwise the CRs can escape within a timescale of τdiff ≈ ℓ2/4D ≈ 3 Myr. Their heating cannot be sustained at a level required for quenching beyond a few Myr. If the CRs heat interstellar gas to above the virial temperature of the host, with an ISM density of around 10 cm−3, from Eq. (4) we estimate that the gas should cool and star-formation would resume in less than a Myr. Given that the quiescent periods for all the post-starburst systems in Table 4 are substantially more than a Myr (generally of order 100 Myr), the cooling of interstellar gas would not enable the reinstatement of star-forming activity. Indeed, the extreme processes at play during the starburst episode may have driven away much of the gas reservoir into the galactic halo by for example the action of starburst-driven outflows.

The discrepancy in the timescales indicates that fresh cold gas is required to rejuvenate the galaxy, and this can only be achieved when circumgalactic conditions allow for the reinstatement of cold filamentary inflows which would likely have been dispersed and evaporated by the processes driving strangulation during the starburst episode. This requires the cooling of the heated circumgalactic gas, which can occur after its irradiation by inverse-Compton X-rays from the host galaxy has declined – again, this is driven by CRs in the galaxy, and will reduce a few Myr after the end of the starburst phase, when CRs have had time to diffuse away. Cooling of the circumgalactic gas can arise on a timescale of 100 Myr, if adopting a characteristic circumgalactic medium density of 10−2 cm−3. Such circumgalactic flows may advect recycled circumgalactic gas with them, as well as ambient magnetic fields and cosmic rays. This could have important implications for future star-forming activities, with such components providing enhanced pressure support against gravitational collapse, or higher metallicities in recycled gas promoting the opposite. Fresh primordial gas could return to feeding the galaxy within a comparable timescale, if returning from a distance of 50 kpc with an inflow velocity of around 400 km s−1 as before.

6.3. Discussion

The timescales in Table 4 suggest that CR feedback operates in all these systems via both the direct quenching and strangulation mechanisms, with some competition between them. This is evident from the similarity between the starburst timescales, tSB and the quenching and strangulation timescales, τQ and τS. We consider three scenarios below:

Low , tSB ≤ {τQ,τS}: At low star-formation rates (e.g. those seen in A1689-zD1 and UDF-983-964), both τS and τQ are dominated by the magnetic saturation timescale, τmag. This would arise over 440 Myr for A1689-zD1, and 320 Myr for UDF-983-964 (see Eq. (54)), with quenching only possible at times close to or after this. The low rates of star-formation might be attributed to a lack of cold gas supply to the galaxies (e.g. by inflows) and suggest it would thus be relatively easy to heat the local reservoir of gas, leading to quenching. If suitably cool gas only persisted in pockets throughout these systems, or if gas were of relatively high temperature (near the virial temperature) so as to cool only enough to yield these relatively low star-formation rates, only a moderate amount of CR heating would then be required for quenching. This could plausibly be achieved even before magnetic saturation had completed (with the associated full CR heating power) and would account for the starburst period being shorter even than the estimated quenching timescale.

Moderate , tSB ≃ τQ: At moderate star-formation rates, τmag is less important to tSB (according to Eq. (54)). In these systems, the quenching timescale τQ and starburst time tSB are reasonably consistent, suggesting they are more dependent on their cold gas reservoirs to form stars, with the quenching mechanism described in this paper operating. Before star-formation can be stopped, much of the cold gas reservoir must be heated by CRs after magnetic saturation. Strangulation is less important here but operates over longer timescales to maintain quenching during the quiescent period.

High , tSB ≃ τS: For systems with high star-formation rates, the collapse of interstellar gas reservoirs may be insufficient to maintain star-forming activity at the estimated levels for a prolonged starburst period (with the quenching timescales being very short). Systems like GNS-zD1 and GNS-zD5 would likely instead rely on the injection of cold gas by filamentary inflows to sustain their extreme starburst. But the continued supply of cold gas allows star-formation to persist under strong CR heating of the ISM, until strangulation can cut off this external gas supply too. This would account for the similarity between tSB and τS in these cases.

After the starburst phase, quiescence would persist until the processes driving star-formation resume. The reasonable consistency between the cooling timescale τcool (the timescale over which circumgalactic gas cools and collapses) and the quiescence timescale16tqui ≈ τ* indicates that star-formation resumes only when cold filamentary inflows are reinstated (after the X-ray heating driven by CRs has subsided).

The quiescent timescale, tqui can be regarded as a measure of the degree of star-formation feedback in the system, rather than simply the time elapsed since the end of the starburst episode (see Appendix D). This may be the case if post-starburst galaxies only tend to exceed detection thresholds when they resume star-formation, so the high-redshift post-starburst galaxies we would observe are preferentially those which have emerged from their quiescent phase. Indeed, this would be supported by the two stellar populations found in MACS1149-JD1, where the youngest population is attributed to the observed lower star-formation rate, and is only a few Myr old (see Hashimoto et al. 2018 for details). Consider the relation of tqui on system parameters: in Fig. 9 we show how it relates to physical galaxy parameters (size rgal and stellar mass M*) as well as the dynamical timescale of each system τdyn, the quenching timescale τQ and the strangulation timescale τS. Here (and subsequently), we have highlighted the high-redshift system MACS1149-JD1, but we note that it appears to be unremarkable compared to its lower-redshift companions in our analysis.

thumbnail Fig. 9.

Quiescence timescales against generic physical (upper row) and derived (lower row) galactic parameters. The data of MACS1149-JD1 are highlighted in red (this galaxy is at much higher redshift than the others in the sample, but does not otherwise appear to be unusual). The value of tqui averaged over all data points is indicated by the horizontal dashed grey line. Upper left: plot against estimated galaxy radius. This indicates that the quiescent time scale is also independent of the size of the system, with a weighted Pearson correlation coefficient of 0.15. Error bars in rgal are as quoted in the literature (A1689-zD1 at 50% and MACS1149-JD1 at 71%), otherwise estimated as 30%. Upper right: plot against stellar mass, which shows a positive correlation (coefficient of 0.58). Lower left: plot against galaxy dynamical timescales, which tqui appears to be independent of, with a correlation of 0.06. Lower centre: plot against quenching timescale, τQ, which shows a negative correlation (coefficient of −0.81). This is consistent with tqui being a marker for the strength of feedback – the longer taken to quench star-formation, the weaker the CR driven feedback and the shorter the quiescence timescale. Uncertainties in τQ are unclear and so are not plotted. Lower right: as per the centre panel, but against the strangulation timescale, τS. Again the negative correlation (of coefficient −0.78) is consistent with tqui being influenced by feedback with shorter values for longer strangulation times. As with τQ, uncertainties in τS are not clear so are omitted.

Open with DEXTER

These panels show a positive correlation between stellar mass and tqui (Pearson coefficient of 0.58), as well as negative correlations with τQ (−0.81) and τS (−0.78). The latter two negative correlations are relatively strong and show that, for longer quenching or strangulation timescales (i.e. when the feedback effect would be weakest), tqui is shorter – this is consistent with tqui being an effective tracer for the level of feedback experienced by a system. The positive correlation between tqui and stellar mass M* is also consistent with feedback, given that the systems with higher stellar mass appear to harbour greater star-formation rates during their starburst episode when normalised by the galaxy size: Fig. 10 shows such a dependence of on stellar mass, while the quiescent star-formation rate appears to simply scale with that during the burst (see Fig. 11). There is no notable correlation between tqui and galaxy size rgal or dynamical timescale τdyn but, in a feedback model, no such correlation would be expected between these parameters. This shows there is no strong dependency of tqui on extensive system quantities, only intensive ones relating to the level of feedback, being τQ and τS. Such results are consistent with the quiescence being strongly influenced by feedback from star-formation.

thumbnail Fig. 10.

Star-formation rate surface density plotted against galaxy stellar mass, indicating these quantities are moderately correlated (weighted Pearson coefficient of 0.71 for the post-starburst systems). Upper panel: for the four galaxies observed during their starburst phase (see Table 2) are plotted in grey, with the post-starbursts are shown in black (see Table 1). We plot the MACS1149-JD1 in red (this system is observed at a much higher redshift than the others in the sample, but its behaviour appears to be relatively consistent). Lower panel: as above, with the removal of the starbursts to more clearly show the distribution of the post-starburst galaxies.

Open with DEXTER

thumbnail Fig. 11.

Observed star-formation rate surface densities for the post-starburst sample plotted against stellar mass. This shows a moderate correlation (weighted Pearson coefficient of 0.55). Again, the high redshift system MACS1149-JD1 is highlighted in red, but does not appear to be inconsistent with the rest of the sample.

Open with DEXTER

The additional four systems observed with ongoing starburst activity are also shown in grey in the upper panel of Fig. 10. These weakly indicate some inverse trend with stellar mass, as would result from the depletion of gas supplies as strangulation takes hold – however it is difficult to draw strong conclusions with this small data set and large uncertainties.

7. Summary

In this study, we have shown the power that can be attained by CR heating in high-redshift starburst galaxies, and how the underlying energy-transfer process in CR heating is split according to two distinctive channels: one in which the secondary electrons arising from the hadronic interactions of CR protons undergo Coulomb scattering to thermalise into their medium over a length-scale of around 0.2 kpc, and the other in which the radiative emission from these secondary electrons (namely inverse-Compton scattering, predominantly with CMB photons) drives a strong X-ray glow which itself can cause a heating effect of the host galaxy and its circumgalactic environment. The second of these channels is particularly important in high-redshift environments, when the energy density of the CMB was substantially higher than it is in the present Universe.

We have demonstrated how these two CR heating channels can result in two distinct forms of CR-induced feedback on star-formation. The direct Coulomb thermalisation process can readily heat the ISM of the host galaxy, up to levels around 10−25 erg cm−3 s−1 in a system exhibiting SN events at a rate of ℛSN = 0.1 yr−1 – this could exceed more conventional feedback mechanisms, for example stellar photons and diffuse X-rays from the hot ISM gases (see Owen et al. 2018). We find that, in a sample of high-redshift starburst and post-starburst galaxies, such heating would be able to increase the ISM temperature to above its virial temperature (and hence decisively quench star-formation in the host) within a few tens of Myr (when star-formation rates during the starburst phase of the system are above a few tens M yr−1). However, this Coulomb thermalisation process is strongly confined by the magnetic fields of the host galaxy, and is relatively ineffective in the wider circumgalactic environment. By contrast, the heating due to inverse-Compton X-rays (referred to as indirect X-ray heating) can operate both within and outside the host system. At z = 12, corresponding to the earliest times at which star-formation has been inferred (see Hashimoto et al. 2018), the CR heating via this channel can attain comparable levels to the direct Coulomb thermalisation process within the ISM. After this (i.e. at lower redshifts), it becomes less important within the galaxy, but remains important outside, where X-ray heating can halt the cold inflowing filaments required to feed star-formation in the most active systems. This can lead to strangulation (e.g. Peng et al. 2015) – by effectively cutting off the supply of cold gas in inflows, this indirect CR X-ray heating can bring star-formation to a halt within just a few hundred Myr. Both of these processes are likely to operate in any given star-forming galaxy at high-redshift. We show that, for more intensely star-forming galaxies, timescales would indicate that CR feedback arises preferentially through strangulation. In less active systems, quenching arises first (with strangulation operating over a longer period of time). We note that, at lower redshifts, the indirect X-ray mechanism does not operate as effectively as it does in the early Universe. This is because it becomes dependent on the up-scattering of stellar photons rather than the CMB, and this yields a substantially diminished X-ray glow compared to early-Universe counterparts. This promotes the idea that high-redshift galaxies are very different environments to their low-redshift descendents and, on galactic scales, gives some insight into how these systems may have evolved.


1

The dynamical mass may be estimated from observations using the velocity dispersion of the line profile. By the virial theorem, Mdyn = f(FWHM)2rgal/G (e.g. Gnerucci et al. 2011), where f ≈ 2.25 (Binney & Tremaine 2008; Förster Schreiber et al. 2009; Epinat et al. 2009), rgal is the approximate radius of the galaxy (inferred from the half-light radius) and G is the Newtonian gravitational constant. The FWHM is the full width of half-maximum of the line profile, being equal to the velocity dispersion of the system. This accounts for the difference in relative velocity due to the circular orbiting motion of the emitting gases around the galaxy. The dynamical timescale is (e.g. Binney & Tremaine 2008).

2

See Fields et al. (2001), Strong et al. (2010), Lemoine-Goumard et al. (2012), Caprioli (2012), Dermer & Powale (2013), Morlino & Caprioli (2012), Wang & Fields (2018) for studies regarding this parameter and a range of possible, reasonable values.

3

This value depends on the SN type as well as the way it interacts with its specific environment (see e.g. Iwamoto & Kunugise 2006). In the case of core-collapse SNe, most of the energy is carried away by neutrinos with perhaps a fraction as low as 0.001, but more likely around 0.01, being retained by the system to, for example, accelerate CRs and drive shocks (see reviews, e.g. Smartt 2009; Janka 2012).

4

Previous studies (see, e.g. Owen et al. 2019) indicate that advection of CRs in outflows does not substantially reduce their abundance within a protogalaxy ISM (by only around 10%). Internal interstellar winds and flows will have a more local impact, in that they push CRs from one part of the galaxy to another – but this arises on a sub-galactic scale, with the galactic-scale picture being unaffected.

5

While the ISM of high redshift protogalaxies cannot be directly observed, we argue that it is not unreasonable to expect that processes driving turbulence within them are not greatly different to those seen in the local Universe. As such, we argue that the diffusion coefficient would likely be similar to that in the Milky Way and, in the absence of more detailed studies or observations, adopting an alternative prescription would not imply more correct physics.

6

Other relevant cooling processes are also applied in calculating these spectra, but their effect is not of great importance.

7

When adopting ℛSN = 1.0 yr−1.

8

This is the energy range which the Chandra X-ray observatory is sensitive to (see http://chandra.harvard.edu) and was chosen in the Schober et al. (2015) study to make observational predictions for the X-ray flux from starbursts to diagnose the presence of CRs and magnetic fields.

9

We note that, even in central galaxy conditions of np = 10 cm−3, the mean free path of a keV photon is around 100 kpc, resulting in the attenuation term being very close to unity with its impact being negligible on our analysis.

10

The results shown in this plot are slightly different to those for CR protons presented in Fig. 5 of Owen et al. (2018), but the two results are not directly comparable. This is because (1) here, we include the attenuation of CRs in our treatment, which was not accounted for in the earlier result (instead this was accounted for later in their simulation pipeline and applied to their final CR heating calculation); (2) the result here integrates over the full energy spectrum of the CR protons and electrons, whereas the result in the previous study shows the CR proton profile when setting the proton distribution to its mean energy of 2.5 GeV; (3) this result is computed using ℛSN = 0.1 yr−1, whereas the previous study showed the case for ℛSN = 1.0 yr−1.

11

As the galactic magnetic field grows in strength, charged CRs can no longer freely stream from their origin. Instead they predominantly diffuse, with their effective diffusion velocity being substantially lower than their free-streaming velocity (see Owen et al. 2018 for details).

12

Since the interfacing region between the ISM and circumgalactic medium (CGM) would be entrained by complicated magnetic field structures – the modelling of which remains beyond the scope of the current paper – we only show the DC channel heating effect up to 0.6 kpc in our results, well within the ISM of the host.

13

We note that regions of increased density would likely harbour stronger magnetic fields than the ISM average. This is seen in Milky Way observations (e.g. Crutcher et al. 2010), and would presumably also arise in gravitationally collapsed clumps and clouds in high-redshift galaxies. While increased magnetic energy densities would yield a larger fraction of CR secondary energy being lost to synchrotron processes, magnetic field strengths are unlikely to be sufficiently strong to compete with inverse-Compton losses until cloud densities exceed nb ≥ 106 cm−3 (Crutcher et al. 2010; Owen et al. 2018). Densities of this level only persist in a very small volume fraction of clouds on sub-pc scales (Draine 2011), and would not be likely to bear a significant influence on the global astrophysics of the host galaxy considered in this study. Moreover, under such conditions, Coulomb scattering would easily dominate over radiative cooling such that the effect of stronger magnetic fields on CR electron losses would be inconsequential.

14

See Sun & Furlanetto (2016), which finds star-formation efficiencies of tens of per cent for halo masses of above 1011 M, corresponding to stellar masses of around 109 M, comparable to the systems used in the sample in this work – see also Behroozi & Silk (2015), and Meier et al. (2002), Turner et al. (2015) which find similar efficiencies in the young starburst NGC 5253 that may also be indicative of the levels to be expected in the types of systems we are modelling.

15

Given the continuous density profile along the filament and into the ISM, a scaling of the internal IX heating parametrisation was considered more meaningful than a re-scaling of the external IX heating power, which would be substantially lower due to the foreground density profile of the galaxy no longer being important at distances above 20 kpc.

16

This lower limit for tqui is comparable to the stellar population age τ* as the age of the stars would correspond to their age by the end of the starburst plus the quiescent time. A more appropriate value for this quantity could be estimated by modelling the star-formation rate within the starburst episode, but this is beyond the scope of this paper and deserving of a dedicated study.

17

These processes arise via the formation and decay of a μ±.

18

– i.e. the energy required for the production of the lowest energy pion (the π0).

19

Dermer & Menon (2009) suggests a lower value of kγe = 2, but we note that this was in the case of ϵr ≳ 40.

20

We note that Eq. (A.11) requires a multiplicative factor if the radiation field is diluted, for example as would be the case with light emitted from a distribution of stars. For the CMB, the field is not diluted and the equation holds without modification.

22

This would not account for the similar behaviour seen in Fig. D.1 where star-formation rates are determined independently. However it is unclear whether this is an indication of a physical process or simply due to low statistics, with some data points aligning to give a false correlation.

Acknowledgments

ERO is supported by a UK Science and Technology Facilities Council PhD Studentship and the International Exchange Studentship of National Tsing Hua University (NTHU), where his visit was kindly hosted by Prof. Shih-Ping Lai. ERO also acknowledges support from a Royal Astronomical Society travel grant and the Institute of Physics C R Barber trust fund and Research Student Conference Fund which supported travel to conferences and meetings where discussion helped to inspire and inform this work. NK’s visit to the Mullard Space Science Laboratory (MSSL) was supported by a Kyoto Sangyo University visiting summer studentship. We thank B. P. Brian Yu (MSSL) and Qin Han (Nanjing University) for critical review of the manuscript, and for constructive comments. We also thank the anonymous referee for their comments and feedback.

References

  1. Abdo, A. A., Ackermann, M., Ajello, M., et al. 2010a, ApJ, 709, L152 [NASA ADS] [CrossRef] [Google Scholar]
  2. Abdo, A. A., Ackermann, M., Ajello, M., et al. 2010b, A&A, 523, L2 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  3. Abdo, A. A., Ackermann, M., Ajello, M., et al. 2010c, A&A, 523, A46 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  4. Abdo, A. A., Ackermann, M., Ajello, M., et al. 2010d, A&A, 512, A7 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  5. Abe, F., Amidei, D., Apollinari, G., et al. 1990, Phys. Rev. D, 41, 2330 [NASA ADS] [CrossRef] [Google Scholar]
  6. Acero, F., Aharonian, F., Akhperjanian, A. G., et al. 2009, Science, 326, 1080 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  7. Ackermann, M., Ajello, M., Allafort, A., et al. 2013, Science, 339, 807 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  8. Aharonian, F. A., & Atoyan, A. M. 2000, A&A, 362, 937 [NASA ADS] [Google Scholar]
  9. Aharonian, F., Bykov, A., Parizot, E., Ptuskin, V., & Watson, A. 2012, Space Sci. Rev., 166, 97 [NASA ADS] [CrossRef] [Google Scholar]
  10. Alatalo, K., Lisenfeld, U., Lanz, L., et al. 2016, ApJ, 827, 106 [NASA ADS] [CrossRef] [Google Scholar]
  11. Albajar, C., Albrow, M. G., Allkofer, O. C., et al. 1990, Nucl. Phys. B, 335, 261 [NASA ADS] [CrossRef] [Google Scholar]
  12. Albini, E., Capiluppi, P., Giacomelli, G., & Rossi, A. M. 1976, Nuovo Cimento A Serie, 32, 101 [NASA ADS] [CrossRef] [Google Scholar]
  13. Alexopoulos, T., Anderson, E. W., Biswas, N. N., et al. 1998, Phys. Lett. B, 435, 453 [NASA ADS] [CrossRef] [Google Scholar]
  14. Allard, D., Parizot, E., & Olinto, A. 2007, Astropart. Phys., 27, 61 [NASA ADS] [CrossRef] [Google Scholar]
  15. Almeida, S. P., Rushbrooke, J. G., Scharenguivel, J. H., et al. 1968, Phys. Rev. D, 174, 1638 [CrossRef] [Google Scholar]
  16. Alner, G. J., Alpgård, K., Ansorge, R. E., et al. 1984, Phys. Lett. B, 138, 304 [NASA ADS] [CrossRef] [Google Scholar]
  17. Alner, G. J., Alpgård, K., Anderer, P., et al. 1985, Phys. Lett. B, 160, 193 [NASA ADS] [CrossRef] [Google Scholar]
  18. Ansorge, R. E., Åsman, B., Burow, L., et al. 1989, Z. Phys. C, 43, 357 [NASA ADS] [CrossRef] [Google Scholar]
  19. Aoyama, S., Hou, K.-C., Shimizu, I., et al. 2017, MNRAS, 466, 105 [NASA ADS] [CrossRef] [Google Scholar]
  20. Arnison, G., Astbury, A., Aubert, B., et al. 1983, Phys. Lett. B, 123, 108 [NASA ADS] [CrossRef] [Google Scholar]
  21. Atoyan, A. M., & Dermer, C. D. 2003, ApJ, 586, 79 [NASA ADS] [CrossRef] [Google Scholar]
  22. Baldini, A., Flaminio, V., Moorhead, W., & Morrison, D. 1987, Total Cross-Sections for Reactions of High Energy Particles, Landolt-Börnstein: Numerical Data and Functional Relationships in Science and Technology - New Series/Elementary Particles, Nuclei and Atoms (Springer) [Google Scholar]
  23. Balsara, D. S., Kim, J., Low, M.-M. M., & Mathews, G. J. 2004, ApJ, 617, 339 [NASA ADS] [CrossRef] [Google Scholar]
  24. Beck, A. M., Lesch, H., Dolag, K., et al. 2012, MNRAS, 422, 2152 [NASA ADS] [CrossRef] [Google Scholar]
  25. Begelman, M. C. 1995, in The Physics of the Interstellar Medium and Intergalactic Medium, eds. A. Ferrara, C. F. C. Heiles, & P. R. Shapiro, ASP Conf. Ser., 80, 545 [NASA ADS] [Google Scholar]
  26. Behroozi, P. S., & Silk, J. 2015, ApJ, 799, 32 [NASA ADS] [CrossRef] [Google Scholar]
  27. Bell, A. R. 1978, MNRAS, 182, 147 [NASA ADS] [CrossRef] [Google Scholar]
  28. Bell, A. R., Schure, K. M., Reville, B., & Giacinti, G. 2013, MNRAS, 431, 415 [NASA ADS] [CrossRef] [Google Scholar]
  29. Benhabiles-Mezhoud, H., Kiener, J., Tatischeff, V., & Strong, A. W. 2013, ApJ, 763, 98 [NASA ADS] [CrossRef] [Google Scholar]
  30. Berezhko, E. G., & Ellison, D. C. 1999, ApJ, 526, 385 [NASA ADS] [CrossRef] [Google Scholar]
  31. Berezinsky, V. S., & Gazizov, A. Z. 1993, Phys. Rev. D, 47, 4206 [NASA ADS] [CrossRef] [Google Scholar]
  32. Berezinsky, V. S., & Grigor’eva, S. I. 1988, A&A, 199, 1 [NASA ADS] [Google Scholar]
  33. Berezinskii, V. S., Bulanov, S. V., Dogiel, V. A., & Ptuskin, V. S. 1990, Astrophysics of Cosmic Rays (Amsterdam: North-Holland) [Google Scholar]
  34. Berezinsky, V., Gazizov, A., & Grigorieva, S. 2006, Phys. Rev. D, 74, 043005 [NASA ADS] [CrossRef] [Google Scholar]
  35. Bernet, M. L., Miniati, F., Lilly, S. J., Kronberg, P. P., & Dessauges-Zavadsky, M. 2008, Nature, 454, 302 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  36. Bethe, H., & Heitler, W. 1934, Proc. R. Soc. A, 146, 83 [NASA ADS] [CrossRef] [Google Scholar]
  37. Bethe, H. A., & Maximon, L. C. 1954, Phys. Rev., 93, 768 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
  38. Bhattacharjee, P., & Sigl, G. 2000, Phys. Rep., 327, 109 [NASA ADS] [CrossRef] [Google Scholar]
  39. Bianchi, S., & Schneider, R. 2007, MNRAS, 378, 973 [NASA ADS] [CrossRef] [Google Scholar]
  40. Binney, J., & Tremaine, S. 2008, Galactic Dynamics, 2nd edn. (Princeton University Press) [Google Scholar]
  41. Birnboim, Y., & Dekel, A. 2003, MNRAS, 345, 349 [NASA ADS] [CrossRef] [Google Scholar]
  42. Blattnig, S. R., Swaminathan, S. R., Kruger, A. T., Ngom, M., & Norbury, J. W. 2000, Phys. Rev. D, 62, 094030 [NASA ADS] [CrossRef] [Google Scholar]
  43. Blumenthal, G. R. 1970, Phys. Rev. D, 1, 1596 [NASA ADS] [CrossRef] [Google Scholar]
  44. Blumenthal, G. R., & Gould, R. J. 1970, Rev. Mod. Phys., 42, 237 [NASA ADS] [CrossRef] [Google Scholar]
  45. Borsellino, A. 1947, Il Nuovo Cimento (1943–1954), 4, 112 [NASA ADS] [CrossRef] [Google Scholar]
  46. Boselli, A., & Gavazzi, G. 2006, PASP, 118, 517 [NASA ADS] [CrossRef] [Google Scholar]
  47. Bournaud, F., & Elmegreen, B. G. 2009, ApJ, 694, L158 [NASA ADS] [CrossRef] [Google Scholar]
  48. Bouwens, R. J., & Illingworth, G. D. 2006, Nature, 443, 189 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  49. Bouwens, R. J., Thompson, R. I., Illingworth, G. D., et al. 2004, ApJ, 616, L79 [NASA ADS] [CrossRef] [Google Scholar]
  50. Bradley, L. D., Bouwens, R. J., Ford, H. C., et al. 2008, ApJ, 678, 647 [NASA ADS] [CrossRef] [Google Scholar]
  51. Brandenburg, A., & Lazarian, A. 2013, Space Sci. Rev., 178, 163 [NASA ADS] [CrossRef] [Google Scholar]
  52. Breakstone, A., Campanini, R., Crawley, H. B., et al. 1984, Phys. Rev. D, 30, 528 [NASA ADS] [CrossRef] [Google Scholar]
  53. Caprioli, D. 2012, J. Cosmol. Astropart. Phys., 7, 038 [NASA ADS] [CrossRef] [Google Scholar]
  54. Ceverino, D., Dekel, A., & Bournaud, F. 2010, MNRAS, 404, 2151 [NASA ADS] [Google Scholar]
  55. Couch, W. J., & Sharples, R. M. 1987, MNRAS, 229, 423 [NASA ADS] [CrossRef] [Google Scholar]
  56. Cowie, L. L., & Songaila, A. 1977, Nature, 266, 501 [NASA ADS] [CrossRef] [Google Scholar]
  57. Cowie, L. L., Songaila, A., Hu, E. M., & Cohen, J. G. 1996, AJ, 112, 839 [NASA ADS] [CrossRef] [Google Scholar]
  58. Crutcher, R. M., Wandelt, B., Heiles, C., Falgarone, E., & Troland, T. H. 2010, ApJ, 725, 466 [NASA ADS] [CrossRef] [Google Scholar]
  59. Dar, A., & de Rújula, A. 2008, Phys. Rep., 466, 179 [NASA ADS] [CrossRef] [Google Scholar]
  60. de Cea del Pozo, E., Torres, D. F., & Rodriguez Marrero, A. Y. 2009, ApJ, 698, 1054 [NASA ADS] [CrossRef] [Google Scholar]
  61. Dehnen, W. 1993, MNRAS, 265, 250 [NASA ADS] [CrossRef] [Google Scholar]
  62. Dekel, A., & Birnboim, Y. 2006, MNRAS, 368, 2 [NASA ADS] [CrossRef] [Google Scholar]
  63. Dekel, A., Sari, R., & Ceverino, D. 2009a, ApJ, 703, 785 [NASA ADS] [CrossRef] [Google Scholar]
  64. Dekel, A., Birnboim, Y., Engel, G., et al. 2009b, Nature, 457, 451 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  65. Dermer, C. D., & Menon, G. 2009, High Energy Radiation from Black Holes: Gamma rays, Cosmic rays, and Neutrinos, Princeton Series in Astrophysics (Princeton, NJ: Princeton Univ. Press) [Google Scholar]
  66. Dermer, C. D., & Powale, G. 2013, A&A, 553, A34 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  67. Di Fazio, A., Occhionero, F., & Vagnetti, F. 1979, A&A, 72, 204 [NASA ADS] [Google Scholar]
  68. Dijkstra, M., Haiman, Z., & Spaans, M. 2006, ApJ, 649, 37 [NASA ADS] [CrossRef] [Google Scholar]
  69. Domingo-Santamaría, E., & Torres, D. F. 2005, A&A, 444, 403 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  70. Draine, B. T. 2011, Physics of the Interstellar and Intergalactic Medium (Princeton University Press) [Google Scholar]
  71. Dressler, A., & Gunn, J. E. 1983, ApJ, 270, 7 [NASA ADS] [CrossRef] [Google Scholar]
  72. Elmegreen, B. G., & Elmegreen, D. M. 2005, ApJ, 627, 632 [NASA ADS] [CrossRef] [Google Scholar]
  73. Elmegreen, D. M., Elmegreen, B. G., & Hirst, A. C. 2004a, ApJ, 604, L21 [NASA ADS] [CrossRef] [Google Scholar]
  74. Elmegreen, D. M., Elmegreen, B. G., & Sheets, C. M. 2004b, ApJ, 603, 74 [NASA ADS] [CrossRef] [Google Scholar]
  75. Elmegreen, D. M., Elmegreen, B. G., Ravindranath, S., & Coe, D. A. 2007, ApJ, 658, 763 [NASA ADS] [CrossRef] [Google Scholar]
  76. Elmegreen, D. M., Elmegreen, B. G., Marcus, M. T., et al. 2009, ApJ, 701, 306 [NASA ADS] [CrossRef] [Google Scholar]
  77. Enßlin, T. A., Pfrommer, C., Springel, V., & Jubelgas, M. 2007, A&A, 473, 41 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  78. Enßlin, T., Pfrommer, C., Miniati, F., & Subramanian, K. 2011, A&A, 527, A99 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  79. Epinat, B., Contini, T., Le Fèvre, O., et al. 2009, A&A, 504, 789 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  80. Falgarone, E., Zwaan, M. A., Godard, B., et al. 2017, Nature, 548, 430 [NASA ADS] [CrossRef] [Google Scholar]
  81. Farber, R., Ruszkowski, M., Yang, H.-Y. K., & Zweibel, E. G. 2018, ApJ, 856, 112 [NASA ADS] [CrossRef] [Google Scholar]
  82. Faucher-Giguère, C.-A., & Kereš, D. 2011, MNRAS, 412, L118 [NASA ADS] [CrossRef] [Google Scholar]
  83. Federrath, C., Chabrier, G., Schober, J., et al. 2011, Phys. Rev. Lett., 107, 114504 [NASA ADS] [CrossRef] [Google Scholar]
  84. Fermi, E. 1949, Phys. Rev., 75, 1169 [NASA ADS] [CrossRef] [Google Scholar]
  85. Ferrara, A., Viti, S., & Ceccarelli, C. 2016, MNRAS, 463, L112 [NASA ADS] [CrossRef] [Google Scholar]
  86. Field, G. B., Goldsmith, D. W., & Habing, H. J. 1969, ApJ, 155, L149 [NASA ADS] [CrossRef] [Google Scholar]
  87. Fields, B. D., Olive, K. A., Cassé, M., & Vangioni-Flam, E. 2001, A&A, 370, 623 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  88. Fiete Grosse-Oetringhaus, J., & Reygers, K. 2010, J. Phys. G Nucl. Phys., 37, 083001 [NASA ADS] [CrossRef] [Google Scholar]
  89. Förster Schreiber, N. M., Genzel, R., Bouché, N., et al. 2009, ApJ, 706, 1364 [NASA ADS] [CrossRef] [Google Scholar]
  90. Förster Schreiber, N. M., Shapley, A. E., Genzel, R., et al. 2011, ApJ, 739, 45 [NASA ADS] [CrossRef] [Google Scholar]
  91. French, K. D., Yang, Y., Zabludoff, A., et al. 2015, ApJ, 801, 1 [NASA ADS] [CrossRef] [Google Scholar]
  92. French, K. D., Zabludoff, A. I., Yoon, I., et al. 2018, ApJ, 861, 123 [NASA ADS] [CrossRef] [Google Scholar]
  93. Fryer, C. L. 1999, ApJ, 522, 413 [NASA ADS] [CrossRef] [Google Scholar]
  94. Fumagalli, M., Prochaska, J. X., Kasen, D., et al. 2011, MNRAS, 418, 1796 [NASA ADS] [CrossRef] [Google Scholar]
  95. Gaggero, D. 2012, Cosmic Ray Diffusion in the Galaxy and Diffuse Gamma Emission, Springer Theses (Springer: Berlin) [CrossRef] [Google Scholar]
  96. Galli, D., & Padovani, M. 2015, ArXiv e-prints [arXiv:1502.03380] [Google Scholar]
  97. Gnerucci, A., Marconi, A., Cresci, G., et al. 2011, A&A, 533, A124 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  98. Goerdt, T., & Ceverino, D. 2015, MNRAS, 450, 3359 [NASA ADS] [CrossRef] [Google Scholar]
  99. Goerdt, T., Dekel, A., Sternberg, A., Gnat, O., & Ceverino, D. 2012, MNRAS, 424, 2292 [NASA ADS] [CrossRef] [Google Scholar]
  100. Goerdt, T., Ceverino, D., Dekel, A., & Teyssier, R. 2015, MNRAS, 454, 637 [NASA ADS] [CrossRef] [Google Scholar]
  101. Goldreich, P., & Sridhar, S. 1995, ApJ, 438, 763 [NASA ADS] [CrossRef] [Google Scholar]
  102. González, V., Labbé, I., Bouwens, R. J., et al. 2010, ApJ, 713, 115 [NASA ADS] [CrossRef] [Google Scholar]
  103. Goto, T. 2006, MNRAS, 369, 1765 [NASA ADS] [CrossRef] [Google Scholar]
  104. Gould, R. J. 1975, ApJ, 196, 689 [NASA ADS] [CrossRef] [Google Scholar]
  105. Gould, R. J., & Burbidge, G. R. 1965, Annu. Astrophys., 28, 171 [NASA ADS] [Google Scholar]
  106. Grazian, A., Castellano, M., Fontana, A., et al. 2012, A&A, 547, A51 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  107. Griffin, R. D., Dai, X., & Thompson, T. A. 2016, ApJ, 823, L17 [NASA ADS] [CrossRef] [Google Scholar]
  108. Gullberg, B., Swinbank, A. M., Smail, I., et al. 2018, ApJ, 859, 12 [NASA ADS] [CrossRef] [Google Scholar]
  109. Gunn, J. E., & Gott, III, J. R. 1972, ApJ, 176, 1 [NASA ADS] [CrossRef] [Google Scholar]
  110. Guo, F., & Oh, S. P. 2008, MNRAS, 384, 251 [NASA ADS] [CrossRef] [Google Scholar]
  111. Hammond, A. M., Robishaw, T., & Gaensler, B. M. 2012, ArXiv eprints [arXiv:1209.1438] [Google Scholar]
  112. Hashimoto, T., Laporte, N., Mawatari, K., et al. 2018, Nature, 557, 392 [NASA ADS] [CrossRef] [Google Scholar]
  113. Haug, E. 1981, Z. Naturforsch, 36, 413 [NASA ADS] [Google Scholar]
  114. Hayashida, M., Stawarz, Ł., Cheung, C. C., et al. 2013, ApJ, 779, 131 [NASA ADS] [CrossRef] [Google Scholar]
  115. Heger, A., Fryer, C. L., Woosley, S. E., Langer, N., & Hartmann, D. H. 2003, ApJ, 591, 288 [NASA ADS] [CrossRef] [Google Scholar]
  116. Hillas, A. M. 1984, ARA&A, 22, 425 [NASA ADS] [CrossRef] [Google Scholar]
  117. Inoue, A. K., Tamura, Y., Matsuo, H., et al. 2016, Science, 352, 1559 [NASA ADS] [CrossRef] [Google Scholar]
  118. Iwamoto, K., & Kunugise, T. 2006, AIP Conf. Proc., 847, 406 [NASA ADS] [CrossRef] [Google Scholar]
  119. Janka, H.-T. 2012, Annu. Rev. Nucl. Part. Sci., 62, 407 [NASA ADS] [CrossRef] [Google Scholar]
  120. Joseph, J., & Rohrlich, F. 1958, Rev. Mod. Phys., 30, 354 [NASA ADS] [CrossRef] [Google Scholar]
  121. Jost, R., Luttinger, J. M., & Slotnick, M. 1950, Phys. Rev., 80, 189 [NASA ADS] [CrossRef] [Google Scholar]
  122. Juvela, M., & Ysard, N. 2011, ApJ, 739, 63 [NASA ADS] [CrossRef] [Google Scholar]
  123. Kafexhiu, E., Aharonian, F., Taylor, A. M., & Vila, G. S. 2014, Phys. Rev. D, 90, 123014 [NASA ADS] [CrossRef] [Google Scholar]
  124. Kamae, T., Karlsson, N., Mizuno, T., Abe, T., & Koi, T. 2006, ApJ, 647, 692 [NASA ADS] [CrossRef] [Google Scholar]
  125. Karlsson, N. 2008, AIP Conf. Proc., 1085, 561 [NASA ADS] [CrossRef] [Google Scholar]
  126. Kaviraj, S., Kirkby, L. A., Silk, J., & Sarzi, M. 2007, MNRAS, 382, 960 [NASA ADS] [CrossRef] [Google Scholar]
  127. Kelner, S. R., Aharonian, F. A., & Bugayov, V. V. 2006, Phys. Rev. D, 74, 034018 [NASA ADS] [CrossRef] [Google Scholar]
  128. Kereš, D., Katz, N., Weinberg, D. H., & Davé, R. 2005, MNRAS, 363, 2 [NASA ADS] [CrossRef] [Google Scholar]
  129. Kereš, D., Katz, N., Fardal, M., Davé, R., & Weinberg, D. H. 2009, MNRAS, 395, 160 [NASA ADS] [CrossRef] [Google Scholar]
  130. Klein, S. R. 2006, Rad. Phys. Chem., 75, 696 [NASA ADS] [CrossRef] [Google Scholar]
  131. Knudsen, K. K., Richard, J., Kneib, J.-P., et al. 2016, MNRAS, 462, L6 [NASA ADS] [CrossRef] [Google Scholar]
  132. Knudsen, K. K., Watson, D., Frayer, D., et al. 2017, MNRAS, 466, 138 [NASA ADS] [CrossRef] [Google Scholar]
  133. Kobayashi, M. A. R., Murata, K. L., Koekemoer, A. M., et al. 2016, ApJ, 819, 25 [NASA ADS] [CrossRef] [Google Scholar]
  134. Kotera, K., & Olinto, A. V. 2011, ARA&A, 49, 119 [NASA ADS] [CrossRef] [Google Scholar]
  135. Kotera, K., Allard, D., & Olinto, A. V. 2010, J. Cosmol. Astropart. Phys., 10, 013 [NASA ADS] [CrossRef] [Google Scholar]
  136. Krumholz, M. R., & Dekel, A. 2010, MNRAS, 406, 112 [NASA ADS] [CrossRef] [Google Scholar]
  137. Lacki, B. C., & Beck, R. 2013, MNRAS, 430, 3171 [NASA ADS] [CrossRef] [Google Scholar]
  138. Lacki, B. C., & Thompson, T. A. 2012, AIP Conf. Proc., 1505, 56 [NASA ADS] [CrossRef] [Google Scholar]
  139. Lacki, B. C., Thompson, T. A., Quataert, E., Loeb, A., & Waxman, E. 2011, ApJ, 734, 107 [NASA ADS] [CrossRef] [Google Scholar]
  140. Larson, R. B., Tinsley, B. M., & Caldwell, C. N. 1980, ApJ, 237, 692 [NASA ADS] [CrossRef] [Google Scholar]
  141. Latif, M. A., Schleicher, D. R. G., Schmidt, W., & Niemeyer, J. 2013, MNRAS, 432, 668 [NASA ADS] [CrossRef] [Google Scholar]
  142. Lazarian, A., & Pogosyan, D. 2000, ApJ, 537, 720 [NASA ADS] [CrossRef] [Google Scholar]
  143. Leite, N., Evoli, C., D’Angelo, M., et al. 2017, MNRAS, 469, 416 [NASA ADS] [CrossRef] [Google Scholar]
  144. Lemoine-Goumard, M., Renaud, M., Vink, J., et al. 2012, A&A, 545, A28 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  145. Lipari, P. 2003, A&A, 545, A28 [Google Scholar]
  146. Loeb, A., & Waxman, E. 2006, J. Cosmol. Astropart. Phys., 2006, 003 [NASA ADS] [CrossRef] [Google Scholar]
  147. Loewenstein, M., Zweibel, E. G., & Begelman, M. C. 1991, ApJ, 377, 392 [NASA ADS] [CrossRef] [Google Scholar]
  148. Mainali, R., Zitrin, A., Stark, D. P., et al. 2018, MNRAS, 479, 1180 [NASA ADS] [Google Scholar]
  149. Mawatari, K., Yamada, T., Fazio, G. G., Huang, J.-S., & Ashby, M. L. N. 2016, PASJ, 68, 46 [NASA ADS] [CrossRef] [Google Scholar]
  150. Meier, D. S., Turner, J. L., & Beck, S. C. 2002, AJ, 124, 877 [NASA ADS] [CrossRef] [Google Scholar]
  151. Moore, B., Lake, G., Quinn, T., & Stadel, J. 1999, MNRAS, 304, 465 [NASA ADS] [CrossRef] [Google Scholar]
  152. Morlino, G. 2017, High-Energy Cosmic Rays from Supernovae, eds. A. W. Murdin, & P. Alsabti, 1711 [Google Scholar]
  153. Morlino, G., & Caprioli, D. 2012, A&A, 538, A81 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  154. Mücke, A., Rachen, J. P., Engel, R., Protheroe, R. J., & Stanev, T. 1999, PASA, 16, 160 [NASA ADS] [CrossRef] [Google Scholar]
  155. Murata, K. L., Kajisawa, M., Taniguchi, Y., et al. 2014, ApJ, 786, 15 [NASA ADS] [CrossRef] [Google Scholar]
  156. Nakamura, K. 2010, J. Phys. G, 37, 075021 [NASA ADS] [CrossRef] [Google Scholar]
  157. Narayanan, D., Cox, T. J., Kelly, B., et al. 2008, ApJS, 176, 331 [NASA ADS] [CrossRef] [Google Scholar]
  158. Nath, B. B., & Biermann, P. L. 1993, MNRAS, 265, 421 [NASA ADS] [CrossRef] [Google Scholar]
  159. Nath, B. B., Laskar, T., & Shull, J. M. 2008, ApJ, 682, 1055 [NASA ADS] [CrossRef] [Google Scholar]
  160. Nozawa, T., Kozasa, T., Habe, A., et al. 2007, ApJ, 666, 955 [NASA ADS] [CrossRef] [Google Scholar]
  161. Nulsen, P. E. J. 1982, MNRAS, 198, 1007 [NASA ADS] [CrossRef] [Google Scholar]
  162. Oesch, P. A., Carollo, C. M., Stiavelli, M., et al. 2009, ApJ, 690, 1350 [NASA ADS] [CrossRef] [Google Scholar]
  163. Oesch, P. A., Bouwens, R. J., Illingworth, G. D., et al. 2014, ApJ, 786, 108 [NASA ADS] [CrossRef] [Google Scholar]
  164. Oesch, P. A., van Dokkum, P. G., Illingworth, G. D., et al. 2015, ApJ, 804, L30 [NASA ADS] [CrossRef] [Google Scholar]
  165. Oesch, P. A., Brammer, G., van Dokkum, P. G., et al. 2016, ApJ, 819, 129 [NASA ADS] [CrossRef] [Google Scholar]
  166. Ono, Y., Ouchi, M., Mobasher, B., et al. 2012, ApJ, 744, 83 [NASA ADS] [CrossRef] [Google Scholar]
  167. Ouchi, M., Ellis, R., Ono, Y., et al. 2013, ApJ, 778, 102 [NASA ADS] [CrossRef] [Google Scholar]
  168. Owen, E. R., Jacobsen, I. B., Wu, K., & Surajbali, P. 2018, MNRAS, 481, 666 [NASA ADS] [CrossRef] [Google Scholar]
  169. Owen, E. R., Jin, X., Wu, K., & Chan, S. 2019, MNRAS, 484, 1645 [NASA ADS] [Google Scholar]
  170. Paglione, T. A. D., Marscher, A. P., Jackson, J. M., & Bertsch, D. L. 1996, ApJ, 460, 295 [NASA ADS] [CrossRef] [Google Scholar]
  171. Papadopoulos, P. P. 2010, ApJ, 720, 226 [NASA ADS] [CrossRef] [Google Scholar]
  172. Papadopoulos, P. P., Thi, W.-F., Miniati, F., & Viti, S. 2011, MNRAS, 414, 1705 [NASA ADS] [CrossRef] [Google Scholar]
  173. Patrignani, C., et al. 2016, Chin. Phys. C, 40, 100001 [NASA ADS] [CrossRef] [Google Scholar]
  174. Peng, Y., Maiolino, R., & Cochrane, R. 2015, Nature, 521, 192 [NASA ADS] [CrossRef] [Google Scholar]
  175. Peng, F.-K., Wang, X.-Y., Liu, R.-Y., Tang, Q.-W., & Wang, J.-F. 2016, ApJ, 821, L20 [NASA ADS] [CrossRef] [Google Scholar]
  176. Persic, M., Rephaeli, Y., & Arieli, Y. 2008, A&A, 486, 143 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  177. Pfrommer, C. 2013, ApJ, 779, 10 [NASA ADS] [CrossRef] [Google Scholar]
  178. Planck Collaboration XIII. 2016, A&A, 594, A13 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  179. Pollack, J. B., & Fazio, G. G. 1963, Phys. Rev., 131, 2684 [NASA ADS] [CrossRef] [Google Scholar]
  180. Protheroe, R. J., & Johnson, P. A. 1996, Astropart. Phys., 4, 253 [NASA ADS] [CrossRef] [Google Scholar]
  181. Pudritz, R., & Fich, M. 2012, Galactic and Extragalactic Star Formation, Nato Science Series C (Springer: Netherlands) [Google Scholar]
  182. Rees, M. J. 1987, QJRAS, 28, 197 [NASA ADS] [Google Scholar]
  183. Rephaeli, Y., & Persic, M. 2014, Nucl. Phys. B–Proc. Suppl., 256, 252 [NASA ADS] [CrossRef] [Google Scholar]
  184. Rephaeli, Y., Arieli, Y., & Persic, M. 2010, MNRAS, 401, 473 [NASA ADS] [CrossRef] [Google Scholar]
  185. Ribaudo, J., Lehner, N., Howk, J. C., et al. 2011, ApJ, 743, 207 [NASA ADS] [CrossRef] [Google Scholar]
  186. Rieder, M., & Teyssier, R. 2016, MNRAS, 457, 1722 [NASA ADS] [CrossRef] [Google Scholar]
  187. Rimondi, F. 1993, 23rd International Symposium on Ultra-High Energy Multiparticle Phenomena Aspen, Colorado, September 12–17, 1993, 0400–404 [Google Scholar]
  188. Robertson, B. E., Ellis, R. S., Dunlop, J. S., McLure, R. J., & Stark, D. P. 2010, Nature, 468, 49 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  189. Rojas-Bravo, C., & Araya, M. 2016, MNRAS, 463, 1068 [NASA ADS] [CrossRef] [Google Scholar]
  190. Roseboom, I. G., Oliver, S., & Farrah, D. 2009, ApJ, 699, L1 [NASA ADS] [CrossRef] [Google Scholar]
  191. Rowlands, K., Wild, V., Nesvadba, N., et al. 2015, MNRAS, 448, 258 [NASA ADS] [CrossRef] [Google Scholar]
  192. Ruszkowski, M., Yang, H.-Y. K., & Reynolds, C. S. 2017, ApJ, 844, 13 [NASA ADS] [CrossRef] [Google Scholar]
  193. Ruszkowski, M., Yang, H.-Y. K., & Reynolds, C. S. 2018, ApJ, 858, 64 [NASA ADS] [CrossRef] [Google Scholar]
  194. Rybicki, G. B., & Lightman, A. P. 1979, Radiative Processes in Astrophysics (Wiley) [Google Scholar]
  195. Salem, M., & Bryan, G. L. 2014, MNRAS, 437, 3312 [NASA ADS] [CrossRef] [Google Scholar]
  196. Samui, S., Subramanian, K., & Srianand, R. 2005, Int. Cosmic Ray Conf., 9, 215 [NASA ADS] [Google Scholar]
  197. Sánchez Almeida, J., Elmegreen, B. G., Muñoz-Tuñón, C., & Elmegreen, D. M. 2014, A&ARv, 22, 71 [NASA ADS] [CrossRef] [Google Scholar]
  198. Sancisi, R., Fraternali, F., Oosterloo, T., & van der Hulst, T. 2008, A&ARv, 15, 189 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  199. Sazonov, S., & Sunyaev, R. 2015, MNRAS, 454, 3464 [NASA ADS] [CrossRef] [Google Scholar]
  200. Schleicher, D. R. G., & Beck, R. 2013, A&A, 556, A142 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  201. Schlickeiser, R. 2002, Cosmic Ray Astrophysics (Berlin: Springer) [CrossRef] [Google Scholar]
  202. Schober, J., Schleicher, D. R. G., & Klessen, R. S. 2013, A&A, 560, A87 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  203. Schober, J., Schleicher, D. R. G., & Klessen, R. S. 2015, MNRAS, 446, 2 [NASA ADS] [CrossRef] [Google Scholar]
  204. Schure, K. M., & Bell, A. R. 2013, MNRAS, 435, 1174 [NASA ADS] [CrossRef] [Google Scholar]
  205. Silvia, D. W., Smith, B. D., & Shull, J. M. 2010, ApJ, 715, 1575 [NASA ADS] [CrossRef] [Google Scholar]
  206. Skorodko, T., Bashkanov, M., Bogoslawsky, D., et al. 2008, Eur. Phys. J. A, 35, 317 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  207. Slattery, P. 1972, Phys. Rev. Lett., 29, 1624 [NASA ADS] [CrossRef] [Google Scholar]
  208. Smartt, S. J. 2009, ARA&A, 47, 63 [NASA ADS] [CrossRef] [Google Scholar]
  209. Smercina, A., Smith, J. D. T., Dale, D. A., et al. 2018, ApJ, 855, 51 [NASA ADS] [CrossRef] [Google Scholar]
  210. Solomon, P. M., & Vanden Bout, P. A. 2005, ARA&A, 43, 677 [NASA ADS] [CrossRef] [Google Scholar]
  211. Stanimirović, S., & Lazarian, A. 2001, ApJ, 551, L53 [NASA ADS] [CrossRef] [Google Scholar]
  212. Stark, D. P., Schenker, M. A., Ellis, R., et al. 2013, ApJ, 763, 129 [NASA ADS] [CrossRef] [Google Scholar]
  213. Stecker, F. W., Tsuruta, S., & Fazio, G. G. 1968, ApJ, 151, 881 [NASA ADS] [CrossRef] [Google Scholar]
  214. Stepney, S., & Guilbert, P. W. 1983, MNRAS, 204, 1269 [NASA ADS] [Google Scholar]
  215. Stewart, K. R., Brooks, A. M., Bullock, J. S., et al. 2013, ApJ, 769, 74 [NASA ADS] [CrossRef] [Google Scholar]
  216. Strong, A. W., Moskalenko, I. V., & Ptuskin, V. S. 2007, Annu. Rev. Nucl. Part. Sci., 57, 285 [NASA ADS] [CrossRef] [Google Scholar]
  217. Strong, A. W., Porter, T. A., Digel, S. W., et al. 2010, ApJ, 722, L58 [NASA ADS] [CrossRef] [Google Scholar]
  218. Sun, G., & Furlanetto, S. R. 2016, MNRAS, 460, 417 [NASA ADS] [CrossRef] [Google Scholar]
  219. Sur, S., Bhat, P., & Subramanian, K. 2018, MNRAS, 475, L72 [NASA ADS] [CrossRef] [Google Scholar]
  220. Tadaki, K.-I., Kodama, T., Tanaka, I., et al. 2014, ApJ, 780, 77 [NASA ADS] [CrossRef] [Google Scholar]
  221. Tang, Q.-W., Wang, X.-Y., & Tam, P.-H. T. 2014, ApJ, 794, 26 [NASA ADS] [CrossRef] [Google Scholar]
  222. Temim, T., Dwek, E., Tchernyshyov, K., et al. 2015, ApJ, 799, 158 [NASA ADS] [CrossRef] [Google Scholar]
  223. Thomas, R., Le Fèvre, O., Scodeggio, M., et al. 2017, A&A, 602, A35 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  224. Thome, W., Eggert, K., Giboni, K., et al. 1977, Nucl. Phys. B, 129, 365 [NASA ADS] [CrossRef] [Google Scholar]
  225. Torres, D. F. 2004, ApJ, 617, 966 [NASA ADS] [CrossRef] [Google Scholar]
  226. Tremaine, S., Richstone, D. O., Byun, Y.-I., et al. 1994, AJ, 107, 634 [NASA ADS] [CrossRef] [Google Scholar]
  227. Tremonti, C. A., Moustakas, J., & Diamond-Stanic, A. M. 2007, ApJ, 663, L77 [NASA ADS] [CrossRef] [Google Scholar]
  228. Turner, J. L., Beck, S. C., Benford, D. J., et al. 2015, Nature, 519, 331 [NASA ADS] [CrossRef] [Google Scholar]
  229. van den Bergh, S., Abraham, R. G., Ellis, R. S., et al. 1996, AJ, 112, 359 [NASA ADS] [CrossRef] [Google Scholar]
  230. VERITAS Collaboration (Acciari, V. A., et al.) 2009, Nature, 462, 770 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  231. Walker, M. A. 2016, ApJ, 818, 23 [NASA ADS] [CrossRef] [Google Scholar]
  232. Wang, C. H. 1991, PhD Thesis, Iowa State University [Google Scholar]
  233. Wang, X., & Fields, B. D. 2014, AIP Conf. Proc., 1595, 231 [NASA ADS] [CrossRef] [Google Scholar]
  234. Wang, X., & Fields, B. D. 2018, MNRAS, 474, 4073 [NASA ADS] [CrossRef] [Google Scholar]
  235. Wardlow, J. L., Cooray, A., Osage, W., et al. 2017, ApJ, 837, 12 [NASA ADS] [CrossRef] [Google Scholar]
  236. Watson, D., Christensen, L., Knudsen, K. K., et al. 2015, Nature, 519, 327 [NASA ADS] [CrossRef] [Google Scholar]
  237. Wentzel, D. G. 1971, ApJ, 163, 503 [NASA ADS] [CrossRef] [Google Scholar]
  238. Whitmore, J. 1974, Phys. Rep., 10, 273 [NASA ADS] [CrossRef] [Google Scholar]
  239. Wiener, J., Zweibel, E. G., & Oh, S. P. 2013, ApJ, 767, 87 [NASA ADS] [CrossRef] [Google Scholar]
  240. Wright, E. L. 2006, PASP, 118, 1711 [NASA ADS] [CrossRef] [Google Scholar]
  241. Yamasawa, D., Habe, A., Kozasa, T., et al. 2011, ApJ, 735, 44 [NASA ADS] [CrossRef] [Google Scholar]
  242. Yan, H., & Lazarian, A. 2004, ApJ, 614, 757 [NASA ADS] [CrossRef] [Google Scholar]
  243. Yan, R., Newman, J. A., Faber, S. M., et al. 2006, ApJ, 648, 281 [NASA ADS] [CrossRef] [Google Scholar]
  244. Yang, Y., Tremonti, C. A., Zabludoff, A. I., & Zaritsky, D. 2006, ApJ, 646, L33 [NASA ADS] [CrossRef] [Google Scholar]
  245. Yoast-Hull, T. M., Gallagher, J. S.III, Aalto, S., & Varenius, E. 2017, MNRAS, 469, L89 [NASA ADS] [CrossRef] [Google Scholar]
  246. Zweibel, E. G. 2013, Phys. Plasmas, 20, 055501 [NASA ADS] [CrossRef] [Google Scholar]
  247. Zweibel, E. 2016, APS Meeting Abstracts, SR1.001 [Google Scholar]

Appendix A: Cosmic ray interactions

Hadronic and leptonic CRs interact through different mechanisms. For instance, their comparatively low mass makes CR electrons considerably more susceptible to scattering and cooling processes compared to their heavier proton counterparts. For protons, absorption processes are more important. We define a “cooling” process as one in which the original CR survives the interaction but loses some fraction of its energy and thus can continue to scatter many times. Conversely, we specify an “absorption” process to be one where either a substantial fraction of the original CR energy is lost in an interaction (such that the particle can no longer be considered high energy afterwards, and is not realistically able to undergo another interaction of the same type), or where the particle is truly destroyed in the process.

A.1. High energy proton interactions

The relevant processes governing hadronic interactions of CRs with radiation and matter fields in astrophysical environments are Photo-pion (pγ) production, Bethe-Heitler pair production and pp pion-production interactions. The pion-producing processes may be considered as absorption mechanisms for hadronic CRs, as they result in substantial energy loss of the particle in a single interaction. The Bethe-Heitler pair production process results in a lower energy loss fraction per interaction and so is better regarded as a continuous cooling mechanism rather than a stochastic absorption process. Hence, we separate our discussion of these mechanisms accordingly.

Each of these interactions produce neutrons, protons, leptons and other hadrons including charged and neutral pions (see Pollack & Fazio 1963; Gould & Burbidge 1965; Stecker et al. 1968; Mücke et al. 1999; Berezinsky et al. 2006; Dermer & Menon 2009). The protons can continue to propagate and undergo further interactions like that of the primary, until the particle energy is below the respective interaction threshold. The neutrons have a half life of 880 s (Nakamura 2010), and so are most likely to undergo a β decay . The charged and neutral pions only have a fleeting lifetime before they undergo decays, so these cannot realistically undergo further interactions. The decay of neutral pions is dominated by the electromagnetic process π0 → 2γ with branching ratio (BR) of 98.8% (Patrignani 2016), arising on a timescale of 8.5 × 10−17 s. Charged pions instead decay to produce leptons (predominantly electrons and positrons) and neutrinos by the weak interaction for positively charged pions, or 17 for negatively charged pions, both with BR = 99.9% (Patrignani 2016) and arising over a longer timescale of 2.6 × 10−8 s. This is the main mechanism by which secondary CR electrons are injected into the environment.

In the following, we denote the energy of the primary CR proton as Ep = γpmpc2, or ϵp = γp(mp/me) if normalised to electron rest mass, where γp is the Lorentz factor. Where interactions involve a target field of photons, the normalised photon energy of frequency ν is written as ϵ = hν/mec2. Alternatively, for CR-baryon interactions, the target baryon (presumed to be a hydrogen nucleus) has a negligible energy compared to the CR primary, and so may be regarded as “at rest”, with γp = 1. In this case, the normalised energy of the target baryon is ϵ ≈ mp/me. It is also useful to define the invariant normalised interaction energy as ϵr = γp(1 − βpcosθ)ϵ, for θ as the angle between the momentum vectors of the incident proton and target field particle (photon or proton). βp is the CR primary velocity in units of c.

A.1.1. Absorption processes

Photo-pion production and pp-pion production are both inelastic processes which ultimately lead to the absorption of the original CR particle and a reprocessing of its energy. While fundamentally similar, the first of these involves the interaction between a CR proton and a low-energy photon, while the second involves a low energy baryon as the target.

Photo-pion production. The dominant interactions in this process are single-pion at lower energies, and multiple-pion production above ϵr ≈ 3500. This occurs both via resonant production and direct production (see Mücke et al. 1999). The resonant production channel (which occurs at around three times the rate of direct production at ϵr ≈ 500) progresses with the formation of the Δ+ resonance which decays by one of two channels, either Δ+ → pπ0 (BR = 2/3), or Δ+ → nπ+ (BR = 1/3) (see Berezinsky & Gazizov 1993), where the pion decay then proceeds as described before. Residual interactions (e.g. pγ → Δ++π, Δ0π+) provide additional charged pions and, when taking all channels into account, the production of each type of pion arises at a similar rate (see Dermer & Menon 2009). Based on experimental data (Baldini et al. 1987), the inelastic cross section for the photo-pion interaction can be reasonably approximated as , where ℋ(...) is the Heaviside step function, with and ϵth = 390. These are the same values as suggested in Dermer & Menon (2009), which are also adopted in Atoyan & Dermer (2003).

The collision rate is defined by the collision timescale as (γp)γπ,coll = [τγπ,coll(γp)]−1 for a CR proton with Lorentz factor γp. Some fraction of these collisions will be elastic, while others will contribute to the pion-producing inelastic pathway. The total collision rate accounts for both of these but, for our purposes, we require the CR proton absorption rate - this is the fraction of collisions that successfully yield pion-production and subsequent CR absorption. The inelastic collision rate follows as:

(A.1)

which is the inverse timescale for the photo-pion interaction. Substituting the dominant CMB radiation field for the photon number density,

(A.2)

where symbols retain their earlier definitions, gives

(A.3)

where the substitution x = ϵ/Θ has been applied. The integral has a standard result of π2/6, yielding the final absorption rate as

(A.4)

pp-pion production. The pp mechanism is similar to the photo-pion process, but with the interaction target as the low energy baryon field of the ISM. The dominant channels proceed via Δ+ and Δ++ resonance production (pp → pΔ+ and pp → nΔ++) which then decay according to

(A.5)

and

(A.6)

with ξ0 and ξ± as the neutral and charged pion multiplicities respectively which themselves are a function of the CR primary energy (see Almeida et al. 1968; Skorodko et al. 2008). Using the empirical parameterisations given by Blattnig et al. (2000), the inelastic cross section for pion production indicates branching ratios for each of the pion species {π+, π, π0} as {0.6, 0.1, 0.3} at 1 GeV, levelling out to {0.3, 0.4, 0.3} at higher energies, above 50 GeV. These suggest that the production rate of each type of pion is similar (as with the photo-pion process). The total inelastic cross section, which is required for the pp absorption rate of CR primaries, is well described by the Kafexhiu et al. (2014) parameterisation:

(A.7)

which we adopt for subsequent calculations in this work. Here, , for Eth as the threshold energy for the process (0.28 GeV18). The rate of CR absorption then follows as:

(A.8)

A.1.2. Cooling processes

Compared to the pion production mechanisms described in the previous section, photo-pair production is better described as a cooling process as each interaction only results in a small change in the CR primary energy. Many cooling processes usually attributed to charged particles (e.g. radiative cooling via inverse-Compton and Synchrotron emission) are not important for protons due to their greater mass. Instead, these are more relevant to the electronic secondaries (see Appendix A.3). Thus, unabsorbed free-streaming protons would more likely to be cooled by the adiabatic cosmological expansion, particularly if they are able to escape from their source galaxy to traverse intergalactic space.

Photo-pair production. The Photo-pair Bethe-Heitler process (Bethe & Heitler 1934) produces predominantly electron-positron pairs over the energies of interest in this study (Blumenthal 1970; Klein 2006) and, for a CR proton, proceeds as pγ → pe+e. When the interaction energy ϵr ≳ 60 (as may be assumed in this work), a reasonable approximation to the fitted inelastic cross section (see Jost et al. 1950; Bethe & Maximon 1954; Blumenthal 1970; Stepney & Guilbert 1983) was found to be:

(A.9)

where αf is the fine structure constant, σT is the Thomson cross-section, and the fitting constant kγe is set to 6.7 (see Owen et al. 2018)19. The timescale for this interaction is

(A.10)

(see Protheroe & Johnson 1996; Dermer & Menon 2009), if assuming that the electron-positron pairs are produced in the zero-momentum frame of the interaction, and the interaction energy is completely dominated by the contribution from the CR primary. The energy loss rate, , is then

(A.11)

where u = mec2/γpkBT and the function

(A.12)

results from evaluating the integral in Eq. (A.10)20, where 𝒞(u) = 0.74 and 𝒟(u) = Γ(5/2)ζ(5/2) when b ≪ 1, 𝒞(u)=u3/2ln(u)eu and 𝒟(u)=u3/2eu when u ≫ 1, and ℰ(u)= − ln[1 − eu] for all values of u (Dermer & Menon 2009).

Adiabatic cosmological expansion. If propagating over cosmological scales, the expansion of the Universe will cool a particle. The rate at which this arises is

(A.13)

(see Gould 1975; Berezinsky & Grigor’eva 1988; Berezinsky et al. 2006), where H0 = (67.8 ± 0.9) km s−1 Mpc−1, Ωm = 0.308 ± 0.012 and ΩΛ = 0.691 ± 0.0062, and negligible curvature and radiation energy densities (Planck Collaboration XIII 2016). This is only important for particles which are not liable to more rapid cooling by other processes (see, e.g. Appendix A.3).

A.2. Electron injection

Secondary high energy (> GeV) CR electrons are injected by the pion-yielding interactions of CR protons with the host environment (e.g. via pp- and photo-pion processes). Such injection of electrons is thought to be even more important than their direct acceleration, particularly in starburst environments where observations indicate that 60–80% of high energy electrons are secondaries (e.g. Loeb & Waxman 2006; Lacki et al. 2011; Lacki & Beck 2013; Schober et al. 2015). This high secondary fraction is also consistent with actively star-forming regions in the more local Universe (see Paglione et al. 1996; Torres 2004; Domingo-Santamaría & Torres 2005; Persic et al. 2008; de Cea del Pozo et al. 2009; Rephaeli et al. 2010) where conditions would be ideal for the hadronic CR interactions thought to produce these electron CR secondaries.

The energy spectrum of the CR secondary particles is governed by the interaction cross sections of the relevant processes driving their formation. These are photo-pion production in radiation dominated systems (i.e. where nph ≫ nb), or pp-pion production in matter dominated systems (nb ≫ nph). In general, we can consider both of these together with their combined effective absorption coefficient as the sum of the contributions from these two processes:

(A.14)

This, coupled with the (energy dependent) multiplicity of the products and the energy spectrum of the primary CRs allows a reasonable model for the differential CR secondary electron production rate (Kelner et al. 2006; Kamae et al. 2006) to be calculated:

(A.15)

This is the injection rate of electrons per unit volume per energy interval (between Ee and Ee + dEe), and introduces the differential CR proton energy density as

(A.16)

This is related to the differential number density of protons (number density per energy interval) by the relation Φp = Ep ⋅ np. In Eq. (A.15), we introduced Fe as the number of electrons produced per interaction per energy interval – effectively the differential multiplicity of electron production. At high energies, above 100 GeV, a parametrisation of the function Fe following Kelner et al. (2006) may be used:

(A.17)

where we have used the substitution x = Ee/Ep. Here,

(A.18)

where, for convenience, we retain the earlier definition for χ as the CR proton energy normalised to the pp interaction threshold energy, . Also

(A.19)

and βe = [0.00042(lnχ)2+0.063lnχ−0.28]−1/4. While valid for high energies, Kelner et al. (2006) shows that this form of Fe to be unsuitable at the lower energies in our range of interest, below 100 GeV. Instead, for this regime they suggest the use of a delta-function (also see Aharonian & Atoyan 2000). In this case, the we break down the process into its composite mechanisms (charged pion production, and their subsequent decay), and model the pion production number by 1–100 GeV CRs per interaction per energy interval (this time between Eπ± and Eπ± + dEπ±) as:

(A.20)

where we have introduced the contraction Ekin = Ep − mpc2. The efficiency term Kπ± encodes the fraction of energy passed to pion production from the CR proton. The branching ratios for the formation of each type of pion are all approximately 1/3 in both the pp-pion and photo-pion interaction channels (see Appendix A.1), and this is not strongly dependent on energy (Aharonian & Atoyan 2000; Blattnig et al. 2000). Accounting also for neutrino production, the total fraction of the CR primary energy passed to charged pions is around 0.6–0.7 (see Dermer & Menon 2009). Here we conservatively adopt the lowest value in this range and set Kπ± = 0.6. Using the delta function approximation (A.20) in Eq. (A.15), we arrive at the injection rate of charged pions per unit volume per energy interval as:

(A.21)

where we have used the contractions α* = α*(Ep, r) and Φp = Φp(Ep, r) and other symbols retain their earlier definitions. Evaluating the integral gives:

(A.22)

(Aharonian & Atoyan 2000), where the energy parameter L is defined as . The charged pions formed in this process rapidly decay into leptons and neutrinos (see Appendix A.1), with the pion energy roughly split equally between the neutrino and lepton products – these form in the ratio of three neutrinos to each lepton. Detailed studies find that the leptons inherit around 26.5% of the pion energy (e.g. Lipari 2003; Lacki & Beck 2013), with the remainder passed to neutrinos and effectively lost from the system.

We may write the rate of energy injection by pions of a specified energy as . By conservation of energy, it follows that this energy is provided by the corresponding decays of electrons such that , where ℳ|Ep is introduced as the electron multiplicity, and Ke = 0.265 is the efficiency of energy transfer from pions to electrons. The multiplicity is determined empirically using a fitting function of the form c1 + c2sc3 as specified in Fiete Grosse-Oetringhaus & Reygers (2010) and Albini et al. (1976) with c1 = 0.0, c2 = 3.102 and c3 = 0.178 (s is the GeV centre-of-mass interaction energy)21. At 1 GeV, the electron production multiplicity is four, with the centre-of-mass energy governed by Ep. This means each injected CR electron inherits around 3–5% of the energy of the CR proton primary initiating the interaction. The electron injection rate (per unit volume) thus follows as:

(A.23)

The pion energy is equal to the fraction of energy passed from the primary proton, as specified by the efficiency parameter Kπ±. Therefore, Eπ± = Kπ±Ep and so we may write

(A.24)

where is the effective overall energy efficiency parameter for the production of electrons from CR protons, accounting for both the efficiencies in the pion production and pion decay steps. In terms of the electron and proton Lorentz factors, Eq. (A.24) may be written as:

(A.25)

A.3. High energy electron interactions

The injected electron secondaries will also interact with their environment. While the hadronic primaries are found to be predominantly absorbed, CR electrons are more susceptible to continuous cooling processes. The main means by which the electrons can lose their energy are: radiatively (by inverse-Compton or synchrotron emission), electro-statically (by Coulomb interactions and free–free processes) and via triplet pair-production (a Bethe-Heitler process).

Coulomb and free–free interactions. The cooling rate of electrons due to Coulomb interactions with a low-density fully-ionised plasma is given by:

(A.26)

(see, e.g. Dermer & Menon 2009; Schleicher & Beck 2013) which we note is independent of the electron energy. Here, lnΛ ≃ 30 is the Coulomb logarithm (which accounts for the ratio between the maximum and minimum impact parameters in a particle-particle interaction). The electron cooling rate due to bremsstrahlung (free–free interactions) is

(A.27)

(e.g. Dermer & Menon 2009; Schleicher & Beck 2013).

Inverse-Compton and synchrotron. In the Thomson limit, cooling rates for inverse-Compton and synchrotron both take the form

(A.28)

(e.g. Blumenthal 1970; Rybicki & Lightman 1979) where Ui is the energy density of the relevant field. For inverse-Compton scattering, the high energy electrons interact with lower energy photons, hence the energy density can be taken as that of the dominant radiation field – in the case of the CMB, UCMB, this can be found by integrating over the black-body spectrum given by Eq. (A.2):

(A.29)

which is a function of redshift (and where symbols retain their earlier meanings). The inverse-Compton cooling function has a more complicated form in the high energy Klein-Nishina limit, as shown in Fig. 1 (see also Blumenthal & Gould 1970; Dermer & Menon 2009).

In synchrotron emission, electrons cool in a magnetic field. The energy density in this case is given by for a (local) magnetic field strength of B.

Triplet pair-production. At energies above a GeV, triplet pair-production can provide some contribution to the cooling energy losses of electrons. This pair-production process operates in a broadly similar manner to the Bethe-Heitler photo-pair production arising with protons. For ϵr ≥ 30, an approximation correct to within 0.1% for the cross section can be found:

(A.30)

(Borsellino 1947; Joseph & Rohrlich 1958; Haug 1981), where ψ = ln(2ϵr). We note that for the energies of interest here (ϵr ≥ 60), we may adopt the same form of approximation as we did with the photo-pair interaction (for the CR protons):

(A.31)

The principal difference is that the fitting constant is set to a lower value of kTPP = 0.05 here. This gives a path length consistent with (Bhattacharjee & Sigl 2000) at 1017 eV and z = 0. The cooling rate then follows as:

(A.32)

where ℱTPP takes the same form as the function ℱγe in Eq. (A.12), but with kγe = 6.7 replaced with kTPP = 0.05 in line with the corresponding cross sections, and u = mec2/γekBT.

Appendix B: Inverse-Compton X-ray luminosity

The inverse-Compton X-ray luminosity of starburst galaxies at high-redshift was calculated analytically by Schober et al. (2015) in the 0.2 to 10 keV band, although this previous study invoked a redshift evolution in the volume of the host galaxy following V ∝  (1 + z)−3. The model used in this paper is otherwise comparable and our benchmark model is of volume 1 kpc3, being similar to the volume of the Schober et al. (2015) model at z = 7. We run our inverse-Compton X-ray luminosity simulation with ℛSN ≈ 1 yr−1 and sum over the galaxy to find the total luminosity of the starburst X-ray glow up to z = 10 to compare with the previous model as a cross-check. We plot our results (black dots) and the previous model when adjusted to account for its additional redshift evolution (grey dashed line) in Fig. B.1, which shows good consistency between the two approaches.

thumbnail Fig. B.1.

Inverse-Compton protogalaxy X-ray luminosity evolution over redshift, as calculated using the approach detailed in this paper with ℛSN ≈ 1 yr−1 (semi-analytical Monte-Carlo simulation) compared to the earlier study, Schober et al. (2015) (dashed grey line), which adopts an analytical model. Our results are consistent with those obtained by Schober et al. (2015).

Open with DEXTER

Appendix C: Heating channel parametrisation

Here we show a cross-check for the two CR heating rate parameterisations introduced in Sect. 5.6. We may check these parameterisations by taking two extreme cases which have been calculated in full, in Sect. 5.4: Consider the ℛSN = 0.1 yr−1 ≡ ℛSF = 16 M yr−1 systems at z = 0 and z = 12 in Fig. 6. The central DC and IX heating levels are shown to be 1.7 × 10−24 erg cm−3 s−1 and 1.4 × 10−29 erg cm−3 s−1 respectively at z = 0, with the IX heating power rising to 4.2 × 10−25 erg cm−3 s−1 by z = 12 (DC being unchanged). Using the parameterisations (50) and (51), our scalings give a central DC and IX heating level of 1.7 × 10−24 erg cm−3 s−1 and 1.6 × 10−29 erg cm−3 s−1 at z = 0 respectively, with the X-ray heating power raising to 4.5 × 10−25 erg cm−3 s−1 by z = 12, which reconcile to within a sufficient degree of accuracy (around 30%). At large distances, the X-ray heating rates at 20 kpc are found to be 6.2 × 10−38 erg cm−3 s−1 and 1.8 × 10−33 erg cm−3 s−1 for z = 0 and z = 12 respectively from our earlier Sect. 5.4 calculation. Instead, from the scaling parametrisation (52), these levels can be estimated as 4.8 × 10−38 erg cm−3 s−1 and 1.4 × 10−33 erg cm−3 s−1 – again, consistent to within similar error bounds. We may thus adopt the parameterisations above to apply our model to the observed galaxies introduced in Sect. 2, avoiding the need to run multiple computationally expensive calculations.

Appendix D: Star-formation rates and feedback

Residual star-formation appears to remain in our sample of post-starburst systems, or at least has re-emerged by the time of observation. This may be indicative of a smoother star-formation history, arising as inflows begin to re-emerge and strengthen over time. In such a scenario, we would expect to see higher rates of star-formation for systems observed at a greater time after the onset of quenching. Indeed, a first view of the data would appear to support this, as seen in the top panel of Fig. D.1 where the star-formation rate ℛSF at the point of observation for each of the high-redshift galaxies is plotted against the time elapsed since the estimated end of their starburst episode, tqui. However, some degeneracies with other effects may be responsible for this trend: in the middle panel we show the star-formation rate surface density, and specific star-formation rate, SSFR = ℛSF/M* in the lower panel. The observed post-starburst shows a more complicated behaviour: while, up to around tqui ≈ 350 M yr−1 kpc−2, increases with tqui. However, for larger tqui values, the star-formation surface density reduces again. This could result from the compactness of some of the galaxies: the behaviour of the SSFR in the lower panel (which does not depend on the physical extent of the system) is different and shows that the longer the time elapsed since the end of star-formation, the lower the observed SSFR. This is intriguing behaviour, and in tension with our first view of the data, and may suggest tqui is better considered as an indicator for the level of feedback experienced by a system (see Sect. 6.3 for further discussion).

thumbnail Fig. D.1.

Star-formation rates in the sample of 12 post-starburst galaxies. MACS1149-JD1 is highlighted in red in all panels as this galaxy is at substantially higher redshift (z = 9.11) than the others in the sample. Upper panel: observed quiescent star-formation rate ℛSF plotted against the length of quiescence period. This period is measured up to the point of observation, so tqui effectively indicates the time elapsed since the end of starburst activity. This shows a strong positive correlation, with a Pearson coefficient (weighted by the error bounds) of 0.93. Middle panel: as above, except the y-axis is star-formation rate surface density (). This removes the correlation seen in the upper panel, with a Pearson coefficient of −0.02 (when removing the clear outlier A1689-zD1, as labelled). Lower panel: y-axis as specific star-formation rate, which is adjusted for galaxy mass where SSFR = ℛSF/M*. If removing the labelled outlier A1689-zD1, this now gives a non-negligible negative correlation of weighted Pearson coefficient −0.64.

Open with DEXTER

One notable outlier in the top panel of Fig. D.1, GNS-zD5, exhibits a remarkably high star-formation rate. This is predominantly due to its large mass – indeed, it falls in line with the rest of the sample in the lower panels when normalised by size or mass. A1689-zD1 is arguably less conformant, with a particularly low star-formation rate for its size and mass. It may be that this system does not truly fit into the burst-mode SFH considered for these galaxies, instead experiencing a more transient and gradual SFH. In this figure (and subsequently), we have highlighted the high-redshift system MACS1149-JD1. Despite having been observed at a considerably higher redshift than the rest of the sample, its star-forming behaviour does not appear to be unusual.

We can also consider the star-formation rate surface density during the starburst episode, (see Fig. D.2). As before, this also correlates with the quiescence timescale, as does the star-formation surface density during the burst. However, the behaviour is again complicated (similar to the middle panel of Fig. D.1), increasing to a peak before falling away again. It is unclear whether this pattern is physical in origin, due to random fluctuations with our limited data points, or whether it is a mathematical feature resulting from the way and associated quantities were calculated: The star-formation rate during the burst is inversely related to tSB which itself is estimated by [t(zobs)−t(zf)] − tqui (see Eq. (1)). This implies an inverse proportionality between (as well as quantities derived from this) and tSB when tqui >  t(zobs)−t(zf), which occurs after around 350 Myr (which is indeed the location of the peaks in Fig. D.2)22. We do not plot the SSFR in Fig. D.2, given that is derived from M* rather than measured independently.

thumbnail Fig. D.2.

As Fig. D.1, but for the starburst phase of each system. The weighted Pearson correlation value for the top panel is found to be 0.71, while the lower panel indicates two behaviours with a less clear trend (overall correlation of 0.39) – see text for details. As SSFR is estimated in the burst phase using tSB and M* rather than being determined independently, it is not physically meaningful to include it here.

Open with DEXTER

All Tables

Table 1.

High-redshift systems with a comparatively suppressed rate of star-formation, but show evidence of a developed stellar population indicative of an earlier starburst episode.

Table 2.

Sample of high-redshift starburst galaxies identified in the literature which appear to be undergoing a starburst event during the epoch at which they are observed.

Table 3.

Peak energy densities of the radiation fields for all 16 of our sample of star-forming galaxies during their starburst phase, compared to the energy density of the CMB radiation field at their observed redshift.

Table 4.

Derived information for high-redshift quenched post-starburst and starburst systems, as introduced in Tables 1 and 2.

All Figures

thumbnail Fig. 1.

Left: CR proton losses in terms of effective interaction path lengths. This results from a protogalaxy model with target density of np = 10 cm−3, a CMB radiation field specified at a redshift of z = 7. The thick black line indicates the total effect that would be experienced by a CR proton at the indicated energy accounting for all energy loss contributions. The lines, as labelled, are (1) pp Pion-production (pp, pπ±π0); (2) CMB Photo-pion (pγ, π±π0); (3) CMB Photo-pair (pγ, e±); (4) adiabatic losses due to cosmological expansion; (5) total. Right: secondary CR electron loss lengths with the same model. The lines, as labelled, are (1) adiabatic losses due to cosmological expansion; (2) inverse-Compton (radiative) cooling due to the CMB with (2a) as the estimated path length in the Thomson limit, and (2b) being the estimated path length if accounting for the transition to Klein-Nishina behaviour at higher energies (we find this does not significantly impact our later results, so the Thomson limit is used in our main calculations); (3) Coulomb losses (responsible for directly heating the ISM plasma); (4) synchrotron (radiative) cooling in a 5 μG magnetic field, being an appropriate strength for this model; (5) Triplet Pair-Production (eγ, e±); (6) free–free (Bremsstrahlung) losses; (7) total.

Open with DEXTER
In the text
thumbnail Fig. 2.

Interaction loss lengths accounting for all relevant processes (i.e. the solid line of Fig. 1), but calculated at redshifts of z = 0 − 10 as shown by the legend. Left: path lengths for protons; right: path lengths for electrons.

Open with DEXTER
In the text
thumbnail Fig. 3.

Electron spectra following their injection from a single discrete source (normalised to the injection level at the minimum energy indicated by the dashed red line, ne, 0), accounting for variations in the energy densities of the radiation fields governing inverse-Compton cooling as well as electron diffusion. In all panels, the minimum energy at which CR secondary electrons may be injected is indicated by the dashed red vertical line – this would correspond to the secondaries produced when a 0.28 GeV CR proton primary undergoes a pp interaction. The injection spectral index would be approximately Γ = −2.1, as indicated in blue in the upper panel. Indices of subsequent steady-state spectra are steeper than this due to the energy dependence of the electron cooling rate. Upper panel: spectral variation with distance: dashed light grey line shows the steady-state spectrum expected at the injection point (we note that this is not the injection spectrum); remaining lines (solid black through to solid grey) denote the spectrum at 0.1–0.5 kpc in 0.1 kpc increments when adopting a CMB energy density at z = 7 and no stellar radiation field. The steepening of the high-end of the spectrum with distance is due to the energy-dependence of the radiative cooling function and of electron diffusion (with higher energy particles diffusing and undergoing radiative cooling more quickly). Middle panel: as above, but spectra calculated at 0.2 kpc (black lines) and 0.5 kpc (blue lines) from their injection location, with CMB at z = 7 and the spectral lines indicating evolution according to: no stellar radiation field (solid lines), stellar radiation field consistent with ℛSF = 3.5 M yr−1 (dashed lines), ℛSF = 35 M yr−1 (dot-dashed lines) and ℛSF = 350 M yr−1 (dotted lines) where star-formation is concentrated in a 1 kpc volume. Lower panel: as above, but with all lines calculated at 0.2 kpc from their injection source location, lines in black corresponding to a CMB-only radiation field and lines in blue corresponding to a strong radiation field with ℛSF = 350 M yr−1, shown at different redshifts over a similar range to our high-redshift galaxy sample: z = 6 (solid line), 7 (dashed line), 8 (dot-dashed line) and 9 (dotted line).

Open with DEXTER
In the text
thumbnail Fig. 4.

Steady-state protogalactic CR proton (primary) and electron (secondary) profile integrated over all energies for a SN rate ℛSN = 0.1 yr−1. The blue line shows the profile for electrons, while the black line corresponds to protons. The solid part of the black line is the proton distribution within the interstellar medium of the protogalaxy, where the magnetic field model is more robust. Outside the protogalaxy, the magnetic field is less well understood and so our result (the dashed part of the black line) should be regarded as an upper limit.10

Open with DEXTER
In the text
thumbnail Fig. 5.

CR heating rates of the protogalaxy ISM by (1) direct Coulomb heating, in black, up to 0.6 kpc, and (2) indirect X-ray heating, in red. The solid lines are the results calculated for ℛSN = 1 yr−1, dotted lines for ℛSN = 0.1 yr−1 and dot-dashed for ℛSN = 0.01 yr−1. ℛSN determines the energy budget of the system, and so the CR heating power in both channels scales directly. The shaded regions show the stellar radiation “floor” for the IX process: any evolution that brings the IX line into this area will instead be dominated by the up-scattering of stellar-radiation by CR electrons, preventing the X-ray heating power from falling into this region. The lower region accounts for a stellar radiation field when ℛSN = 0.01 yr−1, the intermediate case above, for ℛSN = 0.1 yr−1, and the highest region arises for ℛSN = 1.0 yr−1. We note that the DC heating rate falls to negligible levels outside the protogalaxy (not shown) due to the magnetic containment of the CRs.

Open with DEXTER
In the text
thumbnail Fig. 6.

CR heating rates via the direct Coulomb (solid black lines) and indirect X-ray mechanisms (red dashed lines) for a protogalaxy with nb, 0 = 10 cm−3 and ℛSN = 0.1 yr−1 for redshifts z = 0, 2, 4, 6, 8, 10 and 12. The black DC lines (only one line is visible as they overlap) do not show any discernible evolution over this redshift range, but those for the IX effect show significant variation due to the evolving CMB. The central heating rate at z = 12 by both mechanisms is comparable at a level of around 10−25 erg cm−3 s−1. The shaded orange region represents the stellar radiation “floor” for the IX process: any redshift evolution that brings the IX line into this area (i.e. below z ≈ 4) will instead be dominated by the up-scattering of stellar-radiation, preventing the IX power from falling into this region. The DC heating power becomes negligible outside the protogalaxy (not shown).

Open with DEXTER
In the text
thumbnail Fig. 7.

Fraction of 40 MeV CR secondary electron energy passed to the inverse-Compton X-ray heating channel (red solid line) compared to the direct Coulomb heating channel (black solid line) as a function of density. The left hand y-axis corresponds to the relative fraction, and this indicates a strong preference for the IX heating channel in low density media, or the DC channel in higher density media (if ionised). Timescales of respective processes are shown in grey to illustrate the underlying physics: the solid line gives the total energy loss timescale for CR electrons. This is calculated by combining the contributions from inverse-Compton (dashed horizontal line), Coulomb (dotted line) and free–free losses (dot-dashed line). The respective timescale y-axis is shown to the right of the plot.

Open with DEXTER
In the text
thumbnail Fig. 8.

Schematic illustration of the CR heating processes operating in and around a starburst protogalaxy. At A, internal interstellar heating is dominated by the DC channel as CRs and their secondaries are not able to easily escape from the galactic magnetic field of their host system. This increases thermal gas pressure and leads to quenching (French et al. 2015). At B, the IX channel heating mechanism is more important, as this does not rely on CR escape. The X-rays emitted from the galaxy begin to heat and evaporate any cold inflowing filaments, and can do so on a timescale shorter than the inflow time. Such filamentary inflows are required to drive star-formation by providing cold gas which can more easily undergo fragmentation and collapse than ambient interstellar gas (Ribaudo et al. 2011; Sánchez Almeida et al. 2014), so their destruction could strangulate the galaxy (Ribaudo et al. 2011; Sánchez Almeida et al. 2014; Peng et al. 2015).

Open with DEXTER
In the text
thumbnail Fig. 9.

Quiescence timescales against generic physical (upper row) and derived (lower row) galactic parameters. The data of MACS1149-JD1 are highlighted in red (this galaxy is at much higher redshift than the others in the sample, but does not otherwise appear to be unusual). The value of tqui averaged over all data points is indicated by the horizontal dashed grey line. Upper left: plot against estimated galaxy radius. This indicates that the quiescent time scale is also independent of the size of the system, with a weighted Pearson correlation coefficient of 0.15. Error bars in rgal are as quoted in the literature (A1689-zD1 at 50% and MACS1149-JD1 at 71%), otherwise estimated as 30%. Upper right: plot against stellar mass, which shows a positive correlation (coefficient of 0.58). Lower left: plot against galaxy dynamical timescales, which tqui appears to be independent of, with a correlation of 0.06. Lower centre: plot against quenching timescale, τQ, which shows a negative correlation (coefficient of −0.81). This is consistent with tqui being a marker for the strength of feedback – the longer taken to quench star-formation, the weaker the CR driven feedback and the shorter the quiescence timescale. Uncertainties in τQ are unclear and so are not plotted. Lower right: as per the centre panel, but against the strangulation timescale, τS. Again the negative correlation (of coefficient −0.78) is consistent with tqui being influenced by feedback with shorter values for longer strangulation times. As with τQ, uncertainties in τS are not clear so are omitted.

Open with DEXTER
In the text
thumbnail Fig. 10.

Star-formation rate surface density plotted against galaxy stellar mass, indicating these quantities are moderately correlated (weighted Pearson coefficient of 0.71 for the post-starburst systems). Upper panel: for the four galaxies observed during their starburst phase (see Table 2) are plotted in grey, with the post-starbursts are shown in black (see Table 1). We plot the MACS1149-JD1 in red (this system is observed at a much higher redshift than the others in the sample, but its behaviour appears to be relatively consistent). Lower panel: as above, with the removal of the starbursts to more clearly show the distribution of the post-starburst galaxies.

Open with DEXTER
In the text
thumbnail Fig. 11.

Observed star-formation rate surface densities for the post-starburst sample plotted against stellar mass. This shows a moderate correlation (weighted Pearson coefficient of 0.55). Again, the high redshift system MACS1149-JD1 is highlighted in red, but does not appear to be inconsistent with the rest of the sample.

Open with DEXTER
In the text
thumbnail Fig. B.1.

Inverse-Compton protogalaxy X-ray luminosity evolution over redshift, as calculated using the approach detailed in this paper with ℛSN ≈ 1 yr−1 (semi-analytical Monte-Carlo simulation) compared to the earlier study, Schober et al. (2015) (dashed grey line), which adopts an analytical model. Our results are consistent with those obtained by Schober et al. (2015).

Open with DEXTER
In the text
thumbnail Fig. D.1.

Star-formation rates in the sample of 12 post-starburst galaxies. MACS1149-JD1 is highlighted in red in all panels as this galaxy is at substantially higher redshift (z = 9.11) than the others in the sample. Upper panel: observed quiescent star-formation rate ℛSF plotted against the length of quiescence period. This period is measured up to the point of observation, so tqui effectively indicates the time elapsed since the end of starburst activity. This shows a strong positive correlation, with a Pearson coefficient (weighted by the error bounds) of 0.93. Middle panel: as above, except the y-axis is star-formation rate surface density (). This removes the correlation seen in the upper panel, with a Pearson coefficient of −0.02 (when removing the clear outlier A1689-zD1, as labelled). Lower panel: y-axis as specific star-formation rate, which is adjusted for galaxy mass where SSFR = ℛSF/M*. If removing the labelled outlier A1689-zD1, this now gives a non-negligible negative correlation of weighted Pearson coefficient −0.64.

Open with DEXTER
In the text
thumbnail Fig. D.2.

As Fig. D.1, but for the starburst phase of each system. The weighted Pearson correlation value for the top panel is found to be 0.71, while the lower panel indicates two behaviours with a less clear trend (overall correlation of 0.39) – see text for details. As SSFR is estimated in the burst phase using tSB and M* rather than being determined independently, it is not physically meaningful to include it here.

Open with DEXTER
In the text

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

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

Initial download of the metrics may take a while.