Dark-ages Reionization and Galaxy Formation Simulation - XIII. AGN quenching of high-redshift star formation in ZF-COSMOS-20115

TL;DR

This study uses a semi-analytic model to explore the formation and quenching of a massive quiescent galaxy at high redshift, highlighting the role of early mergers, black hole growth, and feedback in halting star formation.

Contribution

It provides a detailed simulation-based explanation for the early quenching of a high-redshift massive galaxy, linking it to merger history and black hole activity.

Findings

ZF-COSMOS-20115 likely hosted by a massive halo (~10^13 M_sun)

Early mergers drove nuclear growth and quenched star formation

Analogues were luminous quasars at z~5, possibly progenitors

Abstract

Massive quiescent galaxies (MQGs) are thought to have formed stars rapidly at early times followed by a long period of quiescence. The recent discovery of a MQG, ZF-COSMOS-20115 at $z\sim4$, only 1.5 Gyr after the big bang, places new constraints on galaxy growth and the role of feedback in early star formation. Spectroscopic follow-up confirmed ZF-COSMOS-20115 as a MQG at $z=3.717$ with an estimated stellar mass of ${\sim}10^{11}\mathrm{M}_\odot$, showing no evidence of recent star formation. We use the Meraxes semi-analytic model to investigate how ZF-COSMOS-20115 analogues build stellar mass, and why they become quiescent. We identify three analogue galaxies with similar properties to ZF-COSMOS-20115. We find that ZF-COSMOS-20115 is likely hosted by a massive halo with virial mass of ${\sim}10^{13}\mathrm{M}_\odot$, having been through significant mergers at early times. These merger events drove intense growth of the nucleus, which later prevented cooling and quenched star formation. ZF-COSMOS-20115 likely remained quiescent at $z<3.7$. We find that the analogues host the most massive black holes in our simulation and were luminous quasars at $z\sim5$, indicating that ZF-COSMOS-20115 and other MQGs may be the descendants of high-redshift quasars. In addition, the model suggests that ZF-COSMOS-20115 formed in a region of intergalactic medium that was reionized early.

Figures

Tables

Table 1. Table 1: A list of the properties of the ZF-COSMOS-20115 analogues at z ∼ 3.7 similar-to 𝑧 3.7 z\sim 3.7 and the time of peak star formation. The estimated properties of ZF-COSMOS-20115 by Glazebrook et al . ( 2017 ) are shown for comparison.

Name	$M_{\ast }\left(z=3.7\right)M_{\odot }$	$\mathrm{SFR}\left(z=3.7\right)/M_{\odot }{\mathrm{yr}}^{-1}$	$M_{\mathrm{vir}}\left(z=3.7\right)/M_{\odot }$	$z_{\mathrm{peak}}$	$\mathrm{SFR}\left(z=z_{\mathrm{peak}}\right)/M_{\odot }{\mathrm{yr}}^{-1}$
ZF-1	$8.92\times {10}^{10}$	0.0	$6.43\times {10}^{12}$	5.1	574.2
ZF-2	$1.02\times {10}^{11}$	0.0	$5.46\times {10}^{12}$	4.9	827.2
ZF-3	$9.38\times {10}^{10}$	0.0	$1.12\times {10}^{13}$	5.7	1000.1
ZF-COSMOS-20115a	$\left(1.46-1.82\right)\times {10}^{11}$		$3.0\times {10}^{12}$	$5-8$	$>990$

Full text

††pagerange: Dark-ages Reionization and Galaxy Formation Simulation - XIII. AGN quenching of high-redshift star formation in ZF-COSMOS-20115–References††pubyear: 2017

1 Introduction

Massive quiescent galaxies (MQGs) are galaxies with stellar masses of the order of $10^{11}\textrm{M}_{\odot}$, and low or null star formation rate (SFR). It is thought that these objects formed rapidly at early times, followed by a long state of quiescence (Ilbert et al., 2013). Due to the absence of recent star-forming events in MQGs, they are observationally challenging to study, especially at high redshift. Using the ultra deep imaging from the FourStar Galaxy Evolution Survey111http://zfourge.tamu.edu/ (ZFOURGE), Straatman et al. (2014) identified 15 MQGs at $z=3.4-4.2$ with an estimated number density of MQGs (hereafter the S14 sample) of $(1.8\pm 0.7)\times 10^{-5}\mathrm{Mpc}^{-3}$. Despite the large uncertainties, this number density is significantly higher than the value expected from extrapolations using observations at lower redshifts (Bell et al., 2003; Muzzin et al., 2013), leading to an unexpectedly high MQG fraction of $34\pm 13$ per cent at $z\sim 4$.

