FULLY COUPLED SIMULATION OF COSMIC REIONIZATION. I. NUMERICAL METHODS AND TESTS

The Astrophysical Journal Supplement Series, 216:16 (24pp), 2015 January C 2015. doi:10.1088/0067-0049/216/1/16 The American Astronomical Society. All rights reserved. FULLY COUPLED SIMULATION OF COSMIC REIONIZATION. I. NUMERICAL METHODS AND TESTS Michael L. Norman 1,2 , Daniel R. Reynolds 3 , Geoffrey C. So 1 , Robert P. Harkness 2,4 , and John H. Wise 5 CASS, University of California, San Diego, 9500 Gilman Drive La Jolla, CA 92093-0424, USA SDSC, University of California, San Diego, 9500 Gilman Drive La Jolla, CA 92093-0505, USA 3 Southern Methodist University, 6425 Boaz Lane, Dallas, TX 75205, USA 4 NICS, Oak Ridge National Laboratory, 1 Bethel Valley Road, Oak Ridge, TN 37831, USA 5 Center for Relativistic Astrophysics, Georgia Institute of Technology, 837 State Street, Atlanta, GA 30332, USA Received 2013 June 8; accepted 2014 November 23; published 2015 January 9 ABSTRACT We describe an extension of the Enzo code to enable fully coupled radiation hydrodynamical simulation of inhomogeneous reionization in large ∼(100 Mpc) 3 cosmological volumes with thousands to millions of point sources. We solve all dynamical, radiative transfer, thermal, and ionization processes self-consistently on the same mesh, as opposed to a postprocessing approach which coarse-grains the radiative transfer. We do, however, employ a simple subgrid model for star formation which we calibrate to observations. The numerical method presented is a modification of an earlier method presented in Reynolds et al. differing principally in the operator splitting algorithm we use to advance the system of equations. Radiation transport is done in the gray flux-limited diffusion (FLD) approximation, which is solved by implicit time integration split off from the gas energy and ionization equations, which are solved separately. This results in a faster and more robust scheme for cosmological applications compared to the earlier method. The FLD equation is solved using the hypre optimally scalable geometric multigrid solver from LLNL. By treating the ionizing radiation as a grid field as opposed to rays, our method is scalable with respect to the number of ionizing sources, limited only by the parallel scaling properties of the radiation solver. We test the speed and accuracy of our approach on a number of standard verification and validation tests. We show by direct comparison with Enzo’s adaptive ray tracing method Moray that the well-known inability of FLD to cast a shadow behind opaque clouds has a minor effect on the evolution of ionized volume and mass fractions in a reionization simulation validation test. We illustrate an application of our method to the problem of inhomogeneous reionization in a 80 Mpc comoving box resolved with 3200 3 Eulerian grid cells and dark matter particles. Key words: cosmology: theory – methods: numerical – radiative transfer redshift interval 6 z 10 using the Hubble Space Telescope support the galaxy reionizer hypothesis, with the caveat that the faint end of the luminosity function which contributes substantially to the number of ionizing photons has not yet been measured (Robertson et al. 2010; Bouwens et al. 2012). Given the paucity of observational information about the pro- cess of cosmic reionization, researchers have resorted to theory and numerical simulation to fill in the blanks. As reviewed by Trac & Gnedin (2011), progress in this area has been dramatic, driven by a synergistic interplay between semi-analytic ap- proaches and numerical simulations. The combination of these two approaches have converged on a qualitative picture of how H reionization proceeds assuming the primary ionizing sources are young, star-forming galaxies. The physics of the reioniza- tion process is determined by the physics of the sources and sinks of ionizing radiation in an expanding universe. Adopting the ΛCDM model of structure formation, galaxies form hierar- chically through the merger of dark matter halos. The structure and evolution of the dark matter density field is now well un- derstood through ultra-high-resolution numerical N-body sim- ulations (Springel et al. 2005; Klypin et al. 2011) and through analytic models based on these simulations (Cooray & Sheth 2002). By making certain assumptions about how ionizing light traces mass and the dynamics of H ii regions, a basic picture of the reionization process has emerged (Furlanetto et al. 2004, 2006; Iliev et al. 2006b; Zahn et al. 2007) that is confirmed by detailed numerical simulations; e.g., Zahn et al. (2011). The basic picture is that galaxies form in the peaks of the dark matter density field and drive expanding H ii regions into their surroundings by virtue of the UV radiation emitted from young, massive stars. These H ii regions are initially isolated, 1. INTRODUCTION The epoch of reionization (EoR) is a current frontier of cos- mological research both observationally and theoretically. Ob- servations constrain the transition from a largely neutral inter- galactic medium (IGM) of primordial gas to a largely ionized one (singly ionized H and He) to the redshift interval z ∼ 11–6, which is a span of roughly 500 Myr. The completion of H reion- ization by z ≈ 6 is firmly established through quasar absorption line studies to luminous, high-redshift quasars which exhibit Lyα Gunn–Peterson absorption troughs (Fan et al. 2006). The precise onset of H reionization (presumably tied to the formation of the first luminous ionizing sources) is presently unknown ob- servationally, however, cosmic microwave background (CMB) measurements of the Thomson optical depth to the surface of last scattering by the Wilkinson Microwave Anisotropy Probe (WMAP) and Planck satellites indicates that the IGM was sub- stantially ionized by z ∼ 10 (Spergel et al. 2003; Komatsu et al. 2009; Jarosik et al. 2011; Planck Collaboration et al. 2014). Since the optical depth measurement is redshift-integrated and averaged over the sky, the CMB observations provide no infor- mation about how reionization proceeded or the nature of the radiation sources that caused it. It is generally believed that reionization begins with the formation of Population III stars at z ∼ 20–30 (Abel et al. 2002; Yoshida et al. 2003; Bromm & Larson 2004; Sokasian et al. 2004), but that soon the ionizing photon budget becomes dominated by young, star forming galaxies (see e.g., Wise et al. 2012; Xu et al. 2013), and to a lesser extent by the first quasars (Madau et al. 1999; Bolton & Haehnelt 2007; Haardt & Madau 2012; Becker & Bolton 2013). Observations of galaxies in the