We present the confusion-limited 1.28 GHz MeerKAT DEEP2 image covering one $\approx 68'$ FWHM primary beam area with $7.6''$ FWHM resolution and $0.55 \pm 0.01$ $\mu$Jy/beam rms noise. Its J2000 center position $\alpha=04^h 13^m 26.4^s$, $\delta=-80^\circ 00' 00''$ was selected to minimize artifacts caused by bright sources. We introduce the new 64-element MeerKAT array and describe commissioning observations to measure the primary beam attenuation pattern, estimate telescope pointing errors, and pinpoint $(u,v)$ coordinate errors caused by offsets in frequency or time. We constructed a 1.4 GHz differential source count by combining a power-law count fit to the DEEP2 confusion $P(D)$ distribution from $0.25$ to $10$ $\mu$Jy with counts of individual DEEP2 sources between $10$ $\mu$Jy and $2.5$ mJy. Most sources fainter than $S \sim 100$ $\mu$Jy are distant star-forming galaxies obeying the FIR/radio correlation, and sources stronger than $0.25$ $\mu$Jy account for $\sim93\%$ of the radio background produced by star-forming galaxies. For the first time, the DEEP2 source count has reached the depth needed to reveal the majority of the star formation history of the universe. A pure luminosity evolution of the 1.4 GHz local luminosity function consistent with the Madau & Dickinson (2014) model for the evolution of star-forming galaxies based on UV and infrared data underpredicts our 1.4 GHz source count in the range $-5 \lesssim \log[S(\mathrm{Jy})] \lesssim -4$.

We examine the properties of the galaxies and dark matter haloes residing in the cluster infall region surrounding the simulated $\Lambda$CDM galaxy cluster studied by Elahi et al. (2016) at z=0. The $1.1\times10^{15}h^{-1}\text{M}_{\odot}$ galaxy cluster has been simulated with eight different hydrodynamical codes containing a variety of hydrodynamic solvers and subgrid schemes. All models completed a dark-matter only, non-radiative and full-physics run from the same initial conditions. The simulations contain dark matter and gas with mass resolution $m_{\text{DM}}=9.01\times 10^8h^{-1}\text{M}_{\odot}$ and $m_{\text{gas}}=1.9\times 10^8h^{-1}\text{M}_{\odot}$ respectively. We find that the synthetic cluster is surrounded by clear filamentary structures that contain ~60% of haloes in the infall region with mass ~$10^{12.5} - 10^{14} h^{-1}\text{M}_{\odot}$, including 2-3 group-sized haloes ($> 10^{13}h^{-1}\text{M}_{\odot}$). However, we find that only ~10% of objects in the infall region are subhaloes residing in haloes, which may suggest that there is not much ongoing preprocessing occurring in the infall region at z=0. By examining the baryonic content contained within the haloes, we also show that the code-to-code scatter in stellar fraction across all halo masses is typically ~2 orders of magnitude between the two most extreme cases, and this is predominantly due to the differences in subgrid schemes and calibration procedures that each model uses. Models that do not include AGN feedback typically produce too high stellar fractions compared to observations by at least ~1 order of magnitude.

Building on the initial results of the nIFTy simulated galaxy cluster comparison, we compare and contrast the impact of baryonic physics with a single massive galaxy cluster, run with 11 state-of-the-art codes, spanning adaptive mesh, moving mesh, classic and modern SPH approaches. For each code represented we have a dark matter only (DM) and non-radiative (NR) version of the cluster, as well as a full physics (FP) version for a subset of the codes. We compare both radial mass and kinematic profiles, as well as global measures of the cluster (e.g. concentration, spin, shape), in the NR and FP runs with that in the DM runs. Our analysis reveals good consistency (<= 20%) between global properties of the cluster predicted by different codes when integrated quantities are measured within the virial radius R200. However, we see larger differences for quantities within R2500, especially in the FP runs. The radial profiles reveal a diversity, especially in the cluster centre, between the NR runs, which can be understood straightforwardly from the division of codes into classic SPH and non-classic SPH (including the modern SPH, adaptive and moving mesh codes); and between the FP runs, which can also be understood broadly from the division of codes into those that include AGN feedback and those that do not. The variation with respect to the median is much larger in the FP runs with different baryonic physics prescriptions than in the NR runs with different hydrodynamics solvers.

