In this work we investigate the properties of the sources that reionized the intergalactic medium (IGM) in the high-redshift Universe. Using a semi-analytical model aimed at reproducing galaxies and black holes in the first 1.5 Gyr of the Universe, we revisit the relative role of star formation and black hole accretion in producing ionizing photons that can escape into the IGM. Both star formation and black hole accretion are regulated by supernova feedback, resulting in black hole accretion being stunted in low-mass halos. We explore a wide range of combinations for the escape fraction of ionizing photons (redshift-dependent, constant and scaling with stellar mass) from both star formation ($\langle f_{\rm esc}^{\rm sf} \rangle$) and AGN ($f_{\rm esc}^{\rm bh}$) to find: (i) the ionizing budget is dominated by stellar radiation from low stellar mass ($M_*<10^9 {\rm M_\odot}$ ) galaxies at $z>6$ with the AGN contribution (driven by $M_{bh}>10^6 {\rm M_\odot}$ black holes in $M_* > 10^9 {\rm M_\odot}$ galaxies) dominating at lower redshifts; (ii) AGN only contribute $10-25\%$ to the cumulative ionizing emissivity by $z=4$ for the models that match the observed reionization constraints; (iii) if the stellar mass dependence of $\langle f_{\rm esc}^{\rm sf} \rangle$ is shallower than $f_{\rm esc}^{\rm bh}$, at $z<7$ a transition stellar mass exists above which AGN dominate the escaping ionizing photon production rate; (iv) the transition stellar mass decreases with decreasing redshift. While AGN dominate the escaping emissivity above the knee of the stellar mass function at $z \sim 6.8$, they take-over at stellar masses that are a tenth of the knee mass by $z=4$.

We explore the connection between the stellar component of galaxies and their host halos during the epoch of reionization ($5 \leq z\leq10$) using the CROC (Cosmic Reionization on Computers) simulations. We compare simulated galaxies with observations and find that CROC underpredicts the abundance of luminous galaxies when compared to observed UV luminosity functions, and analogously the most massive galaxies when compared to observed stellar mass functions. We can trace the deficit of star formation to high redshifts, where the slope of the star formation rate to stellar mass relation is consistent with observations, but the normalization is systematically low. This results in a star formation rate density and stellar mass density that is systematically offset from observations. However, the less luminous or lower stellar mass objects have luminosities and stellar masses that agree fairly well with observational data. We explore the stellar-to-halo mass ratio, a key quantity that is difficult to measure at high redshifts and that models do not consistently predict. In CROC, the stellar-to-halo mass ratio {\it decreases} with redshift, a trend opposite to some abundance matching studies. These discrepancies uncover where future effort should be focused in order to improve the fidelity of modeling cosmic reionization. We also compare the CROC galaxy bias with observational measurements using Lyman-Break Galaxy (LBG) samples. The good agreement of simulation and data shows that the clustering of dark matter halos is properly captured in CROC.

We examine the reionization history of present-day galaxies by explicitly tracing the building blocks of halos from the Cosmic Reionization On Computers project. We track dark matter particles that belong to $z=0$ halos to trace the neutral fractions at corresponding positions during rapid global reionization. The resulting particle reionization histories allow us to explore different definitions of a halo's reionization redshift and to account for the neutral content of the interstellar medium. Consistent with previous work, we find a systematic trend of reionization redshift with mass - present day halos with higher masses have earlier reionization times. Finally, we quantify the spread of reionization times within each halo, which also has a mass dependence.

A few massive ($M_{*} > 10^8 M_{\odot}$), high-redshift ($z = 8-10$) galaxies have recently been discovered to contain stars with ages of several hundred million years, pushing the onset of star formation in these galaxies back to $z\sim15$. The very existence of stars formed so early may serve as a test for cosmological models with little small scale power (and, hence, late formation of cosmic structure). We explore the ages of oldest stars in numerical simulations from the Cosmic Reionization On Computers (CROC) project with Cold Dark Matter (CDM) and two Warm Dark Matter (WDM) cosmologies with 3 keV and 6 keV particles. There are statistically significant differences between average stellar ages of massive galaxies in CDM and 3 keV WDM. However, these differences are much smaller than both the quoted uncertainties in observational data on the ages of galaxies at these redshifts, and the systematic uncertainties in simulation predictions of these ages as assessed by a convergence test. Further theoretical progress will be needed to refine simulation predictions to an accuracy that would enable dark matter particle physics constraints from this probe.

I describe a modification to the original Fast Multipole Method (FMM) of Greengard & Rokhlin that approximates the gravitation field of an FMM cell as a small uniform grid (a "gridlet") of effective masses. The effective masses on a gridlet are set from the requirement that the multipole moments of the FMM cells are reproduced exactly, hence preserving the accuracy of the gravitational field representation. The calculation of the gravitational field from a multipole expansion can then be computed for all multipole orders simultaneously, with a single Fast Fourier Transform, significantly reducing the computational cost at a given value of the required accuracy. The described approach belongs to the class of "kernel independent" variants of the FMM method and works with any Green function.

We investigate Lyman-$\alpha$ transmission spikes at $z > 5$ in synthetic quasar spectra, and discuss their connection to the properties of the intergalactic medium and their ability to constrain reionization models. We use state-of-the-art radiation-hydrodynamic simulations from the Cosmic Reionization On Computers series to predict the number of transmission spikes as a function of redshift, both in the ideal case of infinite spectral resolution and in a realistic observational setting. Transmission spikes are produced in highly-ionized underdense regions located in the vicinity of UV sources. We find that most of the predicted spikes are unresolved by current observations, and show that our mock spectra are consistent with observations of the quasar ULAS J1120+0641 in about 15% of the realizations. The spike height correlates with both the gas density and the ionized fraction, but the former link is erased when synthetic spectra are smoothed to realistically achievable spectral resolutions. There exists a linear relationship between spike width and the extent of the associated underdense region, with a slope that is redshift-dependent. In agreement with observations, the spike transmitted flux is suppressed at small distance from bright galaxies as these reside in overdense regions. We argue that this anti-correlation can be used to constrain large-scale density modes.

The epoch of reionization, when photons from early galaxies ionized the intergalactic medium about a billion years after the Big Bang, is the last major phase transition in the Universe's history. Measuring the characteristics of the transition is important for understanding early galaxies and the cosmic web and for modeling dwarf galaxies in the later Universe. But such measurements require probes of the intergalactic medium itself. Here we describe how the 21-cm line of neutral hydrogen provides a powerful probe of the reionization process and therefore important constraints on both the galaxies and intergalactic absorbers at that time. While existing experiments will make precise statistical measurements over the next decade, we argue that improved 21-cm analysis techniques - allowing imaging of the neutral gas itself - as well as improved theoretical models, are crucial for testing our understanding of this important era.

The early phases of galaxy formation constitute one of the most exciting frontiers in astrophysics. It is during this era that the first luminous sources reionize the intergalactic medium - the moment when structure formation affects every baryon in the Universe. Here we argue that we will obtain a complete picture of this era by combining observations of galaxies with direct measurements of the reionization process: the former will provide a detailed understanding of bright sources, while the latter will constrain the (substantial) faint source population. We further describe how optimizing the comparison of these two measurements requires near-infrared galaxy surveys covering large volumes and retaining redshift information and also improvements in 21-cm analysis, moving those experiments into the imaging regime.

In addition to being a probe of Cosmic Dawn and Epoch of Reionization astrophysics, the 21cm line at $z>6$ is also a powerful way to constrain cosmology. Its power derives from several unique capabilities. First, the 21cm line is sensitive to energy injections into the intergalactic medium at high redshifts. It also increases the number of measurable modes compared to existing cosmological probes by orders of magnitude. Many of these modes are on smaller scales than are accessible via the CMB, and moreover have the advantage of being firmly in the linear regime (making them easy to model theoretically). Finally, the 21cm line provides access to redshifts prior to the formation of luminous objects. Together, these features of 21cm cosmology at $z>6$ provide multiple pathways toward precise cosmological constraints. These include the "marginalizing out" of astrophysical effects, the utilization of redshift space distortions, the breaking of CMB degeneracies, the identification of signatures of relative velocities between baryons and dark matter, and the discovery of unexpected signs of physics beyond the $\Lambda$CDM paradigm at high redshifts.

The "cosmic dawn" refers to the period of the Universe's history when stars and black holes first formed and began heating and ionizing hydrogen in the intergalactic medium (IGM). Though exceedingly difficult to detect directly, the first stars and black holes can be constrained indirectly through measurements of the cosmic 21-cm background, which traces the ionization state and temperature of intergalactic hydrogen gas. In this white paper, we focus on the science case for such observations, in particular those targeting redshifts z $\gtrsim$ 10 when the IGM is expected to be mostly neutral. 21-cm observations provide a unique window into this epoch and are thus critical to advancing first star and black hole science in the next decade.

The surface densities of molecular gas, $\Sigma_{\rm H_2}$, and the star formation rate (SFR), $\dot\Sigma_\star$, correlate almost linearly on kiloparsec scales in observed star-forming (non-starburst) galaxies. We explore the origin of the linear slope of this correlation using a suite of isolated $L_\star$ galaxy simulations. We show that in simulations with efficient feedback, the slope of the $\dot\Sigma_\star$-$\Sigma_{\rm H_2}$ relation on kiloparsec scales is insensitive to the slope of the $\dot\rho_\star$-$\rho$ relation assumed at the resolution scale. We also find that the slope on kiloparsec scales depends on the criteria used to identify star-forming gas, with a linear slope arising in simulations that identify star-forming gas using a virial parameter threshold. This behavior can be understood using a simple theoretical model based on conservation of interstellar gas mass as the gas cycles between atomic, molecular, and star-forming states under the influence of feedback and dynamical processes. In particular, we show that the linear slope emerges when feedback efficiently regulates and stirs the evolution of dense, molecular gas. We show that the model also provides insights into the likely origin of the relation between the SFR and molecular gas in real galaxies on different scales.

We study the mass of quasar-hosting dark matter halos at z $\sim$ 6 and further constrain the fraction of dark matter halos hosting an active quasar $f_{on}$ and the quasar beaming angle $i_{\rm max}$ using observations of CII lines in the literature. We make assumptions that (1) more massive halos host brighter quasars, (2) a fraction of the halos host active quasars with a certain beaming angle, (3) cold gas in galaxies has rotational velocity $V_{\rm circ}=\alpha V_{\rm max}$, and that (4) quasars point randomly on the sky. We find that for a choice of specific $\alpha \gtrsim 1$, the most likely solution has $f_{\rm on} < 0.01$, corresponding to a small duty cycle of quasar activity. However, if we marginalize over $\alpha$, for some choices of a prior a second solution with $f_{\rm on}=1$ appears. Overall, our the constraints are highly sensitive to $\alpha$ and hence inconclusive. Stronger constraints on $f_{\rm on}$ can be made if we better understand the dynamics of cold gas in these galaxies.

The implementation of star formation and stellar feedback in cosmological simulations plays a critical role in shaping galaxy properties. In the first paper of the series, we presented a new method to model star formation as a collection of star clusters. In this paper, we improve the algorithm by eliminating accretion gaps, boosting momentum feedback, and introducing a subgrid initial bound fraction, $f_i$, that distinguishes cluster mass from stellar particle mass. We perform a suite of simulations with different star formation efficiency per freefall time $\epsilon_{\rm ff}$ and supernova momentum feedback intensity $f_{\rm boost}$. We find that the star formation history of a Milky Way-sized galaxy is sensitive to $f_{\rm boost}$, which allows us to constrain its value, $f_{\rm boost}\approx5$, in the current simulation setup. Changing $\epsilon_{\rm ff}$ from a few percent to 200\% has little effect on global galaxy properties. However, on smaller scales, the properties of star clusters are very sensitive to $\epsilon_{\rm ff}$. We find that $f_i$ increases with $\epsilon_{\rm ff}$ and cluster mass. Through the dependence on $f_i$, the shape of the cluster initial mass function varies strongly with $\epsilon_{\rm ff}$. The fraction of clustered star formation and maximum cluster mass increase with the star formation rate surface density, with the normalization of both relations dependent on $\epsilon_{\rm ff}$. The cluster formation timescale systematically decreases with increasing $\epsilon_{\rm ff}$. Local variations in the gas accretion history lead to a 0.25~dex scatter for the integral cluster formation efficiency. Joint constraints from all the observables prefer the runs that produce a median integral efficiency of 16\%.

Using a suite of isolated $L_\star$ galaxy simulations, we show that global depletion times and star-forming gas mass fractions in simulated galaxies exhibit systematic and well-defined trends as a function of the local star formation efficiency per freefall time, $\epsilon_{\rm ff}$, strength of stellar feedback, and star formation threshold. We demonstrate that these trends can be reproduced and explained by a simple physical model of global star formation in galaxies. Our model is based on mass conservation and the idea of gas cycling between star-forming and non-star-forming states on certain characteristic time scales under the influence of dynamical and feedback processes. Both the simulation results and our model predictions exhibit two limiting regimes with rather different dependencies of global galactic properties on the local parameters. When $\epsilon_{\rm ff}$ is small and feedback is inefficient, the total star-forming mass fraction, $f_{\rm sf}$, is independent of $\epsilon_{\rm ff}$ and the global depletion time, $\tau_{\rm dep}$, scales inversely with $\epsilon_{\rm ff}$. When $\epsilon_{\rm ff}$ is large or feedback is very efficient, these trends are reversed: $f_{\rm sf} \propto \epsilon_{\rm ff}^{-1}$ and $\tau_{\rm dep}$ is independent of $\epsilon_{\rm ff}$ but scales linearly with the feedback strength. We also compare our results with the observed depletion times and mass fractions of star-forming and molecular gas and show that they provide complementary constraints on $\epsilon_{\rm ff}$ and the feedback strength. We show that useful constraints on $\epsilon_{\rm ff}$ can also be obtained using measurements of the depletion time and its scatter on different spatial scales.

A deep understanding of the Epoch of Reionization is still missing in our knowledge of the universe. While future probes will allow us to test the precise evolution of the free electron fraction from redshifts between $z\simeq 6$ and $z\simeq 20$, at present one could ask what kind of reionization processes are allowed by present Cosmic Microwave Background temperature and polarization measurements. An early contribution to reionization could imply a departure from the standard picture where star formation determines the reionization onset. BBy considering a broad class of possible reionization parameterizations, we find that current data do not require an early reionization component in our universe and that only one marginal class of models, based on a particular realization of reionization, may point to that. In addition, the frequentist Akaike Information Criterion (AIC) provides strong evidence against alternative reionization histories, favoring the most simple reionization scenario, which describes reionization by means of only one (constant) reionization optical depth $\tau$.

Using an isolated Milky Way-mass galaxy simulation, we compare results from 9 state-of-the-art gravito-hydrodynamics codes widely used in the numerical community. We utilize the infrastructure we have built for the AGORA High-resolution Galaxy Simulations Comparison Project. This includes the common disk initial conditions, common physics models (e.g., radiative cooling and UV background by the standardized package Grackle) and common analysis toolkit yt, all of which are publicly available. Subgrid physics models such as Jeans pressure floor, star formation, supernova feedback energy, and metal production are carefully constrained across code platforms. With numerical accuracy that resolves the disk scale height, we find that the codes overall agree well with one another in many dimensions including: gas and stellar surface densities, rotation curves, velocity dispersions, density and temperature distribution functions, disk vertical heights, stellar clumps, star formation rates, and Kennicutt-Schmidt relations. Quantities such as velocity dispersions are very robust (agreement within a few tens of percent at all radii) while measures like newly-formed stellar clump mass functions show more significant variation (difference by up to a factor of ~3). Systematic differences exist, for example, between mesh-based and particle-based codes in the low density region, and between more diffusive and less diffusive schemes in the high density tail of the density distribution. Yet intrinsic code differences are generally small compared to the variations in numerical implementations of the common subgrid physics such as supernova feedback. Our experiment reassures that, if adequately designed in accordance with our proposed common parameters, results of a modern high-resolution galaxy formation simulation are more sensitive to input physics than to intrinsic differences in numerical schemes.

An optimally efficient explicit numerical scheme for solving fluid dynamics equations, or any other parabolic or hyperbolic system of partial differential equations, should allow local regions to advance in time with their own, locally constrained time steps. However, such a scheme can result in violation of the Courant-Friedrichs-Lewy (CFL) condition, which is manifestly non-local. Although the violations can be considered to be "weak" in a certain sense and the corresponding numerical solution may be stable, such calculation does not guarantee the correct propagation speed for arbitrary waves. We use an experimental fluid dynamics code that allows cubic "patches" of grid cells to step with independent, locally constrained time steps to demonstrate how the CFL condition can be enforced by imposing a condition on the time steps of neighboring patches. We perform several numerical tests that illustrate errors introduced in the numerical solutions by weak CFL condition violations and show how strict enforcement of the CFL condition eliminates these errors. In all our tests the strict enforcement of the CFL condition does not impose a significant performance penalty.