Recently, spectroscopic follow-up (Glazebrook et al., 2017) was performed on the brightest MQG in the S14 sample222The ID of ZF-COSMOS-20115 is 13172 in the S14 sample., ZF-COSMOS-20115. Their analysis revealed that ZF-COSMOS-20115 has a stellar mass of ${M_{*}=1.7^{+0.12}_{-0.24}\times 10^{11}\textrm{M}_{\odot}}$ (see a more recent study of Simpson et al. 2017), and redshift of $3.717\pm 0.001$, but no detectable ongoing star formation based on Balmer absorption lines (i.e. $\textrm{H}\beta$, $\textrm{H}\gamma$ and $\textrm{H}\delta$), implying a current $\mathrm{SFR}<0.2\textrm{M}_{\odot}\mathrm{yr}^{-1}$. Furthermore, modelling the spectral evolution of the stellar population (Fioc & Rocca-Volmerange, 1999) suggests rapid growth of stellar mass with a $\mathrm{SFR}>990\textrm{M}_{\odot}\mathrm{yr}^{-1}$ at the peak of activity and a formation time-scale less than 250 Myr.

In this paper we use the Meraxes semi-analytic model (Mutch et al., 2016a) within the Dark-ages Reionization And Galaxy formation Observables from Numerical Simulations (DRAGONS333http://dragons.ph.unimelb.edu.au) programme to investigate the galaxy formation history of MQGs at high redshift. We note that cosmological simulations have had difficulty producing MQGs at high redshift (Lee & Yi, 2013; Wellons et al., 2015; Davé et al., 2016; Behroozi & Silk, 2016). However, we identify three galaxies with properties similar to ZF-COSMOS-20115.

This paper is organized as follows. We begin with a brief overview of the DRAGONS framework in Section 2 and present the modelled galaxy property in Section 3. We show the ZF-COSMOS-20115 analogues in our model in Section 4 and discuss the history and future of ZF-COSMOS-20115 in Sections 5 and 6. Conclusions are given in Section 7. In this work, we adopt cosmological parameters from the Planck 2015 results ($\Omega_{\mathrm{m}},\Omega_{\mathrm{b}},\Omega_{\mathrm{\Lambda}},h,\sigma_{8},n_{\mathrm{s}}$ = 0.308, 0.0484, 0.692, 0.678, 0.815, 0.968; Planck Collaboration et al. 2015).

2 DRAGONS

The Meraxes semi-analytic model (Mutch et al., 2016a) was specifically designed to study galaxy formation at high redshift and the epoch of reionization (EoR). In this work, we use the updated version of Meraxes (Qin et al., 2017), which includes a detailed prescription of black hole growth and AGN feedback. Our fiducial model was run on the Tiamat-125-HR dark matter halo merger trees. We briefly review the model here and refer the interested reader to the aforementioned references for details.

The underlying Tiamat-125-HR dark matter halo merger trees provide 100 snapshots between $z=35$ and 5 with a time interval of ${\sim}11.1$ Myr and 114 additional snapshots between $z=5$ and 0.56 separated equally in units of Hubble time (Poole et al. in preparation). The particle mass resolution of Tiamat-125-HR is ${\sim}1.2\times 10^{8}h^{-1}\textrm{M}_{\odot}$ and the box size is $125h^{-1}$Mpc. In order to show convergence, Meraxes was also run using the Tiamat halo merger trees (Poole et al., 2016), which were constructed from an N-body simulation sharing identical cosmology with Tiamat-125-HR but at a higher mass resolution ($2.6\times 10^{6}h^{-1}\textrm{M}_{\odot}$) and smaller volume ($67.8h^{-1}$Mpc). These results are presented in Section 5 for comparison.

The Meraxes semi-analytic model consists of a number of important astrophysical processes, including gas infall, cooling, star formation, supernova feedback, AGN feedback, metal enrichment, mergers and reionization (Mutch et al., 2016a; Qin et al., 2017). The model was calibrated against the observed stellar mass function at z\sim 7-{\color[rgb]{0,0,0}0.56} (see the galaxy stellar mass function at $z=2-4$ in the top panels of Fig. 1), black hole mass function at $z\lesssim 0.5$, ionizing emissivity at $z\sim 5-2$ and the Thomson scattering optical depth444The model running with the Tiamat-125-HR merger trees is limited by the resolution and underestimates the ionizing emissivity at higher redshifts (see more in Qin et al. 2017).. The model is in agreement with the observed Magorrian relation at z<{\color[rgb]{0,0,0}0.56} as well as the observed quasar luminosity function across a large redshift range (z\sim 6-{\color[rgb]{0,0,0}0.56}) when an opening angle of 80 deg is chosen to account for obscuration by dust.

Compared to other simulations (e.g. the Millennium simulation; Rong et al. 2017), several features of this model are well suited to study of high-redshift MQGs. We note that the cadence of our model is about 11 Myr at $z\geq 5$ and reaches 30 Myr at $z=3.7$. This high temporal resolution enables us to not only resolve the dynamical time of the disc and simulate the bursty nature of star formation events, but also investigate the evolution of galaxies in more detail. In addition, the agreement between our model and observations, including the stellar mass function at z\sim 7-{\color[rgb]{0,0,0}0.56} and quasar luminosity function at z\sim 6-{\color[rgb]{0,0,0}0.56}, suggests that the model galaxy and AGN catalogues we have constructed are able to represent the relevant observables across cosmic time. This is crucial for investigation of individual galaxy analogues within a cosmological context (e.g. Mutch et al., 2016b; Waters et al., 2016).

3 Galaxy mass function

In the top panels of Fig. 1 we plot the evolution of the model galaxy stellar mass function between $z=2$ and 4. The agreement with observed mass functions demonstrates the ability of the model to reproduce the observed growth of stellar mass across the broad range of cosmic time relevant to this work (Qin et al., 2017). The bottom left panel of Fig. 1 presents the predicted cumulative galaxy stellar mass function at $z=3.7$ for four different samples in our model: i) all galaxies; ii) $\mathrm{SFR}<100\textrm{M}_{\odot}\mathrm{yr}^{-1}$; iii) $\mathrm{SFR}<30\textrm{M}_{\odot}\mathrm{yr}^{-1}$ and iv) $\mathrm{SFR}<1\textrm{M}_{\odot}\mathrm{yr}^{-1}$. We note that there are 15 MQGs detected in the S14 sample, corresponding to a number density of $(1.8\pm 0.7)\times 10^{-5}\mathrm{Mpc}^{-3}$. This is shown in Fig. 1 for comparison. In addition, SFRs derived from stellar population modelling vary from 0 to ${\sim}32\textrm{M}_{\odot}\mathrm{yr}^{-1}$ in the S14 sample. We see that the number density of MQGs predicted by our model with SFRs in this range is in good agreement with the S14 sample. We also show the evolution of the number density of massive galaxies and the MQG fraction in Fig. 1. The model is consistent with observations at lower redshifts ($z<3$; Straatman et al. 2015).