We examine subhaloes and galaxies residing in a simulated LCDM galaxy cluster ($M^{\rm crit}_{200}=1.1\times10^{15}M_\odot/h$) produced by hydrodynamical codes ranging from classic Smooth Particle Hydrodynamics (SPH), newer SPH codes, adaptive and moving mesh codes. These codes use subgrid models to capture galaxy formation physics. We compare how well these codes reproduce the same subhaloes/galaxies in gravity only, non-radiative hydrodynamics and full feedback physics runs by looking at the overall subhalo/galaxy distribution and on an individual objects basis. We find the subhalo population is reproduced to within $\lesssim10\%$ for both dark matter only and non-radiative runs, with individual objects showing code-to-code scatter of $\lesssim0.1$ dex, although the gas in non-radiative simulations shows significant scatter. Including feedback physics significantly increases the diversity. Subhalo mass and $V_{max}$ distributions vary by $\approx20\%$. The galaxy populations also show striking code-to-code variations. Although the Tully-Fisher relation is similar in almost all codes, the number of galaxies with $10^{9}M_\odot/h\lesssim M_*\lesssim 10^{12}M_\odot/h$ can differ by a factor of 4. Individual galaxies show code-to-code scatter of $\sim0.5$ dex in stellar mass. Moreover, strong systematic differences exist, with some codes producing galaxies $70\%$ smaller than others. The diversity partially arises from the inclusion/absence of AGN feedback. Our results combined with our companion papers demonstrate that subgrid physics is not just subject to fine-tuning, but the complexity of building galaxies in all environments remains a challenge. We argue even basic galaxy properties, such as the stellar mass to halo mass, should be treated with errors bars of $\sim0.2-0.4$ dex.

We have simulated the formation of a massive galaxy cluster (M$_{200}^{\rm crit}$ = 1.1$\times$10$^{15}h^{-1}M_{\odot}$) in a $\Lambda$CDM universe using 10 different codes (RAMSES, 2 incarnations of AREPO and 7 of GADGET), modeling hydrodynamics with full radiative subgrid physics. These codes include Smoothed-Particle Hydrodynamics (SPH), spanning traditional and advanced SPH schemes, adaptive mesh and moving mesh codes. Our goal is to study the consistency between simulated clusters modeled with different radiative physical implementations - such as cooling, star formation and AGN feedback. We compare images of the cluster at $z=0$, global properties such as mass, and radial profiles of various dynamical and thermodynamical quantities. We find that, with respect to non-radiative simulations, dark matter is more centrally concentrated, the extent not simply depending on the presence/absence of AGN feedback. The scatter in global quantities is substantially higher than for non-radiative runs. Intriguingly, adding radiative physics seems to have washed away the marked code-based differences present in the entropy profile seen for non-radiative simulations in Sembolini et al. (2015): radiative physics + classic SPH can produce entropy cores. Furthermore, the inclusion/absence of AGN feedback is not the dividing line -as in the case of describing the stellar content- for whether a code produces an unrealistic temperature inversion and a falling central entropy profile. However, AGN feedback does strongly affect the overall stellar distribution, limiting the effect of overcooling and reducing sensibly the stellar fraction.

We have simulated the formation of a galaxy cluster in a $\Lambda$CDM universe using twelve different codes modeling only gravity and non-radiative hydrodynamics (\art, \arepo, \hydra\ and 9 incarnations of GADGET). This range of codes includes particle based, moving and fixed mesh codes as well as both Eulerian and Lagrangian fluid schemes. The various GADGET implementations span traditional and advanced smoothed-particle hydrodynamics (SPH) schemes. The goal of this comparison is to assess the reliability of cosmological hydrodynamical 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 $z=0$, global properties such as mass, and radial profiles of various dynamical and thermodynamical quantities. The underlying gravitational framework can be aligned very accurately for all the codes allowing a detailed investigation of the differences that develop due to the various gas physics implementations employed. As expected, the mesh-based codes ART and AREPO form extended entropy cores in the gas with rising central gas temperatures. Those codes employing traditional SPH schemes show falling entropy profiles all the way into the very centre with correspondingly rising density profiles and central temperature inversions. We show that methods with modern SPH schemes that allow entropy mixing span the range between these two extremes and the latest SPH variants produce gas entropy profiles that are essentially indistinguishable from those obtained with grid based methods.

We present results from a numerical code implementing a new method to solve the master equations describing the evolution of linear perturbations in a spherically symmetric but inhomogeneous background. This method can be used to simulate several configurations of physical interest, such as relativistic corrections to structure formation, the lensing of gravitational waves and the evolution of perturbations in a cosmological void model. This paper focuses on the latter problem, i.e. structure formation in a Hubble scale void in the linear regime. This is considerably more complicated than linear perturbations of a homogeneous and isotropic background because the inhomogeneous background leads to coupling between density perturbations and rotational modes of the spacetime geometry, as well as gravitational waves. Previous analyses of this problem ignored this coupling in the hope that the approximation does not affect the overall dynamics of structure formation in such models. We show that for a giga-parsec void, the evolution of the density contrast is well approximated by the previously studied decoupled evolution only for very large-scale modes. However, the evolution of the gravitational potentials within the void is inaccurate at more than the $10\%$ level, and is even worse on small scales.