In models with dark matter made of particles with keV masses, such as a sterile neutrino, small-scale density perturbations are suppressed, delaying the period at which the lowest mass galaxies are formed and therefore shifting the reionization processes to later epochs. In this study, focusing on Warm Dark Matter (WDM) with masses close to its present lower bound, i.e. around the $3$ keV region, we derive constraints from galaxy luminosy functions, the ionization history and the Gunn-Peterson effect. We show that even if star formation efficiency in the simulations is adjusted to match the observed UV galaxy luminosity functions in both CDM and WDM models, the full distribution of Gunn-Peterson optical depth retains the strong signature of delayed reionization in the WDM model. However, until the star formation and stellar feedback model used in modern galaxy formation simulations is constrained better, any conclusions on the nature of dark matter derived from reionization observables remain model-dependent.

We present a model that explains why galaxies form stars on a time scale significantly longer than the time scales of processes governing the evolution of interstellar gas. We show that gas evolves from a non-star-forming to a star-forming state on a relatively short time scale and thus the rate of this evolution does not limit the star formation rate. Instead, the star formation rate is limited because only a small fraction of star-forming gas is converted into stars before star-forming regions are dispersed by feedback and dynamical processes. Thus, gas cycles into and out of star-forming state multiple times, which results in a long time scale on which galaxies convert gas into stars. Our model does not rely on the assumption of equilibrium and can be used to interpret trends of depletion times with the properties of observed galaxies and the parameters of star formation and feedback recipes in simulations. In particular, the model explains how feedback self-regulates the star formation rate in simulations and makes it insensitive to the local star formation efficiency. We illustrate our model using the results of an isolated $L_*$-sized galaxy simulation that reproduces the observed Kennicutt-Schmidt relation for both molecular and atomic gas. Interestingly, the relation for molecular gas is almost linear on kiloparsec scales, although a nonlinear relation is adopted in simulation cells. We discuss how a linear relation emerges from non-self-similar scaling of the gas density PDF with the average gas surface density.

Baryonic effects are amongst the most severe systematics to the tomographic analysis of weak lensing data which is the principal probe in many future generations of cosmological surveys like LSST, Euclid etc.. Modeling or parameterizing these effects is essential in order to extract valuable constraints on cosmological parameters. In a recent paper, Eifler et al. (2015) suggested a reduction technique for baryonic effects by conducting a principal component analysis (PCA) and removing the largest baryonic eigenmodes from the data. In this article, we conducted the investigation further and addressed two critical aspects. Firstly, we performed the analysis by separating the simulations into training and test sets, computing a minimal set of principle components from the training set and examining the fits on the test set. We found that using only four parameters, corresponding to the four largest eigenmodes of the training set, the test sets can be fitted thoroughly with an RMS $\sim 0.0011$. Secondly, we explored the significance of outliers, the most exotic/extreme baryonic scenarios, in this method. We found that excluding the outliers from the training set results in a relatively bad fit and degraded the RMS by nearly a factor of 3. Therefore, for a direct employment of this method to the tomographic analysis of the weak lensing data, the principle components should be derived from a training set that comprises adequately exotic but reasonable models such that the reality is included inside the parameter domain sampled by the training set. The baryonic effects can be parameterized as the coefficients of these principle components and should be marginalized over the cosmological parameter space.

We present a comparison between several observational tests of the post-reionization IGM and the numerical simulations of reionization completed under the Cosmic Reionization On Computers (CROC) project. The CROC simulations match the gap distribution reasonably well, and also provide a good match for the distribution of peak heights, but there is a notable lack of wide peaks in the simulated spectra and the flux PDFs are poorly matched in the narrow redshift interval 5.5<z<5.7, with the match at other redshifts being significantly better, albeit not exact. Both discrepancies are related: simulations show more opacity than the data.

Using four different suites of cosmological simulations, we generate synthetic spectra for galaxies with different Lyman continuum escape fractions (fesc) at redshifts z=7-9, in the rest-frame wavelength range relevant for the James Webb Space Telescope (JWST) NIRSpec instrument. By investigating the effects of realistic star formation histories and metallicity distributions on the EW(Hb)-beta diagram (previously proposed as a tool for identifying galaxies with very high fesc), we find that neither of these effects are likely to jeopardize the identification of galaxies with extreme Lyman continuum leakage. Based on our models, we expect essentially all z=7-9 galaxies that exhibit rest-frame EW(Hb)< 30 {\AA} to have fesc>0.5. Incorrect assumptions concerning the ionizing fluxes of stellar populations or the dust properties of z>6 galaxies can in principle bias the selection, but substantial model deficiencies of this type will at the same time reveal themselves as an offset between the observed and simulated distribution of z>6 galaxies in the EW(Hb)-beta diagram. Such offsets would thereby allow JWST/NIRSpec measurements of these observables to serve as input for further model refinement.

We present a new implementation of star formation in cosmological simulations, by considering star clusters as a unit of star formation. Cluster particles grow in mass over several million years at the rate determined by local gas properties, with high time resolution. The particle growth is terminated by its own energy and momentum feedback on the interstellar medium. We test this implementation for Milky Way-sized galaxies at high redshift, by comparing the properties of model clusters with observations of young star clusters. We find that the cluster initial mass function is best described by a Schechter function rather than a single power law. In agreement with observations, at low masses the logarithmic slope is $\alpha\approx 1.8-2$, while the cutoff at high mass scales with the star formation rate. A related trend is a positive correlation between the surface density of star formation rate and fraction of stars contained in massive clusters. Both trends indicate that the formation of massive star clusters is preferred during bursts of star formation. These bursts are often associated with major merger events. We also find that the median timescale for cluster formation ranges from 0.5 to 4 Myr and decreases systematically with increasing star formation efficiency. Local variations in the gas density and cluster accretion rate naturally lead to the scatter of the overall formation efficiency by an order of magnitude, even when the instantaneous efficiency is kept constant. Comparison of the formation timescale with the observed age spread of young star clusters provides an additional important constraint on the modeling of star formation and feedback schemes.

The effect of line overlap in the Lyman and Werner bands, often ignored in galactic studies of the atomic-to-molecular transition, greatly enhances molecular hydrogen self-shielding in low metallicity environments, and dominates over dust shielding for metallicities below about 10% solar. We implement that effect in cosmological hydrodynamics simulations with an empirical model, calibrated against the observational data, and provide fitting formulae for the molecular hydrogen fraction as a function of gas density on various spatial scales and in environments with varied dust abundance and interstellar radiation field. We find that line overlap, while important for detailed radiative transfer in the Lyman and Werner bands, has only a minor effect on star formation on galactic scales, which, to a much larger degree, is regulated by stellar feedback.

I show that the Reduced Speed of Light (RSL) approximation, when used properly (i.e. as originally designed - only for the local sources but not for the cosmic background), remains a highly accurate numerical method for modeling cosmic reionization. Simulated ionization and star formation histories from the "Cosmic Reionization On Computers" (CROC) project are insensitive to the adopted value of the reduced speed of light for as long as that value does not fall below about 10% of the true speed of light. A recent claim of the failure of the RSL approximation in the Illustris reionization model appears to be due to the effective speed of light being reduced in the equation for the cosmic background too, and, hence, illustrates the importance of maintaining the correct speed of light in modeling the cosmic background.

We explore the effects of specific assumptions in the subgrid models of star formation and stellar and AGN feedback on intrinsic alignments of galaxies in cosmological simulations of "MassiveBlack-II" family. Using smaller volume simulations, we explored the parameter space of the subgrid star formation and feedback model and found remarkable robustness of the observable statistical measures to the details of subgrid physics. The one observational probe most sensitive to modeling details is the distribution of misalignment angles. We hypothesize that the amount of angular momentum carried away by the galactic wind is the primary physical quantity that controls the orientation of the stellar distribution. Our results are also consistent with a similar study by the EAGLE simulation team.

In this paper, the first in a series on galaxy formation before reionization, we focus on understanding what determines the size and morphology of stellar objects in the first low mass galaxies, using parsec- scale cosmological simulations performed with an adaptive mesh hydrodynamics code. Although the dense gas in which stars are formed tends to have a disk structure, stars are found in spheroids with little rotation. Halos with masses between 10^6 M_sun and 5x10^8 M_sun form stars stochastically, with stellar masses in the range 10^4 M_sun to 2x10^6 M_sun. Nearly independent of stellar mass, we observe a large range of half-light radii for the stars, from a few parsecs to a few hundred parsecs and surface brightnesses and mass-to-light ratios ranging from those typical of globular clusters to ultra-faint dwarfs. In our simulations, stars form in dense stellar clusters with high gas-to-star conversion efficiencies and rather uniform metallicities. A fraction of these clusters remain bound after the gas is removed by feedback, but others are destroyed, and their stars, which typically have velocity dispersions of 20 to 40 km/s, expand until they become bound by the dark matter halo. We thus speculate that the stars in ultra-faint dwarf galaxies may show kinematic and chemical signatures consistent with their origin in a few distinct stellar clusters. On the other hand, some globular clusters may form at the center of primordial dwarf galaxies and may contain dark matter, perhaps detectable in the outer parts.

Using numerical cosmological simulations completed under the "Cosmic Reionization On Computers" (CROC) project, I explore theoretical predictions for the faint end of the galaxy UV luminosity functions at $z\geq 6$. A commonly used Schechter function approximation with the magnitude cut at $M_{\rm CUT}\sim-13$ provides a reasonable fit to the actual luminosity function of simulated galaxies. When the Schechter functional form is forced on the luminosity functions from the simulations, the magnitude cut $M_{\rm CUT}$ is found to vary between $-12$ and $-14$ with a mild redshift dependence. An analytical model of reionization from Madau, Haardt & Rees (1997), as used by Robertson et al. (2015), provides a good description of the simulated results, but only if the redshift dependence of the effective escape fraction (induced by physical processes not captured by the Madau, Haardt & Rees model) is accounted for.

We present a study of a star formation prescription in which star formation efficiency depends on local gas density and turbulent velocity dispersion, as suggested by direct simulations of SF in turbulent giant molecular clouds (GMCs). We test the model using a simulation of an isolated Milky Way-sized galaxy with a self-consistent treatment of turbulence on unresolved scales. We show that this prescription predicts a wide variation of local star formation efficiency per free-fall time, $\epsilon_{\rm ff} \sim 0.1 - 10\%$, and gas depletion time, $t_{\rm dep} \sim 0.1 - 10$ Gyr. In addition, it predicts an effective density threshold for star formation due to suppression of $\epsilon_{\rm ff}$ in warm diffuse gas stabilized by thermal pressure. We show that the model predicts star formation rates in agreement with observations from the scales of individual star-forming regions to the kiloparsec scales. This agreement is non-trivial, as the model was not tuned in any way and the predicted star formation rates on all scales are determined by the distribution of the GMC-scale densities and turbulent velocities $\sigma$ in the cold gas within the galaxy, which is shaped by galactic dynamics. The broad agreement of the star formation prescription calibrated in the GMC-scale simulations with observations, both gives credence to such simulations and promises to put star formation modeling in galaxy formation simulations on a much firmer theoretical footing.

We explore the mean and fluctuating redshifted 21 cm signal in numerical simulations of cosmic reionization from the Cosmic Reionization On Computers (CROC) project. We find that the mean signal varies between about $\pm20\rm{mK}$. Most significantly, we find that the negative pre-reionization dip at $z\sim10-15$ only extends to $\langle\Delta T_B\rangle\sim-20\rm{mK}$, in agreement with prior simulation results and in significant contrast to Pritchard & Loeb analytical model, requiring substantially higher sensitivity from global signal experiments that operate in this redshift range (EDGES-II, LEDA, SCI-HI, and DARE). We also explore the role of dense substructure (filaments and embedded galaxies) in the formation of 21 cm power spectrum. We find that by neglecting the semi-neutral substructure inside ionized bubbles, the power spectrum can be mis-estimated by 25-50\% at scales $k\sim 0.1-1h\rm{Mpc}^{-1}$. This scale range is of a particular interest, because the upcoming 21 cm experiments (MWA, PAPER, HERA) are expected to be most sensitive within it.

This draft report summarizes and details the findings, results, and recommendations derived from the ASCR/HEP Exascale Requirements Review meeting held in June, 2015. The main conclusions are as follows. 1) Larger, more capable computing and data facilities are needed to support HEP science goals in all three frontiers: Energy, Intensity, and Cosmic. The expected scale of the demand at the 2025 timescale is at least two orders of magnitude -- and in some cases greater -- than that available currently. 2) The growth rate of data produced by simulations is overwhelming the current ability, of both facilities and researchers, to store and analyze it. Additional resources and new techniques for data analysis are urgently needed. 3) Data rates and volumes from HEP experimental facilities are also straining the ability to store and analyze large and complex data volumes. Appropriately configured leadership-class facilities can play a transformational role in enabling scientific discovery from these datasets. 4) A close integration of HPC simulation and data analysis will aid greatly in interpreting results from HEP experiments. Such an integration will minimize data movement and facilitate interdependent workflows. 5) Long-range planning between HEP and ASCR will be required to meet HEP's research needs. To best use ASCR HPC resources the experimental HEP program needs a) an established long-term plan for access to ASCR computational and data resources, b) an ability to map workflows onto HPC resources, c) the ability for ASCR facilities to accommodate workflows run by collaborations that can have thousands of individual members, d) to transition codes to the next-generation HPC platforms that will be available at ASCR facilities, e) to build up and train a workforce capable of developing and using simulations and analysis to support HEP scientific research on next-generation systems.

In this paper I show that simulations of reionization performed under the Cosmic Reionization On Computers (CROC) project do converge in space and mass, albeit rather slowly. A fully converged solution (for a given star formation and feedback model) can be determined at a level of precision of about 20%, but such a solution is useless in practice, since achieving it in production-grade simulations would require a large set of runs at various mass and spatial resolutions, and computational resources for such an undertaking are not yet readily available. In order to make progress in the interim, I introduce a weak convergence correction factor in the star formation recipe, which allows one to approximate the fully converged solution with finite resolution simulations. The accuracy of weakly converged simulations approaches a comparable, ~20% level of precision for star formation histories of individual galactic halos and other galactic properties that are directly related to star formation rates, like stellar masses and metallicities. Yet other properties of model galaxies, for example, their HI masses, are recovered in the weakly converged runs only within a factor of two.

We compare the properties of stellar populations of model galaxies from the Cosmic Reionization On Computers (CROC) project with the exiting UV and IR data. Since CROC simulations do not follow cosmic dust directly, we adopt two variants of the dust-follows-metals ansatz to populate model galaxies with dust. Using the dust radiative transfer code Hyperion, we compute synthetic stellar spectra, UV continuum slopes, and IR fluxes for simulated galaxies. We find that the simulation results generally match observational measurements, but, perhaps, not in full detail. The differences seem to indicate that our adopted dust-follows-metals ansatzes are not fully sufficient. While the discrepancies with the exiting data are marginal, the future JWST data will be of much higher precision, rendering highly significant any tentative difference between theory and observations. It is, therefore, likely, that in order to fully utilize the precision of JWST observations, fully dynamical modeling of dust formation, evolution, and destruction may be required.

We revisit the possibility of constraining the properties of dark matter (DM) by studying the epoch of cosmic reionization. Previous studies have shown that DM annihilation was unlikely to have provided a large fraction of the photons that ionized the universe, but instead played a subdominant role relative to stars and quasars. The DM, however, begins to efficiently annihilate with the formation of primordial microhalos at $z\sim100-200$, much earlier than the formation of the first stars. Therefore, if DM annihilation ionized the universe at even the percent level over the interval $z \sim 20-100$, it can leave a significant imprint on the global optical depth, $\tau$. Moreover, we show that cosmic microwave background (CMB) polarization data and future 21 cm measurements will enable us to more directly probe the DM contribution to the optical depth. In order to compute the annihilation rate throughout the epoch of reionization, we adopt the latest results from structure formation studies and explore the impact of various free parameters on our results. We show that future measurements could make it possible to place constraints on the dark matter's annihilation cross section that are at a level comparable to those obtained from the observations of dwarf galaxies, cosmic ray measurements, and studies of recombination.

We use fully self-consistent numerical simulations of cosmic reionization, completed under the Cosmic Reionization On Computers (CROC) project, to explore how well the recombinations in the ionized IGM can be quantified by the effective "clumping factor." The density distribution in the simulations (and, presumably, in a real universe) is highly inhomogeneous and more-or-less smoothly varying in space. However, even in highly complex and dynamic environments the concept of the IGM remains reasonably well-defined; the largest ambiguity comes from the unvirialized regions around galaxies that are over-ionized by the local enhancement in the radiation field ("proximity zones"). That ambiguity precludes computing the IGM clumping factor to better than about 20\%. We also discuss a "local clumping factor," defined over a particular spatial scale, and quantify its scatter on a given scale and its variation as a function of scale.