We note that, despite the large uncertainty, observations may suggest that the MQG fraction stops decreasing at $z\sim 4$ while a continuously decreasing trend towards higher redshifts is usually presented in theoretical models (Glazebrook et al. 2017; see the bottom right panel of Fig. 1). The MQG fraction in the S14 sample is $34\pm 13$ per cent at $z\sim 4$. This was compared to Illustris (Vogelsberger et al., 2014), a suite of hydrodynamical simulations of galaxy formation. Only ${\sim}20$ per cent of massive galaxies in Illustris are MQGs (in agreement with our model, see the $\mathrm{SFR}<30\textrm{M}_{\odot}\mathrm{yr}^{-1}$ line in the bottom right panel of Fig. 1). Based on this, Straatman et al. (2015) suggested that quenching is likely to occur at early times in massive galaxies. However, we note that only one object in the S14 sample, ZF-COSMOS-20115, has been spectroscopically confirmed to be a MQG and, due to the limited volume of the ZFOURGE survey and potential contaminations from dust obscured star-forming galaxies (Straatman et al., 2015), the MQG population at high redshift is not well constrained. Therefore, instead of making attempts to interpret the high MQG fraction at high redshift, in this work, we focus on finding analogues of the spectroscopically confirmed ZF-COSMOS-20115 and discuss the possible evolutionary history and future of this galaxy in order to identify the potential quenching mechanism and possible progenitors of high-redshift MQGs.

4 ZF-COSMOS-20115 analogues

In order to identify the analogues of ZF-COSMOS-20115 in our model, we consider the following two properties.

The stellar mass of ZF-COSMOS-20115 was constrained (Glazebrook et al., 2017) from the equivalent width of the observed spectral energy distribution (SED) and is considered a robust observed property. However, Simpson et al. (2017) recently argued that due to contaminations from a nearby dusty starburst, the stellar mass of ZF-COSMOS-20115 was overestimated by 50 per cent. Therefore, according to Simpson et al. (2017), we set a threshold of M_{*}>{\color[rgb]{0,0,0}0.8}\times 10^{11}\textrm{M}_{\odot}, which returns 139 massive galaxies. Based on the Balmer absorption lines (or the lack of H$\beta$ emission line), the current SFR was reported to be less than 0.2 (or 4)$\textrm{M}_{\odot}\mathrm{yr}^{-1}$. We note that in the alternate case that ZF-COSMOS-20115 is an obscured star-forming galaxy, the SFR would be less than $70-100\textrm{M}_{\odot}\mathrm{yr}^{-1}$, which is inferred from the non-detection threshold of the Herschel/PACS imaging (Straatman et al., 2014). In this work, we focus on the MQG scenario and select ZF-COSMOS-20115 analogues with $\mathrm{SFR}<0.4\ \textrm{M}_{\odot}\mathrm{yr}^{-1}$. This selection returns 10 MQGs, corresponding to 7 per cent of the 139 massive galaxies. 2. 2.

