THESAN-HR: Galaxies in the Epoch of Reionization in warm dark matter, fuzzy dark matter and interacting dark matter

Using high-resolution cosmological radiation-hydrodynamic (RHD) simulations (THESAN-HR), we explore the impact of alternative dark matter (altDM) models on galaxies during the Epoch of Reionization. The simulations adopt the IllustrisTNG galaxy formation model. We focus on altDM models that exhibit small-scale suppression of the matter power spectrum, namely warm dark matter (WDM), fuzzy dark matter (FDM), and interacting dark matter (IDM) with strong dark acoustic oscillations (sDAO). In altDM scenarios, both the halo mass functions and the UV luminosity functions at $z\gtrsim 6$ are suppressed at the low-mass/faint end, leading to delayed global star formation and reionization histories. However, strong non-linear effects enable altDM models to"catch up"with cold dark matter (CDM) in terms of star formation and reionization. The specific star formation rates are enhanced in halos below the half-power mass in altDM models. This enhancement coincides with increased gas abundance, reduced gas depletion times, more compact galaxy sizes, and steeper metallicity gradients at the outskirts of the galaxies. These changes in galaxy properties can help disentangle altDM signatures from a range of astrophysical uncertainties. Meanwhile, it is the first time that altDM models have been studied in RHD simulations of galaxy formation. We uncover significant systematic uncertainties in reionization assumptions on the faint-end luminosity function. This underscores the necessity of accurately modeling the small-scale morphology of reionization in making predictions for the low-mass galaxy population. Upcoming James Webb Space Telescope (JWST) imaging surveys of deep, lensed fields hold potential for uncovering the faint, low-mass galaxy population, which could provide constraints on altDM models.