High-redshift gamma-ray bursts have several advantages for the study of the distant universe, providing unique information about the structure and properties of the galaxies in which they exploded. Spectroscopic identification with large ground-based telescopes has improved our knowledge of the class of such distant events. We present the multi-wavelength analysis of the high-$z$ Swift gamma-ray burst GRB140515A ($z = 6.327$). The best estimate of the neutral hydrogen fraction of the intergalactic medium (IGM) towards the burst is $x_{HI} \leq 0.002$. The spectral absorption lines detected for this event are the weakest lines ever observed in gamma-ray burst afterglows, suggesting that GRB140515A exploded in a very low density environment. Its circum-burst medium is characterised by an average extinction (A$_{\rm V} \sim 0.1$) that seems to be typical of $z \ge 6$ events. The observed multi-band light curves are explained either with a very flat injected spectrum ($p = 1.7$) or with a multi-component emission ($p = 2.1$). In the second case a long-lasting central engine activity is needed in order to explain the late time X-ray emission. The possible origin of GRB140515A from a Pop III (or from a Pop II stars with local environment enriched by Pop III) massive star is unlikely.

What should these lectures be? The subject assigned to us is so broad that many books can be written about it. So, in planning these lectures I had several options. One would be to focus on a narrow subset of topics and to cover them in great detail. Such a subset necessarily would be highly personal and useful to a few read- ers at best. Another option would be to give a very shallow overview of the whole field, but then it won't be very much different from a highly compressed version of a university course (which anyone can take if they wish so). So, I decided to be selfish and to prepare these lectures as if I was teaching my own graduate student. Given my research interests, I selected what the student would need to know to be able to discuss science with me and to work on joint research projects. So, the story presented below is both personal and incomplete, but it does cover several subjects that are poorly represented in the existing textbooks (if at all). Some of topics I focus on below are closely connected, others are disjoint, some are just side detours on specific technical questions. There is an overlapping theme, however. Our goal is to follow the cosmic gas from large scales, low densities, (rel- atively) simple physics to progressively smaller scales, higher densities, closer rela- tion to galaxies, and more complex and uncertain physics. So, we (you - the reader, and me - the author) are going to follow a "yellow brick road" from the gas well be- yond any galaxy confines to the actual sites of star formation and stellar feedback. On the way we will stop at some places for a tour and run without looking back through some others. So, the road will be uneven, but I hope that some readers find it useful.

We compare the results from several sets of cosmological simulations of cosmic reionization, produced under Cosmic Reionization On Computers (CROC) project, with existing observational data on the high-redshift Ly-alpha forest and the abundance of Ly-alpha emitters. We find good consistency with the observational measurements and the previous simulation work. By virtue of having several independent realizations for each set of numerical parameters, we are able to explore the effect of cosmic variance on observable quantities. One unexpected conclusion we are forced into is that cosmic variance is unusually large at z>6, with both our simulations and, most likely, observational measurements are still not fully converged for even such basic quantities as the average Gunn-Peterson optical depth or the volume-weighted neutral fraction. We also find that reionization has little effect on the early galaxies or on global cosmic star formation history, because galaxies whose gas content is affected by photoionization contain no molecular (i.e. star-forming) gas in the first place. In particular, measurements of the faint end of the galaxy luminosity function by JWST are unlikely to provide a useful constraint on reionization.

We use recent proper motion measurements of the tangential velocity of M31, along with its radial velocity and distance, to derive the likelihood of the sum of halo masses of the Milky Way and M31. This is done using a sample halo pairs in the Bolshoi cosmological simulation of $\Lambda$CDM cosmology selected to match properties and environment of the Local Group. The resulting likelihood gives estimate of the sum of masses of $M_{\rm MW,200}+M_{\rm M31,200}=$ $2.40_{-1.05}^{+1.95}\times10^{12}\,M_{\odot}$ ($90\%$ confidence interval). This estimate is consistent with individual mass estimates for the Milky Way and M31 and is consistent, albeit somewhat on the low side, with the mass estimated using the timing argument. We show that although the timing argument is unbiased on average for all pairs, for pairs constrained to have radial and tangential velocities similar to that of the Local Group the argument overestimates the sum of masses by a factor of $1.6$. Using similar technique we estimate the total dark matter mass enclosed within $1$ Mpc from the Local Group barycenter to be $M_{\rm LG}(r<1\, {\rm Mpc})=4.2_{-2.0}^{+3.4}\times10^{12}\,M_{\odot}$ ($90\%$ confidence interval).

Cosmic Reionization On Computers (CROC) is a long-term program of numerical simulations of cosmic reionization. Its goal is to model fully self-consistently (albeit not necessarily from the first principles) all relevant physics, from radiative transfer to gas dynamics and star formation, in simulation volumes of up to 100 comoving Mpc, and with spatial resolution approaching 100 pc in physical units. In this method paper we describe our numerical method, the design of simulations, and the calibration of numerical parameters. Using several sets (ensembles) of simulations in 20 Mpc/h and 40 Mpc/h boxes with spatial resolution reaching 125 pc at z=6, we are able to match the observed galaxy UV luminosity functions at all redshifts between 6 and 10, as well as obtain reasonable agreement with the observational measurements of the Gunn-Peterson optical depth at z<6.

Systematic uncertainties that have been subdominant in past large-scale structure (LSS) surveys are likely to exceed statistical uncertainties of current and future LSS data sets, potentially limiting the extraction of cosmological information. Here we present a general framework (PCA marginalization) to consistently incorporate systematic effects into a likelihood analysis. This technique naturally accounts for degeneracies between nuisance parameters and can substantially reduce the dimension of the parameter space that needs to be sampled. As a practical application, we apply PCA marginalization to account for baryonic physics as an uncertainty in cosmic shear tomography. Specifically, we use CosmoLike to run simulated likelihood analyses on three independent sets of numerical simulations, each covering a wide range of baryonic scenarios differing in cooling, star formation, and feedback mechanisms. We simulate a Stage III (Dark Energy Survey) and Stage IV (Large Synoptic Survey Telescope/Euclid) survey and find a substantial bias in cosmological constraints if baryonic physics is not accounted for. We then show that PCA marginalization (employing at most 3 to 4 nuisance parameters) removes this bias. Our study demonstrates that it is possible to obtain robust, precise constraints on the dark energy equation of state even in the presence of large levels of systematic uncertainty in astrophysical processes. We conclude that the PCA marginalization technique is a powerful, general tool for addressing many of the challenges facing the precision cosmology program.

We discuss the role of recombinations in the IGM, and the related concept of the clumping factor, during cosmic reionization. The clumping factor is, in general, a local quantity that depends on the local over-density and the scale below which the baryon density field can be assumed smooth. That scale, called the filtering scale, is itself depended on over-density and local thermal history. We present a method for building a self-consistent analytical model of inhomogeneous reionization assuming the linear growth rate of the density fluctuation, which accounts for these effects simultaneously. We show that taking into account the local clumping factor introduces significant corrections to the total recombination rate, comparing to the model with a globally uniform clumping factor.

We use simulations of isolated galaxies with a few parsec resolution to explore the connection between the small-scale star formation rate - gas density relation and the induced large-scale correlation between the star formation rate surface density and the surface density of the molecular gas (the Kennicutt-Schmidt relation). We find that, in the simulations, a power-law small-scale "star formation law" directly translates into an identical power-law Kennicutt-Schmidt relation. If this conclusion holds in the reality as well, it implies that the observed approximately linear Kennicutt-Schmidt relation must reflect the approximately linear small-scale "star formation law".

We introduce the AGORA project, a comprehensive numerical study of well-resolved galaxies within the LCDM cosmology. Cosmological hydrodynamic simulations with force resolutions of ~100 proper pc or better will be run with a variety of code platforms to follow the hierarchical growth, star formation history, morphological transformation, and the cycle of baryons in and out of 8 galaxies with halo masses M_vir ~= 1e10, 1e11, 1e12, and 1e13 Msun at z=0 and two different ("violent" and "quiescent") assembly histories. The numerical techniques and implementations used in this project include the smoothed particle hydrodynamics codes GADGET and GASOLINE, and the adaptive mesh refinement codes ART, ENZO, and RAMSES. The codes will share common initial conditions and common astrophysics packages including UV background, metal-dependent radiative cooling, metal and energy yields of supernovae, and stellar initial mass function. These are described in detail in the present paper. Subgrid star formation and feedback prescriptions will be tuned to provide a realistic interstellar and circumgalactic medium using a non-cosmological disk galaxy simulation. Cosmological runs will be systematically compared with each other using a common analysis toolkit, and validated against observations to verify that the solutions are robust - i.e., that the astrophysical assumptions are responsible for any success, rather than artifacts of particular implementations. The goals of the AGORA project are, broadly speaking, to raise the realism and predictive power of galaxy simulations and the understanding of the feedback processes that regulate galaxy "metabolism." The proof-of-concept dark matter-only test of the formation of a galactic halo with a z=0 mass of M_vir ~= 1.7e11 Msun by 9 different versions of the participating codes is also presented to validate the infrastructure of the project.

This document presents (off-line) computing requrements and challenges for Cosmic Frontier science, covering the areas of data management, analysis, and simulations. We invite contributions to extend the range of covered topics and to enhance the current descriptions.

We present results from high-resolution hydrodynamic simulations of isolated SMC- and Milky Way-sized galaxies that include a model for feedback from galactic cosmic rays (CRs). We find that CRs are naturally able to drive winds with mass loading factors of up to ~10 in dwarf systems. The scaling of the mass loading factor with circular velocity between the two simulated systems is consistent with \propto v_c^{1-2} required to reproduce the faint end of the galaxy luminosity function. In addition, simulations with CR feedback reproduce both the normalization and the slope of the observed trend of wind velocity with galaxy circular velocity. We find that winds in simulations with CR feedback exhibit qualitatively different properties compared to SN driven winds, where most of the acceleration happens violently in situ near star forming sites. In contrast, the CR-driven winds are accelerated gently by the large-scale pressure gradient established by CRs diffusing from the star-forming galaxy disk out into the halo. The CR-driven winds also exhibit much cooler temperatures and, in the SMC-sized system, warm (T~10^4 K) gas dominates the outflow. The prevalence of warm gas in such outflows may provide a clue as to the origin of ubiquitous warm gas in the gaseous halos of galaxies detected via absorption lines in quasar spectra.

We investigate the transition from primordial Population III (Pop III) star formation to normal Pop II star formation in the first galaxies using new cosmological hydrodynamic simulations. We find that while the first stars seed their host galaxies with metals, they cannot sustain significant outflows to enrich the intergalactic medium, even assuming a top-heavy initial mass function. This means that Pop III star formation could potentially continue until z~6 in different unenriched regions of the universe, before being ultimately shut off by cosmic reionization. Within an individual galaxy, the metal production and stellar feedback from Pop II stars overtake Pop III stars in 20-200 Myr, depending on galaxy mass.

We revisit the formation and evolution of the first galaxies using new hydrodynamic cosmological simulations with the ART code. Our simulations feature a recently developed model for H2 formation and dissociation, and a star formation recipe that is based on molecular rather than atomic gas. Here, we develop and implement a recipe for the formation of metal-free Population III stars in galaxy-scale simulations that resolve primordial clouds with sufficiently high density. We base our recipe on the results of prior zoom-in simulations that resolved the protostellar collapse in pre-galactic objects. We find the epoch during which Pop III stars dominated the energy and metal budget of the first galaxies to be short-lived. Galaxies which host Pop III stars do not retain dynamical signatures of their thermal and radiative feedback for more than 10^8 yr after the lives of the stars end in pair-instability supernovae, even when we consider the maximum reasonable efficiency of the feedback. Though metals ejected by the supernovae can travel well beyond the virial radius of the host galaxy, they typically begin to fall back quickly, and do not enrich a large fraction of the intergalactic medium. Galaxies with total mass in excess of 3 x 10^6 Msun re-accrete most of their baryons and transition to metal-enriched Pop II star formation.

Previous studies showed that an estimate of the likelihood distribution of the Milky Way halo mass can be derived using the properties of the satellites similar to the Large and Small Magellanic Clouds (LMC and SMC). However, it would be straightforward to interpret such an estimate only if the properties of the Magellanic Clouds (MCs) are fairly typical and are not biased by the environment. In this study we explore whether the environment of the Milky Way affects the properties of the SMC and LMC such as their velocities. To test for the effect of the environment, we compare velocity distributions for MC-sized subhalos around Milky Way hosts in a sample selected simply by mass and in the second sample of such halos selected with additional restrictions on the distance to the nearest cluster and the local galaxy density, designed to mimic the environment of the Local Group (LG). We find that satellites in halos in the LG-like environments do have somewhat larger velocities, as compared to the halos of similar mass in the sample without environmental constraints. We derive the host halo likelihood distribution for the samples in the LG-like envirionment and in the control sample and find that the environment does not significantly affect the derived likelihood. We use the updated properties of the SMC and LMC to derive the constraint on the MW halo mass $\log{({\rm M}_{200} /\msol)}=12.06^{+0.31}_{-0.19}$ (90% confidence interval). We also explore the incidence of close pairs with relative velocities and separations similar to those of the LMC and SMC and find that such pairs are quite rare among $\Lambda$CDM halos. Taking into account the close separation of the MCs in the Busha et al.\ 2011 method results in the shift of the MW halo mass estimate to smaller masses, with the peak shifting approximately by a factor of two.[Abridged]

We study the high column density regime of the HI column density distribution function and argue that there are two distinct features: a turnover at NHI ~ 10^21 cm^-2 which is present at both z=0 and z ~ 3, and a lack of systems above NHI ~ 10^22 cm^-2 at z=0. Using observations of the column density distribution, we argue that the HI-H2 transition does not cause the turnover at NHI ~ 10^21 cm^-2, but can plausibly explain the turnover at NHI > 10^22 cm^-2. We compute the HI column density distribution of individual galaxies in the THINGS sample and show that the turnover column density depends only weakly on metallicity. Furthermore, we show that the column density distribution of galaxies, corrected for inclination, is insensitive to the resolution of the HI map or to averaging in radial shells. Our results indicate that the similarity of HI column density distributions at z=3 and z=0 is due to the similarity of the maximum HI surface densities of high-z and low-z disks, set presumably by universal processes that shape properties of the gaseous disks of galaxies. Using fully cosmological simulations, we explore other candidate physical mechanisms that could produce a turnover in the column density distribution. We show that while turbulence within GMCs cannot affect the DLA column density distribution, stellar feedback can affect it significantly if the feedback is sufficiently effective in removing gas from the central 2-3 kpc of high-redshift galaxies. Finally, we argue that it is meaningful to compare column densities averaged over ~ kpc scales with those estimated from quasar spectra which probe sub-pc scales due to the steep power spectrum of HI column density fluctuations observed in nearby galaxies.

Understanding the gas content of high redshift halos is crucial for studying the formation of the first generation of galaxies and reionization. Recently, Tseliakhovich & Hirata showed that the relative "stream" velocity between the dark matter and baryons at the time of recombination - formally a second order effect, but an unusually large one - can influence the later structure formation history of the Universe. We quantify the effect of the stream velocity on the so-called "characteristic mass" - the minimum mass of a dark matter halo capable of retaining most of its baryons throughout its formation epoch - using three different high-resolution sets of cosmological simulations (with separate transfer functions for baryons and dark matter) that vary in box size, particle number, and the value of the relative velocity between the dark matter and baryons. In order to understand this effect theoretically, we generalize the linear theory filtering mass to properly account for the difference between the dark matter and baryonic density fluctuation evolution induced by the stream velocity. We show that the new filtering mass provides an accurate estimate for the characteristic mass, while other theoretical ansatzes for the characteristic mass are substantially less precise.

We extend the existing analytical model of reionization by Furlanetto et al. (2004) to include the biasing of reionization sources and additional absorption by Lyman Limit systems. Our model is, by construction, consistent with the observed evolution of the galaxy luminosity function at z<8 and with the observed evolution of Ly-{\alpha} forest at z<6. We also find that, for a wide range of values for the relative escape fraction that we consider reasonable, and which are consistent with the observational constraints on the relative escape fraction from lower redshifts, our reionization model is consistent with the WMAP constraint on the Thompson optical depth and with the SPT and EDGES constraints on the duration of reionization. We, therefore, conclude that it is possible to develop physically realistic models of reionization that are consistent with all existing observational constraints.