Using spectral evolution models of stellar population (Fioc & Rocca-Volmerange, 1999), and the assumption that the stellar mass of ZF-COSMOS-20115 increases with a constant SFR for a certain period followed by a long state of quiescence, Glazebrook et al. (2017) found that ZF-COSMOS-20115 has a stellar age of $500-1050$ Myr with a formation time-scale of less than 250 Myr. This implies that ZF-COSMOS-20115 formed at $z\sim 5.8^{+2.3}_{-0.8}$ with a star-forming period of $\Delta z<1.5^{+2.1}_{-0.4}$. We note that strong star formation at $z<4.5$ creates significantly excess flux in the spectrum at wavelengths shorter than $\mathrm{H}_{\alpha}$ compared to the observed SED of ZF-COSMOS-20115. Therefore, we exclude galaxies with $\mathrm{SFR}>100{\textrm{M}_{\odot}}\mathrm{yr}^{-1}$ at $z<4.5$ from the sample.

Based on these specific criteria, we identify three analogues which we will henceforth refer to as ZF-1 to ZF-3. However, we note that the choice of the stellar mass and SFR thresholds does not have a significant impact to our conclusion.

The ZF-COSMOS-20115 analogues in our model have stellar masses of around $10^{11}\mathrm{M}_{\odot}$ and negligible star formation at $z\sim 3.7$ with peak SFRs $\geq 500\mathrm{M}_{\odot}\mathrm{yr}^{-1}$ at $z\sim 5-6$. We summarize their properties at $z\sim 3.7$ and at the time of peak star formation in Table 1, for comparison to the ZF-COSMOS-20115 properties estimated by Glazebrook et al. (2017). We see that, due to the limited volume of our parent N-body simulation, the candidates are slightly less massive (by up to 30 per cent) than ZF-COSMOS-20115. However, their halo masses (${\sim}10^{13}\textrm{M}_{\odot}$) are larger than the value ($10^{12.5}\textrm{M}_{\odot}$) estimated from the galaxy number density and the halo mass function at $z=4$. This implies that the three analogues have larger stellar mass to halo mass ratios by factors of $2-3$, suggesting that these MQGs are hosted in relatively large haloes.

5 The history of ZF-COSMOS-20115

As discussed above, owing to the mass of ZF-COSMOS-20115 and the absence of recent star formation, there must have been a period of intense star formation at early times (Glazebrook et al., 2017). In order to identify why this star formation was subsequently quenched, we investigate the history of the three analogues by tracking their most massive progenitors along the dark matter halo merger trees.

5.1 The role of AGN quenching in ZF-COSMOS-20115

Fig. 2 shows the histories of the three analogues at $z\sim 15-3.5$ (corresponding to a lookback time of $t\sim 13.5-12.0$ Gyr). We show, from top to bottom with thick black lines and star symbols representing the fiducial model, the SFR, stellar mass, cold gas mass, cooling mass, virial mass, central black hole mass, heating mass due to AGN feedback, intrinsic quasar UV magnitude555Note that quasar UV magnitudes induced by accretion events that cannot support the quasar activity for the entire snapshot are not shown here. and baryonic merger ratio666The baryonic merger ratio is defined as the ratio of the total mass of cold gas and stars of the merging companion to that of the target galaxy.. In the background, with thin grey lines and squares we also plot the histories of all massive star-forming galaxies with M_{*}>{\color[rgb]{0,0,0}0.8}\times 10^{11}\textrm{M}_{\odot} and $\mathrm{SFR}>70\textrm{M}_{\odot}\mathrm{yr}^{-1}$. The derived properties and limits of ZF-COSMOS-20115 are shown as circles for comparison. There are several points to note.

The analogues and the massive star-forming galaxies shown in the background of Fig. 2 are selected at $z=3.7$. As shown in the histories of virial mass, the analogues are hosted by relatively massive haloes compared to star-forming galaxies with similar stellar masses, with masses of $M_{\rm vir}\sim 10^{13}{\textrm{M}_{\odot}}$. 2. 2.

Mergers trigger the most intense star formation event in the history of each of ZF-1, ZF-2 and ZF-3 at $z=5.1$, 4.9 and 5.7 with $\mathrm{SFRs}\sim 600$, 800 and 1000 $\textrm{M}_{\odot}\mathrm{yr}^{-1}$, and merger ratios ${\sim}0.18$, 0.34 and 0.14, respectively. However, unlike the star formation history modelling performed in Glazebrook et al. (2017), where the galaxy is assumed to form stars with a constant SFR of $>990{\textrm{M}_{\odot}}\mathrm{yr}^{-1}$ for only $<250$ Myr, the three analogues have longer star formation time-scales of ${\sim}500$ Myr to 1 Gyr with milder evolutions of SFR. 3. 3.