INTRODUCTION
The nature of dark matter (DM) is one of the most outstanding questions in modern cosmology (see Bertone & Hooper 2018, for a review).In the standard cosmological model, DM is assumed to be cold after its decoupling from the primordial plasma and is effectively collisionless (known as CDM).This model is successful at describing cosmic structure formation on large scales ( 1 Mpc; e.g.Blumenthal et al. 1984;Davis et al. 1985) and acts as the foundation for the theory of galaxy formation (e.g.White & Rees 1978;Dekel & Silk 1986;Kauffmann et al. 1993;Springel et al. 2005; Behroozi ★ E-mail: xshen@caltech.eduet al. 2013;Vogelsberger et al. 2014a;Schaye et al. 2015).Many beyond-standard-model theories do predict DM candidates with such properties (see Bertone et al. 2005 for a review), for example, the weakly-interacting massive particles (WIMPs) motivated by the hierarchy problem as well as the light axion(-like) particles (ALPs) motivated by the strong CP problem and string theories.The success of CDM on astrophysical scales and the natural production mechanisms of these candidates have motivated decades of experimental searches (e.g.Ahmed et al. 2009;CDMS II Collaboration et al. 2010;Aprile et al. 2012;Akerib et al. 2014Akerib et al. , 2017;;Aprile et al. 2018).However, as the parameter space and model configuration of these CDM candidates become increasingly constrained given the null results in direct detection (e.g.Aprile et al. 2018;Bertone & Tait 2018), there are strong motivations to explore DM candidates beyond the CDM paradigm.Meanwhile, on small scales (Local dwarf galaxies), several discrepancies between CDM predictions and astrophysical observations have come to light (see Bullock & Boylan-Kolchin 2017, for a review), including but not limited to the missing satellites problem (e.g.Klypin et al. 1999;Moore et al. 1999), the too-big-to-fail problem (TBTF; Boylan-Kolchin et al. 2011;Boylan-Kolchin et al. 2012;Tollerud et al. 2014), the core-cusp problem (e.g.Flores & Primack 1994;Moore 1994;de Blok et al. 2001;Kuzio de Naray et al. 2006;Oh et al. 2015), and recently the diversity problem (e.g.Oman et al. 2015;Kaplinghat et al. 2019;Zavala et al. 2019).DM models alternative to the classical CDM has been proposed as the solution to these anomalies in the Local Universe.However, despite the application to solve the local puzzles, a category of alternative DM (altDM 1 hereafter) models can modify the initial linear matter power spectrum at small, poorly-constrained scales ( 10 ℎ/Mpc) which determines the initial conditions of non-linear structure formation.Such modifications can affect the assembly of haloes/galaxies in the early Universe as well as cosmic reionization.In the following, we will introduce three DM models that fall into this category.
• Warm dark matter (WDM) is a class of DM candidates with large free-streaming velocities.The free-streaming of WDM can suppress the density fluctuations on Mpc scales and below.WDM can be produced through the freeze-out mechanism after reaching thermal equilibrium with the initial plasma (e.g.gravitino; Steffen 2006).In addition to the thermal relic WDM, the sterile neutrino is another well-motivated candidate, which can be non-thermally produced through scattering processes due to their mixing with active neutrinos with the Dodelson-Widrow mechanism (Dodelson & Widrow 1994), through resonant production (Shi & Fuller 1999), or through coupling with other fields (Kusenko 2006;Shaposhnikov & Tkachev 2006).WDM has been empirically motivated as a solution to some of the above-mentioned small-scale astrophysical challenges (e.g.Bode et al. 2001;Lovell et al. 2012Lovell et al. , 2014;;Polisensky & Ricotti 2014;Lovell et al. 2017).The latest constraints on WDM (reported at the 2 level) are e.g.(equivalent) thermal relic mass  WDM > 5.3 keV (3.5 keV, if allowing sudden temperature changes) from Lyman- forest observations (Iršič et al. 2017a), > 5.2 keV from substructures of strong gravitational lenses (Gilman et al. 2020), > 6.5 keV from the ultra-faint dwarf (UFD) satellites of the Milky Way (Nadler et al. 2021) ( 2 keV in Newton et al. 2021).In this work, we choose to study the benchmark model with  WDM = 3 keV.It lies at the edge of the existing constraints bearing in mind the inevitable modeldependence and systematical uncertainties of all reported observational constraints.This value gives prominent signals at the halo mass scale that can be well-resolved by our simulations (see Section 2 for more discussions) and is of interest in explaining the 3.55 keV X-ray emission from the Galatic center as the decay product of sterile neutrinos (e.g.Bulbul et al. 2014;Boyarsky et al. 2014;Schneider 2016;Adhikari et al. 2017).
• Dark Acoustic Oscillations (DAOs) can appear for a generic class of interacting dark matter (IDM) models 2 .The acoustic oscillations can be caused by interactions between DM and Standard Model particles in the early Universe (e.g.Boehm et al. 2002;Sigurdson et al. 2004;Boehm & Schaeffer 2005;Boehm et al. 2014;Schewtschenko et al. 2016;Boddy et al. 2018;Gluscevic & Boddy 2018) and has been considered for thermal WIMP particles (e.g.Loeb & Zaldarriaga   1 Not to be confused with physics models alternative to DM.
2 Not all IDM models feature DAOs.We are studying an empirical benchmark model with strong DAOs in this work.
2005 ;Bertschinger 2006;Bringmann 2009).It was then also generalized to hidden-sector DM models with e.g.DM-dark radiation interactions (e.g.van den Aarssen et al. 2012;Buckley et al. 2014;Cyr-Racine et al. 2014), which has been included in the Effective Theory of Structure Formation (ETHOS) framework (Cyr-Racine et al. 2016;Vogelsberger et al. 2016).DAOs take place on small scales and are accompanied by collisional damping, which can suppress the abundance of small-scale structures as in WDM but with a different physical mechanism.As a continuation of the ETHOS work, Bohr et al. (2020) proposed effective parameters of DAO models that fully characterize structure formation in the non-linear regime at high redshifts ( ≥ 5) in a simple and clear way.We pick the strong DAO (sDAO) model studied in Bose et al. (2019) and Bohr et al. (2021) with the relative strength of the first DAO peak set to unity.This represents a model with strong small-scale residual fluctuations in addition to the major suppression feature at a slightly larger scale.
Structure formation at high redshift serves as a promising complementary channel to constrain the particle nature of DM.The population of faint galaxies during the Epoch of Reionization (EoR) will be sensitive to the altDM physics that suppresses the small-scale power spectrum.The population of high-redshift galaxies has been studied in the context of WDM (e.g.Schultz et al. 2014;Bose et al. 2016b;Dayal et al. 2017;Lopez-Honorez et al. 2017;Lovell et al. 2018;Menci et al. 2018;Bozek et al. 2019), FDM (e.g.Bozek et al. 2015;Menci et al. 2017;Ni et al. 2019;Mocz et al. 2020) and DAO models (Lovell et al. 2018;Bohr et al. 2021;Kurmus et al. 2022).The methods adopted include DM-only simulations with abundance matching (e.g.Schultz et al. 2014;Corasaniti et al. 2017), semi-analytical models of galaxy formation (e.g.Bose et al. 2016b;Dayal et al. 2017;Khimey et al. 2021) as well as cosmological hydrodynamic simulations (e.g.Lovell et al. 2018;Ni et al. 2019).
The suppression of the rest-frame UV luminosity function was identified at faint magnitudes  UV −16 in these models.The UV luminosity functions probed by deep Hubble Space Telescope (HST) programs (e.g.Livermore et al. 2017;Atek et al. 2018;Ishigaki et al. 2018;Bouwens et al. 2022b) have been used to constrain altDM models that suppress small-scale structure formation.In the near future, deep imaging surveys of strongly lensed, inherently faint galaxies concluded with the James Webb Space Telescope (JWST) will push the detection limit by at least three magnitudes and extend the redshift frontier of galaxy searches with its infrared frequency coverage, providing compelling tests on the nature of DM.For example, as part of the JWST Early Release Observations (ERO; Pontoppidan et al. 2022), the galaxy cluster SMACS0723 was imaged using NIRCam, yielding the discovery of galaxies out to  ∼ 16 (e.g.Atek et al. 2023).Another Early Release Science (ERS) program (GLASS-JWST; Treu et al. 2022) pointed to one foreground Hubble Frontier Field (HFF) cluster, Abell 2744.Forthcoming observational programs will extend the measurements of faint-end luminosity functions in these lensed fields.Besides galaxy abundances, many theoretical studies (e.g.Bose et al. 2016bBose et al. , 2017;;Corasaniti et al. 2017;Dayal et al. 2017;Lovell et al. 2018Lovell et al. , 2019;;Ni et al. 2019) found higher occupation fraction of luminous galaxies, higher UV luminosities at a given halo mass in altDM models than in the CDM case.Despite the initial delay in the assembly of DM haloes, the stellar content of galaxies undergoes faster build-up in altDM models, which closes the gap in the global star formation and reionization history between models at lower redshift.The strong, rapid starbursts in low-mass galaxies in these models could be another interesting feature for observations of high-redshift galaxies as well as the Local Group UFDs (e.g.Lovell et al. (2019) and Bozek et al. (2019) found that the simulated dwarfs have delayed but more rapid and diverse star formation histories in WDM).
In addition to the galaxy properties, the neutral gas distribution in the intergalactic medium (IGM) will be sensitive to these altDM physics as well.The mass distribution of baryonic gas is a tracer of the DM density field.The phase of the gas has additional dependencies on the intergalactic ionizing radiation from star formation in galaxies.The neutral phase of the gas can be mapped by the Lyman- forest observations at  5 − 6 (e.g.Cen et al. 1994;Hernquist et al. 1996;Seljak et al. 2005;Viel et al. 2005Viel et al. , 2013;;Delubac et al. 2015;Yang et al. 2020;Bosman et al. 2022) and the intensity mapping of the 21 cm spin-flip transition of the hydrogen atom (Furlanetto et al. 2006;Mellema et al. 2006;Parsons et al. 2010;Pritchard & Loeb 2012).The global signal of the ionized phase can be constrained by the optical depth of the Cosmic Microwave Background (CMB).The impact of altDM models on the IGM properties and various detection signals has been studied in e.g.Viel et al. (2013); Das et al. (2018); Lovell et al. (2018); Bose et al. (2019); Muñoz et al. (2021Muñoz et al. ( , 2022)).In the near future, a more complete picture of the neutral gas distribution in the Universe will be provided by current and upcoming radio instruments, including the Low-Frequency Array (LOFAR; van Haarlem et al. 2013), Square Kilometer Array (SKA; Dewdney et al. 2009), and Hydrogen Epoch of Reionization Array (HERA;DeBoer et al. 2017).These instruments will constrain the timing and morphology of reionization, the properties of the first galaxies, the evolution of large-scale structure, and the early sources of heating throughout the EoR and cosmic dawn.For example, as highlighted in Muñoz et al. (2021), the 21 cm observation can constrain the smallscale fluctuations at  300 ℎ/ Mpc, which has never been reached by any other constraints.
All the EoR constraints involve non-negligible astrophysical uncertainties.This includes the escape fraction of ionizing photons, the star formation efficiency regulated by feedback from star formation as well as the physics that sets the atomic cooling limit for galaxy formation.For example, radiative and supernovae feedback from both Population III (PopIII) and metal-enriched PopII stars can regulate or even prevent star formation in low-mass galaxies and flatten the faint-end UV luminosity function (e.g.Jaacks et al. 2013;Wise et al. 2014;O'Shea et al. 2015;Xu et al. 2016;Ocvirk et al. 2016;Dayal & Ferrara 2018).This can contaminate the signal of altDM physics.The delay in reionization in altDM models can also be compensated by a steeper redshift evolution of the ionizing photon escape fraction (e.g.Dayal et al. 2017).Therefore, it is essential to have theoretical forecasts based on robust galaxy formation and reionization models that have been tested against a wide range of observational constraints.The project (Kannan et al. 2022b;Garaldi et al. 2022;Smith et al. 2022) is a suite of radiation-magneto-hydrodynamic simulations featuring an efficient radiation hydrodynamics solver (Kannan et al. 2019) that precisely captures the interaction between ionizing photons and gas.It is coupled to the well-tested IllustrisTNG galaxy formation model (Pillepich et al. 2018a;Nelson et al. 2019).simulations have been used to make a wide range of predictions for the Universe in the EoR (e.g.Smith et al. 2022;Garaldi et al. 2022;Kannan et al. 2022b,c,a;Borrow et al. 2022;Xu et al. 2022;Yeh et al. 2023).In this work, we use the small-volume, high-resolution variant of the simulations (namely the -project introduced in Borrow et al. 2022) to study galaxy formation in the EoR in altDM models with self-consistent modeling of the reionization process at small scales.
This paper is organized as follows.In Section 2, we introduce the simulations and the altDM models considered.A visual overview of the DM density field and IGM properties is given in Section 3. In Sections 4 and 5, we present predictions of galaxy abundance, including the halo/stellar mass function and rest-frame UV luminosity functions.In Section 6, we take a deeper dive into galaxy properties in different DM models, discussing scaling relations and spatially resolved properties of galaxies.In Section 7, we discuss potential ways to disentangle altDM signatures with various uncertainties in astrophysical processes.The conclusions of the paper are given in Section 8. Throughout this paper, we employ the Planck Collaboration et al. ( 2016) cosmology with ℎ = 0.6774, Ω 0 = 0.3089, and Ω  = 0.0486,  8 = 0.8159,  s = 0.9667.

SIMULATIONS
The project (Kannan et al. 2022b;Garaldi et al. 2022;Smith et al. 2022) is a suite of radiation-magneto-hydrodynamic simulations utilizing the moving-mesh hydrodynamics code AREPO (Springel 2010;Weinberger et al. 2020).Gravity is solved using the hybrid Tree-PM method (Barnes & Hut 1986).The hydrodynamics is solved using the quasi-Lagrangian Godunov method (Godunov & Bohachevsky 1959) on an unstructured Voronoi mesh grid (see Vogelsberger et al. 2020a, for a review).For self-consistent treatment of ionizing radiation, the project employs the radiative transfer (RT) extension A - (Kannan et al. 2019), which solves the first two moments of the RT equation assuming the M1 closure relation (Levermore 1984).The simulation includes the sourcing (from stars and AGNs) and propagation of ionizing photons (in three energy bins relevant for hydrogen and helium photoionization between energy intervals of [13.6, 24.6, 54.4,∞) eV) as well as a non-equilibrium thermochemistry solver to model the coupling of radiation fields to gas.The luminosity and spectral energy density of stars in as a complex function of age and metallicity are calculated using the Binary Population and Spectral Synthesis models (BPASS v2.2.1; Eldridge et al. 2017).The sub-grid escape fraction of stars was set to be  esc = 0.37 in the simulations to match the global reionization history of the Universe in CDM.For details of the simulation methods, we refer to Kannan et al. (2019Kannan et al. ( , 2022b)).
In terms of the galaxy formation model, the simulations employ the IllustrisTNG model (Pillepich et al. 2018a;Nelson et al. 2019), which is an update of the Illustris model (Vogelsberger et al. 2014a,b).The simulations include (1) density-, temperature-, metallicity-and redshift-dependent cooling of metal-enriched gas, (2) a two-phase, --is the run with identical physics and numeric resolution as -but doubling the box side-length.
(2)  box : Side-length of the periodic simulation box.The unit is comoving Mpc (cMpc).(3)  part : Number of particles (cells) in the simulation.In the initial conditions, there are an equal number of DM particles and gas cells.(4)  dm : Mass of DM particles, which is conserved over time.(5)  baryon : Mass of gas cells in the initial conditions as a reference for the baryonic mass resolution.The gas cells are (de-)refined so that the gas mass in each cell is within a factor of two of this target gas mass.Stellar particles stochastically generated out of gas cells can take arbitrary initial masses and are subject to mass loss via stellar evolution (Vogelsberger et al. 2013).( 6)  dm, * : The gravitational softening length for the DM and stellar particles.( 7) ℎ b : The minimum gravitational softening length of gas cells, which are adaptively softened.( 8) Reionization model used.The fiducial model employs on-the-fly radiative transfer of ionizing photons while the uniform UVB model assumes a spatially uniform, time-varying background radiation of ionizing photons.(9) DM model employed.(10) The characteristic wavenumber of the small-scale suppression on the linear matter power spectrum, defined as effective equation of state model for the interstellar medium (ISM) at the sub-resolution level (Springel & Hernquist 2003), (3) star formation in dense gas following the empirically defined Kennicutt-Schmidt relation, (4) thermal and mechanical feedback from supernovae and stellar winds, (5) metal enrichment from stellar evolution and supernovae, and (6) SMBH formation, growth, and feedback (in two regimes described in Weinberger et al. 2017, although SMBH physics have little impact on the low-mass galaxies studied in -).The model has been extensively tested in large-scale simulations and is able to produce realistic galaxies that match a wide range of observations (e.g.Nelson et al. 2018;Springel et al. 2018;Genel et al. 2018;Pillepich et al. 2018b;Naiman et al. 2018;Vogelsberger et al. 2018;Marinacci et al. 2018;Vogelsberger et al. 2020b;Shen et al. 2020).Borrow et al. (2022) introduced a subset of high-resolution smallvolume simulations ( -) using the same numeric setup and physics inputs as the main suite, aiming to explore the formation and evolution of low-mass galaxies in the early Universe.The mass resolution of DM (gas) is set to 4.82×10 5 M (9.04×10 4 M ), which is about 6.4 times better than the flagship -1 simulation and allows atomic cooling haloes ( halo 10 7.5−8 M ; Wise et al. 2014) to be properly resolved ( 100 DM particles).The gravitational softening lengths of DM and stars are set to 0.85 kpc.The gas cells are adaptively softened according to the cell radius with a minimum value set to 0.85 kpc.To explore the impact of reionization on lowmass galaxies, we have included simulations with physics replaced by a uniform UV Background (uniform UVB) model.In this model, we assume reionization is nearly instantaneous at redshift  i = 10 (following the original IllustrisTNG model 3 ) with the strength of the radiation field given by Faucher-Giguère et al. (2009). 3We note that the original Illustris and IllustrisTNG simulations adopted a fixed reionization redshift of  i = 6.However, similar models like EAGLE initiate their UVB at a higher redshift of  i = 11.5 (Schaye et al. 2015).Here we choose  i = 10 since our fully coupled RHD simulations end at  = 5 to ensure efficient use of computational resources.
Simulations with this UVB are performed in exactly the same way as the models.The ionizing radiation flux from the UVB is then passed directly to the thermochemistry module.In most cases, it fully ionizes the gas, though gas can be self-shielded from such radiation if it reaches a high enough density.We employ the model from Rahmati & Schaye (2014) to determine if such gas is self-shielded against the external UVB.
The fiducial -runs can be viewed as putting the small volume in a void of galaxies since no external ionizing photon background is included.On the contrary, the uniform UVB model is equivalent to putting the small volume in a bath of ionizing photons contributed by luminous galaxies nearby.So we effectively bracket the uncertainties from reionization history/morphology in this experiment.This fact is also illustrated in the top panel of Figure 4.The two types of runs have systematically faster/lower reionization history compared to the main suite which has been calibrated to reionization on large scales.A subtle point is that the instantaneous switch-on of the uniform UVB can generate an artificial bump in the star formation history (SFH) of galaxies at  10.This population of old stars should not exist in the physical case when the volume is affected by ionizing feedback from massive galaxies.

Initial conditions and DM models
The initial conditions are generated with the G -4 code (Springel et al. 2021) using the second-order Lagrangian method at the initial redshift  = 49.For the CDM model, we adopt the linear matter power spectrum generated using C (Lewis et al. 2000;Howlett et al. 2012), which is identical to the ones used in the fiducial simulations (Kannan et al. 2022b;Garaldi et al. 2022;Smith et al. 2022) and the IllustrisTNG simulation (Nelson et al. 2019;Pillepich et al. 2018a).Notably, in addition to the vanilla CDM model, the suite includes simulations in WDM, FDM, and sDAO models.These models introduce cutoffs to the linear matter power spectrum at small scales and could affect early structure formation in non-trivial ways.The difference between altDM models and Linear matter power spectrum for the CDM (blue curve), WDM ( WDM = 3 keV; red curve), sDAO (green curve), and FDM ( a = 2 × 10 −21 eV; purple curve) models.Δ 2 () ≡  3  ()/2  2 is the dimensionless power spectrum.For the CDM and FDM models, the linear matter power spectrums are generated using C and A C .For the WDM and sDAO models, we employ analytical transfer functions (see Section 2.1 for details) on top of the CDM power spectrum.The altDM models examined in this study all have small-scale damping of the power spectrum.They differ in the location and shape of the damping wing as well as the oscillatory patterns.
CDM is completely characterized by the transfer function, defined as4

𝑇 (𝑘) ≡
altDM ()  CDM () For the WDM model, the free-streaming of DM particles damps the fluctuations of the density field at small scales.The transfer function for thermal relic WDM can be described as where  = 2,  = −5/ and  takes the value 1.12 (e.g.Viel et al. 2005). controls the characteristic scale of the damping5 .The halfpower wavenumber (where the linear WDM suppression reaches 1/2 in terms of matter power w.r.t. the ΛCDM case) can be approximated as (Viel et al. 2005(Viel et al. , 2013) where  WDM is the mass of the thermal relic WDM particle and Ω DM ≡ Ω 0 − Ω b is the cosmological abundance of DM.We choose  WDM = 3 keV, which is of interest in explaining the 3.55 keV X-ray line from Galactic center (e.g.Bulbul et al. 2014;Boyarsky et al. 2014;Adhikari et al. 2017) and lies at the edge of existing observational constraints (e.g.Iršič et al. 2017a;Gilman et al. 2020;Nadler et al. 2021).The half-power wavenumber defined here is different from the half-mode wavenumber quoted in some literature.The half-mode wavenumber is defined as where the transfer function becomes 1/2 w.r.t. the ΛCDM case and is larger than the half-power wavenumber.
For the IDM model featuring DAOs, we choose the model in the ETHOS framework (Cyr-Racine et al. 2014) studied in Bohr et al. (2021).The small-scale damping of the linear matter power spectrum is caused by interactions between DM and relativistic particles in the early Universe and is accompanied by DAOs.We use the parametrization outlined in Equation (3) of Bohr et al. (2020), which is an extension to Equation (2).The transfer function in the DAO model is mainly controlled by two parameters,  peak (the wavenumber of the first acoustic peak) and ℎ peak (the amplitude of the first acoustic peak).We take the parameter choice,  peak = 40 ℎ Mpc −1 , ℎ peak = 1 in Bohr et al. (2020).In this case, the half-power wavenumber of the first damping wing in the power spectrum is For the FDM model, the quantum mechanical nature of ultralight ALPs moving with the Hubble flow will prevent clustering below the de Broglie wavelength and result in the damping feature at small scales.The transfer function can be approximated as (Hu et al. 2000) where  a is the ALP mass,  eq J 9 ( a /10 −22 eV) 1/2 Mpc −1 is the effective "Jeans" length at matter-radiation equality.The half-power wavenumber is (Hu et al. 2000) To be more accurate, we compute the linear matter power spectrum of FDM with the Boltzmann code A (Hlozek et al. 2015) (developed based on C ; Lewis et al. 2000;Howlett et al. 2012).We choose the ALP particle mass  a = 2 × 10 −21 eV, which yields  FDM 1/2 25 ℎ Mpc −1 similar to the WDM model we study.The isocurvature fluctuations of FDM are ignored, which usually appear at sub-parsec scales (e.g.Hogan & Rees 1988;Kolb & Tkachev 1993;Zurek et al. 2007).Following previous work on high-redshift structure formation in FDM (e.g.Schive et al. 2016;Armengaud et al. 2017), we do not include the quantum pressure in the dynamic of our simulations.This should not alter our predictions for galaxies significantly as found in e.g.Zhang et al. (2018); Nori et al. (2019); Mocz et al. (2020).
In Figure 1, we show the linear matter power spectrum of all the models considered in this paper.The details and numerical parameters of the simulations are summarized in Table 1.In the following parts of the paper, we will discuss the results of these simulations.

DM DENSITY FIELD AND IGM PROPERTIES
The suppression of the linear matter spectrum is directly reflected in the DM density field.In Figure 2, we show images of the DM surface density distributions at  = 6 in different models listed in Table 1.The density map is projected through a slice of thickness ∼ 3 cMpc (roughly half of the box size).Due to the suppression of the power spectrum, the altDM models give smoother DM distributions on small scales and a lower number of filaments and low-mass (sub)haloes.The smoothing of the density field is most noticeable for The suppression of small-scale DM structures will delay the formation of low-mass galaxies, which are the dominant sources of ionizing photons in the early phase of reionization (e.g.Jaacks et al. 2012;Madau & Dickinson 2014;Robertson et al. 2015;Rosdahl et al. 2022;Yeh et al. 2023).Therefore, the morphology of reionization and the thermal properties of the IGM can be affected.In Figure 3  The CDM model has more small clumps of neutral gas retained in subhaloes even in the highly-ionized regions.These low-mass subhaloes and associated cold gas clumps are missing in altDM models.At large scales, WDM and sDAO runs have a higher neutral gas abundance on average and smaller ionizing bubbles at  = 6.The global reionization processes in these two models have been delayed systematically due to the deficiency of ionizing photons contributed by low-mass galaxies.Corresponding features can be found in the temperature distribution as well.The WDM and sDAO runs show larger voids of low-temperature, neutral gas while a lower number of small hot gas bubbles.In the FDM run, although the statistics of high-density gas clumps, filaments, and hot gas bubbles are similar, we find rather different phenomena on large scale.The global ionization fraction at  = 6 has overtaken the CDM value due to the late-time starburst (shown in Figure 4 as well).These strong non-linear effects create larger ionizing bubbles and lower neutral gas column densities in FDM at  6.
(often in the form of dense self-shielded clumps) represents baryons falling into gravitationally bound DM structures and thus traces DM overdensities well.Some less dense neutral gas is retained in filaments around DM haloes, which will be fed into DM haloes fueling subsequent star formation.These two types of neutral gas structures are more abundant with increasing clustering power of DM at small scales and are suppressed in altDM models.On the other hand, the low-density neutral gas filling the IGM is more sensitive to the phase of global reionization that ties to the global SFH.In the WDM and sDAO models, the star formation and reionization history are sys-tematically delayed compared to the CDM case (as will be shown in Figure 4).Therefore, the overall opacity of neutral hydrogen on large scales is enhanced, contrary to the trend at high densities and small scales.In the FDM model, it is the same situation at  7.But between  = 6 and  = 7, the neutral gas preserved by the delayed early-phase reionization leads to enhanced star formation in low-mass galaxies.It elevates the total star formation in the simulation volume (despite lowered number counts of faint galaxies) and accelerates the late-time reionization to even overtake the CDM counterpart.This will be shown in Figure 4 and discussed below.The decreased abundance of dense neutral clumps and filaments at small scales allows the hot ionizing bubbles to penetrate through the IGM in a more violent fashion.This is why the IGM at  = 6 in the FDM model is less opaque than the CDM case.Overall, the IGM properties are more sensitive to the location of the first damping in the matter power spectrum than the small-scale acoustic oscillations.The sDAO model with strong small-scale oscillations (ℎ peak = 1) and steeply rising halo mass function at the low-mass end (compared to WDM and FDM) still has the most delayed star formation and reionization history out of the three altDM models.
The gas temperature map is also shown in Figure 3. Bubbles of gas are heated by stellar feedback and are spatially correlated with DM haloes that are massive enough to host galaxy formation.In the WDM and sDAO runs, we find fewer small-scale clusters of these hot ionized bubbles and larger voids of cold neutral gas on large scales compared to CDM.In the FDM model, due to the strong late-time star formation we mentioned above, the volumefilling fraction of hot ionized bubbles is larger than in the CDM case.The radiative heating of the IGM is dominated by large bubbles from massive galaxies.In general, the temperature distribution of IGM gas is consistent with the trends observed in the neutral gas distribution.The IGM properties in the uniform UVB model are shown in the Appendix A. In short, averaging the ionizing photon background of luminous sources over the entire simulation volume ionizes all of the gas.This produces radically different morphologies of reionization and thermal properties of the IGM.The entire volume is heated to above ∼ 10 4 K.
In the top panel of Figure 4, we show the volume-weighted neutral fraction of hydrogen in simulations.Regardless of the DM model, runs with the uniform UVB have the entire simulation volume ionized right after the activation of the UVB.The -1 result is shown for reference with the dotted line.It emphasizes that the small-volume -runs are tracing a biased volume of the Universe.This bias is mainly caused by the absence of massive, bright sources of ionizing photons that would be present in a real environment.Therefore, we will be cautious about interpreting our results quantitatively as what would be found in real observations but focus on the net differences driven by various of physics effects.
All altDM models cause a delay in cosmic reionization during the early phase when low-mass galaxies primarily contribute to ionizing photons.The extent of the delay in reionization is positively correlated with the damping length/mass scale of the power spectrum.At late stages of reionization ( 7), the FDM model exhibits a period of enhanced star formation and accelerated reionization which enables it to surpass the CDM at  6.However, it is important to note that the phenomena observed in small-volume simulations may not necessarily be representative of the global reionization signal.For example, at  ∼ 6, the SFRD and the ionizing photon budget on ∼ 100 Mpc scales are predominantly influenced by massive galaxies ( * 10 9 M , e.g.Kannan et al. 2022b;Yeh et al. 2023), which are not sampled in -.Therefore, the signal in the global phase of reionization will be much weaker.The solid lines show the results from the simulations with RT, while the dashed lines show the results in the uniform UVB model.In altDM models, reionization histories display systematic delays that correlate with the damping scale of the power spectrum.Interestingly, in the FDM model with the sharpest power spectrum damping, late-time starburst in low-mass galaxies can lead to a global ionized fraction larger than the CDM counterpart at  6. Middle: Evolution of the star formation rate density (SFRD) in the simulation volume.Star formation in the altDM models is suppressed at high redshift due to delayed structure formation at small scales until the late-time starburst closes the gap.In contrast, the uniform UVB model exhibits a bump in SFRD before the activation of the uniform UVB.After the UVB is turned on of the UVB, star formation significantly decreases.Bottom: Evolution of ionizing photon production rate with respect to the value in CDM.Only the results of RT simulations are shown.The late-time starburst in altDM models quickly closes the gap in ionizing photon production.The  ion in FDM surpasses the value in CDM at  ∼ 7 which drives the acceleration of reionization afterward.
In the middle panel of Figure 4, we show the redshift evolution of the star formation rate density (SFRD) in the entire simulation volume.This is derived from the initial-mass-weighted age distribution of all stellar particles collected at  = 5.Again, it is important to note that the results should not be compared to a large-volume-averaged "global" cosmic SFRD due to the limited simulation box size.In the WDM and FDM models, star formation is delayed by approximately 100 Myr, while in the sDAO model, the delay is around 200 Myr.Star formation in these altDM models can quickly "catch up" to the CDM case through a strong late-time starburst.This phenomenon will be discussed in more detail in Section 5 and Section 6.The uniform UVB model exhibits a bump in star formation at  ∼ 12 before the uniform UVB is activated.Star formation at  10 is suppressed as the neutral gas in the diffuse IGM is rapidly ionized by the background radiation.Artificial signals, such as WDM/FDM overtaking CDM at  ∼ 8, can show up since they consume less gas in the early bump of star formation.
In the bottom panel of Figure 4, we present the redshift evolution of the intrinsic ionizing photon production rates, denoted as  ion , normalized by the values in CDM.Within altDM models, the production of ionizing photons is suppressed at high redshifts compared to CDM.At lower redshifts, this disparity narrows due to the strong late-time star formation in altDM models.This trend resembles our findings for the SFRD.Notably, the FDM model surpasses CDM in ionizing photon production rate at  7, which marks the onset of a period of accelerated reionization as shown in the top panel.

HALO AND STELLAR MASS FUNCTIONS
The low-mass end of the (sub)halo mass function is directly affected by the suppression of the small-scale power spectrum in altDM models.In Figure 5, we show the halo mass function in different models at  = 6, 8, 10.The 1 confidence intervals of number densities are estimated using the Gehrels (1986) formulae and are shown as the shaded regions.Throughout this paper, the "halo" is defined as gravitational bound structures identified by the S algorithm (Springel et al. 2001), including both central and satellite galaxies.Halo mass is defined as the total mass of all particles (or cells) gravitationally bound to the halo.The high resolution achieved by the -suite allows us to make reliable predictions of halo properties down to about 5 × 10 7 M (as indicated by the vertical dashed line in the figure).This is close to the atomic cooling limit, ∼ 10 7.5−8 M (Wise et al. 2014), although we do not explicitly model the low-temperature molecular phase of gas and the Lyman-Werner (LW) radiation that dissociates molecular hydrogen.
All altDM models show suppression of the halo mass function in the low-mass end.In Appendix B, we estimate the characteristic mass scale of this suppression according to the half-power wavenumber  1/2 in the linear matter power spectrum.The half-power halo mass log ( halo ( 1/2 )/ M ) is 8.9, 8.8, and 9.4 for the WDM, FDM, and sDAO models studied in this paper, respectively.They accurately predict where the suppression in halo mass function shows up.In addition, different from the cut-off-like feature in WDM and FDM, the sDAO model predicts a steep rise of halo mass function in the lowmass end similar to CDM but with lower normalization.As discussed in Bohr et al. (2021), this is due to the strong DAOs in this model (with the first acoustic peak retaining 100% of the CDM power), which preserve some small-scale fluctuations of the primordial density field.It is known for the WDM and FDM models (or in general any model with steep suppression of power spectrum), spurious structures can form at the limiting mass  lim = 10.1 ρm ( box / White 2007), where  p is where the dimensionless power spectrum (Δ 2 ()) reaches its maximum which is about 10 ℎ/ Mpc for the WDM and FDM models studied here.This limiting mass is calculated to be log ( halo / M ) 7.8, which is consistent with the scale of the uptick in the FDM run in Figure 5.This limiting mass is comparable to the mass of the smallest resolved haloes in the simulation (defined to be haloes of mass equivalent to 100 DM particles).Therefore, in the rest of the paper, we will not consider the problems caused by spurious haloes.
Cosmic variance is important for these small-volume simulations.We estimate this effect by subsampling  = 4 Mpc/ℎ volumes in the --run (listed in Table 1) with the same resolution as -runs but larger box size,  box = 8 Mpc/ℎ.The variance in halo mass function obtained through this subsampling is then multiplied with the variance inflation factor computed in Appendix C. The cosmic variance in CDM is shown with the blue dotted dashed lines in Figure 5.The observed differences between DM models are statistically robust against the cosmic variance.
In Figure 6, we show the stellar mass functions in different models at  = 6, 8, 10.The stellar mass of a galaxy is defined as the total mass of stellar particles within  max6 .In low-mass galaxies, due to the steeply decreasing stellar-to-dark-matter-mass ratio, the statistics of stellar particles are much poorer than the DM counterparts.We mark five times the baryonic mass resolution of simulations as the typical values of 4 (insensitive to halo mass at high redshift, e.g.Bullock et al. 2001).The definition we used here is different from what most literature used.
The main reasons for picking  max as the aperture are: (1) It is less biased.The internal DM distribution in haloes is less affected by altDM physics despite slightly lowered concentration (e.g.Lovell et al. 2014;Bose et al. 2016a;Bohr et al. 2021).We avoid potential indirect influence on physical quantities coming from aperture definitions.(2)  max is a metric directly available in the S catalog.We have compared it to classical definitions e.g.stellar mass within twice the stellar-half-mass radius of the galaxy, and they agree reasonably well with no signs of systematical differences between DM models.vertical dashed line in the figure.The Poisson sampling noise is shown with the shaded region.We estimate the cosmic variance (in CDM) following the same subsampling method used for the halo mass function, and it is illustrated using blue dotted-dashed lines.They are larger compared to those in halo mass functions due to lowered number of halos that host star formation and additional variances from galaxy formation physics.
For altDM models, the suppression of the abundance of low-mass haloes is reflected in the stellar mass function, although to a less extent.For instance, a reduction of approximately 0.5 dex is observed at  * 10 6 M at  = 6 − 10 for the sDAO model.The arrows indicate the estimated half-power stellar masses of altDM models, derived using the methodology detailed in Appendix B. The mass scale of the suppression is consistent with these estimations.However, the signatures of altDM in the stellar mass functions are less prominent compared to those in the halo mass function.They are also less pronounced at lower redshift.Part of the reason is limited statistics (Poisson noise) and additional variance introduced by galaxy formation physics.Another important factor is that galaxies in altDM models are forming stars more rapidly at late times to "catch up" with CDM, as will be discussed in more detail in the following sections.
The uniform UVB model results in the suppression of the stellar mass function at  * 10 6 M at  = 6.This coincides with where altDM signatures become apparent.Once the uniform UVB is activated, reionization occurs quickly throughout the entire simulation volume, with the exception of a few dense, self-shielded clumps in massive galaxies.This rapid reionization suppresses star formation in low-mass galaxies.However, at redshift  = 10, the stellar mass function displays an increase at low masses.Prior to the activation of the UVB, the ignorance of local ionizing sources leads to enhanced star formation.This corresponds to the artificial bump of star formation at  10 in the global SFH, as shown in Figure 4.

REST-FRAME UV LUMINOSITY FUNCTION
To better illustrate the detectability of altDM signatures, in Figure 7, we show the rest-frame UV luminosity functions of galaxies in the simulations.The rest-frame UV (at 1500 Å) luminosities of galaxies are calculated using the Binary Population and Spectral Synthesis models (BPASS v2.2.1; Eldridge et al. 2017), which is consistent with the choice of the model in computing the ionizing photon spectrum.A correction is then applied for galaxies with poor sampling of SFHs (see Section 3 of Smith et al. 2022).Dust attenuation is approximated using the formula in Gnedin (2014)  For reference, the UV luminosity function from the -1 simulation is shown as the black dotted line.On the relatively bright end ( 1500 −13), -1 prediction as well as the small-volume simulations at  = 6 are in good agreement with the observations.Comparing different DM models, the suppression at the faint end due to delayed structure formation at small scales is apparent and similar to what we found in the stellar mass functions.They could be revealed by JWST observations of lense-magnified systems. as which roughly gives  rest UV =  band − 46.7 at  = 6,  rest UV =  band − 47.2 at  = 6 and  rest UV =  band − 47.5 at  = 10.For lensed fields with JWST, we assume it is 5 mag deeper than the unlensed field (corresponding to a typical lens-magnification factor ∼ 100, e.g.Ishigaki et al. 2018;Bouwens et al. 2022b).
As shown in Figure 7, the UV luminosity functions predicted by the fiducial CDM model show no signs of flattening or bending at  1500 −16 at  ≥ 6, which is consistent with the latest HFF results in Bouwens et al. (2022b).A cut-off in the UV luminosity function appears at  1500 ∼ −10, but is likely caused by the limited mass resolution of stellar particles (see the similar feature in the stellar mass function in Figure 6) and we do not explicitly model the physics that lead to a cut-off (e.g.dissociation of molecular hydrogen).
The UV luminosity functions in altDM models show suppression at the faint end (above the resolution limit).The location of the suppression can be reasonably approximated by the half-power halo mass and the scaling relation described in Appendix B. The simulation results are consistent with existing observational constraints.JWST observations conducted in a lensed field will have the potential to reveal the faint-end suppression arising from altDM physics.However, the distinct signature of alternative DM in the UV luminosity function is less prominent than in the halo mass functions.Although the abundance of low-mass galaxies is suppressed in altDM models, the remaining ones are UV brighter.There are three primary reasons for this.(1) Delayed structure formation in altDM models leads to younger stellar populations and higher UV light-to-mass ratios.(2) Owing to the lack of progenitors below the half-power mass scale, haloes in altDM models, WDM and FDM in particular, predominantly assemble through roughly equal-mass (major), gas-rich mergers of haloes near the half-power mass.(3) The delay of reionization (relative to the CDM case under the same numerical configuration) results in a larger supply of cold, neutral, self-shielded gas available for future star formation once accreted into a galaxy.This effect is partly illustrated in Figure 3).These factors collectively contribute to the accelerated formation of high-redshift galaxies in altDM models.Further evidence of this aspect will be shown in the following sections, discussing star formation efficiency and spatially resolved properties of galaxies.Similar phenomena have been found in many previous studies (e.g.Bose et al. 2016b;Corasaniti et al. 2017;Dayal et al. 2017;Lovell et al. 2018Lovell et al. , 2019;;Ni et al. 2019) with widely different theoretical approaches and galaxy formation models.
At  1500 −14 at  ≥ 8, all -simulations fail to sample rare, UV bright sources and tend to underpredict the luminosity function compared to the -1 results and observational data.This is a systematic bias due to the limited volume of -.Towards the faint end, this bias will eventually diminish and has little effect where the altDM signals appear.However, in addition to the bias, the cosmic variances here matter due to not only the small simulation volume but also the small lensed field in observations to push the detection limit down to  UV −12.We estimate this following the subsampling procedure described in Section 4 and the variance inflation factor calculated in Appendix C (only the variance of halo number density is considered, so this is just a "normalization" variance).The effect of cosmic variance is illustrated with the dotted-dashed lines in the figure .The altDM signatures at the faint end could potentially be smeared at  = 6 but are robust against cosmic variances at  8. On the observational side, the typical field of view of plain field surveys is 10 arcmin (e.g.HUDF, Beckwith et al. 2006;XDF, Illingworth et al. 2013), which corresponds to a survey volume of ∼ (60 cMpc) 3 at  ∼ 7 assuming a survey depth of Δ = 1.For lensed fields, the effective survey area will be much smaller on the source plane.For the HFF campaign (Lotz et al. 2017), the typical survey area is 1 arcmin in the source plane at  ∼ 7 (Ishigaki et al. 2018), which corresponds to a survey volume of ∼ (13 cMpc) 3 .Neglecting other sources of observational uncertainties (e.g.completeness, lensing corrections), the cosmic variance will be less problematic compared to the small simulated volume.In the low-mass and intermediate-mass bins, we find patterns that are consistent with the evolution of the overall SFRD shown in Figure 4. Compared to CDM, star formation in the altDM models is suppressed at early times but eventually "catches up" and surpasses the CDM value in the later stages.This late-time starburst is particularly noticeable in the FDM model, where the SFR exceeds the CDM value by a factor of four times at  = 6.The uniform UVB runs display a bump in star formation at  10 regardless of the DM model employed.This is followed by a decline immediately after the UVB activation.Since the differences between DM models primarily arise from structure formation on small scales, they are more prominent in low-mass bins.In the high-mass bin, which is above the half-power mass of all three altDM models, the SFHs converge to the CDM results at  10.
In line with our findings for the stellar mass functions, the uniform UVB model suppresses galaxy abundance at the faint end following the activation of the UVB.By comparing the uniform UVB and the fiducial -run, we effectively bracket the uncertainties arising from the patchy reionization process.In regions close to strong ionizing photon sources, the abundance of faint galaxies will be suppressed similarly to what we observe in the uniform UVB model.In contrast, for regions in voids, only local sources contribute to reionization, and the fiducial -run represents this scenario.The uncertainties associated with the morphology of reionization are larger than the "normalization" uncertainties stemming from the variance of the matter density field (halo number count) as discussed in the paragraph above.Therefore, reliable modeling of the reionization process is essential for predictions of low-mass galaxy populations and quantitative comparisons with observational results.

Star formation history
In Section 3 (see Figure 4), we presented the history of the SFRD in the entire simulation volume.In Figure 8, we calculate the average SFH of haloes in three mass bins selected at  = 5.The bin edges are chosen as 10 8.5 M and 10 9.5 M .For reference, the half-power masses of the three altDM models studied lie between 10 8.5 M and 10 9.5 M , as shown in Figure 5.In the high-mass bin, which is above the half-power masses of all three altDM models, there is a decrease in SFR at  10 in altDM models, which can be attributed to the scarcity of low-mass progenitors at that time.However, at later times ( 10), the average SFHs of haloes appear indistinguishable between the various DM models.In both the low-mass and intermediate-mass bins, we observe a systematic suppression of star formation in the early stages in altDM models.The strength of this suppression is positively correlated with the half-power mass, and the sDAO model, which has the largest half-power mass, demonstrates the most substantial suppression of SFR.However, at later times, starbursts occur in low-mass haloes in altDM models, allowing them to "catch up" with CDM and eventually surpass it.This observation aligns with our findings and discussions about the UV luminosity functions in Section 5.In the FDM model, which features the steepest power spectrum damping, the average SFR in low-mass haloes can be approximately four times higher than that in CDM counterparts at  = 6.The strength of this late-time starburst depends on the shape of the power spectrum damping, rather than simply scaling with the half-power mass.The enhanced star formation rates in altDM models could weaken the distinct signature of altDM in luminosity functions and global star formation/reionization history constraints (e.g.Lovell et al. 2019;Khimey et al. 2021).The dashed lines in Figure 8 show the SFH in the uniform UVB model.Similar to the global evolution of SFRD in this model, there is an increase (decrease) in SFR before (after) the activation of the UVB.This bump-to-suppression feature is amplified in the low-mass halo bin, as is the diversity of SFHs across different DM models.These artificial features in the SFH can contaminate the signature of altDM significantly, which emphasizes the importance of accurately modeling reionization at small scales.
To better illustrate the impact of altDM on properties of the stellar populations, we focus on galaxies at  = 8 where the altDM models exhibit intense starbursts in low-mass haloes and clearly deviate from CDM, as demonstrated in Figure 8.We group them in three host halo mass bins log ( =8 halo / M ) ∈ [7.5, 8.5], [8.5, 9.5] or [9.5, ∞].After stacking the SFH of the galaxies in each group, we determine three key metrics: (1) the fraction of stellar mass formed before  = 10 denoted by  M(z>10); (2) the lookback time that 50% of the stellar mass is formed (mass-weighted median age) represented by  50 ; (3) the lookback time that 90% of the stellar mass is formed represented by  90 .In Figure 9, we show  M ( > 10) versus  90 and  50 in different models.For galaxies with log ( =8 halo / M ) < 9.5, altDM models generate younger stellar populations.For galaxies with log ( =8 halo / M ) ∈ [7.5, 8.5], less than 20% of the stellar mass are formed before  = 10 and the median stellar age is about 80 Myr in the three altDM models, which is about 30% lower than that of the CDM counterparts.These differences are weaker in relatively massive haloes with log ( =8 halo / M ) > 9.5, which is above the halfpower mass of these models.The change in stellar age is systematic, with  90 displaying similar trends as found for  50 .In the uniform UVB model, an artificial bump of star formation at  10 causes  M ( > 10) and  50 to increase up to approximately 0.6 and 200 Myr respectively, regardless of the halo mass.The shaded regions in the figure show the relationship assuming a constant (lower boundary) and an exponentially-rising SFH (upper boundary).The fact that all the points (except for one in the uniform UVB model) lie between the relationship indicates galaxies in general have rising SFHs.The relatively massive galaxies tend to have more steeply-rising SFHs.
The SFH of galaxies will be constrained by JWST observations.Some attempts have been made on bright galaxies using SED fitting on legacy survey data (e.g.Laporte et al. 2021;Tacchella et al. 2022b;Whitler et al. 2023b) or early JWST data (e.g.Tacchella et al. 2022a;Furtak et al. 2023;Whitler et al. 2023a).They reveal the steeply rising SFHs of galaxies in the EoR that is qualitatively consistent with our simulation results.The exact age distribution of samples varies from studies and shows non-trivial dependence on the SFH priors (e.g.Tacchella et al. 2022b).

Star formation efficiency
In the top panel of Figure 10, we show the specific star formation rate (sSFR≡ SFR/ * ) as a function of halo mass.The SFR (instantaneous in gas cells) and  * are all measured using  max as the aperture, following our definition of galaxy stellar mass in Section 4. This choice is less affected by baryonic processes and their non-linear interplay with altDM physics.In general, altDM models predict higher sSFR in haloes below the half-power mass at both  = 6 and  = 10, which implies younger stellar populations.The enhancement in sSFR in low-mass haloes is more pronounced in models with steeper dampings of the linear matter power spectrum (i.e.FDM followed by WDM).As discussed in Section 5, similar phenomena have been found in previous studies of high-redshift structure formation in altDM (e.g.Bose et al. 2016b;Corasaniti et al. 2017;Lovell et al. 2018Lovell et al. , 2019)).But it is the first time that all relevant physics for cosmic reionization are included in a self-consistent fashion with a minimum level of fine-tuning.For example, in semi-analytical studies (e.g.Corasaniti et al. 2017), the enhanced star formation efficiency in low-mass haloes is in fact a consequence of tuning the galaxy formation model in altDM to yield the correct timing of reionization, which is realized in our simulation out-of-the-box.In previous cosmological hydrodynamic simulations (e.g.The specific star formation rate (sSFR) versus halo mass at  = 6 (solid lines) and  = 10 (dotted lines).In altDM models, the sSFRs are higher for halo masses below the half-power mass, ∼ 10 9 M , suggesting younger stellar populations.This enhancement is already noticeable at  = 10.Middle: The gas depletion time  depl.versus halo mass.The depletion time is defined as  gas /SFR.In altDM models, the depletion times are systematically lower than in CDM below the half-power mass.This indicates that a larger fraction of gas in the halo is in a phase suitable for star formation.However, this trend is less evident at  = 10.Bottom: The gas mass ratio  gas versus halo mass. gas is normalized over the universal baryon fraction  b ∼ 0.16.Gas is more abundant in altDM models likely due to more gas-rich major mergers.The uniform UVB model results in extremely small gas fractions in galaxies with  halo < 10 8.5 M , leading to the rapid decline of depletion time observed in the middle panel.

2019)
, where the uniform UVB model was employed, differences in the phase and morphology of reionization as part of the cause of the late-time starburst were not captured.In fact, Ni et al. (2019) found that the luminous fraction of haloes at  ∼ 6 in FDM is lower than the CDM cases, which drives larger discrepancies in the UV luminosity function instead of closing the gap as we have found.
In the middle panel of Figure 10, we present the gas depletion time scale versus halo mass.The gas depletion time is defined as  depl =  gas /SFR, where  gas is defined as the total gas mass within  max .In low-mass galaxies, the depletion time is of order 10 − 100 times the Hubble time  H (), comparable to the dynamical time of the haloes.It is consistent with the 1 − 10% star formation efficiency per free-fall time from the local scaling relations (e.g.Silk 1997; Kennicutt 1998).In altDM models, the gas depletion times are shorter, consistent with their higher star formation efficiency, suggesting not only a larger gas supply but also a higher fraction of the gas in the star-forming phase.For the IllustrisTNG galaxy formation model employed, it means that the neutral gas density is closer to the star formation threshold.The trend is less obvious at  = 10.In the uniform UVB model, the depletion time rapidly declines in low-mass haloes, but it is mainly due to the shortage of gas as will be shown in the bottom panel.
In the bottom panel of Figure 10, we show the gas fraction (normalized using the universal baryon fraction,  b ≡ Ω b /Ω 0 ∼ 0.16) versus halo mass. gas /  b is much below unity since we only measure gas mass within  max rather than the entire halo.The gas reservoir is more abundant in altDM models.The halo mass scale where they deviate from the CDM case is also consistent with the half-power mass.It supports the hypothesis that the haloes below the half-power mass in altDM are assembled mainly through gas-rich major mergers.In the uniform UVB model, the gas abundance is significantly suppressed at  halo 10 8.5 M due to the strong ionizing photon bath, leading to the rapid decline of depletion time seen in the middle panel.