We study how uncertainty in the reionization history of the universe affects estimates of other cosmological parameters from the Cosmic Microwave Background. We analyze WMAP7 data and synthetic Planck-quality data generated using a realistic scenario for the reionization history of the universe obtained from high-resolution numerical simulation. We perform parameter estimation using a simple sudden reionization approximation, and using the Principal Component Analysis (PCA) technique proposed by Mortonson and Hu. We reach two main conclusions: (1) Adopting a simple sudden reionization model does not introduce measurable bias into values for other parameters, indicating that detailed modeling of reionization is not necessary for the purpose of parameter estimation from future CMB data sets such as Planck. (2) PCA analysis does not allow accurate reconstruction of the actual reionization history of the universe in a realistic case.

Stellar feedback plays a key role in galaxy formation by regulating star formation, driving interstellar turbulence and generating galactic scale outflows. Although modern simulations of galaxy formation can resolve scales of 10-100 pc, star formation and feedback operate on smaller, "subgrid" scales. Great care should therefore be taken in order to properly account for the effect of feedback on global galaxy evolution. We investigate the momentum and energy budget of feedback during different stages of stellar evolution, and study its impact on the interstellar medium using simulations of local star forming regions and galactic disks at the resolution affordable in modern cosmological zoom-in simulations. In particular, we present a novel subgrid model for the momentum injection due to radiation pressure and stellar winds from massive stars during early, pre-supernova evolutionary stages of young star clusters. Early injection of momentum acts to clear out dense gas in star forming regions, hence limiting star formation. The reduced gas density mitigates radiative losses of thermal feedback energy from subsequent supernova explosions, leading to an increased overall efficiency of stellar feedback. The detailed impact of stellar feedback depends sensitively on the implementation and choice of parameters. Somewhat encouragingly, we find that implementations in which feedback is efficient lead to approximate self-regulation of global star formation efficiency. We compare simulation results using our feedback implementation to other phenomenological feedback methods, where thermal feedback energy is allowed to dissipate over time scales longer than the formal gas cooling time. We find that simulations with maximal momentum injection suppress star formation to a similar degree as is found in simulations adopting adiabatic thermal feedback.

Molecular hydrogen (H2) is the primary component of the reservoirs of cold, dense gas that fuel star formation in our galaxy. While the H2 abundance is ultimately regulated by physical processes operating on small scales in the interstellar medium (ISM), observations have revealed a tight correlation between the ratio of molecular to atomic hydrogen in nearby spiral galaxies and the pressure in the mid-plane of their disks. This empirical relation has been used to predict H2 abundances in galaxies with potentially very different ISM conditions, such as metal-deficient galaxies at high redshifts. Here, we test the validity of this approach by studying the dependence of the pressure -- H2 relation on environmental parameters of the ISM. To this end, we follow the formation and destruction of H2 explicitly in a suite of hydrodynamical simulations of galaxies with different ISM parameters. We find that a pressure -- H2 relation arises naturally in our simulations for a variety of dust-to-gas ratios or strengths of the interstellar radiation field in the ISM. Fixing the dust-to-gas ratio and the UV radiation field to values measured in the solar neighborhood results in fair agreement with the relation observed in nearby galaxies with roughly solar metallicity. However, the parameters (slope and normalization) of the pressure -- H2 relation vary in a systematical way with ISM properties. A particularly strong trend is the decrease of the normalization of the relation with a lowering of the dust-to-gas ratio of the ISM. We show that this trend and other properties of the pressure -- H2 relation are natural consequences of the transition from atomic to molecular hydrogen with gas surface density.

Cooling and heating functions of cosmic gas are a crucial ingredient for any study of gas dynamics and thermodynamics in the interstellar and intergalactic medium. As such, they have been studied extensively in the past under the assumption of collisional ionization equilibrium. However, for a wide range of applications, the local radiation field introduces a non-negligible, often dominant, modification to the cooling and heating functions. In the most general case, these modifications cannot be described in simple terms, and would require a detailed calculation with a large set of chemical species using a radiative transfer code (the well-known code Cloudy, for example). We show, however, that for a sufficiently general variation in the spectral shape and intensity of the incident radiation field, the cooling and heating functions can be approximated as depending only on several photoionization rates, which can be thought of as representative samples of the overall radiation field. This dependence is easy to tabulate and implement in cosmological or galactic-scale simulations, thus economically accounting for an important but rarely-included factor in the evolution of cosmic gas. We also show a few examples where the radiation environment has a large effect, the most spectacular of which is a quasar that suppresses gas cooling in its host halo without any mechanical or non-radiative thermal feedback.

The baryonic Tully-Fisher relation is an important observational constraint on cosmological and galactic models. However, it is critical to keep in mind that in observations only stars, molecular, and atomic gas are counted, while the contribution of the ionized gas is almost universally missed. The ionized gas is, however, expected to be present in the gaseous disks of dwarf galaxies simply because they are exposed to the cosmic ionizing background and to the stellar radiation that manages to escape from the central regions of the galactic disks into their outer layers. Such an expectation is, indeed, born out both by cosmological numerical simulations and by simple analytical models.

Interactions of cosmic rays with the interstellar gas and radiation fields of the Milky Way provide the majority of the gamma rays observed by the Fermi Gamma Ray Space Telescope. In addition to the gas which is densely concentrated along the Galactic Disk, hydrodynamical simulations and observational evidence favor the presence of a halo of hot (T~10^6 K) ionized hydrogen (H_II), extending with non-negligible densities out to the virial radius of the Milky Way. We show that cosmic ray collisions with this circum-galactic gas should be expected to provide a significant flux of gamma rays, on the order of 10% of the observed isotopic gamma ray background at energies above 1 GeV. In addition, gamma rays originating from the extended H_II halos of other galaxies along a given line-of-sight should contribute to this background at a similar level.

There is ample observational evidence that the star formation rate (SFR) surface density, Sigma_SFR, is closely correlated with the surface density of molecular hydrogen, Sigma_H2. This empirical relation holds both for galaxy-wide averages and for individual >=kpc sized patches of the interstellar medium (ISM), but appears to degrade substantially at a sub-kpc scale. Identifying the physical mechanisms that determine the scale-dependent properties of the observed Sigma_H2-Sigma_SFR relation remains a challenge from a theoretical perspective. To address this question, we analyze the slope and scatter of the Sigma_H2-Sigma_SFR relation using a set of cosmological, galaxy formation simulations with a peak resolution of ~100 pc. These simulations include a chemical network for molecular hydrogen, a model for the CO emission, and a simple, stochastic prescription for star formation that operates on ~100 pc scales. Specifically, star formation is modeled as a Poisson process in which the average SFR is directly proportional to the present mass of H2. The predictions of our numerical model are in good agreement with the observed Kennicutt-Schmidt and Sigma_H2-Sigma_SFR relations. We show that observations based on CO emission are ill suited to reliably measure the slope of the latter relation at low (<20 M_sun pc^-2) H2 surface densities on sub-kpc scales. Our models also predict that the inferred Sigma_H2-Sigma_SFR relation steepens at high H2 surface densities as a result of the surface density dependence of the CO/H2 conversion factor. Finally, we show that on sub-kpc scales most of the scatter in the relation is a consequence of discreteness effects in the star formation process. In contrast, variations of the CO/H2 conversion factor are responsible for most of the scatter measured on super-kpc scales.

It has been recently shown that the relative velocity between the dark matter and the baryons (vbc) at the time of recombination can affect the structure formation in the early universe (Tseliakhovich & Hirata 2010). We statistically quantify this effect using large cosmological simulations. We use three different high resolution sets of simulations (with separate transfer functions for baryons and dark matter) that vary in box size, particle number, and the value of the relative velocity between dark matter and baryons. We show that the total number density of halos is suppressed by ~ 20% at z = 25 for vbc = 1{\sigma}(vbc), where {\sigma}(vbc) is the variance of the relative velocity, while for vbc = 3.4{\sigma}(vbc) the relative suppression at the same redshift reaches 50%, remaining at or above the 30% level all the way to z = 11. We also find high abundance of "empty halos", i.e., halos that have gas fraction below half of the cosmic mean baryonic fraction fb. Specifically we find that for vbc = 1{\sigma}(vbc) all halos below 10^5M\odot are empty at z \geq 19. The high abundance of empty halos results in significant delay in the formation of gas rich mini-halos and the first galaxies.

Characterizing the conversion factor between CO emission and column density of molecular hydrogen, X_CO, is crucial in studying the gaseous content of galaxies, its evolution, and relation to star formation. In most cases the conversion factor is assumed to be close to that of giant molecular clouds (GMCs) in the Milky Way, except possibly for mergers and star-bursting galaxies. However, there are physical grounds to expect that it should also depend on the gas metallicity, surface density, and strength of the interstellar radiation field. The XCO factor may also depend on the scale on which CO emission is averaged due to effects of limited resolution. We study the dependence of X_CO on gas properties and averaging scale using a model that is based on a combination of results of sub-pc scale magneto-hydrodynamic simulations and on the gas distribution from self-consistent cosmological simulations of galaxy formation. Our model predicts a value of X_CO that is consistent with the Galactic value for interstellar medium conditions typical for the Milky Way. For such conditions the predicted X_CO varies by only a factor of two for gas surfaced densities in the range \sim 50 - 500 M_sun / pc^2. However, the model also predicts that more generally on the scale of GMCs, X_CO is a strong function of metallicity, and depends on the column density and the interstellar UV flux. We show explicitly that neglecting these dependencies in observational estimates can strongly bias the inferred distribution of H2 column densities of molecular clouds to have a narrower and offset range compared to the true distribution. We find that when averaged on \sim kpc scales the X-factor depends only weakly on radiation field and column density, but is still a strong function of metallicity. The predicted metallicity dependence can be approximated as X_CO \sim Z^{-{\gamma}} with {\gamma} ~ 0.5 - 0.8.

The South Pole Telescope (SPT), Atacama Cosmology Telescope (ACT), and Wilkinson Microwave Anisotropy Probe (WMAP) have each reported measurements of the cosmic microwave background's (CMB) angular power spectrum which favor the existence of roughly one additional neutrino species, in addition to the three contained in the standard model of particle physics. Neutrinos influence the CMB by contributing to the radiation density, which alters the expansion rate of the universe during the epoch leading up to recombination. In this paper, we consider an alternative possibility that the excess kinetic energy implied by these measurements was possessed by dark matter particles that were produced through a non-thermal mechanism, such as late-time decays. In particular, we find that if a small fraction (<1%) of the dark matter in the universe today were produced through the decays of a heavy and relatively long-lived state, the expansion history of the universe can be indistinguishable from that predicted in the standard cosmological model with an additional neutrino. Furthermore, if these decays take place after the completion of big bang nucleosynthesis, this scenario can avoid tension with the value of three neutrino species preferred by measurements of the light element abundances.

The tight relation of star formation with molecular gas indicated by observations and assumed in recent models implies that the efficiency with which galaxies convert their gas into stars depends on gas metallicity. This is because the abundance of molecular hydrogen is sensitive to the abundance of dust, which catalyzes the formation of H_2 and helps to shield it from dissociating radiation. In this study we point out that in the absence of significant pre-enrichment by Population III stars forming out of zero metallicity gas, such H_2-based star formation is expected to leave an imprint in the form of bi-modality in the metallicity distribution among dwarf galaxies and in the metallicity distribution of stars within individual galaxies. The bi-modality arises because when gas metallicity (and dust abundance) is low, formation of molecular gas is inefficient, the gas consumption time scale is long, and star formation and metal enrichment proceed slowly. When metallicity reaches a critical threshold value star formation and enrichment accelerate, which leads to rapid increase in both stellar mass and metallicity of galaxies. We demonstrate this process both using a simple analytical model and full cosmological simulations. In contrast, observed metallicity distributions of dwarf galaxies or stars within them are not bi-modal. We argue that this discrepancy points to substantial early stochastic pre-enrichment by population III stars to levels Z ~ 0.01 Z_sun in dense, star forming regions of early galaxies.

The condensation of gas and stars in the inner regions of dark matter halos leads to a more concentrated dark matter distribution. While this effect is based on simple gravitational physics, the question of its validity in hierarchical galaxy formation has led to an active debate in the literature. We use a collection of several state-of-the-art cosmological hydrodynamic simulations to study the halo contraction effect in systems ranging from dwarf galaxies to clusters of galaxies, at high and low redshift. The simulations are run by different groups with different codes and include hierarchical merging, gas cooling, star formation, and stellar feedback. We show that in all our cases the inner dark matter density increases relative to the matching simulation without baryon dissipation, at least by a factor of several. The strength of the contraction effect varies from system to system and cannot be reduced to a simple prescription. We present a revised analytical model that describes the contracted mass profile to an rms accuracy of about 10%. The model can be used to effectively bracket the response of the dark matter halo to baryon dissipation. The halo contraction effect is real and must be included in modeling of the mass distribution of galaxies and galaxy clusters.

We study the detailed structure of galaxies at redshifts z > 2 using cosmological simulations with improved modeling of the interstellar medium and star formation. The simulations follow the formation and dissociation of molecular hydrogen, and include star formation only in cold molecular gas. The molecular gas is more concentrated towards the center of galaxies than the atomic gas, and as a consequence, the resulting stellar distribution is very compact. For halos with total mass above 10^{11} Mo, the median half-mass radius of the stellar disks is 0.8 kpc at z = 3. The vertical structure of the molecular disk is much thinner than that of the atomic neutral gas. Relative to the non-radiative run, the inner regions of the dark matter halo change shape from prolate to mildly oblate and align with the stellar disk. However, we do not find evidence for a significant dark disk of dark matter around the stellar disk. The outer halo regions retain the orientation acquired during accretion and mergers, and are significantly misaligned with the inner regions. The radial profile of the dark matter halo contracts in response to baryon dissipation, establishing an approximately isothermal profile throughout most of the halo. This effect can be accurately described by a modified model of halo contraction. The angular momentum of a fixed amount of inner dark matter is approximately conserved over time, while in the dissipationless case most of it is transferred outward during mergers. The conservation of the dark matter angular momentum provides supporting evidence for the validity of the halo contraction model in a hierarchical galaxy formation process.

A basic property of objects, like galaxies and halos that form in cosmological structure formation simulations, is their shape. Here, we critically investigate shape determination methods that are commonly used in the literature. It is found that using an enclosed integration volume and weight factors r^{-2} and r_{ell}^{-2} (elliptical radius) for the contribution of each particle or volume element in the shape tensor leads to biased axis ratios and smoothing of details when calculating the local shape as a function of distance from the center. To determine the local shape of matter distributions as a function of distance for well resolved objects (typically more than O(10^4) particles), we advocate a method that (1) uses an ellipsoidal shell (homoeoid) as an integration volume without any weight factors in the shape tensor and (2) removes subhalos.

As emphasized by previous studies, proper treatment of the density fluctuation on the fundamental scale of a cosmological simulation volume - the "DC mode" - is critical for accurate modeling of spatial correlations on scales > 10% of simulation box size. We provide further illustration of the effects of the DC mode on the abundance of halos in small boxes and show that it is straightforward to incorporate this mode in cosmological codes that use the "supercomoving" variables. The equations governing evolution of dark matter and baryons recast with these variables are particularly simple and include the expansion factor, and hence the effect of the DC mode, explicitly only in the Poisson equation.

It is a well established empirical fact that the surface density of the star formation rate, Sigma_SFR, strongly correlates with the surface density of molecular hydrogen, Sigma_H2, at least when averaged over large (~kpc) scales. Much less is known, however, if (and how) the Sigma_SFR-Sigma_H2 relation depends on environmental parameters, such as the metallicity or the UV radiation field in the interstellar medium (ISM). Furthermore, observations indicate that the scatter in the Sigma_SFR-Sigma_H2 relation increases rapidly with decreasing averaging scale. How the scale-dependent scatter is generated and how one recovers a tight ~ kpc scale Sigma_SFR-Sigma_H2 relation in the first place is still largely debated. Here, these questions are explored with hydrodynamical simulations that follow the formation and destruction of H2, include radiative transfer of UV radiation, and resolve the ISM on ~60 pc scales. We find that within the considered range of H2 surface densities (10-100 Msun/pc^2) the Sigma_SFR-Sigma_H2 relation is steeper in environments of low metallicity and/or high radiation fields (compared to the Galaxy), that the star formation rate at a given H2 surface density is larger, and the scatter is increased. Deviations from a "universal" Sigma_SFR-Sigma_H2 relation should be particularly relevant for high redshift galaxies or for low-metallicity dwarfs at z~0. We also find that the use of time-averaged SFRs produces a large, scale dependent scatter in the Sigma_SFR-Sigma_H2 relation. Given the plethora of observational data expected from upcoming surveys such as ALMA the scale-scatter relation may indeed become a valuable tool for determining the physical mechanisms connecting star formation and H2 formation.