The goal of this work is to compute the number density of far-IR selected galaxies in the comoving frame and along the past lightcone of observationally constrained LTB models and to compare those results with their standard model counterparts. We derived integral number densities and differential number densities using different cosmological distance definitions in the LTB models. Then, we computed selection functions and consistency functions for the luminosity functions in the combined fields of the Herschel/PACS evolutionary probe (PEP) survey in both standard and void cosmologies, from which we derived the observed values of the above-mentioned densities. We used the Kolmogorov-Smirnov statistics to study both the evolution of the consistency functions and its connection to the evolution of the comoving density of sources. Finally, we fitted the power-law behaviour of the densities along the observer's past lightcone. The analysis of the comoving number density shows that the increased flexibility of the LTB models is not enough to fit the observed redshift evolution of the number counts, if it is specialised to a recent best-fit giant void parametrisation. The results for the power-law fits of the densities along the observer's past lightcone show general agreement across both cosmological models studied here around a slope of -2.5 $\pm$ 0.1 for the integral number density on the luminosity-distance volumes. The differential number densities show much bigger slope discrepancies. We conclude that the differential number densities on the observer's past lightcone were still rendered dependent on the cosmological model by the flux limits of the PEP survey. In addition, we show that an intrinsic evolution of the sources must be assumed to fit the comoving number-density redshift evolution in the giant void parametrisation for the LTB models used in this work.

This is the first paper of a series aiming at investigating galaxy formation and evolution in the giant-void class of the Lemaitre-Tolman-Bondi (LTB) models that best fits current cosmological observations. Here we investigate the Luminosity Function (LF) methodology, and how its estimates would be affected by a change on the cosmological model assumed in its computation. Are the current observational constraints on the allowed Cosmology enough to yield robust LF results? We use the far-infrared source catalogues built on the observations performed with the Herschel/PACS instrument, and selected as part of the PACS evolutionary probe (PEP) survey. Schechter profiles are obtained in redshift bins up to z approximately 4, assuming comoving volumes in both the standard model, that is, Friedmann-Lemaitre-Robertson-Walker metric with a perfect fluid energy-momentum tensor, and non-homogeneous LTB dust models, parametrized to fit the current combination of results stemming from the observations of supernovae Ia, the cosmic microwave background, and baryonic acoustic oscillations. We find that the luminosity functions computed assuming both the standard model and LTB void models show in general good agreement. However, the faint-end slope in the void models shows a significant departure from the standard model up to redshift 0.4. We demonstrate that this result is not artificially caused by the used LF estimator which turns out to be robust under the differences in matter-energy density profiles of the models. The differences found in the LF slopes at the faint end are due to variation in the luminosities of the sources, which depend on the geometrical part of the model. It follows that either the standard model is over-estimating the number density of faint sources or the void models are under-estimating it.

A suggested solution to the dark energy problem is the void model, where accelerated expansion is replaced by Hubble-scale inhomogeneity. In these models, density perturbations grow on a radially inhomogeneous background. This large scale inhomogeneity distorts the spherical Baryon Acoustic Oscillation feature into an ellipsoid which implies that the bump in the galaxy correlation function occurs at different scales in the radial and transverse correlation functions. We compute these for the first time, under the approximation that curvature gradients do not couple the scalar modes to vector and tensor modes. The radial and transverse correlation functions are very different from those of the concordance model, even when the models have the same average BAO scale. This implies that if void models are fine-tuned to satisfy average BAO data, there is enough extra information in the correlation functions to distinguish a void model from the concordance model. We expect these new features to remain when the full perturbation equations are solved, which means that the radial and transverse galaxy correlation functions can be used as a powerful test of the Copernican Principle.

Dark energy observations may be explained within general relativity using an inhomogeneous Hubble-scale depression in the matter density and accompanying curvature, which evolves naturally out of an Einstein-de Sitter (EdS) model. We present a simple parameterization of a void which can reproduce concordance model distances to arbitrary accuracy, but can parameterize away from this to give a smooth density profile everywhere. We show how the Hubble constant is not just a nuisance parameter in inhomogeneous models because it affects the shape of the distance-redshift relation. Independent Hubble-rate data from age estimates can in principle serve to break the degeneracy between concordance and void models, but the data is not yet able to achieve this. Using the latest Constitution supernova dataset we show that robust limits can be placed on the size of a void which is roughly independent of its shape. However, the sharpness of the profile at the origin cannot be well constrained due to supernova being dominated by peculiar velocities in the local Universe. We illustrate our results using some recently proposed diagnostics for the Friedmann models.

The Lemaitre-Tolman-Bondi solution has received much attention as a possible alternative to Dark Energy, as it is able to account for the apparent acceleration of the Universe without any exotic matter content. However, in order to make rigorous comparisons between these models and cosmological observations, such as the integrated Sachs-Wolfe effect, baryon acoustic oscillations and the observed matter power spectrum, it is absolutely necessary to have a proper understanding of the linear perturbation theory about them. Here we present this theory in a fully general, and gauge-invariant form. It is shown that scalar, vector and tensor perturbations interact, and that the natural gauge invariant variables in Lemaitre-Tolman-Bondi cosmology do not correspond straightforwardly to the usual Bardeen variables, in the limit of spatial homogeneity. We therefore construct new variables that reduce to pure scalar, vector and tensor modes in this limit.