Mass-metallicity relation
The delayed star formation in altDM models will leave imprints on the metal enrichment of galaxies during or by the end of reionization.Direct observational constraints on galaxy stellar-phase metallicity can be inferred through SED fitting (e.g.Tacchella et al. 2022a;Curtis-Lake et al. 2022;Furtak et al. 2023) and constraints on the gas-phase metallicity can be derived by emission line measurements (e.g.Curti et al. 2023;Fujimoto et al. 2023).
In the left panel of Figure 11, we show galaxy stellar metallicity versus stellar mass (known as the stellar-phase mass-metallicity relation, MZR) of simulated galaxies at  = 6.Results at higher redshifts are similar.The sDAO model results in lower stellar metallicities than CDM at the 1 level for galaxies with  * ∼ 10 6−6.5 M , due to the delayed structure formation.However, the WDM model yields an MZR comparable to that in CDM.The FDM model generates a similar MZR with slightly higher metal enrichment than CDM at  * 10 6 M .The MZRs in these models are affected by the latetime starbursts, which are not as pronounced in the case of sDAO.On the other hand, the uniform UVB model predicts a systematically higher MZR than the CDM run with RT.This discrepancy can be attributed to the artificial star formation that occurs at redshifts  10 before the UVB is activated.
For comparison, we show the early JWST results from Furtak et al. ( 2023) using SED fitting.The galaxies in this sample are in the redshift range 10  15 and are lensed behind the cluster SMACS J0723.3-7327(Atek et al. 2023) with the broadband photometric data from NIRCam and NIRISS fitted by BEAGLE7 (Chevallard & Charlot 2016).The stellar metallicities of galaxies in simulations are in general consistent with these high-redshift sources.The lowmass galaxies at high redshift are quenched by the ionizing radiation background and their metal enrichment patterns will be almost frozen afterward.These galaxies are considered progenitors of the (ultra-)faint dwarf satellite galaxies observed in the Local Group (the Milky Way and M31, e.g.Simon 2019).For comparison, we show the metallicity measurements of local dwarfs presented in Kirby et al. -1 result is shown as the gray dotted line.The faint blue dash-dotted lines indicate cosmic variance.In addition, we show the early JWST measurements of  10 galaxy candidates using SED fitting (Furtak et al. 2023) and the measurements of the Local Group dwarfs compiled in Kirby et al. (2013).Right: Galaxy gas metallicity versus stellar mass at  = 6.The plotting style and label choices are the same as the left panel.For comparison, we show the recent JWST results from Fujimoto et al. (2023) as well as the local dwarf constraints from Lee et al. (2006).The gas and stellar-phase metallicities resemble each other quite well.At  * 10 7 M , all the models give converged MZR.At  * 10 7 M , galaxies in the sDAO model are less metal enriched by about 0.2 dex at mass scale  * ∼ 10 6−7 M while the ones in the uniform UVB model are more metal enriched by about the same level.Signatures are weaker for the WDM and FDM models since the recent burst of star formation can enrich galaxies to the same level as CDM or even more.The simulation results are in good agreement with the latest observational constraints of high redshift galaxies with JWST as well as the local dwarfs, within the scatter of the observed MZR.
(2013)8 , which were derived from Fe abundances of individual stars.We directly convert [Fe/H] to  * / (assuming  = 0.0126; Asplund et al. 2004).The MZRs of simulated galaxies in all the RT models are within the scatter of the observed local dwarfs at  * = 10 5−7 M .The uniform UVB model predicts higher MZR than the local dwarfs.This discrepancy will become even worse considering the potential late-time enrichment of dwarfs.
At the relatively massive end ( * 10 7 M ), the predictions of all models converge and fall within the scatter of the observed sample.The cosmic variance effects are estimated using the same method as in Figure 5 and Figure 7.Meanwhile, there is a clear bias in the small-volume runs, with the -1 prediction being approximately 0.2 dex below the -results.This is due to the absence of massive bright sources in the small volume, which delays reionization and allows for stronger early-phase star formation and metal enrichment.It is worth noting that the observed MZR at  * 10 5 M flattens, which is in contrast to the prediction of cosmological simulations that suggest continuously decreasing metallicities in the low-mass end (e.g.Wheeler et al. 2019).However, the baryonic mass resolution in our simulations is not sufficient to investigate this phenomenon.In future studies, it would be intriguing to explore how patchy reionization or alternative dark matter models might impact the metal enrichment of these low-mass galaxies.
In the right panel of Figure 11, we show the gas-phase MZR of galaxies in simulations.The hydrogen and oxygen abundance are directly taken from gas cells in simulations.The gas-phase MZRs in simulations closely resemble the stellar-phase ones.The Type-II supernovae rates and the continuous mass-loss rate through Asymptotic Giant Branch (AGB) stellar winds drop significantly by ∼ 100 Myr (e.g.Portinari et al. 1998;Leitherer et al. 1999;Vogelsberger et al. 2013), both being the dominant sources of metal enrichment.Therefore, the metal recycling process happens at a relatively small time scale compared to the lifetime of these galaxies.Since these galaxies are predominantly influenced by the recent period of star formation, it is not surprising that the stellar and gas-phase metallicities show great synergy.The simulation results are compared to the JWST spectroscopically-confirmed sources in Fujimoto et al. (2023) and the Local Group dwarfs in Lee et al. (2006).The sources in Fujimoto et al. (2023) are measured using O[III] and H emission lines as tracers.The simulated massive galaxies exhibit slightly higher gas metallicities compared to the JWST sources but within the scatter.Overall, the influence of altDM models is not statistically significant on the MZR plane.A major reason for this is that the late-time starbursts in altDM models wash out the altDM signatures in metal enrichment levels.The precision and statistics of current observations are not sufficient to differentiate between altDM models.