Recent observations indicate that star formation occurs only in the molecular phase of a galaxy's interstellar medium. A realistic treatment of star formation in simulations and analytic models of galaxies therefore requires that one determine where the transition from the atomic to molecular gas occurs. In this paper we compare two methods for making this determination in cosmological simulations where the internal structures of molecular clouds are unresolved: a complex time-dependent chemistry network coupled to a radiative transfer calculation of the dissociating ultraviolet (UV) radiation field, and a simple time-independent analytic approximation. We show that these two methods produce excellent agreement at all metallicities >~10^-2 of the Milky Way value across a very wide range of UV fields. At lower metallicities the agreement is worse, likely because time-dependent effects become important; however, there are no observational calibrations of molecular gas content at such low metallicities, so it is unclear if either method is accurate. The comparison suggests that, in many but not all applications, the analytic approximation provides a viable and nearly cost-free alternative to full time-dependent chemistry and radiative transfer.

A star formation efficiency per free fall time that evolves over the life time of giant molecular clouds (GMCs) may have important implications for models of supersonic turbulence in molecular clouds or for the relation between star formation rate and H2 surface density. We discuss observational data that could be interpreted as evidence of such a time variability. In particular, we investigate a recent claim based on measurements of H2 and stellar masses in individual GMCs. We show that this claim depends crucially on the assumption that H2 masses do not evolve over the life times of GMCs. We exemplify our findings with a simple toy model that uses a constant star formation efficiency and, yet, is able to explain the observational data.

Using numerical simulations of structure formation, we investigate multiple methods of determining the strength of the proximity effect in the HI Lyalpha forest. We analyze three high resolution (~10kpc) redshift snapshots (z=4,3,2.25) of a Hydro-Particle-Mesh simulation to obtain realistic absorption spectra of the HI Lyalpha forest. We begin our analysis investigating the intrinsic biases thought to arise in the widely adopted standard technique of combining multiple lines of sight when searching for the proximity effect. We confirm the existence of this biases. We then concentrate on the analysis of the proximity effect along individual lines of sight. We construct the proximity effect strength distribution (PESD) and confirm that the PESD inferred from a simple averaging technique accurately recovers the input strength of the proximity effect at all redshifts. Moreover, the PESD closely follows the behaviors found in observed samples of quasar spectra. However, the PESD obtained from our new simulated sight lines presents some differences to that of simple Monte Carlo simulations. After developing three new theoretical methods of recovering the strength of the proximity effect on individual lines of sight, we compare their accuracy to the PESD from the simple averaging technique. All our new approaches are based on the maximization of the likelihood function, albeit invoking some modifications. The new techniques presented here fail to recover the input proximity effect in an un-biased way. Thus, employing complex 3D simulations, we provide strong evidence in favor of the proximity effect strength distribution obtained from the simple averaging technique, as method of estimating the UV background intensity, free of any biases.

We present new analytical estimates of the large scale bias of neutral hydrogen (HI). We use a simple, non-parametric model which monotonically relates the total mass of a halo M_tot with its HI mass M_HI at zero redshift; for earlier times we assume limiting models for the HI density evolution consistent with the data presently available, as well as two main scenarios for the evolution of our M_HI - M_tot relation. We find that both the linear and the first nonlinear bias terms exhibit a strong evolution with redshift, regardless of the specific limiting model assumed for the H I density over time. These analytical predictions are then shown to be consistent with measurements performed on the Millennium Simulation. Additionally, we show that this strong bias evolution does not sensibly affect the measurement of the HI power spectrum.

The effect of the cosmic UV background on star formation in high redshift galaxies is explored by means of high resolutions cosmological simulations. The simulations include star formation, 3D radiative transfer, and a highly detailed ISM model, and reach spatial resolution sufficient to resolve formation sites for molecular clouds. In the simulations the local radiation field in the Lyman-Werner band around star-forming molecular clouds dominates over the cosmic UV background by a factor of 100, similarly to the interstellar radiation field in the Milky Way and in a few high redshift galaxies for which measurements exist. The cosmic UV background, therefore, is essentially irrelevant for star formation in normal galaxies.

Using cosmological simulations with a dynamic range in excess of 10 million, we study the transport of gas mass and angular momentum through the circumnuclear region of a disk galaxy containing a supermassive black hole (SMBH). The simulations follow fueling over relatively quiescent phases of the galaxy's evolution (no mergers) and without feedback from active galactic nuclei (AGNs), as part of the first stage of using state-of-the-art, high-resolution cosmological simulations to model galaxy and black hole co-evolution. We present results from simulations at different redshifts (z=6, 4, and 3) and three different black hole masses (30 million, 90 million, and 300 million solar masses; at z=4), as well as a simulation including a prescription that approximates optically thick cooling in the densest regions. The interior gas mass throughout the circumnuclear disk shows transient and chaotic behavior as a function of time. The Fourier transform of the interior gas mass follows a power law with slope -1 throughout the region, indicating that, in the absence of the effects of galaxy mergers and AGN feedback, mass fluctuations are stochastic with no preferred timescale for accretion over the duration of each simulation (~ 1-2 Myr). The angular momentum of the gas disk changes direction relative to the disk on kiloparsec scales over timescales less than 1 Myr, reflecting the chaotic and transient gas dynamics of the circumnuclear region. Infalling clumps of gas, which are driven inward as a result of the dynamical state of the circumnuclear disk, may play an important role in determining the spin evolution of an SMBH, as has been suggested in stochastic accretion scenarios.

We present a detailed description of a phenomenological H2 formation model and local star formation prescription based on the density of molecular (rather than total) gas. Such approach allows us to avoid the arbitrary density and temperature thresholds typically used in star formation recipes. We present results of the model based on realistic cosmological simulations of high-z galaxy formation for a grid of numerical models with varied dust-to-gas ratios and interstellar far UV (FUV) fluxes. Our results show that both the atomic-to-molecular transition on small, ~10 pc scales and the Kennicutt-Schmidt (KS) relation on ~kpc scales are sensititive to the dust-to-gas ratio and the FUV flux. The atomic-to-molecular transition as a function of gas density or column density has a large scatter but is rather sharp and shifts to higher densities with decreasing dust-to-gas ratio and/or increasing FUV flux. Consequently, star formation is concentrated to higher gas surface density regions, resulting in steeper slope and lower amplitude of the KS relation at a given gas surface density, in less dusty and/or higher FUV flux environments. These trends should have a particularly strong effect on the evolution of low-mass, low surface brightness galaxies which typically have low dust content and anemic star formation, but are also likely to be important for evolution of the Milky Way-sized systems. We parameterize the dependencies observed in our simulations in convenient fitting formulae, which can be used to model the dependence of the KS relation on the dust-to-gas ratio and FUV flux in semi-analytic models and in cosmological simulations that do not include radiative transfer and H2 formation.

We present results of self-consistent, high-resolution cosmological simulations of galaxy formation at z~3. The simulations employ recently developed recipe for star formation based on the local abundance of molecular hydrogen, which is tracked self-consistently during the course of simulation. The phenomenological H2 formation model accounts for the effects of dissociating UV radiation of stars in each galaxy, as well as self-shielding and shielding of H2 by dust, and therefore allows us to explore effects of lower metallicities and higher UV fluxes prevalent in high redshift galaxies on their star formation. We compare stellar masses, metallicities, and star formation rates of the simulated galaxies to available observations of the Lyman Break Galaxies (LBGs) and find a reasonable agreement. We find that the Kennicutt-Schmidt (KS) relation exhibited by our simulated galaxies at z~3 is substantially steeper and has a lower amplitude than the z=0 relation at Sigma_gas < 100 Msun/pc^2. The predicted relation, however, is consistent with existing observational constraints for the z~3 Damped Lyman $\alpha$ (DLA) and LBGs. Our tests show that the main reason for the difference from the local KS relation is lower metallicity of the ISM in high redshift galaxies. We discuss several implications of the metallicity-dependence of the KS relation for galaxy evolution and interpretation of observations.

Eulerian hydrodynamical simulations are a powerful and popular tool for modeling fluids in astrophysical systems. In this work, we critically examine recent claims that these methods violate Galilean invariance of the Euler equations. We demonstrate that Eulerian hydrodynamics methods do converge to a Galilean-invariant solution, provided a well-defined convergent solution exists. Specifically, we show that numerical diffusion, resulting from diffusion-like terms in the discretized hydrodynamical equations solved by Eulerian methods, accounts for the effects previously identified as evidence for the Galilean non-invariance of these methods. These velocity-dependent diffusive terms lead to different results for different bulk velocities when the spatial resolution of the simulation is kept fixed, but their effect becomes negligible as the resolution of the simulation is increased to obtain a converged solution. In particular, we find that Kelvin-Helmholtz instabilities develop properly in realistic Eulerian calculations regardless of the bulk velocity provided the problem is simulated with sufficient resolution (a factor of 2-4 increase compared to the case without bulk flows for realistic velocities). Our results reiterate that high-resolution Eulerian methods can perform well and obtain a convergent solution, even in the presence of highly supersonic bulk flows.

We revisit the classical view that quasar contribution to the reionization of hydrogen is unimportant. When secondary ionization are taken into account, in many plausible scenarios for the formation and growth of supermassive black holes quasars contribute substantially or even dominantly at z>8, although their contribution generally falls below that of star-forming galaxies by z=6. Theoretical models that guide the design of the first generation of redshifted 21 cm experiments must, therefore, substantially account for the quasar contribution in order to be even qualitatively accurate.

The cosmic reionization of hydrogen was the last major phase transition in the evolution of the universe, which drastically changed the ionization and thermal conditions in the cosmic gas. To the best of our knowledge today, this process was driven by the ultra-violet radiation from young, star-forming galaxies and from first quasars. We review the current observational constraints on cosmic reionization, as well as the dominant physical effects that control the ionization of intergalactic gas. We then focus on numerical modeling of this process with computer simulations. Over the past decade, significant progress has been made in solving the radiative transfer of ionizing photons from many sources through the highly inhomogeneous distribution of cosmic gas in the expanding universe. With modern simulations, we have finally converged on a general picture for the reionization process, but many unsolved problems still remain in this young and exciting field of numerical cosmology.

The development of radiation hydrodynamical methods that are able to follow gas dynamics and radiative transfer self-consistently is key to the solution of many problems in numerical astrophysics. Such fluid flows are highly complex, rarely allowing even for approximate analytical solutions against which numerical codes can be tested. An alternative validation procedure is to compare different methods against each other on common problems, in order to assess the robustness of the results and establish a range of validity for the methods. Previously, we presented such a comparison for a set of pure radiative transfer tests (i.e. for fixed, non-evolving density fields). This is the second paper of the Cosmological Radiative Transfer (RT) Comparison Project, in which we compare 9 independent RT codes directly coupled to gasdynamics on 3 relatively simple astrophysical hydrodynamics problems: (5) the expansion of an H II region in a uniform medium; (6) an ionization front (I-front) in a 1/r^2 density profile with a flat core, and (7), the photoevaporation of a uniform dense clump. Results show a broad agreement between the different methods and no big failures, indicating that the participating codes have reached a certain level of maturity and reliability. However, many details still do differ, and virtually every code has showed some shortcomings and has disagreed, in one respect or another, with the majority of the results. This underscores the fact that no method is universal and all require careful testing of the particular features which are most relevant to the specific problem at hand.

We describe a phenomenological model for molecular hydrogen formation suited for applications in galaxy formation simulations, which includes on-equilibrium formation of molecular hydrogen on dust and approximate treatment of both its self-shielding and shielding by dust from the dissociating UV radiation. The model is applicable in simulations in which individual star forming regions - the giant molecular complexes - can be identified (resolution of tens of pc) and their mean internal density estimated reliably, even if internal structure is not resolved. In agreement with previous studies, calculations based on our model show that the transition from atomic to fully molecular phase depends primarily on the metallicity, which we assume is directly related to the dust abundance, and clumpiness of the interstellar medium. The clumpiness simply boosts the formation rate of molecular hydrogen, while dust serves both as a catalyst of molecular hydrogen formation and as an additional shielding from dissociating UV radiation. The upshot is that it is difficult to form fully-shielded giant molecular clouds while gas metallicity is low. However, once the gas is enriched to Z ~ 0.01-0.1 solar, the subsequent star formation and enrichment can proceed at a much faster rate. This may keep star formation efficiency in the low-mass, low-metallicity progenitors of galaxies very low for a certain period of time with the effect similar to a strong "feedback" mechanism. [abridged]

Great strides have been made in the last two decades in determining how galaxies evolve from their initial dark matter seeds to the complex structures we observe at z=0. The role of mergers has been documented through both observations and simulations, numerous satellites that may represent these initial dark matter seeds have been discovered in the Local Group, high redshift galaxies have been revealed with monstrous star formation rates, and the gaseous cosmic web has been mapped through absorption line experiments. Despite these efforts, the dark matter simulations that include baryons are still unable to accurately reproduce galaxies. One of the major problems is our incomplete understanding of how a galaxy accretes its baryons and subsequently forms stars. Galaxy formation simulations have been unable to accurately represent the required gas physics on cosmological timescales, and observations have only just begun to detect the star formation fuel over a range of redshifts and environments. How galaxies obtain gas and subsequently form stars is a major unsolved, yet tractable problem in contemporary extragalactic astrophysics. In this paper we outline how progress can be made in this area in the next decade.

Using the 21 cm line, observed all-sky and across the redshift range from 0 to 5, the large scale structure of the Universe can be mapped in three dimensions. This can be accomplished by studying specific intensity with resolution ~ 10 Mpc, rather than via the usual galaxy redshift survey. The data set can be analyzed to determine Baryon Acoustic Oscillation wavelengths, in order to address the question: 'What is the nature of Dark Energy?' In addition, the study of Large Scale Structure across this range addresses the questions: 'How does Gravity effect very large objects?' and 'What is the composition our Universe?' The same data set can be used to search for and catalog time variable and transient radio sources.

One of the next decade's most exciting prospects is to explore the cosmic "dark ages," during which the first stars in the Universe formed, with the 21 cm line of neutral hydrogen. At z>6, this light redshifts into the low-frequency radio regime and can be observed with new generations of wide-field arrays. These experiments have the potential to observe enormous volumes and to improve both cosmological and astrophysical constraints dramatically. In this white paper, we describe how the next decade will see the first steps toward unlocking this potential, allowing us to answer two key questions about the fundamental physics of our Universe and the intergalactic medium: Does the standard cosmological model describe the Universe during the "dark ages?" How did the intergalactic medium evolve during this important time, ending with the reionization of hydrogen?

The cosmic dark ages and the epoch of reionization, during which the first generations of stars and galaxies formed, are among the most compelling frontiers of extragalactic astrophysics and cosmology. Here we describe an exciting new probe of these eras: the 21 cm line of neutral hydrogen, which will allow us to map the neutral intergalactic medium throughout the era of first galaxy formation. In this white paper, we describe how these studies can answer two key questions about galaxy formation: What were the properties of high-redshift galaxies? How did they affect the Universe around them?

(Brief Summary) What is the total radiative content of the Universe since the epoch of recombination? The extragalactic background light (EBL) spectrum captures the redshifted energy released from the first stellar objects, protogalaxies, and galaxies throughout cosmic history. Yet, we have not determined the brightness of the extragalactic sky from UV/optical to far-infrared wavelengths with sufficient accuracy to establish the radiative content of the Universe to better than an order of magnitude. Among many science topics, an accurate measurement of the EBL spectrum from optical to far-IR wavelengths, will address: What is the total energy released by stellar nucleosynthesis over cosmic history? Was significant energy released by non-stellar processes? Is there a diffuse component to the EBL anywhere from optical to sub-millimeter? When did first stars appear and how luminous was the reionization epoch? Absolute optical to mid-IR EBL spectrum to an astrophysically interesting accuracy can be established by wide field imagingat a distance of 5 AU or above the ecliptic plane where the zodiacal foreground is reduced by more than two orders of magnitude.

Using a high resolution cosmological simulation of reionization we have examined the differing structures formed by gas and dark matter at a redshift of 5.1. Baryon-rich regions form a small number of filaments, which connect the largest galaxies in the simulation. More detailed examination of the ten largest galaxies reveals long, slender gaseous filaments about 5 proper kpc in width radiating from the galaxy centers. Extending out from each filament are a few smooth, thin, nearly planar gaseous sheets. By contrast, the dark matter concentrates into quasi-spherical bodies. The results have implications for our understanding of structure formation in the early universe and of the Lyman alpha forest.