After the mergers, cooling is significantly suppressed and the galaxies consume the available cold gas on short time-scales ($\sim 100-300$ Myr). When there is insufficient cold gas in the galaxies ($m_{\mathrm{cold}}<m_{\mathrm{crit}}$, see equation 6 in Mutch et al. 2016a), star formation is quenched and the analogues become fainter in the UV band (the UV magnitude of the host galaxy is shown by the thin black line in the $M_{1450}$ panels of Fig. 2). 4. 4.

Mergers drive black hole growth in our model. The central massive black holes of the three galaxies therefore grow significantly after the merger event. At the time when cooling stops, their black hole masses become $1-1.5$ orders of magnitude larger than in star-forming galaxies with the same stellar mass. During the accretion phase, the central black hole radiates energy and heats the surrounding interstellar medium (ISM). In our model, AGN feedback limits gas condensation by suppressing the cooling flow (Croton et al., 2006). When the black hole becomes massive enough, AGN feedback is able to completely neutralize the cooling mass, leading to quenching of star formation in the three analogues. 5. 5.

Because of the dramatic accretion following a significant merger event, the nucleus is triggered into a quasar phase (Qin et al., 2017). However, the finite accreted mass can only support central massive black hole activity for a finite period of time. In the second last row of Fig. 2, we show the intrinsic quasar UV magnitude (solid thick black lines) in comparison with the magnitude of the host galaxy (solid thin black lines). We see that, before the merger, the total UV magnitude is dominated by stellar light, but that the accretion briefly dominates UV flux during the quasar phase following the merger. However later, at $z=3.7$, the nucleus has become inactive (but still provides AGN feedback through the radio mode, see more details in Qin et al. 2017). This explains the lack of significant UV to optical radiation from the central massive black hole in the observed SED of ZF-COSMOS-20115, which can be well fit with a stellar population model (Glazebrook et al., 2017) without AGN contributions.

5.1.1 Effect of the Eddington ratio

In our fiducial model, black holes are assumed to either accrete and radiate at the Eddington rate ($\epsilon=1$, where $\epsilon$ is the Eddington ratio) or stay quiescent if there is no accretion. This assumption has been shown to provide a good description of black hole growth for the majority of black holes at the relevant redshifts777In our modelling we find that assuming $\epsilon=1$ results in an underestimate of the number of faint AGN at later times ($z<2$) relative to observations. However, Eddington accretion with a value of $\epsilon=1$ is necessary in our model in order to produce a realistic AGN catalogue for the relevant mass range ($M_{\mathrm{BH}}\lesssim 10^{8}\textrm{M}_{\odot}$) and redshifts ($z>3$) (see more in Qin et al. 2017). (Bonoli et al., 2009). In order to directly investigate the importance of the Eddington ratio in determining whether our three analogues are quenched by AGN feedback, we set the Eddington ratio to \epsilon={\color[rgb]{0,0,0}0.3} (dashed thick black lines in Fig. 2) and $\epsilon=0.1$ (dotted thick black lines). We note that without AGN feedback, in the Tiamat-125-HR volume the model can still reproduce the observed galaxy stellar mass function at $z>4$ but overestimates the number of more massive galaxies at $z<2$ (Qin et al., 2017). With $\epsilon=0.1$, we find that while the number of massive galaxies with M_{*}>{\color[rgb]{0,0,0}0.8}\times 10^{11}\textrm{M}_{\odot} increases to 198, the fraction of MQGs decreases significantly at high redshift (see the thin dotted line in the bottom right panel of Fig. 1). Moreover, amongst these galaxies, only five are quiescent galaxies. However, they cannot be considered analogues due to intense recent star formation events at $z<4.5$ that violate the selection criterion 2. Whilst no ZF-COSMOS-20115 analogues are found in the case of sub-Eddington accretion, it is illuminating to look at the variation in the evolution of ZF-1, ZF-2 and ZF-3 with lower black hole growth efficiencies.

We see that when the central massive black hole is sub-Eddington, the black hole growth becomes slower, providing weaker feedback to galaxy formation and having a longer quasar phase with a fainter luminosity. With \epsilon={\color[rgb]{0,0,0}0.3}, the central black holes are inactive at $z=3.7$ but the analogues form stars for a longer period down to $z<4.5$. This violates selection criterion (ii) and results in a stronger UV excess in the spectrum at wavelengths shorter than $\mathrm{H}_{\alpha}$ compared to ZF-COSMOS-20115. With a lower Eddington ratio of $\epsilon=0.1$, black hole growth becomes significantly slower compared to the model with Eddington accretion, and the heating mass due to AGN becomes insignificant. We find that with $\epsilon=0.1$, the analogues are star-forming galaxies at $z=3.7$. Our model therefore predicts that the intense growth of a massive black hole ($\epsilon\sim 1$) in the centre of ZF-COSMOS-20115, likely induced by mergers, resulted in persistent AGN feedback which quenched subsequent star formation.

5.1.2 Evidence for AGN feedback