Spatially-resolved properties
The late-time starburst found in altDM models can influence the stellar distribution within galaxies and the relative contribution of in-situ and ex-situ star formation.In the following, we investigate the spatially resolved properties of galaxies.To enhance the statistics of stellar particles and gas cells, we stack galaxies at  = 8 in three halo mass bins, log ( =8 halo / M ) ∈ [7.5, 8.5], [8.5, 9.5] [9.5, ∞].To quantify galaxy sizes, we measure the half-SFR radius ( SFR 1/2 ), which encloses half of the total SFR within the halo.We combine two proxies for SFR, the instantaneous SFR in gas cells and the ar- Top: Surface density of SFR versus the half-SFR radius for galaxies at  = 8 stacked in three halo mass bins.The half-SFR radius is a proxy for the half-light radius in rest-frame UV which traces recent star formation.We combine the SFR calculated using the instantaneous SFR in gas cells and the one using masses of young stars with age 100 Myr.The FDM model results in smaller galaxy sizes and enhanced SFR surface density.The signature is less prominent in the WDM and sDAO models.Bottom: Fraction of stellar mass with metallicity below 10 −6 versus the metallicity gradient of the same groups of galaxies in the top panel.The metallicity gradient is represented by the difference in the mean metallicity measured within and outside  * 1/2 .Galaxies in the WDM and FDM models have larger metallicity gradients.In low-mass haloes, the fractions of low-metallicity stars are appreciably smaller in WDM and FDM.These findings are consistent with the late-time starbursts, which predominantly occur in the central regions of galaxies.chaeological SFR from young stellar particles (with age 100 Myr).The half-SFR radius offers a better analogy to the half-light radius observed in rest-frame UV, although projection effects could still cause some order-unity differences.In the top panel of Figure 12, we show the SFR surface density (calculated within  SFR 1/2) as a function of  SFR 1/2.The  SFR 1/2 is translated to angular distance on the sky using the angular diameter distance based on the assumed cosmology.In the FDM model where the late-time starburst is most prominent, we find more compact galaxy sizes in host haloes with log ( halo / ) 8.5 compared to their CDM counterparts.In addition, we find enhanced SFR surface densities across all three halo mass bins, including the one above the half-power mass.This increased compactness is likely a result of rapid ongoing star formation in the central regions of galaxies fueled by a series of gas-rich major mergers, as well as a reduced mass of ex-situ stars formed in progenitors.We find similar signatures in the WDM model for relatively low-mass haloes (log ( halo / ) 8.5).However, the sDAO model exhibits less apparent signatures, as it does not show comparably strong late-time starbursts at  = 8 (see Figure 8).
Compact galaxy sizes (roughly 100 pc to 1 kpc) are found for the low-mass ( * 10 7 M ) galaxies in our simulations, regardless of the DM models.These results are consistent with extrapolations (in terms of stellar mass and redshift) of observational findings (e.g.van der Wel et al. 2014;Shibuya et al. 2015) and theoretical expectations (e.g.Mo et al. 1998).Galaxy size is an essential factor in forward modeling lense-magnified sources and determining the luminosity function on the image plane (e.g.Bouwens et al. 2017;Kawamata et al. 2018;Ishigaki et al. 2018;Bouwens et al. 2022a).The upcoming JWST observations will offer more robust morphological constraints on these compact, low-luminosity galaxies.Attempts have been made for the relatively luminous sources (e.g.Ono et al. 2022;Robertson et al. 2022;Yang et al. 2022;Tacchella et al. 2023).These early JWST sources at  UV −18 are all compact with effective radii between 100 pc and 1 kpc.The face values are consistent with our results, but the extrapolation to faint luminosities and potential selection biases are still uncertain.
In the bottom panel of Figure 12, we show the fraction of lowmetallicity stars (log  * < −7) versus the stellar metallicity gradient at the outskirt of galaxies.The metallicity gradient is estimated by calculating the difference between the mass-weighted stellar metallicity measured within and outside the stellar-half-mass radius  * 1/2 .In all three halo mass bins, galaxies in the FDM model exhibit larger stellar metallicity gradients.Similar to what we found for galaxy sizes, the difference is evident above the half-power mass, where the halo mass function has already converged to the CDM results.The metallicity gradient is also driven by centrally dominated late-time star formation as well as lowered metallicities of the ex-situ population of stars at the outskirt of galaxies.Similar trends are found for gas-phase metallicity gradients as well.These signatures are less noticeable in the WDM model and are absent in the sDAO model, where the late-time starburst is not comparably strong at  = 8.Furthermore, in low-mass haloes with log ( halo / M ) ∈ [7.5, 8.5], the WDM and FDM models show a smaller fraction of extremely low-metallicity stars compared to CDM, which are inherited from the earliest generation of star formation.

