We present the Grackle chemistry and cooling library for astrophysical simulations and models. Grackle provides a treatment of non-equilibrium primordial chemistry and cooling for H, D, and He species, including H 2 formation on dust grains; tabulated primordial and metal cooling; multiple UV background models; and support for radiation transfer and arbitrary heat sources. The library has an easily implementable interface for simulation codes written in C, C++, and Fortran as well as a Python interface with added convenience functions for semi-analytical models. As an open-source project, Grackle provides a community resource for accessing and disseminating astrochemical data and numerical methods. We present the full details of the core functionality, the simulation and Python interfaces, testing infrastructure, performance, and range of applicability. Grackle is a fully open-source project and new contributions are welcome.

We use a set of seminumerical simulations based on the Zel'dovich approximation, the friends‐of‐friends algorithm and excursion set formalism to generate reionization maps of high dynamic range with a range of assumptions regarding the distribution and luminosity of ionizing sources and the spatial distribution of sinks for the ionizing radiation. We find that ignoring the inhomogeneous spatial distribution of regions of high gas density where recombinations are important – as is often done in studies of this kind – can lead to misleading conclusions regarding the topology of reionization, especially if reionization occurs in the photon‐starved regime suggested by Lyα forest data. The inhomogeneous spatial distribution of recombinations significantly reduces the mean‐free path of ionizing photons and the typical size of coherently ionized regions. Reionization proceeds then much more as an outside‐in process. Low‐density regions far from ionizing sources become ionized before regions of high gas density not hosting sources of ionizing radiation. The spatial distribution of sinks of ionization radiation also significantly affects the shape and amplitude of the power spectrum of fluctuations of 21 cm emission. The slope of the 21 cm power spectrum as measured by upcoming 21 cm experiments should be able to distinguish to what extent the topology of reionization proceeds outside‐in or inside‐out, while the evolution of the amplitude of the power spectrum with increasing ionized mass fraction should be sensitive to the spatial distribution and the luminosity of ionizing sources.

We introduce a new set of large scale, high resolution hydrodynamical simulations of the intergalactic medium: the Sherwood simulation suite. These are performed in volumes 10 3 -160 3 h −3 comoving Mpc 3 , span almost four orders of magnitude in mass resolution with up to 17.2 billion particles, and employ a variety of physics variations including warm dark matter and galactic outflows. We undertake a detailed comparison of the simulations to high resolution, high signal-to-noise observations of the Lyα forest over the redshift range 2 ≤ z ≤ 5. The simulations are in very good agreement with the observational data, lending further support to the paradigm that the Lyα forest is a natural consequence of the web-like distribution of matter arising in ΛCDM cosmological models. Only a small number of minor discrepancies remain with respect to the observational data. Saturated Lyα absorption lines with column densities N HI > 10 14.5 cm −2 at 2 < z < 2.5 are underpredicted in the models. An uncertain correction for continuum placement bias is required to match the distribution and power spectrum of the transmitted flux, particularly at z > 4. Finally, the temperature of intergalactic gas in the simulations may be slightly too low at z = 2.7 and a flatter temperature-density relation is required at z = 2.4, consistent with the expected effects of non-equilibrium ionisation during He II reionisation.

Large dynamic range numerical simulations of atomic cooling driven collapse of gas in pregalactic dark matter haloes with T vir ∼ 10 000 K show that the gas loses 90 per cent and more of its angular momentum before rotational support sets in. In a fraction of these haloes where the metallicity is low and ultraviolet (UV) radiation suppresses H 2 cooling, conditions are thus very favourable for the rapid build-up of massive black holes. Depending on the progression of metal enrichment, the continued suppression of H 2 cooling by external and internal UV radiation and the ability to trap the entropy produced by the release of gravitational energy, the gas at the centre of the halo is expected to form a supermassive star, a stellar-mass black hole accreting at super-Eddington accretion rates or a compact star-cluster undergoing collisional run-away of massive stars at its centre. In all three cases, a massive black hole of initially modest mass finds itself at the centre of a rapid inflow of gas with inflow rates of 1 M yr −1 . The massive black hole will thus grow quickly to a mass of 10 5 -10 6 M until further inflow is halted either by consumption of gas by star formation or by the increasing energy and momentum feedback from the growing massive black hole. Conditions for the formation of massive seed black holes in this way are most favourable in haloes with T vir ∼ 15 000 K and V vir ∼ 20 km s −1 with less massive haloes not allowing collapse of gas by atomic cooling and more massive haloes being more prone to fragmentation. This should imprint a characteristic mass on the mass spectrum of an early population of massive black hole seeds in pre-galactic haloes which will later grow into the observed population of supermassive black holes in galactic bulges.