In two previous papers, we presented simulations of the first galaxies in a representative volume of the Universe. The simulations are unique because we model feedback-regulated galaxy formation, using time-dependent, spatially-inhomogeneous radiative transfer coupled to hydrodynamics. Here, we study the properties of simulated primordial dwarf galaxies with masses <2x10^8 Msolar and investigate their impact on the intergalactic medium. While many primordial galaxies are dark, about 100--500 per comoving Mpc^3 are luminous but relatively faint. They form preferentially in chain structures, and have low surface brightness stellar spheroids extending to 20% of the virial radius. Their interstellar medium has mean density n_H~10--100 cm^-3, metallicity Z~ 0.01--0.1 Zsolar and can sustain a multi-phase structure. With large scatter, the mean efficiency of star formation scales with halo mass, <f_*> \propto M_dm^2, independent of redshift. Because of feedback, halos smaller than a critical mass, M_crit(z), are devoid of most of their baryons. More interestingly, we find that dark halos have always a smaller M_crit(z) than luminous ones. Metal enrichment of the intergalactic medium is inhomogeneous, with only a 1%--10% volume filling factor of enriched gas with [Z/H]>-3.0 and 10%--50% with [Z/H]>-5.0. At z=10, the fraction of stars with metallicity Z<10^-3 Zsolar is 10^-6 of the total stellar mass. Although detections of high-redshift dwarf galaxies with the James Webb Space Telescope will be a challenge, studies of their fossil records in the local Universe are promising because of their large spatial density.

Using a hydrodynamic adaptive mesh refinement code, we simulate the growth and evolution of a galaxy, which could potentially host a supermassive black hole, within a cosmological volume. Reaching a dynamical range in excess of 10 million, the simulation follows the evolution of the gas structure from super-galactic scales all the way down to the outer edge of the accretion disk. Here, we focus on global instabilities in the self-gravitating, cold, turbulence-supported, molecular gas disk at the center of the model galaxy, which provide a natural mechanism for angular momentum transport down to sub-pc scales. The gas density profile follows a power-law scaling as r^-8/3, consistent with an analytic description of turbulence in a quasi-stationary circumnuclear disk. We analyze the properties of the disk which contribute to the instabilities, and investigate the significance of instability for the galaxy's evolution and the growth of a supermassive black hole at the center.

An estimate for the number of ionizing photons per baryon as a function of redshift is computed based on the plausible extrapolation of the observed galaxy UV luminosity function and the latest results on the properties of the escape fraction of ionizing radiation. It is found that, if the escape fraction for low mass galaxies (Mtot<10^{11}Msun) is assumed to be negligibly small, as indicated by numerical simulations, then there are not enough ionizing photons to reionize the universe by z=6 for the cosmology favored by the WMAP 3rd year results, while the WMAP 1st year cosmology is marginally consistent with the reionization requirement. The escape fraction as a function of galaxy mass would have to be constant to within a factor of two for the whole mass range of galaxies for reionization to be possible within the WMAP 3rd year cosmology.

Nearby dwarf galaxies exhibit tight correlations between their global stellar and dynamical properties, such as circular velocity, mass-to-light ratio, stellar mass, surface brightness, and metallicity. Such correlations have often been attributed to gas or metal-rich outflows driven by supernova energy feedback to the interstellar medium. We use high-resolution cosmological simulations of high-redshift galaxies with and without energy feedback, as well as analytic modeling, to investigate whether the observed correlations can arise without supernova-driven outflows. We find that the simulated dwarf galaxies exhibit correlations similar to those observed as early as z~10, regardless of whether supernova feedback is included. We also show that the correlations can be well reproduced by our analytic model that accounts for realistic gas inflow but assumes no outflows, and star formation rate obeying the Kennicutt-Schmidt law with a critical density threshold. We argue that correlations in simulated galaxies arise due to the increasingly inefficient conversion of gas into stars in low-mass dwarf galaxies rather than supernova-driven outflows. We also show that the decrease of the observed effective yield in low-mass objects, often used as an indicator of gas and metal outflows, can be reasonably reproduced in our simulations without outflows. We show that this trend can arise if a significant fraction of metals in small galaxies is spread to the outer regions of the halo outside the stellar extent via mixing. In this case the effective yield can be significantly underestimated if only metals within the stellar radius are taken into account. Measurements of gas metallicity in the outskirts of gaseous disks of dwarfs would thus provide a key test of such explanation.

We describe a new method to measure the escape fraction fesc of ionizing radiation from distant star-forming galaxies using the afterglow spectra of long-duration gamma-ray bursts (GRBs). Optical spectra of GRB afterglows allow us to evaluate the optical depth of the host ISM, according to the neutral hydrogen column density N(HI) observed along the sightlines toward the star-forming regions where the GRBs are found. Different from previous effort in searching for faint, transmitted Lyman continuum photons, our method is not subject to background subtraction uncertainties and does not require prior knowledge of either the spectral shape of the host galaxy population or the IGM Lya forest absorption along these GRB sightlines. Because most GRBs occur in sub-L_* galaxies, our study also offers the first constraint on fesc for distant low-mass galaxies that dominate the cosmic luminosity density. We have compiled a sample of 27 GRBs at redshift z>2 for which the underlying N(HI) in the host ISM are known. These GRBs together offer a statistical sampling of the integrated optical depth to ionizing photons along random sightlines from star-forming regions in the host galaxies, and allow us to estimate the mean escape fraction <fesc> averaged over different viewing angles. We find <fesc>=0.02\pm 0.02 and place a 95% c.l. upper limit <fesc> <= 0.075 for these hosts. We discuss possible biases of our approach and implications of the result. Finally, we propose to extend this technique for measuring <fesc> at z~0.2 using spectra of core-collapse supernovae.

We model the escape of ionizing radiation from high-redshift galaxies using high-resolution Adaptive Mesh Refinement N-body + hydrodynamics simulations. Our simulations include time-dependent and spatially-resolved transfer of ionizing radiation in three dimensions, including effects of dust absorption. For galaxies of total mass M > 10^11 Msun and star formation rates SFR ~ 1-5 Msun/yr, we find angular averaged escape fractions of 0.01-0.03 over the entire redshift interval studied (3<z<9). In addition, we find that the escape fraction varies by more than an order of magnitude along different lines-of-sight within individual galaxies, from the largest values near galactic poles to the smallest along the galactic disk. The escape fraction declines steeply at lower masses and SFR. We show that the low values of escape fractions are due to a small fraction of young stars located just outside the edge of HI disk. We compare our predicted escape fraction of ionizing photons with previous results, and find a general agreement with both other simulation results and available direct detection measurements at z ~ 3. We also compare our simulations with a novel method to estimate the escape fraction in galaxies from the observed distribution of neutral hydrogen column densities along the lines of sights to long duration gamma-ray bursts. Using this method we find escape fractions of the GRB host galaxies of 2-3%, consistent with our theoretical predictions. [abridged]

We study the luminosity functions of high-redshift galaxies in detailed hydrodynamic simulations of cosmic reionization, which are designed to reproduce the evolution of the Lyman-alpha forest between z=5 and z=6. We find that the luminosity functions and total stellar mass densities are in agreement with observations when plausible assumptions about reddenning at z=6 are made. Our simulations support the conclusion that stars alone reionized the universe.

We used cosmological simulation with self-consistent radiative transfer to investigate the physical nature of Lyman Limit systems at z=4. In agreement with previous studies, we find that most of Lyman Limit systems are ionized by the cosmological background, while higher column density systems seem to be illuminated by the local sources of radiation. In addition, we find that most of Lyman limit systems in our simulations are located within the virial radii of galaxies with a wide range of masses, and are physically associated with them (``bits and pieces'' of galaxy formation). While the finite resolution of our simulations cannot exclude an existence of a second population of self-shielded, neutral gas clouds located in low mass dark matter halos (``minihalos''), our simulations are not consistent with ``minihalos'' dominating the total abundance of Lyman limit systems.

We have investigated the effects of AGN outflows on the amplitude of the matter power spectrum in a simple model of spherically symmetric outflows around realistically clustered AGN population. We find that two competing effects influence the matter power spectrum in two opposite directions. First, AGN outflows move baryons from high to low density regions, decreasing the amplitude of the matter power spectrum by up to 20%. Second, high clustering of the AGN transfers the power from small to larger scales. The exact balance between these two effects depends on the details of outflows on small scales, and quantitative estimates will require much more sophisticated modeling than presented here.

We show that numerical simulations of reionization that resolve the Lyman Limit systems (and, thus, correctly count absorptions of ionizing photons) have converged to about 10% level for 5<z<6.2 and are in reasonable agreement (within 10%) with the SDSS data in this redshift interval. The SDSS data thus constraint the redshift of overlap of cosmic HII regions to z_{OVL} = 6.1+-0.15. At higher redshifts, the simulations are far from convergence on the mean Gunn-Peterson optical depth, but achieve good convergence for the mean neutral hydrogen fraction. The simulations that fit the SDSS data, however, do not have nearly enough resolution to resolve the earliest episodes of star formation, and are very far from converging on the precise value of the optical depth to Thompson scattering - any value between 6 and 10% is possible, depending on the convergence rate of the simulations and the fractional contribution of PopIII stars. This is generally consistent with the third-year WMAP results, but much higher resolution simulation are required to come up with the sufficiently precise value for the Thompson optical depth that can be statistically compared with the WMAP data.

Radiative transfer simulations are now at the forefront of numerical astrophysics. They are becoming crucial for an increasing number of astrophysical and cosmological problems; at the same time their computational cost has come to the reach of currently available computational power. Further progress is retarded by the considerable number of different algorithms (including various flavours of ray-tracing and moment schemes) developed, which makes the selection of the most suitable technique for a given problem a non-trivial task. Assessing the validity ranges, accuracy and performances of these schemes is the main aim of this paper, for which we have compared 11 independent RT codes on 5 test problems: (0) basic physics, (1) isothermal H II region expansion and (2) H II region expansion with evolving temperature, (3) I-front trapping and shadowing by a dense clump, (4) multiple sources in a cosmological density field. The outputs of these tests have been compared and differences analyzed. The agreement between the various codes is satisfactory although not perfect. The main source of discrepancy appears to reside in the multi-frequency treatment approach, resulting in different thicknesses of the ionized-neutral transition regions and different temperature structure. The present results and tests represent the most complete benchmark available for the development of new codes and improvement of existing ones. To this aim all test inputs and outputs are made publicly available in digital form.

We use a combination of high-resolution gasdynamics simulations of high-redshift dwarf galaxies and dissipationless simulations of a Milky Way sized halo to estimate the expected abundance and spatial distribution of the dwarf satellite galaxies that formed most of their stars around z~8 and evolved only little since then. Such galaxies can be considered as fossils of the reionization era, and studying their properties could provide a direct window into the early, pre-reionization stages of galaxy formation. We show that 5-15% of the objects existing at z~8 do indeed survive until the present in the MW like environment without significant evolution. This implies that it is plausible that the fossil dwarf galaxies do exist in the Local Group. Because such galaxies form their stellar systems early during the period of active merging and accretion, they should have spheroidal morphology regardless of their current distance from the host galaxy. We show that both the expected luminosity function and spatial distribution of dark matter halos which are likely to host fossil galaxies agree reasonably well with the observed distributions of the luminous (L_V>10^6 Lsun) Local Group fossil candidates near the host galaxy (d<200 kpc). However, the predicted abundance is substantially larger (by a factor of 2-3) for fainter galaxies (L_V<10^6 Lsun) at larger distances (d>300 kpc). We discuss several possible explanations for this discrepancy.

We use cosmological simulations to explore the large-scale effects of reionization. Since reionization is a process that involves a large dynamic range - from galaxies to rare bright quasars - we need to be able to cover a significant volume of the universe in our simulation without losing hte important small scale effects from galaxies. Here we have taken an approach that uses clumping factors derived from small scale simulations to approximate the radiative transfer on the sub-cell scales. Using this technique, we can cover a simulation size up to $1280 h^{-1} Mpc$ with $10 h^{-1} Mpc$ cells. This allows us to construct synthetic spectra of quasars similar to observed spectra of SDSS quasars at high reshifts and compare them to the observational data. These spectra can then be analyzed for HII region sizes, the presence of the Gunn-Peterson trough and the Lyman-$\alpha$ forest.

We combine cosmological simulations with an AGN luminosity function constrained by optical surveys in order to create a realistic AGN distribution in which to model AGN outflows. The outflows are modeled with assumptions of spherical symmetry and energy conservation. We study the importance of the kinetic luminosity of AGN outflows in determining the fractional volume of the IGM filled with outflows as a function of redshift. We find that kinetic luminosities of < 10% of the bolometric luminosities are required or the entire IGM would be filled with outflows at z=0. We also examine the effects of varying the AGN lifetime and the bias parameter, which describes how tightly correlated AGN are to high density regions. We find that longer lifetimes (tau_{AGN} ~ 1 Gyr), as well as small bias parameters increase the filling fraction.

We present a new, galaxy-halo model of large-scale structure, in which the galaxies entering a given sample are the fundamental objects. Haloes attach to galaxies, in contrast to the standard halo model, in which galaxies attach to haloes. The galaxy-halo model pertains mainly to the relationships between the power spectra of galaxies and mass, and their cross-power spectrum. With surprisingly little input, an intuition-aiding approximation to the galaxy-matter cross-correlation coefficient R(k) emerges, in terms of the halo mass dispersion. This approximation seems valid to mildly non-linear scales (k < ~3 h/Mpc), allowing measurement of the bias and the matter power spectrum from measurements of the galaxy and galaxy-matter power spectra (or correlation functions). This is especially relevant given the recent advances in precision in measurements of the galaxy-matter correlation function from weak gravitational lensing. The galaxy-halo model also addresses the issue of interpreting the galaxy-matter correlation function as an average halo density profile, and provides a simple description of galaxy bias as a function of scale.

We compare the properties of dwarf galaxies in the Local Group with the simulated galaxies formed before reionization in a cosmological simulation of unprecedented spatial and mass resolution,including radiative feedback effects. We find that a subset of the Local Group dwarfs are remarkably similar to the simulated dwarf galaxies in all their properties already before reionization. Based on this similarity, we propose the hypothesis that Local Group dwarfs form in a variety of ways: some of them are ``true fossils'' of the pre-reionization era, some of them form most of their stars later, after reionization (we call them ``survivors'' of the reionization era), and the rest of them form an intermediate group of ``polluted fossils''. We also identify a simple observational test that is able to falsify our hypothesis.

We use cosmological simulations of reionization to predict the effect of large HII regions around individual high-redshift quasars on the possible signal from the redshifted 21cm line of neutral hydrogen in the pre-reionization era. We show that these HII regions appear as ``spectral dips'' in frequency space with equivalent widths in excess of 3 mK MHz or depths in excess of about 1.5 mK, and that they are by far the most prominent cosmological signals in the redshifted HI distribution. These spectral dips are expected to be present in almost every line of sight. If the spectral dips of a large enough sample of HII regions are well resolved in frequency space, the distribution of line depth and equivalent width in frequency with a known observing beamsize can be used to infer the HII region size distribution and the mean difference in neutral hydrogen density between the \HII regions (which may contain self-shielded neutral gas clumps) and the surrounding medium, providing a powerful test for models of reionization.

We use cosmological simulations to study the X-ray ionisation and heating of the intergalactic medium by an early population of accreting black holes. By considering observational constraints from the X-ray background, we find an upper limit for the optical depth to Thompson scattering tau_e~0.17. The redshifted soft X-ray background from these early sources produces:(i) fully ionised atomic hydrogen in the low density intergalactic medium before redshift z~7 (consequently stellar reionisation is characterised by an instantaneous overlap phase of HII regions),(ii) a second HeII reionisation at z~3 and (iii) heats the intergalactic medium to near 10000 K at low redshifts. The typical luminosity in the soft X-ray band of the galaxies hosting the black holes is about one order of magnitude below the sensitivity limit of the Chandra deep field. About a thousand of these sources may be present per square arcmin of the sky, producing detectable fluctuations. Few rarer objects could be luminous enough to be visible in the Chandra deep field. XEUS and Con-X satellites will be able to detect more of these sources that, if radio loud, could be used to study the 21cm forest in absorption.A signature of an early X-ray preionisation is the production of secondary CMB anisotropies on small angular scales. We find that in these models the power spectrum of temperature fluctuations increases with decreasing angular scale (dT~16 muK at 1arcsec scales), while for stellar reionisation scenarios the power decreases on smaller scales. We also show that the redshifted 21 cm radiation from neutral hydrogen can be marginally detected in emission at redshifts 7<z<12. At a redshift of about z~30 a stronger and narrower (in redshift space) signal in absorption against the CMB, that is peculiar to these models,could be detectable.[abridged]