CAN WE DISENTANGLE DM PHYSICS WITH THE UNCERTAINTIES IN THE BARYONIC SECTOR?
In addition to the altDM physics explored in this paper, several other astrophysical factors could also influence the faint-end UV luminosity functions as well as properties of low-mass galaxies in the EoR.The leading three are discussed in the following.
(1) Atomic cooling limit.The low-temperature molecular hydrogen (which is responsible for gas cooling at low temperatures,  10 4 K) can be dissociated with a moderate strength of LW radiation background sourced by PopIII stars (e.g.Haiman et al. 1997Haiman et al. , 2000;;Machacek et al. 2001).Therefore, PopII star formation can only take place in relatively massive haloes ( 10 6 M ) that can host a reservoir of self-shielded, molecular gas or with atomic cooling channels enabled (e.g.Wise & Abel 2007;Greif et al. 2008;O'Shea & Norman 2008).Flattening of the faint-end luminosity functions were predicted at  UV −12 in simulations with explicit models of PopIII star formation & feedback and LW radiation background (e.g.Wise et al. 2014;O'Shea et al. 2015;Ocvirk et al. 2016).In -, the molecular phase of gas and low-temperature cooling are not directly modeled but are partially captured by the effective equation of state of the ISM gas (Springel & Hernquist 2003).As a result, the birth halo masses of stellar particles in our simulations are all above ∼ 10 7.5 M .The mass limit is empirically consistent with the atomic cooling limit ( vir ∼ 10 4 K) but larger than the limit found in e.g.Wise & Abel (2007) due to the lack of PopIII star physics.
(2) Feedback from star formation.This includes stellar (presupernovae) and supernovae feedback that can drive galactic-scale winds as well as radiative feedback from young massive stars (e.g.Gnedin 2000;Wise & Cen 2009;Hopkins et al. 2013;Vogelsberger et al. 2013).Feedback is important in regulating star formation within the halo, leading to systematically lower star formation efficiency in low-mass haloes (e.g.Behroozi et al. 2013).Negative feedback changes the kinematics and the phase of the gas, which limits the effectiveness of gas cooling and gradually boils away the star-forming gas reservoir (e.g.Efstathiou 1992;Thoul & Weinberg 1996;Dijkstra et al. 2004;Shapiro et al. 2004).Radiative feedback such as photoionization and photoheating can have effects on large scales, reducing the neutral gas supply of galaxies in the environment and suppressing star formation activities at later times (e.g.Jeon et al. 2014Jeon et al. , 2015)).These mechanisms all tend to reduce the stellar mass of low-mass haloes or equivalently the number density of galaxies of fixed stellar mass.On the other hand, as galaxies grow more massive, the feedback is eventually dominated by local sources and star formation becomes self-regulating (Wise & Abel 2008;Schaye et al. 2010;Ostriker et al. 2010;Hopkins et al. 2014;Rosdahl et al. 2015;Hopkins et al. 2023), removing the signatures of altDM models at late stages of galaxy evolution.
(3) Morphology of reionization.The phase and morphology of reionization have a strong influence on the abundance of faint galaxies, as emphasized in the previous -study (Borrow et al. 2022).In volumes embedded in the ionizing bubble created by massive bright galaxies, ionizing radiation can deplete the reservoir of cold, neutral gas in low-mass galaxies, thereby suppressing star formation.For example, in the uniform UVB model we tested, the faintend UV luminosity function is elevated (suppressed) before (after) the activation of the UVB, at similar luminosities where altDM signatures appear.We found related signatures in the SFRD evolution and stellar mass functions discussed in previous sections.Conversely, if the volume is situated in a void of ionizing sources (as effectively occurs in the fiducial -runs), the phase of reionization will be delayed, creating a more favorable environment for star formation in low-mass galaxies.This introduces environment-dependent bias and an additional source of uncertainty in the galaxy-halo connection, resulting larger cosmic variance of the UV luminosity functions than the halo mass function.
The physical processes outlined above can generate similar suppression features in the faint end of galaxy UV luminosity function.However, they all impose negative feedback on star formation in a coherent and systematic manner.For instance, a strong, uniform bath of ionizing radiation from bright sources in the environment will consistently advance the phase of reionization and suppress star formation at all stages of reionization.Similarly, a model featuring stronger stellar winds or supernovae or radiative feedback from local sources tends to push gas out of the star-forming phase and suppress star formation in galaxies.For both cases, the negative feedback on star formation is independent of the global phase of halo assembly.
However, as we found for altDM models, a late-time starburst takes place in low-mass galaxies and allows altDM universes to rapidly "catch up" the global phase of reionization.It is effectively a positive feedback that is only associated with low-mass galaxies at late times of the EoR.Importantly, this phenomenon is rooted in the hierarchical assembly pattern of DM haloes rather than merely altering the connection between gas/star formation and haloes.Meanwhile, it is linked to various changes in galaxy properties we observed, including higher star formation efficiency (resulting in younger stellar populations), higher SFR surface density, more compact galaxy sizes, and stronger metallicity gradients at the outskirts of galaxies.The signatures of altDM could be disentangled with uncertainties arising from astrophysical processes through joint analysis of luminosity function measurements as well as constraints on galaxy SFH and morphology.Independent constraints on reionization will also be beneficial.The global signal of the ionized phase from e.g.CMB optical depth, Lyman- emitters at high redshift can help constrain the sub-grid model parameters of the reionization model.In addition, for luminosity function measurements in gravitationally lensed fields, it is important to constrain the local phase of reionization.It will be particularly helpful to reduce the uncertainties from the morphology of reionization.Such measurements can be achieved through Lyman- forest tomography (e.g.Delubac et al. 2015;Yang et al. 2020;Bosman et al. 2022) and the intensity mapping of the 21 cm spinflip transition of the hydrogen atom (Furlanetto et al. 2006;Mellema et al. 2006;Parsons et al. 2010;Pritchard & Loeb 2012).