Our model suggests the existence of ${>}10^{8}{\textrm{M}_{\odot}}$ supermassive black holes and a significant role of AGN feedback in galaxy formation at $z\sim 5-6$ when the Universe was less than ${\sim}1$ Gyr old. In particular, we find that the three ZF-COSMOS-20115 analogues in our model are quenched due to persistent AGN feedback. This could be potentially examined using far-infrared (FIR) and radio observations. For instance, Gobat et al. (2017) constructed a median SED in the range of $10-10^{6}\mu\mathrm{m}$ using stacked images of ${\sim}1000$ MQGs at $z\sim 2$. They discovered significant emissions in the FIR and radio bands of MQGs, which suggests a large content of dust (${\sim}10^{8}\textrm{M}_{\odot}$) and gas (${\sim}10^{10}\textrm{M}_{\odot}$). This is in contrast to local early-type galaxies that are usually found to be gas poor. If present, this gas must be consumed with a low efficiency in order to keep MQGs quiescent across cosmic time. This is likely due to AGN feedback and is supported by a significant excess of radio emission (${\sim}5\times 10^{22}\mathrm{W/Hz}$) in the stacked SED compared to a dust only spectral model.

Conclusive evidence for this becomes observationally more challenging at $z\sim 4$. We note that the near-infrared (NIR) emission of ZF-COSMOS-20115 is less than $3\mu\mathrm{Jy}$ (Spitzer/MIPS888http://www.spitzer.caltech.edu/ 24 $\mu\mathrm{m}$) and no mid-infrared detection was reported by Herschel/PACS999http://sci.esa.int/herschel/ (MIR, $100-160\mu\mathrm{m}$, Straatman et al. 2015). In addition, within a 20 arcsec box around ZF-COSMOS-20115 there is no radio source in the COSMOS VLA Deep catalogue101010http://irsa.ipac.caltech.edu/Missions/cosmos.html (Schinnerer et al., 2010), the typical rms and resolution of which are 7.5 $\mu\mathrm{Jy}$ and 1.7 arcsec. Next-generation telescopes are expected to provide enormously improved sensitivity and resolution covering a wide range of wavelengths.

We summarize the detection capabilities of JWST111111https://jwst.nasa.gov/, WFIRST121212https://wfirst.gsfc.nasa.gov/, ASKAP131313http://www.atnf.csiro.au/projects/askap/index.html and SKA141414http://skatelescope.org/ in Fig. 3 compared to the SED of ZF-COSMOS-20115 in the observed frame (Glazebrook et al., 2017), the sensitivities of the NIR and MIR instruments employed by ZFOURGE (Dickinson & Giavalisco, 2003; Straatman et al., 2016) and the VLA-COSMOS Deep Survey. Within the context of high-redshift quiescent galaxies, JWST will provide broad-band imaging from 0.7 to 28.5 $\mu\mathrm{m}$ with sensitivities of ${\sim}10^{-2}\mu\mathrm{Jy}$ to $1-10\mu\mathrm{Jy}$ (3h exposure) from NIR to MIR (Greene et al., 2007; Rieke, 2011; Bouchet et al., 2015). With a sufficient high-redshift MQG sample becoming available from JWST, radio flux can be estimated using the stacking technique with telescopes such as ASKAP (Norris et al., 2011). In addition, the high latitude survey of WFIRST (Spergel et al., 2013; Gehrels & Spergel, 2015) will cover 6 per cent of the entire sky, a region that is 20000 times larger than ZFROUGE with AB magnitude limits of $\sim 27$ in the NIR ($0.7-2\mu\mathrm{m}$). Together with the planned SKA-SUR all sky survey (Kapinska et al., 2015), this will enable direct evaluation of AGN activity in MQGs.

5.2 The possible connection between MQGs and high-redshift quasars

Fig. 2 shows that the UV fluxes of the three analogues for most of their histories are dominated by stellar light, and that at $z\sim 3.7$ their central nuclei have become inactive. However, the central massive black holes of the three analogues have previously grown significantly (${\sim}10^{8}\textrm{M}_{\odot}$) during the quasar phase at $z\sim 5$ when accretion becomes important to or even dominates the total UV luminosity. Fig. 4 shows the correlation between black hole mass and stellar mass of the three analogues at $z=4.9$, 4.5 and 5.3 when their AGN activities reached a peak. Note that the model has been shown to successfully reproduce the observed Magorrian relation at $z\lesssim 0.5$ (Qin et al., 2017). We see that at $z\sim 5$, the three analogues have larger black holes than other galaxies in our model, and in particular, that ZF-3 hosts the most massive black hole in the model at $z=5.3$, and has a luminosity corresponding to that of a SDSS quasar (Glikman et al., 2011; Shen & Kelly, 2012; McGreer et al., 2013).

Indeed at the time of peak AGN activity, the quasar UV luminosities of ZF-1, ZF-2 and ZF-3 increase to $M_{1450}=-23.0,-24.7$ and $-24.8$, up to ${\sim}2.5$ magnitudes brighter than the host galaxies. We show the predicted quasar luminosity function at $z=5$ in the top panel of Fig. 5. The result using the Tiamat halo merger trees is shown for comparison and to indicate the convergence of the quasar luminosity function at $z\sim 5$. The model is in agreement with observations at the bright end (Shen & Kelly, 2012; McGreer et al., 2013). We indicate the quasar luminosities of the three analogues in the top panel of Fig. 5 and find that the progenitors of ZF-1, ZF-2 and ZF-3 are bright quasars in the range that optical surveys such as SDSS are able to detect. This suggests that ZF-COSMOS-20115 is likely the descendent of a high-redshift quasar and illustrates the possible connection between MQGs and high-redshift quasars (Maiolino et al., 2012).

In order to better demonstrate this, for each massive galaxy selected at $z=5.3$, we present the correlation between black hole mass and subsequent SFR at $z=3.7$ in the bottom panel of Fig. 5. ZF-3 is indicated with a filled circle. We see that galaxies with less massive black holes ($M_{\mathrm{BH}}<10^{7.5}\textrm{M}_{\odot}$) can have a range of subsequent star formation levels. This is due to self-regulated stellar feedback being the dominant mechanism in these galaxies. However, due to persistent AGN feedback, galaxies with more massive black holes do not have intense subsequent star formation. This indicates that a high-redshift bright quasar will likely become a MQG at later times, when its accretion disc has been completely consumed.

5.3 The reionization history

The three ZF-COSMOS-20115 analogues are hosted by massive haloes, which provide deep gravitational potentials that efficiently accrete baryons. The high baryon accretion efficiency induces intense star formation at high redshift, providing a significant number of ionizing photons that reionize surrounding hydrogen at a very early time. In this section, we discuss the reionization history of ZF-COSMOS-20115.

In Fig. 6, we show the evolution of the average neutral hydrogen fraction within boxes with different sizes around the three analogues. The neutral hydrogen fraction is calculated using 21cmFAST (Mesinger et al., 2011) with the galaxy catalogue provided by Meraxes (see more details in Mutch et al. 2016a). We see that the analogues start ionizing the intergalactic medium (IGM) at $z\gtrsim 14$ with inner regions becoming ionized first. For example, ZF-1, ZF-2 and ZF-3 (and their nearby galaxies) reionize 3.6Mpc regions at $z\sim 9.5$, 9 and 10.5, respectively, and 9.4Mpc regions by $z\sim 7$, 6.5 and 8. We also present the average neutral hydrogen fraction of the entire Tiamat-125-HR box in Fig. 6.

Due to the simulation not resolving dwarf galaxies, reionization only finishes at $z<5$ in the Tiamat-125-HR volume. We therefore also show the result calculated using the higher resolution Tiamat halo merger trees for comparison. This higher resolution model has been shown (Qin et al., 2017) to match the observed ionizing emissivity at $z\sim 5-2$ (Becker & Bolton, 2013) and the Thomson scattering optical depth (Planck Collaboration et al., 2016). Comparing the average evolution in ionized fraction between the two simulations illustrates that reionization predicted are delayed by $\Delta z\sim 1$ in the low-resolution simulation. The model suggests that ZF-COSMOS-20115 is located in a region of IGM that was ionized early in the reionization era.

6 the subsequent evolution of ZF-COSMOS-20115

In this section, we use the Meraxes semi-analytic model to explore possible scenarios for the subsequent evolution of ZF-COSMOS-20115 at $z<3.7$. Fig. 7 shows the future of the three analogues from $z\sim 4$ to z\sim{\color[rgb]{0,0,0}0.56}. We see that AGN feedback keeps gas in the analogues from cooling, and that without a major merger it is unlikely that a galaxy like ZF-COSMOS-20115 will have a further strong star formation event at $z<3.7$. However, depending on the halo properties or merger ratios, starbursts can in some cases be triggered. Because ZF-COSMOS-20115 is likely to be cold gas poor due to AGN heating, subsequent starburst events would only persist for a short time period. We also see the following trends.

When there is a merger, the three analogues increase their stellar masses by absorbing baryons from their merging companions. The UV flux changes due to different star-forming histories of the merging companion (most likely the companion is a star-forming galaxy and therefore the UV flux increases after the merger), as well as the merger-triggered starburst. However, depending on the merger ratio, only a small fraction of mergers can trigger additional large starbursts. For example, ZF-1 stays quiescent until $z\sim 1.1$ when several galaxies merge with it, inducing a minor star-forming event with a $\mathrm{SFR}\sim 25{\textrm{M}_{\odot}}\mathrm{yr}^{-1}$. 2. 2.

With the central supermassive black hole heating the ISM, continuous star formation is unlikely to occur. However, there are three short ($\Delta z<0.1$) star formation events in ZF-3 at $z\sim 2-3$, when the cold gas mass, which gradually increases from stellar mass recycling (Mutch et al., 2016a), exceeds the critical value. We note that the star formation law in our model depends on a critical mass, above which galaxies are able to form stars. This critical mass is inferred from the host halo properties and increases with virial mass and halo angular momentum (based on Kennicutt 1998 and Kauffmann 1996). While the halo mass of ZF-3 grows smoothly, the spin parameter, which is a measure of specific angular momentum, decreases by a factor of 2 during the star-forming events. Assuming full conservation of specific angular momentum for newly cooling gas (Mutch et al., 2016a), the cold gas disc contracts, resulting in higher gas densities. Under circumstances where the cold gas density becomes larger than the critical value, star formation can occur. 3. 3.

ZF-2 merges151515Note that using a central weighting algorithm, the Tiamat-125-HR halo merger trees were constructed by tracking halo cores instead of the majority of halo particles. Since the core is where a galaxy forms and evolves, this strategy improves the semi-analytic galaxy formation modelling. In this work, subfind (Springel et al., 2008), a standard halo finder, is adopted, and in some circumstances subfind redistributes the dark matter particles of a central halo to its satellites and leads to a central halo not being the most massive halo (central-satellite switching, see more details in Poole et al. in preparation). This results in unrealistic changes in halo properties. For clarity, we do not show the properties when this occurs and we note that this does not have any impact on our conclusions or predictions of the future of ZF-COSMOS-20115. with a larger galaxy at $z\sim 2.5-1.2$. The host halo of ZF-2 initially becomes smaller due to tidal stripping from the merging companion, and eventually merges into the more massive halo with a baryonic merger ratio of ${\sim}1.4$. When the halo is being stripped, the host galaxy loses its gas component which becomes unbound. Consequently, without gas falling into the centre, the central massive black hole stops accreting and there is no AGN feedback or heating in ZF-2 during the stripping period. On the other hand, the cold gas disc is considered to be gravitationally bound in our model and remains intact during galaxy stripping. When the halo mass becomes smaller, the disc size and hence the critical mass required for star formation decreases. This results in star formation activity which quickly consumes the available gas.

7 conclusions

ZF-COSMOS-20115 is a MQG at $z\sim 3.7$. We model MQGs using the Meraxes semi-analytic model and identify three ZF-COSMOS-20115 analogues in the fiducial model presented in Qin et al. (2017). These analogues have properties that are in agreement with the observed constraints on ZF-COSMOS-20115 inferred from the recent spectroscopic follow-up of Glazebrook et al. (2017).

We find that the three analogues are hosted by more massive haloes ($10^{13}\textrm{M}_{\odot}$) compared to other star-forming galaxies with similar stellar masses. Following the most massive progenitor, we track the history of the three analogues and identify significant merger events at $z\sim 5$. We find that when the mergers drive intense Eddington-limited growth of the central massive black hole, the cooling flow is significantly suppressed by the resulting AGN feedback. We further investigate scenarios when black holes accrete at sub-Eddington rates and find insufficient heating energy from AGN feedback to suppress star formation. In particular, with black holes accreting at a third of the Eddington rate, the ZF-COSMOS-20115 analogues still have high SFRs at $z\sim 4.5$, while, with the Eddington ratio of 0.1, they are star-forming galaxies at $z\sim 3.7$. Our model therefore suggests that there was a period when the central massive black hole grew rapidly in ZF-COSMOS-20115, probably triggered by mergers, and that the persistent feedback from AGN quenched the subsequent star formation.

In addition, we find the three analogues host the most massive black holes in our simulation, and that they were luminous quasars ($M_{1450}\lesssim-23$) at $z\sim 5$. Moreover, all galaxies with massive black holes at $z\sim 5$ have quenched star formation by $z\sim 3.7$. This suggests that ZF-COSMOS-20115 is the descendent of a high-redshift quasars similar to those detected in SDSS. By investigating the ionizing regions around the three analogues, we find that ZF-COSMOS-20115 formed in a region that was reionized early. We also follow the future of our ZF-COSMOS-20115 analogues down to z={\color[rgb]{0,0,0}0.56} and find that further strong or continuous star formation events are unlikely to occur in ZF-COSMOS-20115 at $z<3.7$.

Acknowledgements

This research was supported by the Victorian Life Sciences Computation Initiative (VLSCI), grant ref. UOM0005, on its Peak Computing Facility hosted at the University of Melbourne, an initiative of the Victorian Government, Australia. Part of this work was performed on the gSTAR national facility at Swinburne University of Technology. gSTAR is funded by Swinburne and the Australian Governments Education Investment Fund. This work was supported by the Flagship Allocation Scheme of the NCI National Facility at the ANU, generous allocations of time through the iVEC Partner Share and Australian Supercomputer Time Allocation Committee. AM acknowledges support from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation program (Grant No. 638809 – AIDA).