We explore the possibility of detecting HI 21cm absorption by the neutral intergalactic medium (IGM) toward very high redshift radio sources, and by gas associated with the first collapsed structures, using the Square Kilometer Array at low frequency (100 to 200 MHz). The epoch considered is between the time when the first ionizing sources form and when the bulk of the neutral IGM becomes ionized. Expected IGM absorption signal includes $\sim 1%$ absorption by the mean neutral IGM (the radio 'Gunn-Peterson' effect'), plus deeper, narrow lines ($\ge 5%$, a few km s$^{-1}$) arising in mild density inhomogeneities with typical values of cosmic overdensity $\delta\sim 10$, precisely the structures that at later times give rise to the Ly$\alpha$ forest (the `21cm forest'). Absorption can also arise in gas associated with collapsed structures ($\delta\ge 100$), including 'minihalos' ($\le 10^7$ M$_\odot$) and protodisks ($\ge 10^8$ M$_\odot$). We consider SKA sensitivity limits and the evolution of radio source populations, and conclude that it is reasonable to hypothesize the existence of an adequate number of high-$z$ radio sources against which such absorption studies could be performed, provided that reionization occurs at $z < 10$. Lastly, we discuss the possibility of `line confusion' due to radio recombination lines arising in the ionized IGM.

We have developed an algorithm to find haloes in an N-body dark matter simulation, called VOBOZ (VOronoi BOund Zones), which has as little dependence on free parameters as we can manage. By using the Voronoi diagram, we achieve nonparametric, `natural' measurements of each particle's density and set of neighbors. We then eliminate much of the ambiguity in merging sets of particles together by identifying every possible density peak, and measuring the probability that each does not arise from Poisson noise. This probability could be interpreted as a measure of a halo's sub-haloness, and might be of use to researchers of the halo model. The first parameter in VOBOZ controls the subtlety of particle unbinding, and may be eliminated if one is cavalier with processor time; even if one is not, the results saturate to the parameter-free answer when the parameter is sufficiently small. The only parameter which remains, an outer density cut-off, does not influence whether or not haloes are identified, nor does it have any effect on subhaloes; it only affects the masses returned for supercluster haloes.

I show that advanced simulations of cosmological reionization are able to fit the observed data on the mean transmitted flux in the hydrogen Lyman-alpha line at z~6. At the same time, posteriori models can be constructed that also produce a large value (20%) for the Thompson scattering optical depth, consistent with the WMAP measurements. Thus, it appears that a consistent picture emerges in which early reionization (as suggested by WMAP) is complete by z~6 in accord with the SLOAN data.

We use a numerical simulation of cosmological reionization to estimate the likelihood of detecting Lyman-alpha emitting galaxies during the pre-reionization era. We show that it is possible to find galaxies even at z~9 that are barely affected by the dumping wing of the Gunn-Peterson absorption from the neutral IGM outside of their HII regions. The damping wing becomes rapidly more significant at z>9, but even at z>10 is it not inconceivable (although quite hard) to see a Lyman-alpha emission line from a star-forming galaxy.

We use cosmological simulations of reionization to predict the possible signal from the redshifted 21 cm line of neutral hydrogen in the pre-reionization era and to investigate the observability of this signal. We show that the separation of the mean (global) signal over the whole sky from the known foreground contamination may be feasible, but very challenging. In agreement with previous studies, we demonstrate that measuring angular fluctuations in the HI signal is likely to be extremely difficult if not impossible because of the overwhelming contamination from the galactic and extragalactic foregrounds. However,we show that the sharp HI fluctuations in the frequency domain should be easily separable from the relatively smooth spectra of the foregrounds, and that these fluctuations should be detectable even at moderate angular resolution (10-20 arcmin).

By comparing the PSCz galaxy power spectrum with the results of nested pure dark matter N-body simulations, we try to understand how infrared-selected galaxies populate dark-matter haloes, paying special attention to the method of halo identification in the simulations. We thus test the hypothesis that baryonic physics negligibly affects the distribution of galaxies down to the smallest scales yet observed. We are successful in reproducing the PSCz power spectrum on scales < ~40 h/Mpc, near our resolution limit, by imposing a central density cut-off on simulated haloes, which gives a rough minimum mass and circular velocity of haloes in which PSCz galaxies formed.

We present the star formation histories, luminosities, colors, mass to light ratios, and halo masses of "galaxies" formed in a simulation of cosmological reionization. We compare these galaxies with Lyman Break Galaxies observed at high redshift. While the simulation is severely limited by the small box size, the simulated galaxies do appear as fainter cousins of the observed ones. They have similar colors, and both simulated and observed galaxies can be fitted with the same luminosity function. This is significant to the process of reionization since in the simulation the stars alone are responsible for reionization at a redshift consistent with that of the Gunn-Peterson trough seen in the absorption spectra of QSO's. The brightest galaxies at z = 4 are indeed quite young in accord with observational studies. But these brightest galaxies are free riders - they contributed only about 25% to the reionization of the universe at z > 6. Instead, the bulk of the work was done by dimmer galaxies, those that fall within the 28 < I < 30 magnitude range at z = 4. These dimmer galaxies are not necessarily less massive than the brightest ones.

By comparing the PSCz galaxy power spectrum with haloes from nested and phased N-body simulations, we try to understand how IRAS infrared-selected galaxies populate dark-matter haloes. We pay special attention to the way we identify haloes in the simulations.

We report on the linear matter power spectrum reconstructed from the Lyman alpha forest, on baryonic wiggles in the galaxy power spectrum, and on how to measure real space (as opposed to redshift space) power.

We investigate the relationship between the dark matter and baryons in the linear regime. This relation is quantified by the so-called ``filtering scale''. We show that a simple gaussian ansatz which uses the filtering scale provides a good approximation to the exact solution.

We show how the morphological analysis of the maps of the secondary CMB anisotropies can detect an extended period of ``smoldering'' reionization, during which the universe remains partially ionized. Neither radio observations of the redshifted 21cm line nor IR observations of the redshifted Lyman-alpha forest will be able to detect such a period. The most sensitive to this kind of non-gaussianity parameters are the number of regions in the excursion set, the perimeter of the excursion set, and the genus of the largest (by area) region. For example, if the universe reionized fully at z=6, but maintained about 1/3 mean ionized fraction since z=20, then a 2 arcmin map with 500x500 pixel resolution and a signal-to-noise ratio S/N=1/2 allows to detect the non-gaussianity due to reionization with better than 99% confidence level.

We explore the possibility of detecting HI 21cm absorption by the neutral intergalactic medium (IGM) toward very high redshift radio sources. The epoch considered is between the time when the first ionizing sources form and when the bulk of the neutral IGM becomes ionized. Due to the extreme Ly-alpha opacities of the neutral IGM, objects within this 'gray age' can only be observed at wavelengths longer than about 1micron. We use the latest simulations of the evolution of the IGM in the context of LambdaCDM structure formation models constrained by observations of the highest redshift QSOs to predict the optical depth as a function of frequency of the neutral IGM due to the HI 21cm line. We then simulate radio spectra assuming observational parameters for future large area radio telescopes. These spectra show that HI 21cm absorption studies could be a powerful probe of the rich structure of the neutral IGM prior to the epoch of reionization, including 1% absorption by the mean neutral IGM, plus deeper, narrow lines (5% and a few km/s). Most of the variations in transmissivity are due to the mild density inhomogeneities with typical values of the cosmic overdensity of about 10, precisely the structures that at later times give rise to the Ly-alpha forest. We also consider sensitivity limits and the evolution of radio source populations, and conclude that it is reasonable to hypothesize the existence of an adequate number of high-z radio sources against which such absorption studies could be performed.

[abridged]We use 3D cosmological simulations with radiative transfer to study the formation and evolution of the first galaxies in a LCDM cosmology. We find that the first luminous objects ("small-halos") are characterized by a bursting star formation (SF) that is self-regulated by a feedback process acting on cosmological instead of galactic scales. The global star formation history is regulated by the mean number of ionizing photons that escape from each source, \epsilon_{UV}\fesc. It is almost independent of the assumed star formation efficiency parameter, \epsilon_*, and the intensity of the dissociating background. The main feedback process that regulates the SF is the re-formation of H_2 in front of HII regions and inside relic HII regions. The HII regions remain confined inside filaments, maximizing the production of H_2 in overdense regions through cyclic destruction/reformation of H_2. If \epsilon_{UV}\fesc > 10^{-8}/\epsilon_* the SF is self-regulated, photo-evaporation of "small-halo" objects dominate the metal pollution of the low density IGM, and the mass of produced metals depends only on \fesc. If \epsilon_{UV}\fesc \simlt 10^{-8}/\epsilon_*, positive feedback dominates, and "small-halo" objects constitute the bulk of the mass in stars and metals at least until redshift z \sim 10. "Small-halo" objects cannot reionize the universe because the feedback mechanism confines the HII regions inside the large scale structure filaments. In contrast to massive objects, which can reionize voids, "small-halo" objects partially ionize only the dense filaments while leaving the voids mostly neutral.

In cold dark matter (CDM) cosmogonies, low-mass objects play an important role in the evolution of the universe. Not only are they the first luminous objects to shed light in a previously dark universe, but, if their formation is not inhibited by their own feedback, they dominate the galaxy mass function until redshift z \sim 5. In this paper we present and discuss the implementation of a 3D cosmological code that includes most of the needed physics to simulate the formation and evolution of the first galaxies with a self-consistent treatment of radiative feedback. The simulation includes continuum radiative transfer using the ``Optically Thin Variable Eddington Tensor'' (OTVET) approximation and line-radiative transfer in the H_2 Lyman-Werner bands of the background radiation. We include detailed chemistry for H_2 formation/destruction, molecular and atomic cooling/heating processes, ionization by secondary electrons, and heating by Ly\alpha resonant scattering. We find that the first galaxies ("small-halos") are characterized by a bursting star formation, self-regulated by a feedback process that acts on cosmological scales. Their formation is not suppressed by feedback processes; therefore, their impact on cosmic evolution cannot be neglected. The main focus of this paper is on the methodology of the simulations, and we only briefly introduce some of the results. An extensive discussion of the results and the nature of the feedback mechanism are the focus of a companion paper.

We investigate possible systematic errors in the recent measurement of the matter power spectrum from the Lyman-alpha forest by Croft et al. (2001). We find that for a large set of prior cosmological models the Croft et al. result holds quite well, with systematic errors being comparable to random ones, when a dependence of the recovered matter power spectrum on the cosmological parameters at z~3 is taken into account. We find that peculiar velocities cause the flux power spectrum to be smoothed over about 100-300 km/s, dependng on scale. Consequently, the recovered matter power spectrum is a smoothed version of the underlying true power spectrum. Uncertainties in the recovered power spectrum are thus correlated over about 100-300 km/s. As a side effect, we find that residual fluctuations in the ionizing background, while having almost no effect on the recovered matter power spectrum, significantly bias estimates of the baryon density from the Lyman-alpha forest data. We therefore conclude that the Croft et al. result provides a powerful new constraint on cosmological parameters and models of structure formation.

Comparison of the observed evolution of the Ly-alpha transmitted flux in the spectra of four highest redshift quasars discovered by SLOAN survey with the theoretical prediction for this evolution based on the state-of-the-art numerical simulations of cosmological reionization already allows one to constrain the redshift of reionization to z = 6.2 +- 0.1s +- 0.2r, where systematic and random errors are given respectively.

We use one-dimensional radiative transfer simulations to study the evolution of H_2 gas-phase (H^- catalyzed) formation and photo-dissociation regions in the primordial universe. We find a new positive feedback mechanism capable of producing shells of H_2 in the intergalactic medium, which are optically thick in some Lyman-Werner bands. While these shells exist, this feedback effect is important in reducing the H_2 dissociating background flux and the size of photo-dissociation spheres around each luminous object. The maximum background opacity of the IGM in the H_2 Lyman-Werner bands is \tau_{H_2} ~ 1-2 for a relic molecular fraction x_{H_2}=2 x 10^{-6}, about 6 times greater than found by Haiman, Abel & Rees. Therefore, the relic molecular hydrogen can decrease the photo-dissociation rate by about an order of magnitude. The problem is relevant to the formation of small primordial galaxies with masses M_{DM} < 10^8 M_\odot, that rely on molecular hydrogen cooling to collapse. Alternatively, the universe may have remained dark for several hundred million years after the birth of the first stars, until galaxies with virial temperature T_{vir} > 10^4 K formed.

We present a new approach to numerically model continuum radiative transfer based on the Optically Thin Variable Eddington Tensor (OTVET) approximation. Our method insures the exact conservation of the photon number and flux (in the explicit formulation) and automatically switches from the optically thick to the optically thin regime. It scales as N logN with the number of hydrodynamic resolution elements and is independent of the number of sources of ionizing radiation (i.e. works equally fast for an arbitrary source function). We also describe an implementation of the algorithm in a Soften Lagrangian Hydrodynamic code (SLH) and a multi--frequency approach appropriate for hydrogen and helium continuum opacities. We present extensive tests of our method for single and multiple sources in homogeneous and inhomogeneous density distributions, as well as a realistic simulation of cosmological reionization.

Future observations with ALMA will be able to detect star-forming primordial gal axies, and perhaps even their dwarf spheroidal satellites, in CO emission lines at redshifts approaching 10. These observations will compliment other tools designed to study the dawn of galaxy formation, such as NGST and FIRST.

It is well known that the mass function for_halos_ in CDM cosmology is a relatively steep power law for low masses, possibly too steep to be consistent with observations. But how steep is the_galaxy_ mass function? We have analyzed the stellar and gas mass functions of the first massive luminous objects formed in a \Lambda CDM universe, as calculated in the numerical simulation described in Gnedin (2000ab). We found that while the dark matter mass function is steep, the stellar and gas mass functions are flatter for low mass objects. The stellar mass function is consistently flat at the low mass end. Moreover, while the gas mass function follows the dark matter mass function until reionization at z~7, between z=7 and z=4, the gas mass function also flattens considerably at the low mass end. At z=4, the gas and stellar mass functions are fit by a Schechter function with \alpha ~ -1.2 +/- 0.1, significantly shallower than the dark matter halo mass function and consistent with some recent observations. The baryonic mass functions are shallower because (a) the dark matter halo mass function is consistent with the Press-Schechter formulation at low masses n(M) M^-2 and (b) heating/cooling and ionization processes appear to cause baryons to collect in halos with the relationship M_b M_d^4 at low masses. Combining (a) and (b) gives n(M_b) M_b^-5/4, comparable to the simulation results. Thus, the well known observational fact that low mass galaxies are underabundant as compared to expectations from numerical dark matter simulations or Press-Schechter modeling of CDM universes emerges naturally from these results, implying that perhaps no ``new physics'' beyond the standard model is needed.

We use numerical simulation of cosmological reionization to calculate the secondary CMB anisotropies in a representative flat low density cosmological model. We show that the kinetic Sunyaev-Zel'dovich effect (scattering off of moving electrons in the ionized intergalactic medium) is dominated by the nonlinear hydrodynamic and gravitational evolution of the density and velocity fields, rather than the detailed distribution of the ionization fraction (``patchy reionization'') on all angular scales. Combining our results with the recent calculation of secondary CMB anisotropies by Springel et al., we are able to accurately predict the power spectrum of the kinetic SZ effect on almost all angular scales.

I use simulations of cosmological reionization to quantify the effect of photoionization on the gas fraction in low mass objects, in particular the characteristic mass scale below which the gas fraction is reduced compared to the universal value. I show that this characteristic scale can be up to an order of magnitude lower than the linear theory Jeans mass, and that even if one defines the Jeans mass at a higher overdensity, it does not track the evolution of this characteristic suppression mass. Instead, the filtering mass, which corresponds directly to the scale over which baryonic perturbations are smoothed in linear perturbation theory, provides a remarkably good fit to the characteristic mass scale. Thus, it appears that the effect of reionization on structure formation in both the linear and nonlinear regimes is described by a single characteristic scale, the filtering scale of baryonic perturbations. In contrast to the Jeans mass, the filtering mass depends on the full thermal history of the gas instead of the instantaneous value of the sound speed, so it accounts for the finite time required for pressure to influence the gas distribution in the expanding universe. In addition to the characteristic suppression mass, I study the full shape of the probability distribution to find an object with a given gas mass among all the objects with the same total mass, and I show that the numerical results can be described by a simple fitting formula that again depends only on the filtering mass. This simple description of the probability distribution may be useful for semi-analytical modeling of structure formation in the early universe.

I show how the existing observational data on Local Group dwarf galaxies can be used to estimate the average star formation law during the first 3 Gyr of the history of the universe. I find that the observational data are consistent with the orthodox Schmidt law with a star formation efficiency of about 4 percent if the star formation is continuous (during the first 3 Gyr). The efficiency is proportionally higher if most of the gas in the dwarfs was consumed (and never replenished) in a short time interval well before the universe turned 3 Gyr.