CONCLUSIONS
In this paper, we utilized -, a suite of high-resolution, cosmological, hydrodynamic simulations of galaxies in the EoR, to explore the impact of alternative DM models that suppress the linear matter power spectrum.The inclusion of the radiative transfer method alongside the galaxy formation model and the high numeric resolution featured in these simulations provide a solid foundation for us to study various galaxy properties as well as the morphology of reionization in altDM models.In this paper, we focus on the properties of galaxies that can be revealed by the JWST.Our major findings can be summarized as follows.
• Galaxy abundance: In altDM models, the imprint of the suppression in the linear matter power spectrum is found in the halo mass functions, the stellar mass functions, and the UV luminosity functions at the low-mass/faint end.The mass/luminosity scale of the suppression correlates well with the half-power wavenumber  1/2 (see Appendix B for details).In addition, the shape of the suppression in galaxy abundance reflects the shape of the damping wing of the matter power spectrum.For example, a steeper damping of the power spectrum in FDM compared to WDM results in a steeper suppression in the mass functions as well.The residual acoustic oscillation features result in continuously rising mass functions in the low-mass end.JWST observations of lense-magnified fields are promising in detecting the suppression features due to altDM models, but one needs to be aware of various astrophysical effects influencing low-mass galaxy abundance.
• Complexities from reionization modeling: In terms of the abundance of faint galaxies, one important astrophysical uncertainty is the morphology of the reionization process.In simulations using a uniform UVB model, a strong bath of ionizing photons contributed by massive bright non-local sources is applied to low-mass haloes.As a result, the abundance of low-mass, faint galaxies at  * 10 6 M or  1500 −14 are suppressed.The location of this suppression is similar to where altDM signatures show up.Therefore, it is important to correctly model the morphology of reionization to make robust predictions for the faint-end luminosity functions and use that to constrain the nature of DM.
• Complexities from non-linear effects: In all altDM models, in particular the ones with steep power spectrum dampings, we find positive feedback on the late-time star formation in galaxies with host halo mass below the half-power mass.This manifests as brighter UV luminosities, higher sSFR, and gas abundance as well as shorter gas depletion time in haloes below the half-power mass scale.The latetime starbursts allow the SFRD and ionizing photon production rates in altDM to "catch up" with or even surpass the CDM value at lower redshifts.Notably, the volume-weighted ionized fraction in the FDM model is accelerated to exceed the CDM counterpart at  ∼ 6.
• Metal enrichment: The delayed star formation in altDM models can lead to suppression of metal abundances in low-mass galaxies until the late-time starbursts take over.For the sDAO model, the signal is at the 1 level at the mass range  * ∼ 10 6−6.5 at  = 6.The WDM and FDM models on the other hand experience a boost in late-time star formation and metal enrichment, and can thus "catch up" with CDM or even exceed it.However, the differences are still encoded in the spatial distribution of metals in galaxies.In altDM models, especially the FDM model, stellar/gas metallicities show stronger gradients at the outskirts of galaxies, and the population of "old" stars with extremely-low metallicities is suppressed.
• Galaxy sizes: In altDM models, we find more compact galaxy sizes in low-mass haloes below the half-power mass.This is associated with enhanced SFR surface density.The differences are driven by the late-time starburst, which takes place in the central regions of galaxies, as a consequence of more frequent gas-rich major mergers in altDM models.
• Ways to disentangle DM physics and astrophysical uncertainties: Many astrophysical processes could generate suppression features in the faint-end luminosity function, which can contaminate the altDM signatures.Two examples are a strong ionizing photon background and stronger stellar and supernovae feedback.Most of these processes will cast negative feedback on star formation in a coherent and systematic fashion.However, for altDM models, we find positive feedback on star formation that is only associated with low-mass galaxies at late times of the EoR.It is a unique feature of altDM physics that delays early structure formation and is linked to various changes in galaxy SFH and morphology.A joint analysis of the luminosity functions and galaxy properties could help disentangle these effects.
In this study, we investigated the properties of faint galaxies during the EoR with a focus on the influence of altDM models.Our findings corroborate several previous assertions, such as the reduced abundance of faint galaxies and the promise of detecting this with JWST observations.By combining the radiation-hydrodynamic solver and the IllustrisTNG galaxy formation model, we can study the effects of reionization on low-mass galaxies, as well as their reciprocal influence on reionization, within altDM models.We identify enhanced late-time star formation in altDM models, consistent with earlier studies, which arises from both different patterns of hierarchical assembly and delayed reionization process in altDM.This can introduce non-linear effects on the reionization history at small scales.We also draw attention to the complexities inherent in both reionization and galaxy formation models and propose ways to disentangle them with altDM physics.
Looking ahead, it is crucial to further explore the impact of baryonic physics in a more controlled setting.We plan to conduct zoom-in RHD simulations with multi-phase ISM physics, allowing us to investigate the impact of diverse star formation and feedback prescriptions on low-mass galaxies, providing valuable insights into the interplay between the "luminous" and DM physics.In addition, simulations with larger and more representative volumes will allow us to make quantitative predictions for the properties of the IGM.This will enable us to examine metrics such as the HI clumping factor and make predictions for the Lyman- forest and 21 cm line intensity mapping.The DM distribution on large scales is nearly identical to the fiducial CDM run using physics.Bottom: Neutral hydrogen column density and gas temperature distribution at  = 6.The IGM is almost completely ionized after the uniform UVB is activated, with the exception of some dense, selfshielded gas in the most massive halo and the surrounding filaments.The IGM temperature has been elevated to 10 4 K by the strong ionizing radiation bath.
Empirically, we find  0 = 1.9 gives the best prediction on where the suppression on UV luminosity functions shows up.We use this relation to estimate the characteristic UV luminosity where differences between DM models show up.

APPENDIX C: ESTIMATING COSMIC VARIANCE
In Section 4, we describe the subsampling approach to estimate cosmic variances.Here we provide details on estimating the variance inflation factor.The halo mass function in a biased density field is known as the conditional halo mass function.As shown in e.g.Mo & White (1996); Trapp & Furlanetto (2020) This paper has been typeset from a T E X/L A T E X file prepared by the author.
Figure1.Linear matter power spectrum for the CDM (blue curve), WDM ( WDM = 3 keV; red curve), sDAO (green curve), and FDM ( a = 2 × 10 −21 eV; purple curve) models.Δ 2 () ≡  3  ()/2  2 is the dimensionless power spectrum.For the CDM and FDM models, the linear matter power spectrums are generated using C and A C .For the WDM and sDAO models, we employ analytical transfer functions (see Section 2.1 for details) on top of the CDM power spectrum.The altDM models examined in this study all have small-scale damping of the power spectrum.They differ in the location and shape of the damping wing as well as the oscillatory patterns.

Figure 2 .
Figure 2. DM surface density map at  = 6 centered on the most massive halo in the simulation volume.Here we show the runs with the fiducial physics (with RT).The DM density fields in runs with uniform UVB are identical to their RT counterparts.In the top left corner of each panel, we show a zoom-in image of the structures around the massive halo.Due to the suppressed small-scale power spectrum in altDM models, fewer low-mass haloes and filamentary structures formed.The sDAO model has residual fluctuations at small scales due to the DAOs.It has subtle differences from the WDM and FDM models with a single damping of the power spectrum.The overall density field in sDAO is smoothed out while a limited number of small-scale structures can still form.
, we show the column density map of neutral hydrogen at  = 6 in the simulations.The camera position and the slice of volume visualized are identical to the DM image (Figure 2).The neutral gas in the IGM is a highly-biased tracer of the underlying DM distribution (e.g.Furlanetto et al. 2006; Pritchard & Loeb 2012; Rahmati et al. 2013; Villaescusa-Navarro et al. 2018).The high-density neutral gas