We use cosmological hydrodynamic simulations to investigate formation of galactic bulges within the framework of hierarchical clustering in a representative CDM cosmological model. We show that largest objects forming at cosmological redshifts z~4 resemble observed bulges of spiral galaxies or moderate size ellipticals in their general properties like sizes, shapes, and density profiles. This is consistent with observational data indicating the existence of ``old'' bulges and ellipticals at more moderate redshifts. These bulges are gas dominated at redshift z=3, with high rates of star formation and would appear to be good candidates for small blue galaxies seen in the Hubble Deep Field.

We investigate the generation of magnetic field by the Biermann battery in cosmological ionization fronts, using new simulations of the reionization of the universe by stars in protogalaxies. Two mechanisms are primarily responsible for magnetogenesis: i) the breakout of I-fronts from protogalaxies, and ii) the propagation of I-fronts through the high density neutral filaments which are part of the cosmic web. The first mechanism is dominant prior to overlapping of ionized regions (z ~ 7), whereas the second continues to operate even after that epoch. However, after overlap the field strength increase is largely due to the gas compression occurring as cosmic structures form. As a consequence, the magnetic field at z ~ 5 closely traces the gas density, and it is highly ordered on megaparsec scales. The mean mass-weighted field strength is B_0 ~ 10^{-19} G in the simulation box. There is a relatively well-defined, nearly linear correlation between B_0 and the baryonic mass of virialized objects, with B_0 ~ 10^{-18} G in the most massive objects (M ~ 10^9 M_sun) in our simulations. This is a lower limit, as lack of numerical resolution prevents us from following small scale dynamical processes which could amplify the field in protogalaxies. Although the field strengths we compute are probably adequate as seed fields for a galactic dynamo, the field is too small to have had significant effects on galaxy formation, on thermal conduction, or on cosmic ray transport in the intergalactic medium. It could, however, be observed in the intergalactic medium through innovative methods based on analysis of gamma-ray burst photon arrival times.

We develop a method to extract the "effective equation of state" of the intergalactic medium from the doppler b parameter distribution of the low-density Lyman-alpha forest. We test the method on numerical simulations and then apply it to published observations of the Lyman-alpha forest at redshifts z from 0 to 4. We find that the effective equation of state is close to isothermal at redshift z=3, indicating that a second reheating of the IGM took place at z=3. This reheating can plausibly be identified with the reionization of HeII observed to occur at z about 3.

I use cosmological simulations that incorporate a physically motivated approximation to three-dimensional radiative transfer that recovers correct asymptotic ionization front propagation speeds for some cosmologically relevant density distributions transfer to investigate the process of the reionization of the universe by ionizing radiation from proto-galaxies. Reionization proceeds in three stages and occupies a large redshift range from z~15 until z~5. During the first, ``pre-overlap'' stage, HII regions gradually expand into the low density IGM, leaving behind neutral high density protrusions. During the second, ``overlap'' stage, that occurs in about 10% of the Hubble time, HII regions merge and the ionizing background rises by a large factor. During the third, ``post-overlap'' stage, remaining high density regions are being gradually ionized as the required ionizing photons are being produced. Residual fluctuations in the ionizing background reach significant (more than 10%) levels for the Lyman-alpha forest absorption systems with column densities above 10^14 - 10^15 cm^-2 at z=3 to 4.

We solve the Boltzmann equation for cosmological neutrinos around the epoch of the electron-positron annihilation in order to verify the freeze-out approximation and to compute accurately the cosmological neutrino distribution function. We find the radiation energy density to be about 0.3% higher than the one predicted by the freeze-out approximation. As a result, the spectrum of the Cosmic Microwave Background anisotropies changes by 0.3-05%, depending on the angular scale, and the amplitude of the mass fluctuations on scales below about 100 h^{-1} Mpc decreases by about 0.2-0.3%.

We have simulated the formation of an X-ray cluster in a cold dark matter universe using 12 different codes. The codes span the range of numerical techniques and implementations currently in use, including SPH and grid methods with fixed, deformable or multilevel meshes. The goal of this comparison is to assess the reliability of cosmological gas dynamical simulations of clusters in the simplest astrophysically relevant case, that in which the gas is assumed to be non-radiative. We compare images of the cluster at different epochs, global properties such as mass, temperature and X-ray luminosity, and radial profiles of various dynamical and thermodynamical quantities. On the whole, the agreement among the various simulations is gratifying although a number of discrepancies exist. Agreement is best for properties of the dark matter and worst for the total X-ray luminosity. Even in this case, simulations that adequately resolve the core radius of the gas distribution predict total X-ray luminosities that agree to within a factor of two. Other quantities are reproduced to much higher accuracy. For example, the temperature and gas mass fraction within the virial radius agree to about 10%, and the ratio of specific kinetic to thermal energies of the gas agree to about 5%. Various factors contribute to the spread in calculated cluster properties, including differences in the internal timing of the simulations. Based on the overall consistency of results, we discuss a number of general properties of the cluster we have modelled.

We use numerical simulations of structure formation in a Cold Dark Matter model to predict the absorption lines in the soft X-rays produced by heavy elements in the shock-heated intergalactic medium at low redshift. The simulation incorporates a model for metal production in galaxies and the subsequent dispersion of the metals to the intergalactic medium. We analyze in particular absorption lines produced by oxygen, and calculate the ionization stage taking into account the observed X-ray background at the present time. We find that oxygen is fully ionized by the X-ray background in low-density voids, and is mostly in the form of OVII and OVIII in the sheets and filamentary regions. Strong absorption lines of OVII and OVIII with equivalent widths W ~ 100 km/s are produced in filamentary regions of overdensities ~ 100 and temperatures ~ 10^6 K, located in the outskirts of groups and clusters of galaxies. The OVII line at E = 574 eV is generally the strongest one in these systems. Our model predicts that any X-ray source (such as a quasar) should typically show about one OVII absorption line with W > 100 km/s in the interval from z=0 to z=0.3. These lines could be detected with the upcoming generation of X-ray telescopes, and their origin in intervening systems could be confirmed by the association with groups of galaxies and X-ray emitting halos near the line-of-sight at the same redshift. The hot intergalactic medium may be one of the main reservoirs of baryons in the present universe, and the heavy element X-ray absorption lines offer a promising possibility of detecting this new component in the near future.

Emission lines in hydrogen can be used to measure the approximate redshift of the reionization of the universe. This is an important measurement given the lack of a convincing theoretical prediction of this epoch. There is a rapid change in absorption at the reionization event, which can leave sharp features in the spectrum of recombining hydrogen. These features in the near infrared background permit us to measure the redshift at which reionization occurs. We have calculated these signals based on numerical simulations of the Lambda-CDM cosmogony.

The paper is withdrawn.

I apply the well controlled Hydro-PM approximation of Gnedin & Hui to model the column density distribution of the Lyman-alpha forest for 25 different flat cosmological scenarios, including variants of the standard CDM, tilted CDM, CDM with a cosmological constant, and CHDM models. I show that within the accuracy of the HPM approximation the slope of the column density distribution reflects the degree of nonlinearity of the cosmic gas distribution and is a function of the rms linear density fluctuation at the characteristic filtering scale only. The amplitude of the column density distribution, expressed as the value for the ionizing intensity, is derived as a function of the cosmological parameters (to about 40% accuracy). The observational data are currently consistent with the value for the ionizing intensity being constant in the redshift interval z~2-4.

I demonstrate by means of high resolution cosmological simulations, which include modelling of a two-phase interstellar medium, that the dominant mechanism for transporting heavy elements from the proto-galaxies into the IGM is the merger mechanism as discovered by Gnedin & Ostriker. Direct ejection of the interstellar gas by supernovae plays only a minor role in transporting metals into the IGM: for a realistic cosmological scenario only a small fraction of all metals in the IGM is delivered by the supernova-driven winds, while most of all metals in the IGM are transported by the merger mechanism. As the result, the metallicity distribution in the IGM is highly inhomogeneous, in agreement with studies of the QSO metal absorption systems, and the predicted metallicity distribution of Lyman-alpha absorbers as a function of their column density is in excellent agreement with the observational data.

We develop a method to calculate the column density distribution of the Lyman-alpha forest for column densities in the range $10^{12.5} - 10^{14.5} cm^{-2}$. The Zel'dovich approximation, with appropriate smoothing, is used to compute the density and peculiar velocity fields. The effect of the latter on absorption profiles is discussed and it is shown to have little effect on the column density distribution. An approximation is introduced in which the column density distribution is related to a statistic of density peaks (involving its height and first and second derivatives along the line of sight) in real space. We show that the slope of the column density distribution is determined by the temperature-density relation as well as the power spectrum on scales $2 h Mpc^{-1} < k < 20 h Mpc^{-1}$. An expression relating the three is given. We find very good agreement between the column density distribution obtained by applying the Voigt-profile-fitting technique to the output of a full hydrodynamic simulation and that obtained using our approximate method for a test model. This formalism then is applied to study a group of CDM as well as CHDM models. We show that the amplitude of the column density distribution depends on the combination of parameters $(\Omega_b h^2)^2 T_0^{-0.7} J_{HI}^{-1}$, which is not well-constrained by independent observations. The slope of the distribution, on the other hand, can be used to distinguish between different models: those with a smaller amplitude and a steeper slope of the power spectrum on small scales give rise to steeper distributions, for the range of column densities we study. Comparison with high resolution Keck data is made.

We introduce an efficient and accurate alternative to full hydrodynamic simulations, Hydro-PM (HPM), for the study of the low column density Lyman-alpha forest. It consists of a Particle-Mesh solver, modified to compute, in addition to the gravitational potential, an effective potential due to the gas pressure. Such an effective potential can be computed from the density field because of a tight correlation between density and pressure in the low density limit, which can be calculated for any photo-reionization history by a method outlined in Hui & Gnedin (1997). Such a correlation exists, in part, because of minimal shock-heating in the low density limit. We compare carefully the density and velocity fields as well as absorption spectra, computed using HPM versus hydrodynamic simulations, and find good agreement. We show that HPM is capable of reproducing measurable quantities, such as the column density distribution, computed from full hydrodynamic simulations, to a precision comparable to that of observations. We discuss how, by virtue of its speed and accuracy, HPM can enable us to use the Lyman-alpha forest as a cosmological probe. We also discuss in detail the smoothing of the gas (or baryon) fluctuation relative to that of the dark matter on small scales due to finite gas pressure. It is shown the conventional wisdom that the linear gas fluctuation is smoothed on the Jeans scale is incorrect for general reionization (or reheating) history; the correct linear filtering scale is in general smaller than the Jeans scale after reheating, but larger prior to it. (abridged)

A possible inconsistency arising when a radiative cooling term is incorporated in a finite resolution self-gravitating hydrodynamic code is discussed. The inconsistency appears when the heating-cooling balance within the cooling and collapsing gas cloud is broken near the resolution limit of a numerical code. As the result, the cooling time of a fluid element increases enormously leading to the unphysical conclusion that the fluid element does not cool and is therefore stable against the collapse. A special cooling consistency condition is introduced which approximately restores the heating-cooling balance and leads to a numerical solution that closely mimics the exact (infinite resolution) solution.

Absorption lines in the Lyman-alpha forest can be thought of as peaks in neutral hydrogen density along lines of sight. The column density distribution (the number density of absorption lines as a function of column density) is then a statistic of density peaks, which contains information about the underlying power spectrum. In particular, we show that the slope of the distribution provides a measure of power on scales smaller than those probed by studies of present-day large scale structure.

We develop an efficient method to study the effects of reionization history on the temperature-density relation of the intergalactic medium in the low density limit (overdensity less than 5). It is applied to the study of photo-reionization models in which the amplitude, spectrum and onset epoch of the ionizing flux, as well as the cosmology, are systematically varied. We find that the mean temperature-density relation at z=2-4 is well approximated by a power-law equation of state for uniform reionization models. We derive analytical expressions for its evolution and exhibit its asymptotic behavior: it is found that for sufficiently early reionization, imprints of reionization history prior to z=10 on the temperature-density relation are washed out. In this limit the temperature at cosmic mean density is proportional to (\Omega_b h/\sqrt\Omega_0)^{1/1.7}. While the amplitude of the radiation flux at the ionizing frequency of HI is found to have a negligible effect on the temperature-density relation as long as the universe reionizes before z=5, the spectrum can change the overall temperature by about 20%, through variations in the abundances of helium species. However the slope of the mean equation of state is found to lie within a narrow range for all reionization models we study, where reionization takes place before z=5. We discuss the implications of these findings for the observational properties of the Lyman-alpha forest. In particular, uncertainties in the temperature of the intergalactic medium, due to the uncertain reionization history of our universe, introduces a 30% scaling in the amplitude of the column density distribution while the the slope of the distribution is only affected by about 5%. Finally, we discuss how a fluctuating ionizing field affects the above results. We argue that under

We simulate a plausible cosmological model in considerable physical and numerical detail through the successive phases of reheating (at 10<z<20), formation of Pop III stars at z=15 (due to molecular hydrogen cooling), with subsequent reionization at z=7. We assume an efficiency of high mass star formation appropriate to leave the universe, after it becomes transparent, with an ionizing background J_21=0.4 (at z=4), near (and perhaps slightly below) the observed value. Since the same stars produce the ionizing radiation and the first generation of heavy elements, a mean metallicity of Z/Z_\sun=1/200 is produced in this early phase, but there is a large variation about this mean, with the high density regions having Z/Z_\sun=1/30 and low density regions essentially no metals. Reionization, when it occurs, is very rapid, which will leave a signature which may be detectable by very large area meter-wavelength radio instruments. Also, the background UV radiation field will show a sharp drop from 1Ryd to 4Ryd due to absorption edges. The simulated volume is too small to form L_* galaxies, but the smaller objects which are found in the simulation obey the Faber-Jackson relation. In order to explore theoretically this domain of "the end of the dark ages" quantitatively, numerical simulations must have a mass resolution of the order of 10^{4.5} M_\sun} in baryons, high spatial resolution (1 kpc) to resolve strong clumping, and allow for detailed and accurate treatment of both the radiation field and atomic/molecular physics.

We develop an analytical method based on the lognormal approximation to compute the column density distribution of the Lyman-alpha forest in the low column density limit. We compute the column density distributions for six different cosmological models and found that the standard, COBE-normalized CDM model cannot fit the observations of the Lyman-alpha forest at z=3. The amplitude of the fluctuations in that model has to be lowered by a factor of almost 3 to match observations. However, the currently viable cosmological models like the lightly tilted COBE-normalized CDM+Lambda model, the CHDM model with 20% neutrinos, and the low-amplitude Standard CDM model are all in agreement with observations, to within the accuracy of our approximation, for the value of the cosmological baryon density at or higher than the old Standard Bing Bang Nucleosynthesis value of 0.0125 for the currently favored value of the ionizing radiation intensity. With the low value for the baryon density inferred by Hogan & Rugers (1996), the models can only marginally match observations.

We note that current observational evidence strongly favors a conventional recombination of ionized matter subsequent to redshift z=1200, followed by reionization prior to redshift z=5 and compute how this would have occurred in a standard scenario for the growth of structure. Extending prior semi-analytic work, we show by direct, high-resolution numerical simulations (of a COBE normalized CDM+Lambda model) that reheating, will occur in the interval 15>z>7, followed by reionization and accompanied by a significant increase in the Jeans mass. However, the evolution of the Jeans mass does not significantly affect star formation in dense, self-shielded clumps of gas, which are detached from the thermal evolution of the rest of the universe. On average, the growth of the Jeans mass tracks the growth of the nonlinear mass scale, a result we suspect is due to nonlinear feedback effects. Cooling on molecular hydrogen leads to a burst of star formation prior to reheating which produces Population III stars with Omega_* reaching 10^{-5.5} and Z/Z_sun reaching 10^{-3.7} by z=14. Star formation subsequently slows down as molecular hydrogen is depleted by photo-destruction and the rise of the temperature. At later times, z<10, when the characteristic virial temperature of gas clumps reach 10,000 degrees, star formation increases again as hydrogen line cooling become efficient. Objects containing Pop III stars accrete mass with time and, as soon as they reach 10,000 K virial temperature, they engage in renewed star formation and turn into normal Pop II objects having an old Pop III metal poor component.

Building a self-gravitating hydrodynamic code as a combination of a hydrodynamic solver and a gravity solver is discussed. We show that straightforward combining those two solvers generally leads to a code that does not conserve energy locally, and instead a special Consistency Condition ought to be satisfied. A particular example of combining Softened Lagrangian Hydrodynamics (SLH) with a P3M gravity solver is used to demonstrate the effect of the Consistency Condition for a self-gravitating hydrodynamic code. The need to supplement the SLH method with the P3M gravity solver arose because the Moving Mesh Gravity solver, used in conjunction with the SLH method previously, was found to produce inaccurated results. We also show that most existing cosmological hydrodynamic codes implicitly satisfy the Consistency Condition.