Figure 3 .
Figure 3. Overview of the properties of the IGM in simulations at  = 6.Each panel corresponds to one DM model studied.In each panel, the top left part shows the column density distribution of neutral hydrogen.The lower right part shows the (mass-weighted) gas temperature of the same volume as the top left part mirrored.The CDM model has more small clumps of neutral gas retained in subhaloes even in the highly-ionized regions.These low-mass subhaloes and associated cold gas clumps are missing in altDM models.At large scales, WDM and sDAO runs have a higher neutral gas abundance on average and smaller ionizing bubbles at  = 6.The global reionization processes in these two models have been delayed systematically due to the deficiency of ionizing photons contributed by low-mass galaxies.Corresponding features can be found in the temperature distribution as well.The WDM and sDAO runs show larger voids of low-temperature, neutral gas while a lower number of small hot gas bubbles.In the FDM run, although the statistics of high-density gas clumps, filaments, and hot gas bubbles are similar, we find rather different phenomena on large scale.The global ionization fraction at  = 6 has overtaken the CDM value due to the late-time starburst (shown in Figure4as well).These strong non-linear effects create larger ionizing bubbles and lower neutral gas column densities in FDM at  6.

Figure 4 .
Figure 4. Top: Evolution of the global volume-weighted neutral hydrogen fraction in simulations.The solid lines show the results from the simulations with RT, while the dashed lines show the results in the uniform UVB model.In altDM models, reionization histories display systematic delays that correlate with the damping scale of the power spectrum.Interestingly, in the FDM model with the sharpest power spectrum damping, late-time starburst in low-mass galaxies can lead to a global ionized fraction larger than the CDM counterpart at  6. Middle: Evolution of the star formation rate density (SFRD) in the simulation volume.Star formation in the altDM models is suppressed at high redshift due to delayed structure formation at small scales until the late-time starburst closes the gap.In contrast, the uniform UVB model exhibits a bump in SFRD before the activation of the uniform UVB.After the UVB is turned on of the UVB, star formation significantly decreases.Bottom: Evolution of ionizing photon production rate with respect to the value in CDM.Only the results of RT simulations are shown.The late-time starburst in altDM models quickly closes the gap in ionizing photon production.The  ion in FDM surpasses the value in CDM at  ∼ 7 which drives the acceleration of reionization afterward.

Figure 5 .
Figure 5. Halo mass function at  =6, 8 and 10 from top to bottom.We show the halo mass functions from the simulations in different DM models.For each model, results from the model (the uniform UVB model) are shown as solid (dashed) lines.The faint blue dotted-dashed lines show the cosmic variance estimated in CDM.The halo mass function in the -1 simulation is shown as the black dotted line.The arrows at the bottom of each panel indicate the half-power mass of altDM models.Below the half-power mass, we find suppression to the halo mass function due to altDM physics, of which the shape depends on the shape of the damping feature in the power spectrum.

Figure 6 .
Figure 6.Galaxy stellar mass function at  =6, 8 and 10 in different models.The labeling follows Figure 5.The vertical lines indicate the mass roughly corresponding to 5 stellar particles in the galaxy.Similar to the halo mass function, we find suppression at the low-mass end of the stellar mass functions in altDM models below the half-power stellar mass.The signature is less prominent due to increased statistical noises and cosmic variance.The uniform UVB model can suppress the stellar mass function in a similar fashion to altDM after the activation of the UVB.

Figure 7 .
Figure 7. Rest-frame UV luminosity function at  = 6, 8, 10.The UV luminosities are calculated using the population synthesis code and an empirical dust correction.The cosmic variance effect is illustrated with the dotteddashed line.The observational constraints are taken from the compilation in Vogelsberger et al. (2020b) and Bouwens et al. (2022b) (see the main text for details).The vertical lines show the detection limits of HST and JWST.For reference, the UV luminosity function from the -1 simulation is shown as the black dotted line.On the relatively bright end ( 1500 −13),-1 prediction as well as the small-volume simulations at  = 6 are in good agreement with the observations.Comparing different DM models, the suppression at the faint end due to delayed structure formation at small scales is apparent and similar to what we found in the stellar mass functions.They could be revealed by JWST observations of lense-magnified systems.

Figure 8 .
Figure 8.Average SFH of haloes selected at  = 5 binned by halo mass.We choose three mass bins log (  =5 halo / M ) ∈ [7.5, 8.5], [8.5, 9.5] and [9.5, ∞].In the low-mass and intermediate-mass bins, we find patterns that are consistent with the evolution of the overall SFRD shown in Figure4.Compared to CDM, star formation in the altDM models is suppressed at early times but eventually "catches up" and surpasses the CDM value in the later stages.This late-time starburst is particularly noticeable in the FDM model, where the SFR exceeds the CDM value by a factor of four times at  = 6.The uniform UVB runs display a bump in star formation at  10 regardless of the DM model employed.This is followed by a decline immediately after the UVB activation.Since the differences between DM models primarily arise from structure formation on small scales, they are more prominent in low-mass bins.In the high-mass bin, which is above the half-power mass of all three altDM models, the SFHs converge to the CDM results at  10.

Figure 9 .
Figure9.Fraction of stellar mass formed before  = 10,  M ( > 10), versus the mass-weighted age of the stellar population ( x is defined as the lookback time when x% of the stellar mass is formed) for galaxies at  = 8.Different DM or reionization models are distinguished by different colors, while haloes in different mass bins are distinguished by different markers (as labelled).In terms of the shaded band, the upper boundary represents the relationship obtained by assuming an exponentially rising SFH, while the lower boundary represents the results assuming a constant SFH.In general, points from RT simulations are above the lower boundary of this band, which indicates rising SFH.In relatively low-mass haloes, the stellar populations are younger and more dominated by the recent burst of star formation.In relatively massive haloes, the differences in stellar ages are smaller and the SFHs are close to or even steeper than exponentially rising.Conversely, the uniform UVB model (activated at  ∼ 10) creates an artificial bump in the SFH at  10 and a population of old stars.

Figure 10 .
Figure10.Top: The specific star formation rate (sSFR) versus halo mass at  = 6 (solid lines) and  = 10 (dotted lines).In altDM models, the sSFRs are higher for halo masses below the half-power mass, ∼ 10 9 M , suggesting younger stellar populations.This enhancement is already noticeable at  = 10.Middle: The gas depletion time  depl.versus halo mass.The depletion time is defined as  gas /SFR.In altDM models, the depletion times are systematically lower than in CDM below the half-power mass.This indicates that a larger fraction of gas in the halo is in a phase suitable for star formation.However, this trend is less evident at  = 10.Bottom: The gas mass ratio  gas versus halo mass. gas is normalized over the universal baryon fraction  b ∼ 0.16.Gas is more abundant in altDM models likely due to more gas-rich major mergers.The uniform UVB model results in extremely small gas fractions in galaxies with  halo < 10 8.5 M , leading to the rapid decline of depletion time observed in the middle panel.

Figure 11 .
Figure 11.Left: Galaxy stellar metallicity versus stellar mass at  = 6.The median metallicities are shown as solid lines.The 1 dispersion is shown as the shaded region, for CDM only for clarity.-1result is shown as the gray dotted line.The faint blue dash-dotted lines indicate cosmic variance.In addition, we show the early JWST measurements of  10 galaxy candidates using SED fitting(Furtak et al. 2023) and the measurements of the Local Group dwarfs compiled inKirby et al. (2013).Right: Galaxy gas metallicity versus stellar mass at  = 6.The plotting style and label choices are the same as the left panel.For comparison, we show the recent JWST results fromFujimoto et al. (2023) as well as the local dwarf constraints fromLee et al. (2006).The gas and stellar-phase metallicities resemble each other quite well.At  * 10 7 M , all the models give converged MZR.At  * 10 7 M , galaxies in the sDAO model are less metal enriched by about 0.2 dex at mass scale  * ∼ 10 6−7 M while the ones in the uniform UVB model are more metal enriched by about the same level.Signatures are weaker for the WDM and FDM models since the recent burst of star formation can enrich galaxies to the same level as CDM or even more.The simulation results are in good agreement with the latest observational constraints of high redshift galaxies with JWST as well as the local dwarfs, within the scatter of the observed MZR.
Figure12.Top: Surface density of SFR versus the half-SFR radius for galaxies at  = 8 stacked in three halo mass bins.The half-SFR radius is a proxy for the half-light radius in rest-frame UV which traces recent star formation.We combine the SFR calculated using the instantaneous SFR in gas cells and the one using masses of young stars with age 100 Myr.The FDM model results in smaller galaxy sizes and enhanced SFR surface density.The signature is less prominent in the WDM and sDAO models.Bottom: Fraction of stellar mass with metallicity below 10 −6 versus the metallicity gradient of the same groups of galaxies in the top panel.The metallicity gradient is represented by the difference in the mean metallicity measured within and outside  * 1/2 .Galaxies in the WDM and FDM models have larger metallicity gradients.In low-mass haloes, the fractions of low-metallicity stars are appreciably smaller in WDM and FDM.These findings are consistent with the late-time starbursts, which predominantly occur in the central regions of galaxies.

Figure A1 .
Figure A1.Top: DM surface density in the CDM-uniform UVB run at  = 6.The DM distribution on large scales is nearly identical to the fiducial CDM run using physics.Bottom: Neutral hydrogen column density and gas temperature distribution at  = 6.The IGM is almost completely ionized after the uniform UVB is activated, with the exception of some dense, selfshielded gas in the most massive halo and the surrounding filaments.The IGM temperature has been elevated to 10 4 K by the strong ionizing radiation bath.
, the normalization of the conditional halo mass function at (, ) scales linearly (when  is much smaller than the total mass in the volume) with the mass overdensity of the volume as (, )/ n(, ) = 1 + (, )  m , where (, ) is the bias function and n is the mean halo mass function.Therefore, the variance of the mass function log (, ) equals the variance of the overdensity field multiplied by (, ).The ratio between the total variance of the overdensity field over what we measure from subsampling  = √︄  2 ( = 4 Mpc/ℎ)  2 ( = 4 Mpc/ℎ) −  2 ( = 8 Mpc/ℎ) not vary in altDM models (since 4 Mpc/ℎ is much above the suppression scale).The variance of the halo mass function linearly scales with the variance of the overdensity field and therefore can be estimated as the variance calculated from subsampling times this inflation factor,√︃Var sample [log (, )] × .

Table 1 . Simulations of the -suite. Each column corresponds to the following information: (
1) Name of the simulation.