Over recent decades cosmological simulations of galaxy formation have been instrumental in advancing our understanding of structure and galaxy formation in the Universe. These simulations follow the nonlinear evolution of galaxies modelling a variety of physical processes over an enormous range of time and length scales. A better understanding of the relevant physical processes improved numerical methods and increased computing power have led to simulations that can reproduce a large number of the observed galaxy properties. Modern simulations model dark matter dark energy and ordinary matter in an expanding spacetime starting from welldefined initial conditions. The modelling of ordinary matter is most challenging due to the large array of physical processes affecting this component. Cosmological simulations have also proven useful to study alternative cosmological models and their impact on the galaxy population. This Technical Review presents a concise overview of the methodology of cosmological simulations of galaxy formation and their different applications.
['Astrophysics - Astrophysics of Galaxies', 'Astrophysics - Cosmology and Nongalactic Astrophysics']
Cosmological simulations of galaxy formation
{'Cosmological Simulations of Galaxy Formation': 'Mark Vogelsberger 1 , Federico Marinacci 2 , Paul Torrey 3 , and Ewald Puchwein 4  \n- 1 Department of Physics, Kavli Institute for Astrophysics and Space Research, Massachusetts Institute of Technology, Cambridge, MA 02139, USA\n- 2 Department of Physics and Astronomy, University of Bologna, via Gobetti 93/2, 40129 Bologna, Italy\n- 3 Department of Astronomy, University of Florida, 316 Bryant Space Sciences Center, Gainesville, FL 32611 USA\n- 4 Leibniz-Institut fur Astrophysik Potsdam (AIP), An der Sternwarte 16, 14482 Potsdam, Germany', 'ABSTRACT': 'Over the last decades, cosmological simulations of galaxy formation have been instrumental for advancing our understanding of structure and galaxy formation in the Universe. These simulations follow the nonlinear evolution of galaxies modeling a variety of physical processes over an enormous range of scales. A better understanding of the physics relevant for shaping galaxies, improved numerical methods, and increased computing power have led to simulations that can reproduce a large number of observed galaxy properties. Modern simulations model dark matter, dark energy, and ordinary matter in an expanding space-time starting from well-defined initial conditions. The modeling of ordinary matter is most challenging due to the large array of physical processes affecting this matter component. Cosmological simulations have also proven useful to study alternative cosmological models and their impact on the galaxy population. This review presents a concise overview of the methodology of cosmological simulations of galaxy formation and their different applications.', '1 Introduction': 'Modern astronomical surveys provide enormous amounts of observational data confronting our theories of structure and galaxy formation. Interpreting these observations demands accurate theoretical predictions. However, galaxy formation is a challenging problem due to its intrinsic multi-scale and multi-physics character. Cosmological computer simulations are, hence, the method of choice for tackling these complexities when studying the properties, growth and evolution of galaxies. These simulations are important to understand the detailed workings of structure and galaxy formation. Dark matter builds the backbone for structure formation and is therefore a key ingredient of these simulations. In addition, dark energy is responsible for the accelerated expansion of the Universe and must also be considered. Despite the fact that the nature of dark matter and dark energy are not known, simulations can make detailed and reliable predictions for these dark components based on their general characteristics. Ordinary matter, e.g. stars and gas, contribute only about five percent to the energy budget of the Universe. Nevertheless, simulating this matter component is essential to study galaxies, but, unfortunately, it is also the most challenging aspect of galaxy formation. Recent simulations follow the formation of individual galaxies and galaxy populations from well-defined initial conditions and yield realistic galaxy properties 1 . Visual representations of the predictions of some of these simulations are shown in Figure 1. At the heart of these simulations are detailed galaxy formation models. Among others, these models describe the cooling of gas, the formation of stars, and the energy and momentum injection caused by supermassive black holes and massive stars 2 . More recent simulations also model the impact of radiation fields, relativistic particles, and magnetic fields leading to a more and more complex description of the galactic ecosystem and the detailed evolution of galaxies in the cosmological context. Galaxy formation simulations have also  \nbecome important for cosmological studies since they can, for example, explore the impact of alternative cosmological models on the galaxy population. Cosmological simulations of galaxy formation therefore provide important insights into a wide range of problems in astrophysics and cosmology. The most important components of cosmological galaxy formation simulations are discussed in this review. A schematic overview of the different ingredients of cosmological simulations is presented in Figure 2.', '2 Cosmological Framework': 'Cosmological simulations of galaxy formation are performed within a cosmological model and start from specific initial conditions. Both of these ingredients are now believed to be known to high precision.', '2.1 Cosmological Model': 'Various observations revealed that our Universe is geometrically flat and dominated by dark matter and dark energy accounting for about ∼ 95% of the energy density. Standard model particles make up for the remaining ∼ 5% and are collectively referred to as baryons. The leading model for structure formation assumes that dark matter is cold, with negligible random motions when decoupled from other matter, and collisionless, so-called cold dark matter, and dark energy is represented by a cosmological constant L , which drives the accelerated expansion of the Universe. This leads to the concordance L CDMmodel, which builds the framework for galaxy formation. Measurements of the cosmic microwave background combined with other observations such as the distance-redshift relation from Type Ia supernovae, abundances of galaxy clusters, and galaxy clustering constrain the fundamental parameters of the L CDMmodel 3 .', '2.2 Initial Conditions': 'Initial conditions for cosmological simulations specify the perturbations imposed on top of a homogeneous expanding background. The background model is generally taken to be a spatially flat FriedmannLemaˆıtre-Robertson-Walker space-time with a defined composition of dark matter, dark energy and baryons. Inflation predicts Gaussian perturbations, where the joint probability distribution of density fluctuations is a multidimensional Gaussian completely specified by its matter power spectrum P ( | k | ) . The post-recombination density field is the linear convolution of the primordial fluctuation field as predicted by inflation with a transfer function T ( k ) 4-7 . Therefore, the power spectrum used to initialize simulations generally takes the form P ( k ) = Ak n | T ( k ) | 2 with n ≈ 1. Once the linear density fluctuation field has been specified at some initial time, typically at redshift z ∼ 100, dark matter particle positions and velocities are assigned along with baryon density, velocity, and temperature fields. The standard approach for dark matter is to displace simulation particles from a uniform Cartesian lattice or glass-like 8, 9 particle configuration using a linear theory approximation 10 or low-order perturbation theory 11-14 . A gravitational glass is made by advancing particles from random positions using the opposite sign of gravity until they freeze in comoving coordinates. Baryon positions and velocities are set in a similar way, and the baryon temperature is often roughly initialized to the redshift-dependent microwave background temperature. Two types of initial conditions are commonly employed: uniformly sampled periodic large volumes or zoom initial conditions, where a low resolution background realization of the density fields surrounds a high resolution region of interest. The computational cost of these zoom simulations increases with the mass of the object that is studied for a given mass resolution. Zoom simulations of dwarf galaxies are therefore computationally less expensive than zoom simulations of large galaxy clusters given the larger number of resolution elements. Some simulations also employ constrained initial conditions to mimic, for example, the local Universe, e.g. nearby dark matter overdensities 15,16 .', 'Generating initial conditions': 'initial positions: x = q + D ( t ) Y ( q )  \ninitial velocities: a ( t ) ˙ x = a ( t ) d D ( t ) d t Y ( q ) = a ( t ) H ( t ) dln D dln a D ( t ) Y ( q )  \nComoving initial positions, x , are assigned based on the unperturbed particle position, q , the linear growth factor, D ( t ) , and the scale factor, a , which is related to the initial redshift, z = 1 / a -1. The curl-free displacement field Y is computed by solving the linearized continuity equation GLYPH<209> · Y = -d / D ( t ) , where d is the relative density fluctuation.', '3 Simulating Dark Matter': "Dark matter builds the backbone for the formation of galaxies, which are expected to form at the centers of dark matter overdensities, so-called halos. The continuum limit of non-interacting dark matter particles is described by the collisionless Boltzmann equation coupled to Poisson's equation. This pair of equations has to be solved in an expanding background Universe dictated by the Friedmann equations, which are derived from the field equations of general relativity. Most cosmological simulations employ Newtonian rather than relativistic gravity, which provides a good approximation since linear structure growth is identical in the matter dominated regime in the two theories, and non-linear large-scale structure induces velocities far below the speed of light. Cosmological simulations are also typically performed with periodic boundary conditions to mimic the large-scale homogeneity and isotropy of the matter distribution of the Universe, i.e. the cosmological principle.", 'Modeling dark matter': "collisionless Boltzmann equation: d f = ¶ f + v ¶ f -¶ F ¶ f = 0  \nd t ¶ t ¶ r ¶ r ¶ v  \nPoisson's equation: GLYPH<209> 2 F = 4 p G ∫ f d v  \nThe collisionless Boltzmann equation describes the evolution of the phase-space density or distribution function of dark matter, f = f ( r , v , t ) , under the influence of the collective gravitational potential, F , given by Poisson's equation. The collisionless Boltzmann equation states the conservation of the local phase-space density; i.e. Liouville's theorem.", '3.1 Numerical Techniques': "The high dimensionality of the collisionless Boltzmann equation prohibits efficient numerical solution methods based on standard discretization techniques for partial differential equations. Therefore, over the past decades, other numerical techniques have been developed to solve this problem more efficiently. An overview of some selected simulation codes and the employed dark matter simulation techniques is presented in Table 1.  \nThe N-Body method: N -body methods are often employed to follow the collisionless dynamics of dark matter, where the phase-space density is sampled by an ensemble of N phase-space points r i , ˙ r i , i = 1 . . . N with masses m i . The conservation of f along the flow implies that the masses m i remain unchanged along each trajectory. N -body methods therefore solve the collisionless Boltzmann equation by the method of characteristics. Alternatively, this method can also be interpreted as a Monte Carlo technique since any initial sample of N phase-space points drawn from the same phase-space density at t = 0 results in an N -body model for the time evolution of f ( r , ˙ r , t ) . The ensemble of all N particles together represents the coarse grained phase-space density 〈 f 〉 ≈ GLYPH<229> i m i f ( r i ( t ) , ˙ r i ( t )) . The latter represents a typical Monte Carlo estimate that can be applied also to other quantities, like the configuration space density. This  \nFigure 1. Visual representations of some selected recent structure and galaxy formation simulations. The simulations are divided in large volume simulations providing statistical samples of galaxies, and zoom simulations resolving smaller scales in more detail. Furthermore, they are also divided in dark matter-only, i.e. N -body, and dark matter plus baryons, i.e. hydrodynamical simulations. Dark matter-only simulations have now converged on a wide range of predictions for the large-scale clustering of dark matter and the dark matter distribution within gravitationally bound dark matter halos. Recent hydrodynamical simulations reproduce galaxy populations that agree remarkably well with observational data. However, many detailed predictions of these simulations are still sensitive to the underlying implementation of baryonic physics.  \n<!-- image -->  \nsampling is subject to Poisson noise, and high particle numbers are therefore desirable to reduce noise in these estimates. To avoid unphysical two-body scatterings between nearby particles, gravitational interactions are softened on small scales so that the particle collection represents a smoothed density field.  \nA variety of kernel-based smoothing techniques are implemented, and some simulations also implement adaptive softening schemes to reduce the softening length in high density regions to reach higher spatial force resolution 17 . The main challenge of N -body simulations is to efficiently calculate the gravitational force that governs the motion of the dark matter sample particles. Once the forces have been calculated, the particles are advanced based on symplectic integration schemes commonly implemented through a Leapfrog integrator. Symplectic integrators exactly solve an approximate Hamiltonian such that the numerical time evolution is a canonical map and preserves certain conserved quantities, such as the total angular momentum, and the phase-space volume. Cosmological simulations are further confronted with a large dynamic range in timescales; i.e. in high-density regions orders of magnitude smaller timesteps are required than in low-density regions. Integration schemes with individual timesteps are therefore typically employed. The time integration is no longer symplectic in a formal sense when individual short-range timesteps are chosen for different particles. Methods to calculate gravitational forces of the N -body system can roughly be divided in two groups: approaches to accelerate the direct summation problem through approximations, or mesh-based methods to calculate the forces. The former approaches aim for efficient numerical solutions of the integral form of Poisson's equation. The latter methods aim for efficient techniques to solve the differential form of Poisson's equation.  \nSolving the integral form of Poisson's equation: The integral form of Poisson's equation, F ( r ) = -G ∫ d r ' r ( r ' ) / | r -r ' | , can be translated to a discrete direct summation problem with complexity O ( N 2 ) . Solving this problem directly results in the so-called particle-particle scheme, and the earliest simulations employed this brute-force approach. The most common method to accelerate the direct summation through approximations is the so-called tree approach 18 . Here, contributions to the gravitational potential from distant particles are approximated by the lowest order terms in a multipole expansion of the mass distribution at a coarse level of the tree reducing the computational cost to O ( N log N ) . The approximation used in the tree method is formally obtained by Taylor expanding the force around some expansion center of the particle group. Often an octree is implemented in cosmological simulations, where each cubic cell is split into up to eight child cells resulting in a tree-like hierarchy of cubic nodes with the root node containing all particles at its bottom. The particles within each of the tree nodes constitute a well-defined and localized group that build the basis for the tree force calculation. A further improvement to O ( N ) complexity is possible through the use of the fast multipole method, where forces are calculated between two tree nodes rather than between individual particles and nodes. This method is best implemented using a tree structure 19 , although the original proposed method was based on a fixed mesh 20 . Implementing periodic boundary conditions for these direct summation-based schemes typically requires Ewald summation techniques 21 originally developed for solid-state physics 22 .  \nSolving the differential form of Poisson's equation: Mesh-based methods aim to solve the differential form of Poisson's equation, GLYPH<209> 2 F ( r ) = 4 p G r ( r ) . This equation can be solved efficiently through fast Fourier transform-based methods, with Poisson's equation in Fourier space k 2 ˜ F ( k ) = -4 p G ˜ r ( k ) , leading to the so-called particle-mesh method 23 . To obtain forces, the potential is then differentiated using a finite difference approximation and the forces are interpolated to the particle positions. The calculation of the gravitational forces via a fast Fourier transform has only a O ( N log N ) complexity, where N is the number of mesh cells. The computational cost does not depend on the details of the particle distribution, and no explicit force softening is necessary for this scheme since the force is automatically softened on the grid scale. Combining the particle-mesh method with a set of nested grids of increasing resolution enables an efficient force solver for inhomogeneous systems resulting in adaptive-mesh-refinement schemes. Multigrid or multilevel methods, which solve the discretized form of Poisson's equation using relaxation methods, such as Gauss-Seidel iterations, are also commonly employed 24 . The advantage of this technique  \nover the fast Fourier transform approach is that the grid does not need to be equidistant, but can be locally adapted according to the particle density. The structure of such an adaptively refined mesh is identical to that of a shallow octree.  \nHybrid schemes: A variety of schemes combine direct summation-based techniques, for short range forces, with Fourier transform-based methods, for long range forces. The most basic example of this is the particle-particle plus particle-mesh method 25 . A common hybrid scheme is the tree-particle-mesh method 26 where the direct summation for short range interactions is approximated by a tree-like method. Combinations of the multigrid method with the fast Fourier transform are also employed, where the Fourier transform is used as a force solver on the coarsest grid 27 . Most modern simulations implement these hybrid solvers to achieve high efficiency.  \nBeyond N-Body method: Conceptually different methods to solve the collisionless Boltzmann equation have also been developed. However, none of these alternatives have so far been widely used for general structure formation simulations. These different methods are motivated, among others, by the desire to resolve the fine-grained structure of the phase-space density and to avoid numerical inaccuracies of the N -body approach like the artificial clumping of simulation particles for dark matter models with a cut-off in the initial power spectrum 28 . Among the methodological alternatives to the N -body method are, for example, a reformulation of the Boltzmann-Poisson system as a Schrodinger equation 29-31 , the waterbag method 32,33 , geodesic deviation equation-based methods 34,35 , Lagrangian tessellation techniques 36 , and direct integration schemes using finite volume approaches based on positive flux conservation methods of plasma physics 37 .  \nTable 1: Major galaxy formation simulation codes  \n- a PM: particle-mesh; TreePM: tree + PM, FM: fast multipole, P 3 M: particle-particle-particle-mesh; ML: multilevel; MG: multigrid  \n- c data-based: data parallelism focuses on distributing data across different nodes, which operate on the data in parallel; task-based: task parallelism focuses on distributing tasks concurrently performed\n- d private: private code; public: publicly available code (in some cases with limited functionality)\n- e gravity solver is based on PKDGRAV3\n- f based on the GADGET-3 code", '3.2 Some Key Results of N -body Simulations': 'The earliest dark matter simulations studied halo population models 50 , the assembly of massive clusters 51 , and the growth of large-scale structure 52,53 . Since then the resolution of these simulations has grown exponentially starting from a few thousand to multi-trillion particle simulations today 54 . Table 2 presents some selected recent structure and galaxy formation simulations. The findings of these simulations can  \nroughly be divided in two categories: the large-scale distribution of dark matter and the structure of dark matter halos. The interaction between baryons and dark matter does affect the structure of dark matter on smaller scales, which is especially important for the internal structure of dark matter halos 55-61 . Studying these phenomena requires, however, simulations that model both dark matter and baryons.  \nThe large-scale distribution of dark matter: Cold dark matter simulations predict that the large-scale distribution of dark matter is not completely homogeneous, but instead exhibits a web-like structure consisting of voids, walls, filaments, and halos quantified through, among others, the halo mass and matter correlation functions.  \nHalo mass function: The halo mass function quantifies the comoving number density of dark matter halos as a function of their virial mass, M vir , typically defined as the mass, M 200 , enclosed within a radius r 200 containing a mean density 200 times the critical density of the Universe. Recently, also other halo boundary definitions like the splashback radius 62,63 , which corresponds to the outermost caustic originally discussed in symmetric analytic halo formation models 64,65 , have been proposed to avoid, for example, the pseudo-evolution of the halo mass and radius 66 . In simulations, dark matter halos are identified through cluster finding methods like the friend-of-friends algorithm 67 and extensions of this based on gravitational unbinding 68 or phase-space structure finding taking into account also velocity space information 69 . In the cold dark matter cosmogony, structure forms through the hierarchical merging of dark matter halos 67 , and the corresponding evolution of the halo mass function has been studied extensively 50,70-77 . Most importantly these studies revealed that the low-mass end of the halo mass function has a power law slope close to -2. Furthermore, the high mass end of the halo mass function is exponentially suppressed. The halo mass function is also an important probe of the nature of dark matter since many particle candidates predict strong, scale-dependent deviations from the expectations of the cold dark matter model 78,79 . The high mass shape and evolution of the halo mass functions also constrains cosmological parameters 80 . Simulation-based empirical halo mass functions are often expressed as M dn / d M = r 0 dln s -1 / d Mf ( s ( M )) , where r 0 is the mean mass density of the Universe, s ( M ) is the variance of the linear density field within a top-hat filter containing mass M , and f ( s ) is a function that is determined empirically by fitting the simulation results. This functional form of the halo mass function is motivated and also predicted by the analytic Press-Schechter model 50 . However, the shape of f ( s ) found in simulations differs significantly from the analytic model originally proposed, agreeing much better with a version motivated by ellipsoidal rather than spherical collapse 81 . The detailed form of f ( s ) depends, among others, on simulation details and halo mass definitions, and a variety of empirical fitting functions have been published 77,82-85 .  \nDark matter distribution: A major success of cold dark matter simulations is their ability to predict the matter distribution on large scales 86 , which is described through the two-point correlation function x ( r ) . For a set of points this function is defined as x ( r ) = 〈 N p 〉 / N m -1, where 〈 N p 〉 indicates the average number of pairs in a thin shell of radius r centered on one point of the set and N m the expected number of pairs in the same shell given a uniform distribution of points. Although this function can be estimated analytically in the linear regime, dark matter simulations are needed to probe its evolution into the non-linear regime. The dark matter correlation function signal grows with time and develops a characteristic shoulder at small scales 87 . This effect can be explained by the relative contribution of the one-halo term, i.e. pairs composed of particles within the same halo, and the two-halo term, i.e. pairs formed by particles in different halos, to the clustering signal 88 . Finally, the dark matter correlation function has a markedly different shape than the galaxy correlation function. The latter has a power law shape over a significant range of scales and an amplitude nearly constant at all redshifts 86 . Indeed, galaxies trace the highest peaks of the dark matter distribution, and their clustering does not change significantly with time, as more and more dark matter  \nstructures grow. This bias needs to be taken properly into account when estimating the large-scale total matter distribution using galaxy tracers 89 .  \nStructure of dark matter halos: Cold dark matter simulations have also established multiple characteristics of the dark matter distribution within collapsed and virialized dark matter halos. This has most importantly led to the discovery of a nearly universal radial density profile of dark matter halos.  \nInternal halo structure: The dark matter mass distribution within halos is well described by a nearuniversal spherically averaged density profile, the so-called Navarro-Frenk-White profile 90,91 : r ( r ) = r s / [( r / r s )( 1 +( r / r s )) 2 ] with a characteristic density r s , and a transition radius r s . This form of density profiles has been shown to arise also in the absence of hierarchical growth like, for example, in hot dark matter models 28 , or models with truncated initial power spectra 92 . The central slope of dark matter halos has been debated for a long time and is also affected by baryonic physics effects that require hydrodynamical simulations. More recent higher resolution simulations found a central slope shallower than -1, indicating that the density profile is better described by a functional form with a gradually changing slope profile 93 : ln ( r ( r ) / r -2 ) = ( -2 / a )[( r / r -2 ) a -1 ] with slope a and transition radius r -2 . This profile had previously been used to fit star counts in the Milky Way 94 , and is known as the Einasto profile. The adjustable shape parameter, a , shows considerable scatter but increases systematically with halo mass at z = 0. The ratio of the virial radius, r vir , and the transition radius, r s , is called the concentration parameter, c , that correlates with the mass of the halo leading to the mass-concentration relation 95-97 ( c GLYPH<181> M -d , d ≈ 0 . 1). Simulations demonstrated that the dependence of halo concentration on mass, initial fluctuation spectrum and cosmological parameters all reflect a dependence of concentration on the actual halo formation time 91 . Specifically, lower mass halos typically assemble earlier, and thus have higher concentration, due to the higher density of the Universe at the time of their formation. The shapes of halos have also been studied, and those depart from sphericity, with halos typically being prolate and increasingly so towards their centers. Major-to-minor axis ratios of two or greater are not uncommon, and more massive halos tend to be less spherical than lower mass halos 98-100 . The exact shapes of dark matter halos also depend on the dark matter particle physics model. Simulations also provide information on the velocity structure of halos. The averaged velocity anisotropy profile, b ( r ) = 1 -0 . 5 s 2 t / s 2 r , grows from zero, i.e. isotropic, to about 0.5, i.e. mild radial anisotropy, towards the outer regions 101,102 . Here, s t denotes the tangential and s r the radial velocity dispersion, with the total velocity dispersion being s 2 = s 2 t + s 2 r . A b value of 1 and b →-¥ correspond to systems where dark matter particles have purely radial and purely circular orbits, respectively. Simulated dark matter halos therefore turn out to be almost isotropic in their inner regions and to be somewhat radially biased at larger radii. Although both r ( r ) and s ( r ) are not close to a power law, the combination f ( r ) = r ( r ) / s 3 ( r ) , also called pseudo-phase-space density, is remarkably close to a power law, with slope ≈-1 . 875 102 . This power law index is identical to that of solutions for self-similar infall onto a point mass from an otherwise uniform Einstein-de Sitter Universe 65 .  \nHalo substructure: As the resolution of dark matter simulations increased, halos within halos, so-called subhalos, could be resolved 92,103 . Subhalos have cuspy, Navarro-Frenk-White-like density profiles but they tend to be less extended than comparable halos in the field due to tidal stripping 101,104 . Bound subhalos with m = M sub / M vir > 10 -7 contain about 10% of the halo mass within the virial radius 104 . Lower mass halos tend to have fewer subhalos and lower subhalo mass fractions at a given m . This shift is due to the difference in the relative dynamical age of halos; e.g. substructure is more effectively destroyed by tides in older, galactic halos compared to more massive galaxy cluster halos. The cumulative subhalo mass function is a power law N ( > m ) GLYPH<181> m -s for m glyph[lessmuch] 1, with s close to the critical value of unity 104,105 . For s = 1 each logarithmic mass bin contributes equally to the total mass in substructure. This  \nis logarithmically divergent as m approaches zero, and implies that a significant fraction of the mass could, in principle, be locked in halos too small to be resolved by the simulations. This can, for example, have important implications for the prediction of dark matter annihilation signals since these small unresolved halos can boost the overall resolved annihilation emission 106 . The abundance of subhalos also varies systematically with other properties of the parent halo, like, for example, the concentration leading to a lower amount of substructure with increasing halo concentration 107 . The radial distribution of subhalos varies only little with the mass or concentration of the parent halo. It is much less centrally concentrated than the overall dark matter profile 104 .  \nFigure 2. Overview of the key ingredients of cosmological simulations. These simulations are performed within a given cosmological framework, and start from specific initial conditions. This framework includes physical models for gravity, dark matter, dark energy, and the type of initial conditions. Two types of simulations are typically performed: either large volume simulations or zoom simulations. The evolution equations of the main matter components, dark matter and gas, are discretized using different techniques and evolved forward in time. The dark matter component follows the equations of collisionless gravitational dynamics that are in most cases solved through the N -body method using different techniques to calculate the gravitational forces. The gas component of baryons is described through the equations of hydrodynamics that are solved, for example, with Lagrangian or Eulerian methods. Various astrophysical processes must also be considered to achieve a realistic galaxy population. Many of these are implemented through effective sub-resolution models.  \n<!-- image -->', '4 Simulating Baryons': 'Dark matter and dark energy dominate the energy budget of the Universe, but the visible components of galaxies consist of baryons. Simulating baryons is therefore crucial to make predictions for the visible Universe. Initially, the baryon component is solely comprised of gas, mostly hydrogen and helium. Some of this gas eventually turns into stars during structure formation. Astrophysical gases in cosmological simulations are typically described as inviscid ideal gases following the Euler equations, which can be expressed in different forms leading to different numerical discretization schemes. Hydrodynamics in cosmological simulations is numerically demanding due to the large dynamic range, highly supersonic flows, and large Reynolds numbers.', 'Modeling cosmic gas': 'Different forms of the hydrodynamical equations. D / d t ≡ ¶ / ¶ t + v · GLYPH<209> denotes the Lagrangian derivative and e = u + v 2 / 2 the total energy per unit mass. The equations are closed through P =( g -1 ) r u with g = 5 / 3. For the arbitrary Lagrangian-Eulerian formulation the grid moves with velocity w and cell volumes evolve as d V / d t = ∫ V ( GLYPH<209> · w ) d V .', '4.1 Numerical Techniques': 'The hydrodynamical equations can be discretized in different ways employing methods that roughly fall into three classes: Lagrangian, Eulerian or arbitrary Lagrange-Eulerian techniques. The Lagrangian specification of the field assumes an observer that follows an individual fluid parcel, with its own properties like density, as it moves through space and time. The Eulerian specification, on the other hand, focuses on specific locations in space through which the fluid flows as time passes. In addition, numerical approaches can also be distinguished between mesh-free and mesh-based algorithms. Mesh-free methods do not require connections between nodes, but are rather based on interactions of each node with its neighbors. An overview of some selected simulation codes and the employed hydrodynamical simulation techniques is shown in Table 1.  \nEulerian Methods: Eulerian methods are the traditional method to solve the system of hyperbolic partial differential equations that constitute ideal hydrodynamics. The most common approaches include finite volume, finite difference, finite element, spectral or wavelet methods. For example, accurate Godunov finite volume schemes offer high-order spatial accuracy, have negligible post-shock oscillations and low numerical diffusivity. For these methods a Riemann problem is solved across cell faces, which yields the required fluxes at each cell face to update the conserved quantities. If the cell is assumed to have uniform properties, this is called a first-order Godunov solver. Modern implementations employ parabolic interpolation, known as the piecewise parabolic method 108,109 . The large dynamic range of cosmological simulations requires adaptive meshes, where the mesh size can be reduced based on some refinement criterion. This leads to the class of adaptive-mesh-refinement schemes 110-113 , which were first developed for solving general problems involving hyperbolic partial differential equations, and then later were also  \napplied to cosmological simulations. Recently also discontinuous Galerkin methods 114-116 became more popular in computational astrophysics since they offer a framework for discretizing hyperbolic problems at any order of spatial accuracy, together with attractive data locality by combining features of spectral element and finite volume methods.  \nLagrangian Methods: Smoothed particle hydrodynamics is a widely used mesh-free Lagrangian technique for approximating the continuum dynamics of fluids through the use of sampling particles, which may also be viewed as interpolation points, following the equations of motion derived from the hydrodynamical equations 117-120 . Energy, linear momentum, angular momentum, mass, and entropy, assuming no artificial viscosity operates, are all simultaneously conserved. The local resolution follows the mass flow, which is convenient to represent large density contrasts. Over the last years various improved formulations of the smoothed particle hydrodynamics method have been developed and applied to cosmological simulations 121-126 . A few cosmological simulations have also employed Lagrangian mesh-based hydrodynamics schemes, which are based on grid deformation techniques 127,128 . However, mesh-tangling effects are a major problem of such multi-dimensional mesh-based Lagrangian hydrodynamics methods.  \nArbitrary Lagrangian-Eulerian Methods: For arbitrary Lagrangian-Eulerian methods, the grid velocity can be freely chosen. For astrophysical applications, such a scheme has recently been realized through a Voronoi tessellation of a set of discrete mesh-generating points, which are allowed to move freely 40 . A finite volume hydrodynamic scheme with the Voronoi cells as control volumes can then be consistently defined. Most importantly, due to the mathematical properties of the Voronoi tessellation, the mesh continuously deforms and changes its topology as a result of the point motion, without ever leading to problematic mesh-tangling effects. Similar methods have over the last years also been implemented in other simulation codes 129,130 . Most recently new types of arbitrary Lagrangian-Eulerian, mesh-free, finite mass and finite volume methods have been successfully applied to astrophysical and galaxy formation problems 45 .', '4.2 Baryonic Physics': 'The hydrodynamical equations have to be complemented by various astrophysical processes that shape the galaxy population. Most of these processes are implemented through effective, so-called sub-resolution models, which are necessary due to the limited numerical resolution of simulations.  \nGas cooling: Gas dissipates its internal energy through cooling processes, like collisional excitation and ionization, inverse Compton, recombination and free-free emission. Cooling processes are coupled to the energy equation using cooling functions that are either tabulated or extracted from chemical networks. Cosmological simulations often assume that the gas is optically thin and in ionization equilibrium and neglect three-body processes that are typically unimportant. In addition to primordial cooling also cooling due to heavy elements, so-called metals, is important. Metal line cooling dominates for temperatures 10 5 glyph[lessorsimilar] T glyph[lessorsimilar] 10 7 K. Early simulations typically employed cooling rates assuming collisional ionization equilibrium 131 , but most later galaxy formation models account for the photo-ionization of metals by the metagalactic radiation field 132 . For most post-reionization simulations this metagalactic radiation field is assumed to be spatially uniform but time-dependent 133 . Simulations that resolve the cold phase of the interstellar medium also include gas cooling below 10 4 K via fine-structure and molecular cooling. In neutral atomic gas, the efficiency of cooling is sensitive to the residual ionization degree. In molecular gas ( n > ∼ 100cm -3 , T < ∼ 50K), the CO molecule dominates the cooling at low densities while at higher densities CI, O 2 and H 2 O start to contribute 134 . Gas cooling is a direct physical process that is not implemented through a sub-resolution model. However, following all cooling processes in detail requires sufficient numerical resolution to resolve the different gas phases.  \nInterstellar medium: Carefully modeling the interstellar medium is important since its properties directly impact star formation. However, simulating the interstellar medium is challenging due to its complex multiphase structure including magnetic fields and relativistic particles. Especially modeling the cold phase is technically difficult because of the short timescales associated with the dense gas. These timescales require very small time-steps to reliably follow the cold gas evolution. Moreover, the implementation of additional physical processes is needed to accurately model such a phase. To circumvent this problem, this dense gas phase is often not directly modeled but rather described by an effective polytropic equation of state 135-137 ; i.e. T GLYPH<181> r g ( r ) , which naturally emerges from an equilibrium two-phase interstellar medium where a hot, supernova-heated and volume-filling phase co-exists with a colder phase containing the bulk of the mass 135 . More recent modeling efforts started to abandon the effective equation of state approach and instead aimed towards resolving the multi-phase structure directly. Such simulations are starting to be able to resolve the Jeans mass of gas, corresponding to the scale of molecular cloud complexes. Therefore, a more direct modeling of the multi-phase interstellar medium is possible 138 . In such simulations, the gas density and temperature distributions follow a multi-modal distribution 138-140 . Generally, the cold gas phase dominates ( ∼ 90%) the gas mass budget, but occupies a very small volume fraction ( ∼ 1%), which is mostly comprised of hotter gas 141 . Simulating the molecular phase of the interstellar medium is challenging because it requires detailed modeling of the interaction between gas, dust, and radiation, which tends to destroy molecules unless gas is able to effectively self-shield from ionizing radiation 142 . Detailed models of the interstellar medium have also to take into account the various feedback sources that ultimately shape the structure of the interstellar medium. Thus, future simulations have to consider how this complex interplay of such a wide range of physical processes affects the properties of the interstellar medium.  \nStar formation: Cold and dense gas eventually forms stars, and simulations therefore transform a portion of this gas into collisionless star particles, representing co-eval, single-metallicity stellar populations described by an underlying initial stellar mass function. Observations support a nearly universal star formation efficiency in molecular gas, where about 1% of the gas is converted into stars per free fall time 143,144 . Based on a calculated star formation rate, the gas is converted into star particles typically using a probabilistic sampling scheme. The star formation rate is usually computed based on a KennicuttSchmidt type relation as dM glyph[star] / dt = e M g / t ff , where M g is the gas cell/particle mass, t ff is the gravitational free fall time and e is a conversion efficiency typically in the range 0 . 01 -1 135,145,146 . However, not all the gas elements are eligible for star formation. Commonly adopted criteria are based on: a density threshold 124,125,135,145-150 , restricting star formation to gravitationally bound regions identified via the virial parameter - that quantifies the degree of pressure support against gravitational collapse 146,151 , Jeans length-based criteria - that is gas must be prone to gravitational instability 146,148,152,153 , restricting star formation to the molecular gas phase 145,146,154-158 or converging flows ( GLYPH<209> · v < 0 147 ). Alternative to the probabilistic sampling scheme, and to better model the clustered nature of star formation, a few simulations also consider star clusters as the unit of star formation by allowing the growth of star particles through accretion from the ambient medium 159 . Once stellar particles have been formed, modern galaxy formation models also track the stellar evolution and mass return of these stars to the gas component. This leads to an enrichment of the gas with metals. Early models tracked only Type II supernova enrichment, but recent models also follow asymptotic giant branch stars 160 , Type Ia supernovae, which are important for iron enrichment 161 , and neutron star mergers for r-process element enrichment 162 . The actual enrichment is based on metal yield models derived from detailed stellar evolution calculations. These yields are however still rather uncertain, at least by a factor of two, particularly at low metallicities and for more  \nmassive stars. This uncertainty then propagates into predictions for metal abundances in simulations. Future cosmological simulations will still have to implement star formation as sub-resolution models with individual stars as their building blocks.  \nStellar feedback: Stars interact with their surrounding gas through the injection of energy and momentum leading to a feedback loop regulating star formation. To regulate star formation, stellar feedback must be efficient in launching galactic-scale outflows to eject gas from galaxies, and a plethora of sub-resolution schemes exists to achieve an efficient generation of galactic winds. Those differ in the way energy and momentum, most notably in the form of supernova explosions, are coupled to the surrounding gas. Essentially the energy can be deposited thermally or kinetically. In the first case, excessive radiative gas cooling must be avoided. While cooling in dense and cold gas is physically expected, at the comparatively low resolution of cosmological simulations it cannot be modeled reliably. The result is then an artificial excessive cooling of the gas, which leads to the unphysical loss of the supernova feedback energy via radiation and greatly reduces its effectiveness. Some approaches therefore disable the radiative cooling of the affected gas for a prescribed amount of time ( ∼ 10 7 yr) 147 , or heat the gas probabilistically to reach high enough temperatures ( T ∼ 10 6 K) for radiative cooling to become ineffective on time scales of ∼ 10 7 yr 137 . In the second case, kinetic energy cannot be radiated away until it thermalizes. However, the use of hydrodynamically-decoupled wind particles, to realize a non-local injection of momentum in the gas surrounding active star forming regions, can still be necessary to obtain large-scale galactic outflows 135,150,161,163 . Recently, more explicit models for stellar feedback have been developed. In addition to supernova feedback they also take into account other feedback channels, such as energy and momentum injection by stellar winds and photoionization and radiation pressure due to radiation emitted by young, massive stars 139,145,146,164,165 . The combination of these processes then leads to a regulation of star formation to the observed low gas to star conversion efficiency of 1% per free-fall time 143,144 . Stellar feedback must be efficient in launching galactic-scale outflows to eject gas from galaxies, thereby also explaining the low baryon retention fraction in galaxies 166,167 . Recent explicit feedback models can make direct predictions for the outflow rates of these outflows 168 , whereas older models typically prescribe the mass loading of these galactic-scale outflows close to the galaxies. Sub-resolution models of stellar feedback vary widely among different galaxy formation models. More work is required to understand in detail which stellar feedback channels are most important for shaping the different types of galaxies.  \nSupermassive black holes: Supermassive black holes are observed in massive galaxies 169,170 , in small, bulge-less disc galaxies 171,172 as well as in dwarf galaxies 173,174 . Simulations therefore include models for supermassive black holes, and numerically seed them typically in dark matter haloes with masses glyph[greaterorsimilar] 10 10 -10 11 M glyph[circledot] since the true seeds cannot be resolved, and their origin is not yet fully understood. They then accrete mass often based on an Eddington-rate-capped Bondi-Hoyle-like accretion rate: ˙ M BH =( 4 p G 2 M 2 BH r ) / ( c 2 s + v 2 rel ) 3 / 2 , where r and c s are the gas density and gas sound speed, respectively, and v rel denotes the relative velocity between the gas and the black hole. Depending on the numerical resolution this accretion rate is sometimes artificially increased, possibly in a density-dependent fashion, to compensate for the inability of simulations to resolve the multi-phase structure of gas 175 . Many simulations also explored variations of the Bondi-Hoyle model to overcome its limitations. The Bondi model, for example, implicitly assumes that the accreting gas has negligible angular momentum, which is most likely unrealistic. Some models therefore assume that black holes might be primarily fed by gas driven to the centers by gravitational torques from non-axisymmetric perturbations 176 , which have more recently been explored in simulations 177-182 . Black holes also grow through mergers, which are modeled in cosmological simulations as well. Due to resolution limitations, general relativistic effects are not taken into account and it is assumed that the black holes of the two galaxies merge instantly once they come  \nclose enough, i.e. within their numerical accretion radius, which is typically calculated based on a local gas resolution element nearest neighbor search.  \nFeedback from active galactic nuclei: Active galactic nuclei are related to observational phenomena associated with accreting supermassive black holes including electromagnetic radiation, relativistic jets, and less-collimated non-relativistic outflows 183 . The resulting energy and momentum couple with the surrounding gas leading to a regulation of black hole growth and star formation in more massive halos ( M > ∼ 10 12 M glyph[circledot] ). This feedback is commonly divided in two modes that are implemented differently in simulations: quasar and radio mode. However, some galaxy formation models do not make this distinction arguing that cosmological simulations lack the resolution to properly distinguish the two feedback modes, and to limit the number of feedback channels to the minimum required to match the observational data 124 . Quasar mode feedback is associated with the radiatively efficient mode of black hole growth and is often implemented through energy or momentum injection assuming that the bolometric luminosity is proportional to the accretion rate, and a fixed fraction of this luminosity is deposited into the neighboring gas 184,185 . Recent works have also implemented momentum-driven winds via radiation pressure on dust 186-188 and via broad-line-region winds 189 . Radio mode feedback is caused by highly-collimated jets of relativistic particles, which are often associated with X-ray bubbles with enough energy to offset cooling losses. Therefore, this feedback mode is assumed to be important for the regulation of star formation in massive galaxies. Radio mode feedback is often implemented as a second sub-resolution feedback channel once the accretion rate is below a critical value 190,191 . Jets themselves cover an enormous dynamic range, being launched at several Schwarzschild radii, and propagating outwards to tens of kpc. Directly resolving them in detail in cosmological simulations is therefore currently not feasible. The sub-resolution models for supermassive black holes are therefore still very uncertain since they have to bridge a very large scale gap between the actual accretion and feedback, and the scales that can be resolved with simulations.  \nMagnetic fields: Magnetic fields permeate the Universe on all scales and impact the motion of ionized gas. Conversely, gas dynamics affects the topology and strength of magnetic fields. Cosmological simulations typically employ the ideal magnetohydrodynamics approach, which is a good approximation for cosmological magnetic fields. This approach assumes that the plasma is perfectly conducting and that relativistic effects, i.e. terms GLYPH<181> ( v / c ) 2 such as the displacement current c -1 ¶ E / ¶ t , are negligible. However, for other situations the ideal magnetohydrodynamics approximation breaks down and non-ideal terms, such as ohmic resistivity, ambipolar diffusion and the Hall effect, must be taken into account. These effects are important especially at very small spatial scales, e.g. for individual star formation, causing a diffusion of the magnetic field. On large cosmological scales the impact of magnetic fields on the dynamics of gas is rather limited 192 . However, magnetic fields are an essential constituent of the interstellar medium, providing both pressure support against gravity 193 and influencing the propagation of cosmic rays 194 . Cosmological simulations including magnetic fields through the ideal magnetohydrodynamics are typically initialized with a certain magnetic seed field, since the approximations and assumptions of ideal magnetohydrodynamics do not permit the self-consistent generation of magnetic fields. Some simulations also consider source terms like the Biermann battery effect or field injection from stellar winds as the source for initial magnetic fields 195 . In most cases, the initial conditions of such cosmological simulations contain however a small seed field of the order of roughly 10 -10 Gauss at a redshift of around z ∼ 100. The simulation results are not sensitive to this seed field as long as its value is not significantly too large, close to violating observational constraints 192,196 , or vanishingly small. The reason for this insensitivity lies in the strong amplification processes that occur during structure formation. This amplification typically occurs in two phases. At high redshifts, a turbulent dynamo leads to an exponential amplification of the magnetic fields in halos. Once the initial turbulent amplification phase has saturated, a second  \nphase of magnetic field amplification starts leading to a linear growth caused by a galactic dynamo 197 . The numerical discretization of the ideal magnetohydrodynamics equations is challenging because of the solenoidal constraint GLYPH<209> · B = 0. Two main families of discretization techniques exist: divergence cleaning schemes and constrained transport. For the cleaning approach, source terms are added to the underlying magnetohydrodynamics equations to correct for divergence errors 198,199 . Constrained transport discretizations 200 , on the other hand, guarantee that the divergence is zero by construction. However, a more complex implementation is required in that case. For instance, either vector potentials 201 , Euler potentials 202,203 or staggered discretizations of the magnetic field components 204-208 must be employed.', 'Ideal magnetohydrodynamics equations:': "MHDMaxwell equations:  \n¶r ¶ t + GLYPH<209> · ( r v ) = 0 ¶r v ¶ t + GLYPH<209> · ( r v ⊗ v + P 1 ) = J × B c ¶ ( r e + e B ) ¶ t + GLYPH<209> · [ ( r e + P ) v + c E × B 4 p ] =  \n0  \nThe evolution of the magnetic field, B , is given by the induction equation, ¶ B / ¶ t = GLYPH<209> × ( v × B ) . Magnetic fields act on gas through the Lorentz force, J × B / c with the current density, J = c GLYPH<209> × B / ( 4 p ) . The energy equation contains the magnetic energy density, e B = || B || 2 / 8 p , and the Poynting vector, c ( E × B / 4 p ) , in the flux part.  \nCosmic rays: Relativistic nuclei and electrons, known as cosmic rays, are another important component of the galactic ecosystem. They are accelerated through diffusive shock acceleration mostly in supernova remnants and jets of active galactic nuclei (first-order Fermi acceleration) and turbulence (second-order Fermi acceleration). Cosmic rays contribute to the pressure in the interstellar medium 209,210 , provide an important heating channel 211,212 , and potentially play a role in driving galactic gas outflows 213-222 due to their shallow equation of state ( P cr GLYPH<181> r 4 / 3 cr ), their long cooling time, and their ability to transfer energy to outflows outside of star-forming discs 223 . The propagation of cosmic rays is dictated by the strength and topology of the underlying magnetic fields.  \nReliably modeling the propagation of cosmic rays therefore requires a detailed modeling of magnetic fields. To capture all these effects self-consistently, the injection, acceleration and the transport of cosmic rays, through anisotropic diffusion and streaming, must be included in simulations. This requires, in principle, a detailed knowledge of the cosmic ray energy spectrum to accurately estimate energy losses and heating rates. The discretization of the cosmic ray transport terms is difficult. For example, anisotropic diffusion requires discretization techniques that avoid the violation of the entropy condition by limiting the transverse fluxes 217,224-226 . Modeling cosmic ray streaming is particularly challenging because of the discontinuous dependence of the streaming velocity on the sign of the scalar product between the magnetic field and the cosmic ray pressure gradient in the one-moment formulation of cosmic ray hydrodynamics. This leads to unphysical oscillations of the solution and small time steps especially near cosmic ray pressure maxima if not addressed in form of regularization techniques - such as replacing the sign function with the hyperbolic tangent function that ensures a smooth dependence of the streaming velocity on cosmic rays and gas properties 219,226,227 , albeit at the expense of a dependence of the solution on a numerical parameter. An elegant solution of this problem is to replace the equation for the cosmic ray energy by two equations for cosmic ray energy and flux that are coupled to the MHD system of equations 228,229 . This two-moment formulation can be derived from quasi-linear theory of cosmic ray transport and describes  \nGLYPH<209> × B = 4 p c J GLYPH<209> · B = 0 1 c ¶ B ¶ t + GLYPH<209> × E = 0 E = -v × B c  \ncosmic ray streaming and diffusion self-consistently with a hyperbolic set of equations, which also contains the evolution equations for Alfv'en waves that are self-generated by the streaming cosmic rays 229 .", 'MHDMaxwell equations:': "¶r ¶ t + GLYPH<209> · ( r v ) = 0 GLYPH<209> × B = 4 p c J GLYPH<209> · B = 0 ¶r v ¶ t + GLYPH<209> · ( r v ⊗ v + P 1 ) = J × B c -GLYPH<209> P cr 1 c ¶ B ¶ t + GLYPH<209> × E = 0 E = -v × B c ¶ ( r e + e B ) ¶ t + GLYPH<209> · [ ( r e + P ) v + c E × B 4 p ] = -( v + v st ) · GLYPH<209> P cr + L th + G th  \nCosmic rays energy density evolution:  \n¶e cr ¶ t + GLYPH<209> · [ e cr ( v + v st ) -ke b ( b · GLYPH<209> e cr )] = -P cr GLYPH<209> · ( v + v st ) + L cr + G cr .  \nCosmic rays exhibit a force on the gas through GLYPH<209> P cr. Their energy density is influenced by streaming with velocity v st ( e cr [ v + v st ] ), anisotropic diffusion with coefficient ke ( ke b [ b · GLYPH<209> e cr ]) , and adiabatic processes due to the compression of the Alfv'en frame ( P cr GLYPH<209> · [ v + v st ] ). The terms L th , L cr, G th and G cr, represent non-adiabatic source and sink terms.  \nRadiation Hydrodynamics: Radiation alters the thermal, kinetic, and chemical state of the gas. Radiation hydrodynamics simulations are required to capture this self-consistently. In the context of cosmological simulations, radiation hydrodynamics simulations have so far primarily been employed to study the epoch of reionization 153,230,231 . These simulations are aimed at exploring the high redshift Universe and are typically not evolved towards the low redshift regime. Consequently, the employed galaxy formation models within these simulations can also not be tested against low redshift predictions. Only a limited number of simulations have studied the impact of radiation in the context of galaxy formation simulations 140,232 . The main reason for this lack of detailed radiation hydrodynamics studies is that numerical radiative transfer is challenging because of the high dimensionality caused by the frequency and directional dependencies of photon propagation. Even more challenging is the fact that in general the speed of light poses severe constraints on the timesteps of these simulations, which can however be circumvented to some degree through the application of a reduced speed of light approximation 233-236 . The most common numerical methods for radiation hydrodynamics are ray-tracing, Monte Carlo, and moment-based methods. The ray-tracing method discretizes the radiative transfer equation along individual directions from each source. Long characteristic ray-tracing schemes 237-239 , cast rays from the source through the whole simulation domain, and the transport, absorption and emission of radiation is computed along each ray. Long characteristic schemes are accurate but computationally expensive, since they scale as O ( N s × N p c ) , where N s is the number of sources, N c is the number of underling discretization resolution elements, for example cells, and p is a method- and geometry-dependent exponent 240,241 . Short characteristic methods 242-245 , on the other hand, solve the radiative transport only along rays that connect nearby cells allowing an efficient parallelization and merging procedures to break the O ( N s × N p c ) scaling. Monte Carlo methods 246-251 , often only applied in post-processing, emit photon packets and propagate them probing the gas opacity, interaction lengths and scattering angles from underlying probability density functions thus stochastically solving the radiative transfer equation. One drawback of Monte Carlo schemes is that the signal-to-noise ratio improves only as the square root of the number of photon packets due to Poisson noise. Still, Monte Carlo is highly accurate and photon weighting, path-based estimators, and discrete diffusion schemes help overcome the efficiency barriers that inhibit convergence 252-255 . Moment-based methods became popular  \nover the last years due to superior scalability 140,256-259 . They are based on a fluid-like description of radiation fields by taking zeroth, first and second moments of the radiation specific intensity with respect to the angular variable. This defines a radiation energy density E n , flux F n , pressure tensor P n and hyperbolic conservation laws for the energy density and the radiation flux. Similar to the hydrodynamical case, where an equation of state is required to relate gas pressure and density, a non-unique closure relation is required to relate P n to E n and F n . A widely used approach is to define P n ≡ E n D , where D is the Eddington tensor that can be estimated with different methods, for example through flux-limited diffusion 117,260 , the optically thin variable Eddington tensor approach 233,261,262 or the M1 closure 140,256,259,263,264 . For the former methods D is estimated assuming that the gas between sources of radiation is always optically thick or thin. The M1 method, instead, computes the Eddington tensor by using local radiation quantities.", 'Radiative transfer equation:': '¶r ¶ t + GLYPH<209> · ( r v ) = 0 ¶ ( r v ) ¶ t + GLYPH<209> · ( r v ⊗ v + P 1 ) = G p ¶ ( r e ) ¶ t + GLYPH<209> · ( r e + P ) v = -L + G E  \n1 c ¶ I n ¶ t + n · ¶ I n ¶ r = -kn I n + j n  \nThe radiative transfer equation relates the specific radiation intensity, I n , with the absorption coefficient, kn , and the specific emissivity, j n . The radiation direction of propagation is represented by the unit vector n . L is the cooling function, G p and G E are source terms that describe the transfer of momentum and energy from the radiation to the gas.  \nOther physics: Additional physical processes considered in some cosmological simulations of galaxy formation are, for example, dust physics 265-273 , thermal conduction 188,225,274-278 , and viscosity 279-282 . Dust has typically been neglected in galaxy formation simulations since it contributes only about ∼ 1% to the mass budget of the interstellar medium. However, dust plays an important role for the evolution of the interstellar medium affecting the thermochemistry and radiation processing. Therefore, recently galaxy formation began to incorporate first simple dust models to follow its production, growth and destruction in the interstellar medium. Most of these implementations treat dust as a passive scalar and model the processes affecting the dust population through effective rate equations. Thermal conduction is another physical effect that is often neglected in cosmological simulations of galaxy formation. However, in hot plasmas of galaxy clusters, conduction can affect the thermodynamic properties of galaxy clusters as has recently been demonstrated 278,283-285 . Simulating thermal conduction requires a precise numerical magnetohydrodynamics implementation to resolve the strength and topology of the magnetic field, and an efficient anisotropic diffusion solver to model the conduction 225 .  \nCaveats and limitations: Simulations of the dark matter component typically boil down to implementing efficient N -body methods and parallelization schemes. Simulations of the baryonic matter component are however more challenging, since they require reliable hydrodynamics numerical schemes and well-posed sub-resolution models. These additional complications lead to some caveats and limitations of such simulations.  \nCalibration: The numerical implementation of baryonic physics is based on sub-resolution models due to the intrinsic resolution limitations of any simulation. These effective models depend on a certain number of adjustable parameters. Depending on the exact galaxy formation model implementation, these parameters can either be chosen based on physical arguments or they require a certain calibration procedure. The latter  \napproach is often employed in large volume simulations, where the sub-resolution models are less detailed compared to those of zoom simulations. The calibration process consists of a parameter exploration for the effective models through a large number of simulations. These simulations typically cover a smaller volume compared to production simulations. The calibration is then based on a comparison to some key observables of the galaxy population like the star formation rate density as a function of cosmic time, the galaxy stellar mass function at z = 0 and the present-day stellar-to-halo mass relation.  \nNumerical convergence: Cosmological simulations have to cover a wide range of spatial and time scales. This implies that simulations have to aim for the highest possible number of resolution elements. However, even state-of-the-art simulations cannot capture all relevant scales. Simulations are therefore often performed at different resolution levels to understand the exact dependence of the results on the number of resolution elements. A simulation prediction is then said to be converged once this prediction does not significantly change anymore if the numerical resolution is further increased.  \nDiverging results: Various simulations now agree on a wide range of predictions. This is especially the case for predictions of the stellar content of galaxies and related observables. However, there is also a wide range of predictions that diverge among different simulations. For example, the characteristics of gas around galaxies are very sensitive to the feedback implementations used in the different galaxy formation models. This can lead to rather different outcomes for the thermodynamic structure of gas around galaxies. Such difference can then be used to differentiate and test galaxy formation models.', '4.3 Some Key Results of Hydrodynamical Simulations': 'The results of hydrodynamical simulations can directly be confronted with observational data providing important tests for galaxy formation models. This often involves the construction of detailed mock observations based on the simulated data 286,287 . Early simulations successfully reproduced properties of the intergalactic medium such as the column density distribution of the Lymana forest 288 . Many simulations also focused on the formation of individual galaxies 289-293 . However, such simulations suffered for a long time from, for example, inconsistent stellar masses, galaxy sizes, star formation histories and galaxy morphologies 163,294-296 . Only recently simulations began to produce realistic galaxies 87,124,125,146,149,297,298 . However, different sub-resolution implementations of astrophysical processes remain a major source of uncertainties. Results of hydrodynamical simulations can be grouped into those for global properties for the whole galaxy population, and those for the properties of individual galaxies.  \nGlobal Properties: Large volume simulations are ideally suited to explore global properties of the galaxy population due to their large statistical sample size. This enables direct comparisons to astronomical galaxy surveys. Table 2 presents some selected recent structure and galaxy formation simulations.  \nStellar content of galaxies: One of the most fundamental properties of the galaxy population is the galaxy stellar mass function, which quantifies the comoving number density of galaxies as a function of galaxy stellar mass. Stellar mass functions are frequently described by a Schechter function 299 with parameter M ∗ , a characteristic mass scale above which the distribution is exponentially suppressed, a normalization f ∗ , and a ∗ setting the low-mass slope. Observed low redshift parameters are roughly given by log ( M ∗ / M glyph[circledot] ) ≈ 11, log ( F ∗ / Mpc -3 ) ≈-2 . 7, a ∗ ≈-1 . 2 300 . However, double Schechter functions provide an even better description of low redshift galaxy stellar mass functions 301-305 . The halo mass function exhibits a steeper low-mass slope, ≈-2, than the galaxy stellar mass function and the exponential suppression occurs at a lower volume density. Reproducing the observed stellar mass function therefore requires a strong suppression of star formation at both the low and high mass ends. Galaxy formation models assume that supernova feedback flattens out the low-mass (M < ∼ 10 12 M glyph[circledot] ) slope by suppressing star formation 306-308 while the suppression of bright and high-mass (M > ∼ 10 12 M glyph[circledot] ) galaxies  \nis regulated by feedback from active galactic nuclei. Energetically plausible forms of supernova and active galactic nuclei feedback in simulations resulted in galaxy stellar mass functions that are consistent with observational data. Simulation predictions are also often confronted with empirical constraints on the relationship between stellar mass and halo mass, which are derived based on various galaxy-halo mapping techniques 166,167 . This ratio of stellar mass to halo mass peaks around halo masses of roughly ∼ 10 12 M glyph[circledot] , where star formation is most efficient. For higher and lower halo masses, the star formation rates are reduced due to feedback processes. Modern large volume simulations reproduce the stellar to halo mass relationship at low and high redshifts reasonably well 124,309 .  \nGas around galaxies: One of the key advantages of hydrodynamical simulations compared to semianalytic models (see Box 1) is their ability to make detailed predictions for the distribution and properties of gas around galaxies including the circumgalactic medium, the intracluster medium, and the intergalactic medium. The circumgalactic and intergalactic media are quite diffuse ( n ∼ 10 -3 -10 -7 cm -3 ) and cool ( T ∼ 10 4 -6 K) and observations in emission, like Lymana and metal lines, are therefore rather challenging. However, absorption line observations from background quasars can probe the distribution, enrichment, and ionization state of this gas. One of the first successes of hydrodynamical simulations has been the reproduction of the declining trend of the number of absorbing clouds per unit redshift and linear interval of HI column density with column density in the Lymana forest 288 . Reproducing properties of the circumgalactic medium, however, is significantly more challenging. Observations of this gas indicate that it features a rich multi-phase structure where individual lines of sight simultaneously contain highly ionized, warm, and cool atomic species 310,311 . The coolest and densest parts of this gas have spatial scales of 10 -100pc 312 , although the coherence scale can reach up to ∼ 1kpc 313 . These spatial scales are below the typical circumgalactic gas resolution limits of galaxy formation simulations. More recently, cosmological simulations with special circumgalactic gas refinement schemes have been employed to overcome some of the resolution limitations. Such simulations increase the numerical resolution in the circumgalactic gas reaching smaller spatial scales 314-317 . At z = 2 such simulations can reach a spatial resolution below ∼ 100pc 316 , and at z = 0 below ∼ 1kpc within the circumgalactic medium 315 . In addition to resolution concerns, the circumgalactic medium is influenced by feedback-driven outflows from galaxies, whose characteristics are not yet properly understood and modeled. The circumgalactic medium can therefore also be used to constrain feedback mechanisms. The intracluster medium can directly be observed via X-ray observations due to the much higher gas temperatures ( T ∼ 10 7 -8 K). Many properties of the intracluster medium , like X-ray and Sunyaev-Zeldovich scaling relations or the iron distribution, can be accurately modeled in simulations 318-321 . However, significant challenges remain for galaxy formation models to reproduce cluster entropy profiles and, in particular, distinct cool-core, and non cool-core clusters 126,322 .  \nGalaxy clustering: Galaxy clustering varies as a function of galaxy mass and galaxy properties, e.g. formation time, star formation rate, color. Simulations now reproduce a number of features in the galaxy clustering signal including the mass dependent two-point correlation length 87 , which increases with increasing masses 323-325 , the clustering signal for non- and star-forming galaxies 87,326 , and the steepening of the power law slope g of the galaxy correlation function with declining redshift ( g ∼ 1 . 8 at z glyph[similarequal] 0 and g ∼ 1 . 6 at z glyph[similarequal] 1 87,327 ).  \nScaling relations: Galaxies exhibit a wide range of scaling relations linking various observables constituting another important test for galaxy formation models. Modern large volume hydrodynamical simulations broadly reproduce many galaxy scaling relations including the mass-size 328 , the supermassive black hole mass-stellar velocity dispersion relation 329 , and the mass-metallicity 330 relation 331-334 . Also other galaxy characteristics like the color of galaxies as a function of galaxy stellar mass can now be reasonably well reproduced by cosmological simulations 298,335,336 . However, there are still points of tension including, for  \nexample, the magnitude of the scatter, the detailed shape, or the dependence on additional galaxy properties.', 'Box 1: Semi-analytic modeling of galaxy formation': "Studying baryonic physics through hydrodynamical simulations is computationally expensive compared to dark matter-only N -body simulations. An alternative approach is to model baryonic physics on top of N -body dark matter simulations through analytic models. This combination of numerical dark matter-only simulations, and analytic models for the prescription of baryonic physics, is known as semi-analytic modeling 337-341 . These semi-analytic models track, for example, how much gas accretes onto halos, how much hot gas cools and turns into stars, or how feedback processes remove cold gas from the galaxy or heat the halo gas. The models are based on the merger history of dark matter halos extracted from N -body simulations. The result of such a calculation is a predicted galaxy population that can be compared to observational data in a similar way as the output of full hydrodynamical simulations. The key advantage of semi-analytical models is their efficiency. It is therefore possible to perform a wide range of calculations, using different model variations. However, a disadvantage of semi-analytic models is that they are less self-consistent compared to hydrodynamical simulations. Furthermore, studying detailed gas properties, for example, the circumgalactic gas with these models is not directly possible since the gas component is not resolved.  \nGalaxy Properties: The detailed properties of late-type disc-like and early-type spheroid-dominated galaxies have been studied extensively using simulations.  \nProperties of late-type galaxies: Simulating the formation of star-forming, late-type galaxies has been one of the most pressing challenges of computational galaxy formation. For a long time, simulations struggled to form galaxies with extended and rotationally-supported stellar and gaseous discs as observed in the Universe. These discs are expected to form through angular momentum conservation of the cooling gas in dark matter halos 342,343 . However, realizing this mechanism in cosmological simulations turned out to be difficult, and early works produced galaxies dominated by a stellar spheroidal component, with a sub-dominant disc only 294,344 . More efficient stellar feedback schemes were required to offset runaway radiative losses of the star-forming gas, the so-called overcooling catastrophe 345 , and to eject the low-angular momentum material responsible for the creation of the dominant stellar bulge , the so-called angular momentum catastrophe 346 . The success of modern simulations in producing late-type disc galaxies is largely due to the ability of stellar feedback to regulate star formation efficiently 125,145,146,164,297,347-352 . More recently, magnetic fields in late-type galaxies have also been studied to understand their topology and field strengths 197,353-357 . Furthermore, the impact of cosmic rays in galaxies has been studied in more detail over the last years 217,219 . These results indicate that cosmic rays are potentially important for driving galactic outflows.  \nProperties of early-type galaxies: Simulations can also reproduce spheroid-dominated early-type systems, which broadly match the early formation history 358 , scaling relations (e.g. the mass and size or velocity dispersion) 359,360 , and the metallicity distribution 361 of observed early-type galaxies. The assembly of such large objects proceeds in two phases 362-365 . At high redshift ( z > ∼ 1 . 5), galaxies grow predominantly in-situ by efficiently converting gas into stars. At later times, mass is predominantly gained through accretion of smaller substructures, i.e. mergers, which also considerably increases galaxy sizes. Spatiallyresolved spectral observations have shown that spheroid-dominated galaxies have diverse kinematics and shapes. The kinematics is usually described through the so-called spin parameter l R . This quantity is used to split galaxies into fast ( l R > 0 . 1) and slow ( l R < 0 . 1) rotator classes. The ellipticity e , instead, is used to define the spheroid's shape. Simulations have played a major role in building a physical picture to explain the diversity in kinematics and shapes of spheroid-dominating galaxies that are produced based on their formation histories 358,366,367 , in particular through gas dissipation, which builds rotationally-supported  \nstructures, and mergers, which set the late-time spin parameter.  \nTable 2: Recent structure and galaxy formation simulations  \n- c for particle based codes, the minimum softening length is reported; for mesh codes, the minimum cell size is quoted (dark matter/gas)", '5 Simulations of Alternative Cosmological Models': 'Cosmological simulations of galaxy formation have also been used to explore alternative cosmological models. At the most basic level the cosmological model can be altered in three different ways: alternative forms of dark matter, alternative forms of dark energy, or alternative forms of gravity. We note that many simulations of alternative cosmological models typically only consider the dark matter component and do not model baryons. However, these simulations then neglect the important backreaction between baryons and dark matter. Similarly, simulations including baryons are also now important to infer cosmological parameters. For example, DESI, LSST and Euclid will rely on models based on galaxy formation simulations to achieve their forecasted precision. Future explorations of alternative cosmologies have to consider and include these effects by also modeling the baryon component.', 'Box 2: Small-scale problems of cold dark matter': 'The cold dark matter paradigm correctly describes the large-scale distribution of galaxies. On sub-galactic scales however, some problems have been identified over the last decades 392 . Among the most relevant challenges are: the under-abundance of dwarf galaxies in the Milky Way and in the field (the missing satellites problem 92,393-396 ), the inconsistency of inner dark matter density profiles in low surface brightness and dwarf galaxies (the cuspcore problem 397,398 ), the deficit of dark matter in the inner regions of massive dwarf galaxies (the too-big-to-fail problem 399,400 ), and the large variety of shapes of dwarf rotation curves (the diversity problem 401 ). Most of these problems have been found by contrasting dark matter-only simulations with observations, which do not take into account the complex baryonic dark matter interactions. It is therefore possible that these challenges can be solved through the proper modeling of baryonic physics. For instance, the existence of dark matter cores can potentially be explained by the gravitational transfer of energy from supernovae into the orbits of dark matter particles 402-406 . Alternatively, these discrepancies between observations and cold dark matter simulations can also be explored through alternative dark matter models. These small-scale problems have therefore generated significant interest in the exploration of alternative dark matter scenarios.', '5.1 Alternative Forms of Dark Matter': 'A wide range of alternative dark matter models have been proposed over the last decades. However, not all of these models have been studied in detail through simulations. Mostly three main classes of alternative dark matter models have been simulated: warm dark matter, self-interacting dark matter, and fuzzy dark matter. Many of these models have been invoked to address small-scale problems of the cold dark matter paradigm (see Box 2).  \nWarm Dark Matter: Cold dark matter models exhibit a highk cut-off in the initial power spectrum due to free-streaming or collisional damping. For a canonical weakly interactive massive particle this cut-off is of the order of 1 comoving parsec corresponding to a mass scale of 10 -6 M glyph[circledot] 407 . Warm dark matter models, on the other hand, have an effective free-streaming length l fs that scales inversely with particle mass 408 . For recent cosmologies 3 , this relation is approximately l fs = 33 ( m WDM / 1keV ) -1 . 11 kpc and the corresponding free-streaming mass is M fs = 2 × 10 7 ( m WDM / 1keV ) -3 . 33 M glyph[circledot] . The reduction of small-scale power within warm dark matter models has two consequences: first, a reduction of low mass halos, and second a reduction of the central density of halos. Simulations of warm dark matter models are typically carried out with the same numerical methods as cold dark matter simulations, but with modified initial conditions. However, the power spectrum cut-off leads to artificial and numerical discreteness effects in N -body simulations 28 . Special care is then required to avoid a contamination of results in that case. More recently alternative methods based on phase-space tessellation techniques have been employed  \nto study warm dark matter models avoiding these numerical artifacts 78 .  \nSelf-Interacting Dark Matter: Dark matter models that involve dark matter self-interactions 409,410 have also been explored extensively. Self-interactions are commonly quantified in terms of the cross section per unit particle mass, s / m . Models with constant and velocity-dependent cross sections have both been studied with simulations 411 . The high central dark matter densities observed in clusters exclude self-interacting dark matter models with s / m glyph[greaterorsimilar] 0 . 5cm 2 / g for these cluster mass scales. Recently, more general self-interacting dark matter models have been suggested. Those have both truncated power spectra and self-interactions 412,413 . Such models affect the internal structure of dark matter halos through the scattering of particles that cause the formation of density cores. On the other hand, the truncated power spectra also lead, similar to warm dark matter models, to a suppression of halo substructure. Various recent simulations have demonstrated that models with s / m ≈ 0 . 5 -10cm 2 / g produce dark matter cores in dwarf galaxies with sizes ∼ 0 . 3 -1 . 5kpc and central densities 2 -0 . 2 × 10 8 M glyph[circledot] / kpc 3 = 7 . 4 -0 . 74GeV / cm 3 that can alleviate some cold dark matter small-scale problems 414-417 . Simulations of self-interacting dark matter are based on the N -body approach coupled to a local Monte Carlo-based probabilistic scattering scheme to model particle self-interactions.  \nFuzzy Dark Matter: An ultralight bosonic scalar field is a completely different alternative to the cold dark matter paradigm 418 , where a bosonic fluid with a particle mass of m ∼ 10 -22 eV suppresses small-scale structure owing to macroscopic quantum properties 419-421 with a typical de Broglie wavelength of l DB ∼ 1kpc 422,423 . The dark matter fluid forms in this case a cosmological Bose-Einstein condensate 424-426 . Such an ultralight scalar field of spin-0 at zero temperature is described in the non-relativistic limit by the Schrodinger-Poisson equations 419,420,427,428 : i ¯ h ¶y / ¶ t = -¯ h 2 / 2 m GLYPH<209> 2 y + mV y and GLYPH<209> 2 V = 4 p G ( r -¯ r ) , where r = | y | 2 is the fluid density, ¯ r is the mean density, and V is the potential. One consequence of the macroscopic quantum behavior of the fluid is that the fluid admits stable, minimum-energy soliton configurations forming at the centers of self-gravitating halos. These kpc-scale soliton cores offer one possible solution to the cusp-core problem of cold dark matter. Numerically, the Schrodinger-Poisson equations can, for example, be solved through adaptive spectral methods or through a reformulation into a hydrodynamics problem, that can be solved with hydrodynamical discretization techniques, based on the Madelung formulation 429,430 .', '5.2 Alternative Forms of Dark Energy': 'Cosmological simulations must include at least a cosmological constant to account for the accelerated expansion of the Universe. A wide range of alternative dark energy models have, however, been considered in the literature 431 and a number of these have also been studied with simulations 432 .  \nDynamical dark energy: The most simple extension in the dark energy sector is to assume a dark energy density that is time dependent but still spatially homogeneous - at least on sub-horizon scales. This behavior can, for example, be obtained in scalar field models of dark energy 432 . Cosmic structure growth is then only affected via an altered background expansion. The only change required to perform cosmological simulations of such models is then to modify the calculation of the Hubble expansion rate in the numerical integration 433,434 . As the growth function is different than in L CDM, extra care is also required when choosing the amplitude of matter density fluctuations in the initial conditions, i.e. taking into account at what redshift observational constraints on the amount of fluctuations are aimed to be matched. For example, models with a higher dark energy density at early times suppress structure growth and have hence a lower amplitude of fluctuations at redshift zero for the same scalar amplitude in the cosmic microwave background 433,434 . Dynamical dark energy can have a surprisingly large impact on galaxy properties in simulations 435 . In practice, this results in degeneracies between cosmology and the feedback physics that  \nis required to match observations.  \nInhomogeneous dark energy: Models of dark energy that exhibit sizable spatial fluctuations within the horizon represent the next level of complexity. For such models, and even more so for the coupled dark energy models, a clear distinction between dark energy and modified gravity is often not possible as accelerations arising from spatial fluctuations in the dark energy field can also be interpreted as modifications to the laws of gravity. Relatively little simulation work has been done on models in which inhomogeneous dark energy interacts with matter only gravitationally, such as, for example, in the clustering dark energy scenario 436 .  \nCoupled dark energy: In the hope to alleviate the puzzle of the similar energy density of matter and dark energy at the present cosmic epoch, additional non-gravitational couplings between these sectors have been proposed 437 . Such a coupling of dark energy to matter could either be universal, i.e. involving all matter species, or non-universal, with dark energy, for example, coupling only to dark matter but not to baryons. Models with a universal coupling typically require a screening mechanism that hides its effects in dense environments like the solar system, where experimental tests of gravity tightly constrain a direct coupling to baryons. In contrast, models with a coupling only to dark matter are observationally much less constrained. In both cases, growing perturbations in the matter density field can naturally give rise to corresponding fluctuations in the coupled dark energy field. Coupled dark energy scenarios have been widely studied with simulations, either avoiding 438,439 or including 440,441 a treatment of the spatial fluctuations of dark energy. In the former case, the main effects of coupling terms are a time dependence of the gravitating particle mass of the coupled matter species, as well as a velocity dependent friction term. Accounting for the spatial fluctuation additionally results in a fifth force proportional to the gradient of the dark energy field. These effects have, for example, been found to lower the concentrations and baryon fraction of halos 439 , thereby reducing potential tensions compared to a L CDMcosmology.', '5.3 Alternative Forms of Gravity': "While general relativity has been tested to high precision within the solar system, constraints on galactic and intergalactic scales are much weaker. Indeed, additional components that have so far not been directly observed, dark matter and dark energy, need to be added to allow a viable description of cosmology by general relativity. As an alternative, modifications of the laws of gravity have been proposed, which could make at least one of these components obsolete.  \nModified gravity as an alternative to dark matter: Dark matter models successfully explain observations on many different scales, including the cosmic microwave background, the Lymana forest, the clustering of galaxies, and the internal dynamics of galaxies and galaxy clusters. Most work aimed at replacing the role of dark matter by a modification of the laws of gravity has focused only on a subset of these areas. For example, modified Newtonian dynamics 442,443 , a change in Newton's second law at small acceleration values ( glyph[vector] F = m m ( | glyph[vector] a | / a 0 ) glyph[vector] a , with a 0 ∼ 10 -10 ms -2 and m ( x ) → 1 for x glyph[greatermuch] 1 and m ( x ) → x for x glyph[lessmuch] 1), or alternatively a change in Poisson's equation of Newtonian gravity ( glyph[vector] GLYPH<209> · ( m ( | glyph[vector] a | / a 0 ) glyph[vector] a ) = 4 pr ), has been proposed to account for the flat rotation curves of galaxies at large radii. Since this modified Poisson's equation is non-linear, gravity algorithms that are based on the principle of linear force superposition such as direct summation, tree and Fourier transform-based schemes are not suitable to simulate these kinds of models. Simulations have therefore been performed with the multigrid method with the full approximation scheme 24 . The non-linear partial differential equation is then discretized on a grid with a finite difference representation of the differential operator and iteratively solved using Gauss-Seidel relaxation. Since modified Newtonian dynamics is not a relativistic theory, relativistic extensions of it have also been proposed, for example, tensor-vector-scalar gravity (TeVeS) 444,445 . Here gravity is mediated by a tensor  \n(metric), vector and scalar field. However, these models have not been widely studied in full cosmological simulations yet. Models without a dark matter component such as modified Newtonian dynamics also naturally account for the tight observed relation between the gravitational acceleration inferred from galaxy rotation curves and that expected from the observed baryonic mass 446,447 . However, galaxy formation simulations within the L CDMframework can also produce a sufficiently tight relation 448-450 .  \nModified gravity as an alternative to dark energy: Dark energy has only been observed through its impact on the background expansion of the Universe. Replacing dark energy with a modification to the laws of gravity is, compared to replacing dark matter, easier. In fact, a cosmological constant in the Einstein field equations can also be interpreted as modified gravity rather than an unexpectedly small zero-point energy of a quantum field. In the literature, a wide range of much more sophisticated modified gravity theories have been considered. While many of these can account for the observed accelerated expansion of the Universe, Occam's razor would typically disfavor them compared to a cosmological constant in the absence of observational evidence beyond the observed background expansion. Cosmological simulations have been widely used to investigate the observational signatures of such extended gravity models to guide observational searches for potential modifications of gravity over a wide range of scales. Many modified gravity models that exhibit interesting behavior on, for example, galactic and intergalactic scales have been designed such that they approach general relativity in dense environments such as the solar system to avoid violating experimental constraints. Such screening mechanisms typically involve non-linear partial differential equations, which renders them numerically challenging and requires tailored techniques 451 . Most schemes resort to the multigrid method with the full approximation scheme 24 , e.g. employed on an adaptively refining mesh 452-454 . With such methods cosmological simulations have been carried out for a number of screened modified gravity models, including Chameleonf ( R ) , DGP, symmetron, dilaton, and Galileon gravity 455-457 . Most such studies focused on collisionless simulations. Semi-analytical galaxy formation models combined with Chameleonf ( R ) gravity demonstrated that the gravity modification effects on basic properties such as galaxy stellar mass functions and cosmic star formation rate densities are rather small and comparable to the uncertainties of the semi-analytical models 458 . Clustering signals and relative velocities of halo pairs can, however, change notably 458,459 . Post-processing L CDMgalaxy formation simulations with a modified gravity solver suggests that there should be characteristic changes in the internal kinematics of galaxies such as features in their rotation curves near the screening threshold 460 , which can also result in degeneracies with the core/cusp problem 461 . Fully self-consistent simulation studies of galaxy formation in such screened modified gravity models have only started very recently 462 . Such simulations should in principle also take into account effects that modified gravity has on stellar physics 463 .", '6 Conclusions and Outlook': 'Cosmological simulations of galaxy formation play a crucial role for our understanding of galaxy formation. Especially, the last years have seen enormous progress on two fronts: large volume simulations modeling large samples of galaxies, and zoom simulations with refined galaxy formation models that resolve the physical processes in more detail. Modern galaxy formation simulations reproduce now a plethora of observational results and create virtual universes that are to first order nearly identical to the real Universe. At the same time, these simulations are also employed to explore alternative cosmological models with modifications to the nature of dark matter, dark energy and gravity. This progress in the field of galaxy formation simulations has mostly been driven by a better understanding of galaxy formation physics, refined numerical methods, and the steady growth of computing power.  \nCosmological simulations of galaxy formation use a variety of different numerical methods, and different implementations of galaxy formation physics. Despite these differences, such simulations have now converged on a wide range of predictions for the evolution of galaxies. It therefore seems that the basic physical mechanisms that shape the galaxy population have been identified, and that their current modeling is sufficient to produce realistic galaxy populations. However, these physical processes are implemented through still rather crude sub-resolution models. Sub-resolution models aim to capture the relevant physics through an effective description. In fact, cosmological simulations will always have to rely on these sub-resolution models since truly ab initio cosmological simulations of galaxy formation are and will remain impossible. One danger associated with the application of sub-resolution models is the belief that the reproduction of large amounts of observational data automatically implies a correct and physically plausible effective model and therefore detailed understanding of galaxy formation. This is problematic since sub-resolution models contain per construction a certain number of adjustable and degenerate parameters, and at the same time do not really capture the detailed physics at play but only provide an effective coarse description. Caution is therefore required to not over-interpret some of the recent successes generated by these models.  \nOne of the next goals of computational galaxy formation is to understand which detailed physical processes drive the outcomes of effective physical models. For example, many simulations employ rather crude and incomplete models for the generation of galactic outflows without a detailed modeling of the driving process. Future simulations should aim at understanding these processes in more detail to illuminate the true physical processes at work going beyond the crude effective models to gain more fundamental insights. This will also lead to a better understanding of what physics actually drives the overall behavior of currently existing coarse-grained effective sub-resolution models. While constructing new models and simulations, it is important to keep in mind that the major goal of simulations is not primarily to fit observed data, but rather to gain insights into galaxy formation physics. Advances in this direction benefit often more from failures of certain ideas or models, rather than a perfect reproduction of observational data that is to some degree subject to the calibration of free model parameters and the coarse-grained nature of the employed models. Another frontier of cosmological galaxy formation simulations is the desire to provide large volume simulations that match the statistical sample sizes of upcoming large observational surveys. This requires very large volume simulations with well-understood sub-resolution models. The development and better understanding of refined sub-resolution models, the desire to achieve higher numerical resolution, and simulations with larger volumes represent the main challenges of cosmological simulations of the next decade.', 'Acknowledgements': "We thank David Barnes, Mike Boylan-Kolchin, Lars Hernquist, Rahul Kannan, Hui Li, Stephanie O'Neil, Rudiger Pakmor, Christoph Pfrommer, Laura Sales, Aaron Smith, Volker Springel, and Rainer Weinberger for useful comments. We also thank the reviewers for helpful feedback. MV acknowledges support through an MIT RSC award, a Kavli Research Investment Fund, NASA ATP grant NNX17AG29G, and NSF grants AST-1814053 and AST-1814259. FM acknowledges support through the Program 'Rita Levi Montalcini' of the Italian MIUR.", 'References': "- 1. Somerville, R. S. & Dav'e, R. Physical Models of Galaxy Formation in a Cosmological Framework. ARA&A 53 , 51-113 (2015). DOI 10.1146/annurev-astro-082812-140951. 1412.2712 .\n- 2. Naab, T. & Ostriker, J. P. Theoretical Challenges in Galaxy Formation. ARA&A 55 , 59-109 (2017). DOI 10.1146/annurev-astro-081913-040019. 1612.06891 .\n- 3. Planck Collaboration et al. Planck 2015 results. XIII. Cosmological parameters. A&A 594 , A13 (2016). DOI 10.1051/0004-6361/201525830. 1502.01589 .\n- 4. Seljak, U. & Zaldarriaga, M. A Line-of-Sight Integration Approach to Cosmic Microwave Background Anisotropies. ApJ 469 , 437 (1996). DOI 10.1086/177793. astro-ph/ 9603033 .\n- 5. Peacock, J. A. The evolution of galaxy clustering. MNRAS 284 , 885-898 (1997). DOI 10.1093/mnras/284.4.885. astro-ph/9608151 .\n- 6. Eisenstein, D. J. & Hu, W. Baryonic Features in the Matter Transfer Function. ApJ 496 , 605-614 (1998). DOI 10.1086/305424. astro-ph/9709112 .\n- 7. Eisenstein, D. J. & Hu, W. Power Spectra for Cold Dark Matter and Its Variants. ApJ 511 , 5-15 (1999). DOI 10.1086/306640. astro-ph/9710252 .\n- 8. Baugh, C. M., Gaztanaga, E. & Efstathiou, G. A comparison of the evolution of density fields in perturbation theory and numerical simulations - II. Counts-in-cells analysis. MNRAS 274 , 1049-1070 (1995). DOI 10.1093/mnras/274.4.1049. astro-ph/9408057 .\n- 9. White, S. D. M. Cosmology and large scale structure. In Schaeffer, R., Silk, J., Spiro, M. & Zinn-Justin, J. (eds.) Cosmology and Large Scale Structure , 349-430 (1996).\n- 10. Zel'dovich, Y. B. Gravitational instability: An approximate theory for large density perturbations. A&A 5 , 84-89 (1970).\n- 11. Bertschinger, E. Multiscale Gaussian Random Fields and Their Application to Cosmological Simulations. ApJS 137 , 1-20 (2001). DOI 10.1086/322526. astro-ph/0103301 .\n- 12. Jenkins, A. Second-order Lagrangian perturbation theory initial conditions for resimulations. MNRAS 403 , 1859-1872 (2010). DOI 10.1111/j.1365-2966.2010.16259.x. 0910.0258 .\n- 13. Hahn, O. & Abel, T. Multi-scale initial conditions for cosmological simulations. MNRAS 415 , 2101-2121 (2011). DOI 10.1111/j.1365-2966.2011.18820.x. 1103.6031 .\n- 14. Garrison, L. H., Eisenstein, D. J., Ferrer, D., Metchnik, M. V. & Pinto, P. A. Improving initial conditions for cosmological N-body simulations. MNRAS 461 , 4125-4145 (2016). DOI 10.1093/mnras/stw1594. 1605.02333 .\n- 15. Hoffman, Y. & Ribak, E. Constrained realizations of Gaussian fields - A simple algorithm. ApJ 380 , L5-L8 (1991). DOI 10.1086/186160.\n- 16. Salmon, J. Generation of Correlated and Constrained Gaussian Stochastic Processes for NBody Simulations. ApJ 460 , 59 (1996). DOI 10.1086/176952.\n- 17. Price, D. J. & Monaghan, J. J. An energy-conserving formalism for adaptive gravitational force softening in smoothed particle hydrodynamics and N-body codes. MNRAS 374 , 13471358 (2007). DOI 10.1111/j.1365-2966.2006.11241.x. astro-ph/0610872 .\n- 18. Barnes, J. & Hut, P. A hierarchical O(N log N) force-calculation algorithm. Nat. 324 , 446-449 (1986). DOI 10.1038/324446a0.\n- 19. Dehnen, W. A Very Fast and Momentum-conserving Tree Code. ApJ 536 , L39-L42 (2000). DOI 10.1086/312724. astro-ph/0003209 .\n- 20. Greengard, L. & Rokhlin, V. A fast algorithm for particle simulations. J. Comput. Phys. 73 , 325-348 (1987). DOI 10.1016/0021-9991(87)90140-9.\n- 21. Hernquist, L., Bouchet, F. R. & Suto, Y. Application of the Ewald method to cosmological N-body simulations. ApJS 75 , 231-240 (1991). DOI 10.1086/191530.\n- 22. Ewald, P. P. Die Berechnung optischer und elektrostatischer Gitterpotentiale. Annalen der Physik 369 , 253-287 (1921). DOI 10.1002/andp.19213690304.\n- 23. Hockney, R. W. & Eastwood, J. W. Computer Simulation Using Particles (New York: McGrawHill, 1981).\n- 24. Brandt, A. Multi-level adaptive solutions to boundary-value problems. Math. Comp. 31 , 333390 (1977). DOI 10.1090/S0025-5718-1977-0431719-X.\n- 25. Efstathiou, G., Davis, M., White, S. D. M. & Frenk, C. S. Numerical techniques for large cosmological N-body simulations. ApJS 57 , 241-260 (1985). DOI 10.1086/191003.\n- 26. Bode, P. & Ostriker, J. P. Tree Particle-Mesh: An Adaptive, Efficient, and Parallel Code for Collisionless Cosmological Simulation. ApJS 145 , 1-13 (2003). DOI 10.1086/345538. astro-ph/0302065 .\n- 27. Kravtsov, A. V., Klypin, A. A. & Khokhlov, A. M. Adaptive Refinement Tree: A New High-Resolution N-Body Code for Cosmological Simulations. ApJS 111 , 73-94 (1997). DOI 10.1086/313015. astro-ph/9701195 .\n- 28. Wang, J. & White, S. D. M. Discreteness effects in simulations of hot/warm dark matter. MNRAS 380 , 93-103 (2007). DOI 10.1111/j.1365-2966.2007.12053.x. astro-ph/0702575 .\n- 29. Widrow, L. M. & Kaiser, N. Using the Schroedinger Equation to Simulate Collisionless Matter. ApJ 416 , L71 (1993). DOI 10.1086/187073.\n- 30. Schaller, M., Becker, C., Ruchayskiy, O., Boyarsky, A. & Shaposhnikov, M. A new framework for numerical simulations of structure formation. MNRAS 442 , 3073-3095 (2014). DOI 10.1093/mnras/stu1069. 1310.5102 .\n- 31. Uhlemann, C., Kopp, M. & Haugg, T. Schrodinger method as N-body double and UV completion of dust. Phys. Rev. D 90 , 023517 (2014). DOI 10.1103/PhysRevD.90.023517. 1403.5567 .\n- 32. Colombi, S. & Touma, J. Vlasov-Poisson in 1D: waterbags. MNRAS 441 , 2414-2432 (2014). DOI 10.1093/mnras/stu739.\n- 33. Colombi, S. Vlasov-Poisson in 1D for initially cold systems: post-collapse Lagrangian perturbation theory. MNRAS 446 , 2902-2920 (2015). DOI 10.1093/mnras/stu2308. 1411.4165 .\n- 34. Vogelsberger, M., White, S. D. M., Helmi, A. & Springel, V. The fine-grained phase-space structure of cold dark matter haloes. MNRAS 385 , 236-254 (2008). DOI 10.1111/j.13652966.2007.12746.x. 0711.1105 .\n- 35. Vogelsberger, M. & White, S. D. M. Streams and caustics: the fine-grained structure of L cold dark matter haloes. MNRAS 413 , 1419-1438 (2011). DOI 10.1111/j.1365-2966.2011.18224.x. 1002.3162 .\n- 36. Hahn, O., Abel, T. & Kaehler, R. A new approach to simulating collisionless dark matter fluids. MNRAS 434 , 1171-1191 (2013). DOI 10.1093/mnras/stt1061. 1210.6652 .\n- 37. Yoshikawa, K., Yoshida, N. & Umemura, M. Direct Integration of the Collisionless Boltzmann Equation in Six-dimensional Phase Space: Self-gravitating Systems. ApJ 762 , 116 (2013). DOI 10.1088/0004-637X/762/2/116. 1206.6152 .\n- 38. Teyssier, R. Cosmological hydrodynamics with adaptive mesh refinement. A new high resolution code called RAMSES. A&A 385 , 337-364 (2002). DOI 10.1051/0004-6361:20011817. arXiv:astro-ph/0111367 .\n- 39. Springel, V. The cosmological simulation code GADGET-2. MNRAS 364 , 1105-1134 (2005). DOI 10.1111/j.1365-2966.2005.09655.x. arXiv:astro-ph/0505010 .\n- 40. Springel, V. E pur si muove: Galilean-invariant cosmological hydrodynamical simulations on a moving mesh. MNRAS 401 , 791-851 (2010). DOI 10.1111/j.1365-2966.2009.15715.x. 0901.4107 .\n- 41. Bryan, G. L. et al. ENZO: An Adaptive Mesh Refinement Code for Astrophysics. ApJS 211 , 19 (2014). DOI 10.1088/0067-0049/211/2/19. 1307.2265 .\n- 42. Jetley, P., Kale, L. V., Gioachin, F., Quinn, T. & Mendes, C. Massively parallel cosmological simulations with changa. In 2008 IEEE International Parallel & Distributed Processing Symposium , 1-12 (IEEE Computer Society, Los Alamitos, CA, USA, 2008). DOI 10.1109/IPDPS.2008.4536319.\n- 43. Gioachin, F., Kal'e, L. V., Quinn, T. R., Jetley, P. & Wesolowski, L. Scaling hierarchical n-body simulations on gpu clusters. In SC Conference , 1-11 (IEEE Computer Society, Los Alamitos, CA, USA, 2010). DOI 10.1109/SC.2010.49.\n- 44. Menon, H. et al. Adaptive techniques for clustered N-body cosmological simulations. Comput. Astrophys. Cosmol. 2 , 1 (2015). DOI 10.1186/s40668-015-0007-9. 1409.1929 .\n- 45. Hopkins, P. F. A new class of accurate, mesh-free hydrodynamic simulation methods. MNRAS 450 , 53-110 (2015). DOI 10.1093/mnras/stv195. 1409.7395 .\n- 46. Habib, S. et al. HACC: Simulating sky surveys on state-of-the-art supercomputing architectures. New A 42 , 49-65 (2016). DOI 10.1016/j.newast.2015.06.003. 1410.2805 .\n- 47. Potter, D., Stadel, J. & Teyssier, R. PKDGRAV3: beyond trillion particle cosmological simulations for the next era of galaxy surveys. Comput. Astrophys. Cosmol. 4 , 2 (2017). DOI 10.1186/s40668-017-0021-1. 1609.08621 .\n- 48. Wadsley, J. W., Keller, B. W. & Quinn, T. R. Gasoline2: a modern smoothed particle hydrodynamics code. MNRAS 471 , 2357-2369 (2017). DOI 10.1093/mnras/stx1643. 1707.03824 .\n- 49. Schaller, M., Gonnet, P., Chalk, A. B. G. & Draper, P. W. SWIFT: SPH With Inter-dependent Fine-grained Tasking (2018). 1805.020 .\n- 50. Press, W. H. & Schechter, P. Formation of Galaxies and Clusters of Galaxies by Self-Similar Gravitational Condensation. ApJ 187 , 425-438 (1974). DOI 10.1086/152650.\n- 51. White, S. D. M. The dynamics of rich clusters of galaxies. MNRAS 177 , 717-733 (1976). DOI 10.1093/mnras/177.3.717.\n- 52. Aarseth, S. J., Gott, J. R., III & Turner, E. L. N-body simulations of galaxy clustering. I - Initial conditions and galaxy collapse times. ApJ 228 , 664-683 (1979). DOI 10.1086/156892.\n- 53. Efstathiou, G. The clustering of galaxies and its dependence upon Omega. MNRAS 187 , 117127 (1979). DOI 10.1093/mnras/187.2.117.\n- 54. Skillman, S. W. et al. Dark Sky Simulations: Early Data Release. arXiv e-prints (2014). 1407.2600 .\n- 55. Gnedin, O. Y., Kravtsov, A. V., Klypin, A. A. & Nagai, D. Response of Dark Matter Halos to Condensation of Baryons: Cosmological Simulations and Improved Adiabatic Contraction Model. ApJ 616 , 16-26 (2004). DOI 10.1086/424914. astro-ph/0406247 .\n- 56. Di Cintio, A. et al. The dependence of dark matter profiles on the stellar-to-halo mass ratio: a prediction for cusps versus cores. MNRAS 437 , 415-423 (2014). DOI 10.1093/mnras/stt1891. 1306.0898 .\n- 57. Zhu, Q. et al. Baryonic impact on the dark matter distribution in Milky Way-sized galaxies and their satellites. MNRAS 458 , 1559-1580 (2016). DOI 10.1093/mnras/stw374.\n- 58. Ben'ıtez-Llambay, A., Frenk, C. S., Ludlow, A. D. & Navarro, J. F. Baryon-induced dark matter cores in the EAGLE simulations. MNRAS 488 , 2387-2404 (2019). DOI 10.1093/mnras/stz1890. 1810.04186 .\n- 59. Bose, S. et al. No cores in dark matter-dominated dwarf galaxies with bursty star formation histories. MNRAS 486 , 4790-4804 (2019). DOI 10.1093/mnras/stz1168. 1810.03635 .\n- 60. Katz, H. et al. Stellar feedback and the energy budget of late-type Galaxies: missing baryons and core creation. MNRAS 480 , 4287-4301 (2018). DOI 10.1093/mnras/sty2129.\n- 61. Read, J. I., Walker, M. G. & Steger, P. Dark matter heats up in dwarf galaxies. MNRAS (2019). DOI 10.1093/mnras/sty3404. 1808.06634 .\n- 62. Diemer, B. & Kravtsov, A. V. Dependence of the Outer Density Profiles of Halos on Their Mass Accretion Rate. ApJ 789 , 1 (2014). DOI 10.1088/0004-637X/789/1/1. 1401.1216 .  \n- 63. More, S., Diemer, B. & Kravtsov, A. V. The Splashback Radius as a Physical Halo Boundary and the Growth of Halo Mass. ApJ 810 , 36 (2015). DOI 10.1088/0004-637X/810/1/36. 1504.05591 .\n- 64. Fillmore, J. A. & Goldreich, P. Self-similar gravitational collapse in an expanding universe. ApJ 281 , 1-8 (1984). DOI 10.1086/162070.\n- 65. Bertschinger, E. Self-similar secondary infall and accretion in an Einstein-de Sitter universe. ApJS 58 , 39-65 (1985). DOI 10.1086/191028.\n- 66. Diemer, B., More, S. & Kravtsov, A. V. The Pseudo-evolution of Halo Mass. ApJ 766 , 25 (2013). DOI 10.1088/0004-637X/766/1/25. 1207.0816 .\n- 67. Davis, M., Efstathiou, G., Frenk, C. S. & White, S. D. M. The evolution of large-scale structure in a universe dominated by cold dark matter. ApJ 292 , 371-394 (1985). DOI 10.1086/163168.\n- 68. Springel, V., White, S. D. M., Tormen, G. & Kauffmann, G. Populating a cluster of galaxies - I. Results at z=0. MNRAS 328 , 726-750 (2001). DOI 10.1046/j.1365-8711.2001.04912.x. astro-ph/0012055 .\n- 69. Behroozi, P. S., Wechsler, R. H. & Wu, H.-Y. The ROCKSTAR Phase-space Temporal Halo Finder and the Velocity Offsets of Cluster Cores. ApJ 762 , 109 (2013). DOI 10.1088/0004637X/762/2/109. 1110.4372 .\n- 70. Bond, J. R., Cole, S., Efstathiou, G. & Kaiser, N. Excursion set mass functions for hierarchical Gaussian fluctuations. ApJ 379 , 440-460 (1991). DOI 10.1086/170520.\n- 71. Jenkins, A. et al. The mass function of dark matter haloes. MNRAS 321 , 372-384 (2001). DOI 10.1046/j.1365-8711.2001.04029.x. astro-ph/0005260 .\n- 72. Sheth, R. K. & Tormen, G. An excursion set model of hierarchical clustering: ellipsoidal collapse and the moving barrier. MNRAS 329 , 61-75 (2002). DOI 10.1046/j.13658711.2002.04950.x. astro-ph/0105113 .\n- 73. White, M. The Mass Function. ApJS 143 , 241-255 (2002). DOI 10.1086/342752. astro-ph/0207185 .\n- 74. Reed, D. et al. Evolution of the mass function of dark matter haloes. MNRAS 346 , 565-572 (2003). DOI 10.1046/j.1365-2966.2003.07113.x. astro-ph/0301270 .\n- 75. Warren, M. S., Abazajian, K., Holz, D. E. & Teodoro, L. Precision Determination of the Mass Function of Dark Matter Halos. ApJ 646 , 881-885 (2006). DOI 10.1086/504962. astro-ph/0506395 .\n- 76. Reed, D. S., Governato, F., Quinn, T., Stadel, J. & Lake, G. The age dependence of galaxy clustering. MNRAS 378 , 777-784 (2007). DOI 10.1111/j.1365-2966.2007.11826.x. astro-ph/0602003 .\n- 77. Tinker, J. et al. Toward a Halo Mass Function for Precision Cosmology: The Limits of Universality. ApJ 688 , 709-728 (2008). DOI 10.1086/591439. 0803.2706 .\n- 78. Angulo, R. E., Hahn, O. & Abel, T. The warm dark matter halo mass function below the cut-off scale. MNRAS 434 , 3337-3347 (2013). DOI 10.1093/mnras/stt1246. 1304.2406 .\n- 79. Schneider, A. Structure formation with suppressed small-scale perturbations. MNRAS 451 , 3117-3130 (2015). DOI 10.1093/mnras/stv1169. 1412.2133 .\n- 80. Eke, V. R., Cole, S. & Frenk, C. S. Cluster evolution as a diagnostic for Omega. MNRAS 282 (1996). DOI 10.1093/mnras/282.1.263. astro-ph/9601088 .\n- 81. Sheth, R. K., Mo, H. J. & Tormen, G. Ellipsoidal collapse and an improved model for the number and spatial distribution of dark matter haloes. MNRAS 323 , 1-12 (2001). DOI 10.1046/j.1365-8711.2001.04006.x. astro-ph/9907024 .\n- 82. Angulo, R. E. et al. Scaling relations for galaxy clusters in the Millennium-XXL simulation. MNRAS 426 , 2046-2062 (2012). DOI 10.1111/j.1365-2966.2012.21830.x. 1203.3216 .\n- 83. Watson, W. A. et al. The halo mass function through the cosmic ages. MNRAS 433 , 1230-1245 (2013). DOI 10.1093/mnras/stt791. 1212.0095 .\n- 84. Heitmann, K. et al. The Q Continuum Simulation: Harnessing the Power of GPU Accelerated Supercomputers. ApJS 219 , 34 (2015). DOI 10.1088/0067-0049/219/2/34. 1411.3396 .\n- 85. Bocquet, S., Saro, A., Dolag, K. & Mohr, J. J. Halo mass function: baryon impact, fitting formulae, and implications for cluster cosmology. MNRAS 456 , 2361-2373 (2016). DOI 10.1093/mnras/stv2657. 1502.07357 .\n- 86. Springel, V., Frenk, C. S. & White, S. D. M. The large-scale structure of the Universe. Nat. 440 , 1137-1144 (2006). DOI 10.1038/nature04805. arXiv:astro-ph/0604561 .\n- 87. Springel, V. et al. First results from the IllustrisTNG simulations: matter and galaxy clustering. MNRAS 475 , 676-698 (2018). DOI 10.1093/mnras/stx3304. 1707.03397 .\n- 88. Cooray, A. & Sheth, R. Halo models of large scale structure. Phys. Reports 372 , 1-129 (2002). DOI 10.1016/S0370-1573(02)00276-4. astro-ph/0206508 .\n- 89. Benson, A. J., Cole, S., Frenk, C. S., Baugh, C. M. & Lacey, C. G. The nature of galaxy bias and clustering. MNRAS 311 , 793-808 (2000). DOI 10.1046/j.1365-8711.2000.03101.x. astro-ph/9903343 .\n- 90. Navarro, J. F., Frenk, C. S. & White, S. D. M. The Structure of Cold Dark Matter Halos. ApJ 462 , 563 (1996). DOI 10.1086/177173. astro-ph/9508025 .\n- 91. Navarro, J. F., Frenk, C. S. & White, S. D. M. A Universal Density Profile from Hierarchical Clustering. ApJ 490 , 493-508 (1997). DOI 10.1086/304888. astro-ph/9611107 .\n- 92. Moore, B. et al. Dark Matter Substructure within Galactic Halos. ApJ 524 , L19-L22 (1999). DOI 10.1086/312287. astro-ph/9907411 .\n- 93. Navarro, J. F. et al. The inner structure of L CDM haloes - III. Universality and asymptotic slopes. MNRAS 349 , 1039-1051 (2004). DOI 10.1111/j.1365-2966.2004.07586.x. astro-ph/0311231 .\n- 94. Einasto, J. On the Construction of a Composite Model for the Galaxy and on the Determination of the System of Galactic Parameters. Trudy Astrofizicheskogo Instituta Alma-Ata 5 , 87-100 (1965).\n- 95. Ludlow, A. D. et al. The mass-concentration-redshift relation of cold dark matter haloes. MNRAS 441 , 378-388 (2014). DOI 10.1093/mnras/stu483. 1312.0945 .\n- 96. Ludlow, A. D. et al. The mass-concentration-redshift relation of cold and warm dark matter haloes. MNRAS 460 , 1214-1232 (2016). DOI 10.1093/mnras/stw1046. 1601.02624 .\n- 97. Klypin, A., Yepes, G., Gottlober, S., Prada, F. & Heß, S. MultiDark simulations: the story of dark matter halo concentrations and density profiles. MNRAS 457 , 4340-4359 (2016). DOI 10.1093/mnras/stw248. 1411.4001 .\n- 98. Jing, Y. P. Intrinsic correlation of halo ellipticity and its implications for large-scale weak lensing surveys. MNRAS 335 , L89-L93 (2002). DOI 10.1046/j.1365-8711.2002.05899.x. astro-ph/0206098 .\n- 99. Allgood, B. et al. The shape of dark matter haloes: dependence on mass, redshift, radius and formation. MNRAS 367 , 1781-1796 (2006). DOI 10.1111/j.1365-2966.2006.10094.x. astro-ph/0508497 .\n- 100. Bett, P. et al. The spin and shape of dark matter haloes in the Millennium simulation of a L cold dark matter universe. MNRAS 376 , 215-232 (2007). DOI 10.1111/j.1365-2966.2007.11432.x. astro-ph/0608607 .\n- 101. Diemand, J., Kuhlen, M. & Madau, P. Formation and Evolution of Galaxy Dark Matter Halos and Their Substructure. ApJ 667 , 859-877 (2007). DOI 10.1086/520573. astro-ph/ 0703337 .\n- 102. Navarro, J. F. et al. The diversity and similarity of simulated cold dark matter haloes. MNRAS 402 , 21-34 (2010). DOI 10.1111/j.1365-2966.2009.15878.x. 0810.1522 .\n- 103. Ghigna, S. et al. Dark matter haloes within clusters. MNRAS 300 , 146-162 (1998). DOI 10.1046/j.1365-8711.1998.01918.x. astro-ph/9801192 .\n- 104. Springel, V. et al. The Aquarius Project: the subhaloes of galactic haloes. MNRAS 391 , 16851711 (2008). DOI 10.1111/j.1365-2966.2008.14066.x. 0809.0898 .\n- 105. Gao, L. et al. The Phoenix Project: the dark side of rich Galaxy clusters. MNRAS 425 , 21692186 (2012). DOI 10.1111/j.1365-2966.2012.21564.x. 1201.1940 .\n- 106. Springel, V. et al. Prospects for detecting supersymmetric dark matter in the Galactic halo. Nat. 456 , 73-76 (2008). DOI 10.1038/nature07411. 0809.0894 .\n- 107. Gao, L., White, S. D. M., Jenkins, A., Stoehr, F. & Springel, V. The subhalo populations of L CDM dark haloes. MNRAS 355 , 819-834 (2004). DOI 10.1111/j.1365-2966.2004.08360.x. astro-ph/0404589 .\n- 108. Colella, P. & Woodward, P. R. The Piecewise Parabolic Method (PPM) for Gas-Dynamical Simulations. J. Comput. Phys. 54 , 174-201 (1984). DOI 10.1016/0021-9991(84)90143-8.\n- 109. Woodward, P. R. PPM: Piecewise-Parabolic Methods for Astrophysical Fluid Dynamics. In Winkler, K.-H. A. & Norman, M. L. (eds.) NATO Advanced Science Institutes (ASI) Series C , vol. 188 of NATO Advanced Science Institutes (ASI) Series C , 245 (1986).\n- 110. Berger, M. J. & Oliger, J. Adaptive Mesh Refinement for Hyperbolic Partial Differential Equations. J. Comput. Phys. 53 , 484-512 (1984). DOI 10.1016/0021-9991(84)90073-1.\n- 111. Berger, M. J. & Colella, P. Local adaptive mesh refinement for shock hydrodynamics. J. Comput. Phys. 82 , 64-84 (1989). DOI 10.1016/0021-9991(89)90035-1.\n- 112. Klein, R. I., McKee, C. F. & Colella, P. On the hydrodynamic interaction of shock waves with interstellar clouds. 1: Nonradiative shocks in small clouds. ApJ 420 , 213-236 (1994). DOI 10.1086/173554.\n- 113. Bryan, G. L. & Norman, M. L. Simulating X-ray Clusters with Adaptive Mesh Refinement. In American Astronomical Society Meeting Abstracts , vol. 27 of Bulletin of the American Astronomical Society , 1421 (1995).\n- 114. Cockburn, B., Lin, S.-Y. & Shu, C.-W. TVB Runge Kutta Local Projection Discontinuous Galerkin Finite Element Method for Conservation Laws III: One-Dimensional Systems. J. Comput. Phys. 84 , 90-113 (1989). DOI 10.1016/0021-9991(89)90183-6.\n- 115. Mocz, P., Vogelsberger, M., Sijacki, D., Pakmor, R. & Hernquist, L. A discontinuous Galerkin method for solving the fluid and magnetohydrodynamic equations in astrophysical simulations. MNRAS 437 , 397-414 (2014). DOI 10.1093/mnras/stt1890. 1305.5536 .\n- 116. Guillet, T., Pakmor, R., Springel, V., Chandrashekar, P. & Klingenberg, C. High-order Magnetohydrodynamics for Astrophysics with an Adaptive Mesh Refinement Discontinuous Galerkin Scheme. MNRAS (2019). DOI 10.1093/mnras/stz314. 1806.02343 .\n- 117. Lucy, L. B. A numerical approach to the testing of the fission hypothesis. AJ 82 , 1013-1024 (1977). DOI 10.1086/112164.\n- 118. Gingold, R. A. & Monaghan, J. J. Smoothed particle hydrodynamics - Theory and application to non-spherical stars. MNRAS 181 , 375-389 (1977). DOI 10.1093/mnras/181.3.375.\n- 119. Springel, V. Smoothed Particle Hydrodynamics in Astrophysics. ARA&A 48 , 391-430 (2010). DOI 10.1146/annurev-astro-081309-130914. 1109.2219 .\n- 120. Price, D. J. Smoothed particle hydrodynamics and magnetohydrodynamics. J. Comput. Phys. 231 , 759-794 (2012). DOI 10.1016/j.jcp.2010.12.011. 1012.1885 .\n- 121. Read, J. I. & Hayfield, T. SPHS: smoothed particle hydrodynamics with a higher order dissipation switch. MNRAS 422 , 3037-3055 (2012). DOI 10.1111/j.1365-2966.2012.20819.x. 1111.6985 .\n- 122. Hopkins, P. F. A general class of Lagrangian smoothed particle hydrodynamics methods and implications for fluid mixing problems. MNRAS 428 , 2840-2856 (2013). DOI 10.1093/mnras/sts210. 1206.5006 .\n- 123. Keller, B. W., Wadsley, J., Benincasa, S. M. & Couchman, H. M. P. A superbubble feedback model for galaxy simulations. MNRAS 442 , 3013-3025 (2014). DOI 10.1093/mnras/stu1058. 1405.2625 .  \n- 124. Schaye, J. et al. The EAGLE project: simulating the evolution and assembly of galaxies and their environments. MNRAS 446 , 521-554 (2015). DOI 10.1093/mnras/stu2058. 1407.7040 .\n- 125. Wang, L. et al. NIHAO project - I. Reproducing the inefficiency of galaxy formation across cosmic time with a large sample of cosmological hydrodynamical simulations. MNRAS 454 , 83-94 (2015). DOI 10.1093/mnras/stv1937. 1503.04818 .\n- 126. Barnes, D. J. et al. The Cluster-EAGLE project: global properties of simulated clusters with resolved galaxies. MNRAS 471 , 1088-1106 (2017). DOI 10.1093/mnras/stx1647. 1703.10907 .\n- 127. Gnedin, N. Y. Softened Lagrangian hydrodynamics for cosmology. ApJS 97 , 231-257 (1995). DOI 10.1086/192141.\n- 128. Pen, U.-L. A High-Resolution Adaptive Moving Mesh Hydrodynamic Algorithm. ApJS 115 , 19-34 (1998). DOI 10.1086/313074. astro-ph/9704258 .\n- 129. Duffell, P. C. & MacFadyen, A. I. TESS: A Relativistic Hydrodynamics Code on a Moving Voronoi Mesh. ApJS 197 , 15 (2011). DOI 10.1088/0067-0049/197/2/15. 1104.3562 .\n- 130. Vandenbroucke, B. & De Rijcke, S. The moving mesh code SHADOWFAX. Astron. Comput. 16 , 109-130 (2016). DOI 10.1016/j.ascom.2016.05.001. 1605.03576 .\n- 131. Sutherland, R. S. & Dopita, M. A. Cooling functions for low-density astrophysical plasmas. ApJS 88 , 253-327 (1993). DOI 10.1086/191823.\n- 132. Wiersma, R. P. C., Schaye, J. & Smith, B. D. The effect of photoionization on the cooling rates of enriched, astrophysical plasmas. MNRAS 393 , 99-107 (2009). DOI 10.1111/j.13652966.2008.14191.x. 0807.3748 .\n- 133. Haardt, F. & Madau, P. Radiative Transfer in a Clumpy Universe. IV. New Synthesis Models of the Cosmic UV/X-Ray Background. ApJ 746 , 125 (2012). DOI 10.1088/0004-637X/746/2/125. 1105.2039 .\n- 134. Tielens, A. G. G. M. The Physics and Chemistry of the Interstellar Medium (Cambridge, UK: Cambridge University Press, 2010).\n- 135. Springel, V. & Hernquist, L. Cosmological smoothed particle hydrodynamics simulations: a hybrid multiphase model for star formation. MNRAS 339 , 289-311 (2003). DOI 10.1046/j.13658711.2003.06206.x. arXiv:astro-ph/0206393 .\n- 136. Agertz, O., Teyssier, R. & Moore, B. The formation of disc galaxies in a L CDM universe. MNRAS 410 , 1391-1408 (2011). DOI 10.1111/j.1365-2966.2010.17530.x. 1004.0005 .\n- 137. Dalla Vecchia, C. & Schaye, J. Simulating galactic outflows with thermal supernova feedback. MNRAS 426 , 140-158 (2012). DOI 10.1111/j.1365-2966.2012.21704.x. 1203.5667 .\n- 138. Hopkins, P. F., Quataert, E. & Murray, N. The structure of the interstellar medium of starforming galaxies. MNRAS 421 , 3488-3521 (2012). DOI 10.1111/j.1365-2966.2012.20578.x. 1110.4636 .\n- 139. Agertz, O., Kravtsov, A. V., Leitner, S. N. & Gnedin, N. Y. Toward a Complete Accounting of Energy and Momentum from Stellar Feedback in Galaxy Formation Simulations. ApJ 770 , 25 (2013). DOI 10.1088/0004-637X/770/1/25. 1210.4957 .\n- 140. Rosdahl, J., Schaye, J., Teyssier, R. & Agertz, O. Galaxies that shine: radiationhydrodynamical simulations of disc galaxies. MNRAS 451 , 34-58 (2015). DOI 10.1093/mnras/stv937. 1501.04632 .\n- 141. Emerick, A., Bryan, G. L. & Mac Low, M.-M. Simulating an isolated dwarf galaxy with multichannel feedback and chemical yields from individual stars. MNRAS 482 , 1304-1329 (2019). DOI 10.1093/mnras/sty2689. 1807.07182 .\n- 142. Krumholz, M. R. & Gnedin, N. Y. A Comparison of Methods for Determining the Molecular Content of Model Galaxies. ApJ 729 , 36 (2011). DOI 10.1088/0004-637X/729/1/36. 1011.4065 .\n- 143. Bigiel, F. et al. A Constant Molecular Gas Depletion Time in Nearby Disk Galaxies. ApJ 730 , L13 (2011). DOI 10.1088/2041-8205/730/2/L13. 1102.1720 .\n- 144. Krumholz, M. R., Dekel, A. & McKee, C. F. A Universal, Local Star Formation Law in Galactic Clouds, nearby Galaxies, High-redshift Disks, and Starbursts. ApJ 745 , 69 (2012). DOI 10.1088/0004-637X/745/1/69. 1109.4150 .\n- 145. Hopkins, P. F. et al. Galaxies on FIRE (Feedback In Realistic Environments): stellar feedback explains cosmologically inefficient star formation. MNRAS 445 , 581-603 (2014). DOI 10.1093/mnras/stu1738. 1311.2073 .\n- 146. Hopkins, P. F. et al. FIRE-2 simulations: physics versus numerics in galaxy formation. MNRAS 480 , 800-863 (2018). DOI 10.1093/mnras/sty1690. 1702.06148 .\n- 147. Stinson, G. et al. Star formation and feedback in smoothed particle hydrodynamic simulations I. Isolated galaxies. MNRAS 373 , 1074-1090 (2006). DOI 10.1111/j.1365-2966.2006.11097.x. arXiv:astro-ph/0602350 .\n- 148. Teyssier, R., Pontzen, A., Dubois, Y. & Read, J. I. Cusp-core transformations in dwarf galaxies: observational predictions. MNRAS 429 , 3068-3078 (2013). DOI 10.1093/mnras/sts563. 1206.4895 .\n- 149. Vogelsberger, M. et al. Properties of galaxies reproduced by a hydrodynamic simulation. Nat. 509 , 177-182 (2014). DOI 10.1038/nature13316. 1405.1418 .\n- 150. Pillepich, A. et al. Simulating galaxy formation with the IllustrisTNG model. MNRAS 473 , 4077-4106 (2018). DOI 10.1093/mnras/stx2656. 1703.02970 .\n- 151. Semenov, V. A., Kravtsov, A. V. & Gnedin, N. Y. The Physical Origin of Long Gas Depletion Times in Galaxies. ApJ 845 , 133 (2017). DOI 10.3847/1538-4357/aa8096. 1704.04239 .\n- 152. Trebitsch, M., Blaizot, J., Rosdahl, J., Devriendt, J. & Slyz, A. Fluctuating feedback-regulated escape fraction of ionizing radiation in low-mass, high-redshift galaxies. MNRAS 470 , 224-239 (2017). DOI 10.1093/mnras/stx1060. 1705.00941 .  \n- 182. Dav'e, R. et al. SIMBA: Cosmological simulations with black hole growth and feedback. MNRAS 486 , 2827-2849 (2019). DOI 10.1093/mnras/stz937. 1901.10203 .  \n- 183. Krolik, J. H. Active galactic nuclei : from the central black hole to the galactic environment (Princeton, N. J. : Princeton University Press, 1999).\n- 184. Springel, V., Di Matteo, T. & Hernquist, L. Modelling feedback from stars and black holes in galaxy mergers. MNRAS 361 , 776-794 (2005). DOI 10.1111/j.1365-2966.2005.09238.x. astro-ph/0411108 .\n- 185. Di Matteo, T., Springel, V. & Hernquist, L. Energy input from quasars regulates the growth and activity of black holes and their host galaxies. Nat. 433 , 604-607 (2005). DOI 10.1038/nature03335. astro-ph/0502199 .\n- 186. Debuhr, J., Quataert, E., Ma, C.-P. & Hopkins, P. Self-regulated black hole growth via momentum deposition in galaxy merger simulations. MNRAS 406 , L55-L59 (2010). DOI 10.1111/j.1745-3933.2010.00881.x. 0909.2872 .\n- 187. Costa, T., Rosdahl, J., Sijacki, D. & Haehnelt, M. G. Driving gas shells with radiation pressure on dust in radiation-hydrodynamic simulations. MNRAS 473 , 4197-4219 (2018). DOI 10.1093/mnras/stx2598. 1703.05766 .\n- 188. Barnes, D. J., Kannan, R., Vogelsberger, M. & Marinacci, F. Radiative AGN feedback on a moving mesh: the impact of the galactic disc and dust physics on outflow properties. arXiv e-prints (2018). 1812.01611 .\n- 189. Choi, E., Ostriker, J. P., Naab, T. & Johansson, P. H. Radiative and Momentum-based Mechanical Active Galactic Nucleus Feedback in a Three-dimensional Galaxy Evolution Code. ApJ 754 , 125 (2012). DOI 10.1088/0004-637X/754/2/125. 1205.2082 .\n- 190. Sijacki, D., Springel, V., Di Matteo, T. & Hernquist, L. A unified model for AGN feedback in cosmological simulations of structure formation. MNRAS 380 , 877-900 (2007). DOI 10.1111/j.1365-2966.2007.12153.x. 0705.2238 .\n- 191. Weinberger, R. et al. Simulating galaxy formation with black hole driven thermal and kinetic feedback. MNRAS 465 , 3291-3308 (2017). DOI 10.1093/mnras/stw2944. 1607.03486 .\n- 192. Marinacci, F. & Vogelsberger, M. Effects of simulated cosmological magnetic fields on the galaxy population. MNRAS 456 , L69-L73 (2016). DOI 10.1093/mnrasl/slv176. 1508. 06631 .\n- 193. Ferri'ere, K. M. The interstellar environment of our galaxy. Rev. Mod. Phys. 73 , 1031-1066 (2001). DOI 10.1103/RevModPhys.73.1031. astro-ph/0106359 .\n- 194. Kotera, K. & Olinto, A. V. The Astrophysics of Ultrahigh-Energy Cosmic Rays. ARA&A 49 , 119-153 (2011). DOI 10.1146/annurev-astro-081710-102620. 1101.4256 .\n- 195. Donnert, J., Dolag, K., Lesch, H. & Muller, E. Cluster magnetic fields from galactic outflows. MNRAS 392 , 1008-1021 (2009). DOI 10.1111/j.1365-2966.2008.14132.x. 0808.0919 .\n- 196. Marinacci, F. et al. First results from the IllustrisTNG simulations: radio haloes and magnetic fields. MNRAS 480 , 5113-5139 (2018). DOI 10.1093/mnras/sty2206. 1707.03396 .\n- 197. Pakmor, R. et al. Magnetic field formation in the Milky Way like disc galaxies of the Auriga project. MNRAS 469 , 3185-3199 (2017). DOI 10.1093/mnras/stx1074. 1701.07028 .\n- 198. Powell, K. G., Roe, P. L., Linde, T. J., Gombosi, T. I. & De Zeeuw, D. L. A Solution-Adaptive Upwind Scheme for Ideal Magnetohydrodynamics. J. Comput. Phys. 154 , 284-309 (1999). DOI 10.1006/jcph.1999.6299.\n- 199. Dedner, A. et al. Hyperbolic Divergence Cleaning for the MHD Equations. J. Comput. Phys. 175 , 645-673 (2002). DOI 10.1006/jcph.2001.6961.\n- 200. Evans, C. R. & Hawley, J. F. Simulation of magnetohydrodynamic flows - A constrained transport method. ApJ 332 , 659-677 (1988). DOI 10.1086/166684.\n- 201. Mocz, P. et al. A moving mesh unstaggered constrained transport scheme for magnetohydrodynamics. MNRAS 463 , 477-488 (2016). DOI 10.1093/mnras/stw2004. 1606.02310 .\n- 202. Rosswog, S. & Price, D. MAGMA: a three-dimensional, Lagrangian magnetohydrodynamics code for merger applications. MNRAS 379 , 915-931 (2007). DOI 10.1111/j.13652966.2007.11984.x. 0705.1441 .\n- 203. Dolag, K. & Stasyszyn, F. An MHD GADGET for cosmological simulations. MNRAS 398 , 1678-1697 (2009). DOI 10.1111/j.1365-2966.2009.15181.x. 0807.3553 .\n- 204. Stone, J. M. & Norman, M. L. ZEUS-2D: A Radiation Magnetohydrodynamics Code for Astrophysical Flows in Two Space Dimensions. II. The Magnetohydrodynamic Algorithms and Tests. ApJS 80 , 791 (1992). DOI 10.1086/191681.\n- 205. Londrillo, P. & del Zanna, L. On the divergence-free condition in Godunov-type schemes for ideal magnetohydrodynamics: the upwind constrained transport method. J. Comput. Phys. 195 , 17-48 (2004). DOI 10.1016/j.jcp.2003.09.016. astro-ph/0310183 .\n- 206. Fromang, S., Hennebelle, P. & Teyssier, R. A high order Godunov scheme with constrained transport and adaptive mesh refinement for astrophysical magnetohydrodynamics. A&A 457 , 371-384 (2006). DOI 10.1051/0004-6361:20065371. astro-ph/0607230 .\n- 207. Teyssier, R., Fromang, S. & Dormy, E. Kinematic dynamos using constrained transport with high order Godunov schemes and adaptive mesh refinement. J. Comput. Phys. 218 , 44-67 (2006). DOI 10.1016/j.jcp.2006.01.042. astro-ph/0601715 .\n- 208. Balsara, D. S. A two-dimensional HLLC Riemann solver for conservation laws: Application to Euler and magnetohydrodynamic flows. J. Comput. Phys. 231 , 7476-7503 (2012). DOI 10.1016/j.jcp.2011.12.025. 1110.0750 .\n- 209. Ferriere, K. Magnetic Fields in Galaxies: Their Origin and Their Impact on the Interstellar Medium. Acta Astron. Sinica 44 , 115-122 (2003).\n- 210. Cox, D. P. The Three-Phase Interstellar Medium Revisited. ARA&A 43 , 337-385 (2005). DOI 10.1146/annurev.astro.43.072103.150615.\n- 211. Field, G. B., Goldsmith, D. W. & Habing, H. J. Cosmic-Ray Heating of the Interstellar Gas. ApJ 155 , L149 (1969). DOI 10.1086/180324.\n- 212. Wolfire, M. G., Hollenbach, D., McKee, C. F., Tielens, A. G. G. M. & Bakes, E. L. O. The neutral atomic phases of the interstellar medium. ApJ 443 , 152-168 (1995). DOI 10.1086/175510.\n- 213. Uhlig, M. et al. Galactic winds driven by cosmic ray streaming. MNRAS 423 , 2374-2396 (2012). DOI 10.1111/j.1365-2966.2012.21045.x. 1203.1038 .\n- 214. Booth, C. M., Agertz, O., Kravtsov, A. V. & Gnedin, N. Y. Simulations of Disk Galaxies with Cosmic Ray Driven Galactic Winds. ApJ 777 , L16 (2013). DOI 10.1088/2041-8205/777/1/L16. 1308.4974 .\n- 215. Hanasz, M. et al. Cosmic Rays Can Drive Strong Outflows from Gas-rich High-redshift Disk Galaxies. ApJ 777 , L38 (2013). DOI 10.1088/2041-8205/777/2/L38. 1310.3273 .\n- 216. Salem, M., Bryan, G. L. & Hummels, C. Cosmological Simulations of Galaxy Formation with Cosmic Rays. ApJ 797 , L18 (2014). DOI 10.1088/2041-8205/797/2/L18. 1412.0661 .\n- 217. Pakmor, R., Pfrommer, C., Simpson, C. M. & Springel, V. Galactic Winds Driven by Isotropic and Anisotropic Cosmic-Ray Diffusion in Disk Galaxies. ApJ 824 , L30 (2016). DOI 10.3847/2041-8205/824/2/L30. 1605.00643 .\n- 218. Simpson, C. M. et al. The Role of Cosmic-Ray Pressure in Accelerating Galactic Outflows. ApJ 827 , L29 (2016). DOI 10.3847/2041-8205/827/2/L29. 1606.02324 .\n- 219. Ruszkowski, M., Yang, H.-Y. K. & Zweibel, E. Global Simulations of Galactic Winds Including Cosmic-ray Streaming. ApJ 834 , 208 (2017). DOI 10.3847/1538-4357/834/2/208. 1602.04856 .\n- 220. Farber, R., Ruszkowski, M., Yang, H.-Y. K. & Zweibel, E. G. Impact of Cosmic-Ray Transport on Galactic Winds. ApJ 856 , 112 (2018). DOI 10.3847/1538-4357/aab26d. 1707.04579 .\n- 221. Girichidis, P., Naab, T., Hanasz, M. & Walch, S. Cooler and smoother - the impact of cosmic rays on the phase structure of galactic outflows. MNRAS 479 , 3042-3067 (2018). DOI 10.1093/mnras/sty1653. 1805.09333 .\n- 222. Jacob, S., Pakmor, R., Simpson, C. M., Springel, V. & Pfrommer, C. The dependence of cosmic ray-driven galactic winds on halo mass. MNRAS 475 , 570-584 (2018). DOI 10.1093/mnras/stx3221. 1712.04947 .\n- 223. Pfrommer, C., Pakmor, R., Schaal, K., Simpson, C. M. & Springel, V. Simulating cosmic ray physics on a moving mesh. MNRAS 465 , 4500-4529 (2017). DOI 10.1093/mnras/stw2941. 1604.07399 .\n- 224. Sharma, P. & Hammett, G. W. Preserving monotonicity in anisotropic diffusion. J. Comput. Phys. 227 , 123-142 (2007). DOI 10.1016/j.jcp.2007.07.026. 0707.2616 .\n- 225. Kannan, R., Springel, V., Pakmor, R., Marinacci, F. & Vogelsberger, M. Accurately simulating anisotropic thermal conduction on a moving mesh. MNRAS 458 , 410-424 (2016). DOI 10.1093/mnras/stw294. 1512.03053 .\n- 226. Butsky, I. S. & Quinn, T. R. The Role of Cosmic-ray Transport in Shaping the Simulated Circumgalactic Medium. ApJ 868 , 108 (2018). DOI 10.3847/1538-4357/aaeac2. 1803.06345 .\n- 227. Sharma, P., Colella, P. & Martin, D. Numerical implementation of streaming down the gradient: Application to fluid modeling of cosmic rays and saturated conduction. Siam J. on Sci. Comput. 32 , 3476-3494 (2009). DOI 10.1137/100792135.\n- 228. Jiang, Y.-F. & Oh, S. P. A New Numerical Scheme for Cosmic-Ray Transport. ApJ 854 , 5 (2018). DOI 10.3847/1538-4357/aaa6ce. 1712.07117 .\n- 229. Thomas, T. & Pfrommer, C. Cosmic-ray hydrodynamics: Alfv'en-wave regulated transport of cosmic rays. MNRAS 485 , 2977-3008 (2019). DOI 10.1093/mnras/stz263. 1805.11092 .\n- 230. Gnedin, N. Y. & Kaurov, A. A. Cosmic Reionization on Computers. II. Reionization History and Its Back-reaction on Early Galaxies. ApJ 793 , 30 (2014). DOI 10.1088/0004637X/793/1/30. 1403.4251 .\n- 231. Ocvirk, P. et al. Cosmic Dawn II (CoDa II): a new radiation-hydrodynamics simulation of the self-consistent coupling of galaxy formation and reionization. arXiv e-prints (2018). 1811.11192 .\n- 232. Wise, J. H., Abel, T., Turk, M. J., Norman, M. L. & Smith, B. D. The birth of a galaxy - II. The role of radiation pressure. MNRAS 427 , 311-326 (2012). DOI 10.1111/j.13652966.2012.21809.x. 1206.1043 .\n- 233. Gnedin, N. Y. & Abel, T. Multi-dimensional cosmological radiative transfer with a Variable Eddington Tensor formalism. New A 6 , 437-455 (2001). DOI 10.1016/S1384-1076(01)00068-9. astro-ph/0106278 .\n- 234. Gnedin, N. Y. On the Proper Use of the Reduced Speed of Light Approximation. ApJ 833 , 66 (2016). DOI 10.3847/1538-4357/833/1/66. 1607.07869 .\n- 235. Deparis, N., Aubert, D., Ocvirk, P., Chardin, J. & Lewis, J. Impact of the reduced speed of light approximation on ionization front velocities in cosmological simulations of the epoch of reionization. A&A 622 , A142 (2019). DOI 10.1051/0004-6361/201832889. 1803.01634 .\n- 236. Ocvirk, P., Aubert, D., Chardin, J., Deparis, N. & Lewis, J. Impact of the reduced speed of light approximation on the post-overlap neutral hydrogen fraction in numerical simulations of the epoch of reionization. A&A 626 , A77 (2019). DOI 10.1051/0004-6361/201832923. 1803.02434 .\n- 237. Mihalas, D. & Mihalas, B. W. Foundations of radiation hydrodynamics (New York, Oxford University Press, 1984).\n- 238. Abel, T., Norman, M. L. & Madau, P. Photon-conserving Radiative Transfer around Point Sources in Multidimensional Numerical Cosmology. ApJ 523 , 66-71 (1999). DOI 10.1086/307739. astro-ph/9812151 .\n- 239. Jaura, O., Glover, S. C. O., Klessen, R. S. & Paardekooper, J.-P. SPRAI: coupling of radiative feedback and primordial chemistry in moving mesh hydrodynamics. MNRAS 475 , 2822-2834 (2018). DOI 10.1093/mnras/stx3356. 1711.02542 .\n- 240. Rijkhorst, E.-J., Plewa, T., Dubey, A. & Mellema, G. Hybrid characteristics: 3D radiative transfer for parallel adaptive mesh refinement hydrodynamics. A&A 452 , 907-920 (2006). DOI 10.1051/0004-6361:20053401. astro-ph/0505213 .  \n- 241. Rosen, A. L., Krumholz, M. R., Oishi, J. S., Lee, A. T. & Klein, R. I. Hybrid Adaptive RayMoment Method (HARM 2 ): A highly parallel method for radiation hydrodynamics on adaptive grids. J. Comput. Phys. 330 , 924-942 (2017). DOI 10.1016/j.jcp.2016.10.048. 1607.01802 .\n- 242. Whalen, D. & Norman, M. L. A Multistep Algorithm for the Radiation Hydrodynamical Transport of Cosmological Ionization Fronts and Ionized Flows. ApJS 162 , 281-303 (2006). DOI 10.1086/499072. astro-ph/0508214 .\n- 243. Trac, H. & Cen, R. Radiative Transfer Simulations of Cosmic Reionization. I. Methodology and Initial Results. ApJ 671 , 1-13 (2007). DOI 10.1086/522566. astro-ph/0612406 .\n- 244. Pawlik, A. H. & Schaye, J. TRAPHIC - radiative transfer for smoothed particle hydrodynamics simulations. MNRAS 389 , 651-677 (2008). DOI 10.1111/j.1365-2966.2008.13601.x. 0802.1715 .\n- 245. Petkova, M. & Springel, V. A novel approach for accurate radiative transfer in cosmological hydrodynamic simulations. MNRAS 415 , 3731-3749 (2011). DOI 10.1111/j.13652966.2011.18986.x. 1012.1017 .\n- 246. Ciardi, B., Ferrara, A., Marri, S. & Raimondo, G. Cosmological reionization around the first stars: Monte Carlo radiative transfer. MNRAS 324 , 381-388 (2001). DOI 10.1046/j.13658711.2001.04316.x. astro-ph/0005181 .\n- 247. Oxley, S. & Woolfson, M. M. Smoothed particle hydrodynamics with radiation transfer. MNRAS 343 , 900-912 (2003). DOI 10.1046/j.1365-8711.2003.06751.x.\n- 248. Tasitsiomi, A. Ly a Radiative Transfer in Cosmological Simulations and Application to a z ˜= 8 Ly a Emitter. ApJ 645 , 792-813 (2006). DOI 10.1086/504460. astro-ph/0510347 .\n- 249. Semelin, B., Combes, F. & Baek, S. Lyman-alpha radiative transfer during the epoch of reionization: contribution to 21-cm signal fluctuations. A&A 474 , 365-374 (2007). DOI 10.1051/0004-6361:20077965. 0707.2483 .\n- 250. Dullemond, C. P. et al. RADMC-3D: A multi-purpose radiative transfer tool. Astrophysics Source Code Library (2012). 1202.015 .\n- 251. Smith, A. et al. The physics of Lyman a escape from high-redshift galaxies. MNRAS 484 , 39-59 (2019). DOI 10.1093/mnras/sty3483. 1810.08185 .\n- 252. Lucy, L. B. Improved Monte Carlo techniques for the spectral synthesis of supernovae. A&A 345 , 211-220 (1999).\n- 253. Gentile, N. A. Implicit Monte Carlo Diffusion-An Acceleration Method for Monte Carlo Time-Dependent Radiative Transfer Simulations. J. Comput. Phys. 172 , 543-571 (2001). DOI 10.1006/jcph.2001.6836.\n- 254. Densmore, J. D., Urbatsch, T. J., Evans, T. M. & Buksas, M. W. A hybrid transport-diffusion method for Monte Carlo radiative-transfer simulations. J. Comput. Phys. 222 , 485-503 (2007). DOI 10.1016/j.jcp.2006.07.031.\n- 255. Smith, A., Tsang, B. T. H., Bromm, V. & Milosavljevi'c, M. Discrete diffusion Lyman a radiative transfer. MNRAS 479 , 2065-2078 (2018). DOI 10.1093/mnras/sty1509. 1709.10187 .\n- 256. Levermore, C. D. Relating Eddington factors to flux limiters. J. Quant. Spectrosc. Radiat. Transf. 31 , 149-160 (1984). DOI 10.1016/0022-4073(84)90112-2.\n- 257. Gonz'alez, M., Audit, E. & Huynh, P. HERACLES: a three-dimensional radiation hydrodynamics code. A&A 464 , 429-435 (2007). DOI 10.1051/0004-6361:20065486.\n- 258. Rosdahl, J., Blaizot, J., Aubert, D., Stranex, T. & Teyssier, R. RAMSES-RT: radiation hydrodynamics in the cosmological context. MNRAS 436 , 2188-2231 (2013). DOI 10.1093/mnras/stt1722. 1304.7126 .\n- 259. Kannan, R. et al. AREPO-RT: radiation hydrodynamics on a moving mesh. MNRAS 485 , 117-149 (2019). DOI 10.1093/mnras/stz287. 1804.01987 .\n- 260. Krumholz, M. R., Klein, R. I., McKee, C. F. & Bolstad, J. Equations and Algorithms for Mixed-frame Flux-limited Diffusion Radiation Hydrodynamics. ApJ 667 , 626-643 (2007). DOI 10.1086/520791. astro-ph/0611003 .\n- 261. Finlator, K., Ozel, F. & Dav'e, R. A new moment method for continuum radiative transfer in cosmological re-ionization. MNRAS 393 , 1090-1106 (2009). DOI 10.1111/j.13652966.2008.14190.x. 0808.3578 .\n- 262. Petkova, M. & Springel, V. An implementation of radiative transfer in the cosmological simulation code GADGET. MNRAS 396 , 1383-1403 (2009). DOI 10.1111/j.13652966.2009.14843.x. 0812.1801 .\n- 263. Dubroca, B. & Feugeas, J. Etude th'eorique et num'erique d'une hi'erarchie de mod'eles aux moments pour le transfert radiatif. Acad. des Sci. Paris Comptes Rendus Ser. Sci. Math. 329 , 915-920 (1999). DOI 10.1016/S0764-4442(00)87499-6.\n- 264. Ripoll, J.-F., Dubroca, B. & Duffa, G. Modelling radiative mean absorption coefficients. Combust. Theory Model. 5 , 261-274 (2001). DOI 10.1088/1364-7830/5/3/301.\n- 265. McKinnon, R., Torrey, P. & Vogelsberger, M. Dust formation in Milky Way-like galaxies. MNRAS 457 , 3775-3800 (2016). DOI 10.1093/mnras/stw253. 1505.04792 .\n- 266. Aoyama, S. et al. Galaxy simulation with dust formation and destruction. MNRAS 466 , 105121 (2017). DOI 10.1093/mnras/stw3061. 1609.07547 .\n- 267. Hou, K.-C., Hirashita, H., Nagamine, K., Aoyama, S. & Shimizu, I. Evolution of dust extinction curves in galaxy simulation. MNRAS 469 , 870-885 (2017). DOI 10.1093/mnras/stx877. 1704.01769 .\n- 268. McKinnon, R., Torrey, P., Vogelsberger, M., Hayward, C. C. & Marinacci, F. Simulating the dust content of galaxies: successes and failures. MNRAS 468 , 1505-1521 (2017). DOI 10.1093/mnras/stx467. 1606.02714 .\n- 269. Aoyama, S., Hou, K.-C., Hirashita, H., Nagamine, K. & Shimizu, I. Cosmological simulation with dust formation and destruction. MNRAS 478 , 4905-4921 (2018). DOI 10.1093/mnras/sty1431. 1802.04027 .\n- 270. McKinnon, R., Vogelsberger, M., Torrey, P., Marinacci, F. & Kannan, R. Simulating galactic dust grain evolution on a moving mesh. MNRAS 478 , 2851-2886 (2018). DOI 10.1093/mnras/sty1248. 1805.04521 .\n- 271. Vogelsberger, M. et al. Dust in and around galaxies: dust in cluster environments and its impact on gas cooling. MNRAS 487 , 4870-4883 (2019). DOI 10.1093/mnras/stz1644. 1811.05477 .\n- 272. Gjergo, E. et al. Dust evolution in galaxy cluster simulations. MNRAS 479 , 2588-2606 (2018). DOI 10.1093/mnras/sty1564. 1804.06855 .\n- 273. Hou, K.-C., Aoyama, S., Hirashita, H., Nagamine, K. & Shimizu, I. Dust scaling relations in a cosmological simulation. MNRAS 485 , 1727-1744 (2019). DOI 10.1093/mnras/stz121. 1901.02886 .\n- 274. Ruszkowski, M., Lee, D., Bruggen, M., Parrish, I. & Oh, S. P. Cosmological Magnetohydrodynamic Simulations of Cluster Formation with Anisotropic Thermal Conduction. ApJ 740 , 81 (2011). DOI 10.1088/0004-637X/740/2/81. 1010.2277 .\n- 275. Smith, B., O'Shea, B. W., Voit, G. M., Ventimiglia, D. & Skillman, S. W. Cosmological Simulations of Isotropic Conduction in Galaxy Clusters. ApJ 778 , 152 (2013). DOI 10.1088/0004637X/778/2/152. 1306.5748 .\n- 276. Arth, A., Dolag, K., Beck, A. M., Petkova, M. & Lesch, H. Anisotropic thermal conduction in galaxy clusters with MHD in Gadget. arXiv e-prints (2014). 1412.6533 .\n- 277. Yang, H.-Y. K. & Reynolds, C. S. Interplay Among Cooling, AGN Feedback, and Anisotropic Conduction in the Cool Cores of Galaxy Clusters. ApJ 818 , 181 (2016). DOI 10.3847/0004637X/818/2/181. 1512.05796 .\n- 278. Kannan, R. et al. Increasing Black Hole Feedback-induced Quenching with Anisotropic Thermal Conduction. ApJ 837 , L18 (2017). DOI 10.3847/2041-8213/aa624b. 1612.01522 .\n- 279. Parrish, I. J., McCourt, M., Quataert, E. & Sharma, P. The effects of anisotropic viscosity on turbulence and heat transport in the intracluster medium. MNRAS 422 , 704-718 (2012). DOI 10.1111/j.1365-2966.2012.20650.x. 1201.0754 .\n- 280. Suzuki, K., Ogawa, T., Matsumoto, Y. & Matsumoto, R. Magnetohydrodynamic Simulations of the Formation of Cold Fronts in Clusters of Galaxies: Effects of Anisotropic Viscosity. ApJ 768 , 175 (2013). DOI 10.1088/0004-637X/768/2/175.\n- 281. ZuHone, J. A., Kunz, M. W., Markevitch, M., Stone, J. M. & Biffi, V. The Effect of Anisotropic Viscosity on Cold Fronts in Galaxy Clusters. ApJ 798 , 90 (2015). DOI 10.1088/0004637X/798/2/90. 1406.4031 .\n- 282. Su, K.-Y. et al. Feedback first: the surprisingly weak effects of magnetic fields, viscosity, conduction and metal diffusion on sub-L* galaxy formation. MNRAS 471 , 144-166 (2017). DOI 10.1093/mnras/stx1463. 1607.05274 .\n- 283. Balbus, S. A. Stability, Instability, and 'Backward' Transport in Stratified Fluids. ApJ 534 , 420-427 (2000). DOI 10.1086/308732. astro-ph/9906315 .\n- 284. Quataert, E. Buoyancy Instabilities in Weakly Magnetized Low-Collisionality Plasmas. ApJ 673 , 758-762 (2008). DOI 10.1086/525248. 0710.5521 .\n- 285. Barnes, D. J. et al. Enhancing AGN efficiency and cool-core formation with anisotropic thermal conduction. arXiv e-prints (2018). 1805.04109 .\n- 286. Torrey, P. et al. Synthetic galaxy images and spectra from the Illustris simulation. MNRAS 447 , 2753-2771 (2015). DOI 10.1093/mnras/stu2592. 1411.3717 .\n- 287. Trayford, J. W. et al. Optical colours and spectral indices of z = 0.1 eagle galaxies with the 3D dust radiative transfer code skirt. MNRAS 470 , 771-799 (2017). DOI 10.1093/mnras/stx1051. 1705.02331 .\n- 288. Hernquist, L., Katz, N., Weinberg, D. H. & Miralda-Escud'e, J. The Lyman-Alpha Forest in the Cold Dark Matter Model. ApJ 457 , L51 (1996). DOI 10.1086/309899. astro-ph/ 9509105 .\n- 289. Katz, N. & Gunn, J. E. Dissipational galaxy formation. I - Effects of gasdynamics. ApJ 377 , 365-381 (1991). DOI 10.1086/170367.\n- 290. Navarro, J. F. & Benz, W. Dynamics of cooling gas in galactic dark halos. ApJ 380 , 320-329 (1991). DOI 10.1086/170590.\n- 291. Katz, N. Dissipational galaxy formation. II - Effects of star formation. ApJ 391 , 502-517 (1992). DOI 10.1086/171366.\n- 292. Katz, N., Hernquist, L. & Weinberg, D. H. Galaxies and gas in a cold dark matter universe. ApJ 399 , L109-L112 (1992). DOI 10.1086/186619.\n- 293. Navarro, J. F. & Steinmetz, M. The Effects of a Photoionizing Ultraviolet Background on the Formation of Disk Galaxies. ApJ 478 , 13-28 (1997). DOI 10.1086/303763. astro-ph/ 9605043 .\n- 294. Navarro, J. F. & Steinmetz, M. Dark Halo and Disk Galaxy Scaling Laws in Hierarchical Universes. ApJ 538 , 477-488 (2000). DOI 10.1086/309175. arXiv:astro-ph/0001003 .\n- 295. Abadi, M. G., Navarro, J. F., Steinmetz, M. & Eke, V. R. Simulations of Galaxy Formation in a L Cold Dark Matter Universe. I. Dynamical and Photometric Properties of a Simulated Disk Galaxy. ApJ 591 , 499-514 (2003). DOI 10.1086/375512. astro-ph/0211331 .\n- 296. Scannapieco, C., Tissera, P. B., White, S. D. M. & Springel, V. Effects of supernova feedback on the formation of galaxy discs. MNRAS 389 , 1137-1149 (2008). DOI 10.1111/j.13652966.2008.13678.x. 0804.3795 .\n- 297. Grand, R. J. J. et al. The Auriga Project: the properties and formation mechanisms of disc galaxies across cosmic time. MNRAS 467 , 179-207 (2017). DOI 10.1093/mnras/stx071.\n- 298. Kaviraj, S. et al. The Horizon-AGN simulation: evolution of galaxy properties over cosmic time. MNRAS 467 , 4739-4752 (2017). DOI 10.1093/mnras/stx126. 1605.09379 .\n- 299. Schechter, P. An analytic expression for the luminosity function for galaxies. ApJ 203 , 297306 (1976). DOI 10.1086/154079.  \n- 300. Panter, B., Jimenez, R., Heavens, A. F. & Charlot, S. The star formation histories of galaxies in the Sloan Digital Sky Survey. MNRAS 378 , 1550-1564 (2007). DOI 10.1111/j.13652966.2007.11909.x. astro-ph/0608531 .\n- 301. Pozzetti, L. et al. zCOSMOS - 10k-bright spectroscopic sample. The bimodality in the galaxy stellar mass function: exploring its evolution with redshift. A&A 523 , A13 (2010). DOI 10.1051/0004-6361/200913020. 0907.5416 .\n- 302. Baldry, I. K. et al. Galaxy And Mass Assembly (GAMA): the galaxy stellar mass function at z < 0.06. MNRAS 421 , 621-634 (2012). DOI 10.1111/j.1365-2966.2012.20340.x. 1111.5707 .\n- 303. Ilbert, O. et al. Mass assembly in quiescent and star-forming galaxies since z glyph[similarequal] 4 from UltraVISTA. A&A 556 , A55 (2013). DOI 10.1051/0004-6361/201321100. 1301.3157 .\n- 304. Muzzin, A. et al. The Evolution of the Stellar Mass Functions of Star-forming and Quiescent Galaxies to z = 4 from the COSMOS/UltraVISTA Survey. ApJ 777 , 18 (2013). DOI 10.1088/0004-637X/777/1/18. 1303.4409 .\n- 305. Weigel, A. K., Schawinski, K. & Bruderer, C. Stellar mass functions: methods, systematics and results for the local Universe. MNRAS 459 , 2150-2187 (2016). DOI 10.1093/mnras/stw756. 1604.00008 .\n- 306. Dekel, A. & Silk, J. The Origin of Dwarf Galaxies, Cold Dark Matter, and Biased Galaxy Formation. ApJ 303 , 39 (1986). DOI 10.1086/164050.\n- 307. Larson, R. B. Effects of supernovae on the early evolution of galaxies. MNRAS 169 , 229-246 (1974). DOI 10.1093/mnras/169.2.229.\n- 308. White, S. D. M. & Frenk, C. S. Galaxy Formation through Hierarchical Clustering. ApJ 379 , 52 (1991). DOI 10.1086/170483.\n- 309. Pillepich, A. et al. First results from the IllustrisTNG simulations: the stellar mass content of groups and clusters of galaxies. MNRAS 475 , 648-675 (2018). DOI 10.1093/mnras/stx3112. 1707.03406 .\n- 310. Werk, J. K. et al. The COS-Halos Survey: Physical Conditions and Baryonic Mass in the Low-redshift Circumgalactic Medium. ApJ 792 , 8 (2014). DOI 10.1088/0004-637X/792/1/8. 1403.0947 .\n- 311. Werk, J. K. et al. The COS-Halos Survey: Origins of the Highly Ionized Circumgalactic Medium of Star-Forming Galaxies. ApJ 833 , 54 (2016). DOI 10.3847/1538-4357/833/1/54. 1609.00012 .\n- 312. Stern, J., Hennawi, J. F., Prochaska, J. X. & Werk, J. K. A Universal Density Structure for Circumgalactic Gas. ApJ 830 , 87 (2016). DOI 10.3847/0004-637X/830/2/87. 1604.02168 .\n- 313. Rubin, K. H. R., Diamond-Stanic, A. M., Coil, A. L., Crighton, N. H. M. & Stewart, K. R. Galaxies Probing Galaxies in PRIMUS. II. The Coherence Scale of the Cool Circumgalactic Medium. ApJ 868 , 142 (2018). DOI 10.3847/1538-4357/aad566. 1806.08801 .\n- 314. Peeples, M. S. et al. Figuring Out Gas &amp; Galaxies in Enzo (FOGGIE). I. Resolving Simulated Circumgalactic Absorption at 2 ≤ z ≤ 2.5. ApJ 873 , 129 (2019). DOI 10.3847/15384357/ab0654. 1810.06566 .\n- 315. van de Voort, F., Springel, V., Mandelker, N., van den Bosch, F. C. & Pakmor, R. Cosmological simulations of the circumgalactic medium with 1 kpc resolution: enhanced H I column densities. MNRAS 482 , L85-L89 (2019). DOI 10.1093/mnrasl/sly190. 1808.04369 .\n- 316. Suresh, J., Nelson, D., Genel, S., Rubin, K. H. R. & Hernquist, L. Zooming in on accretion - II. Cold circumgalactic gas simulated with a super-Lagrangian refinement scheme. MNRAS 483 , 4040-4059 (2019). DOI 10.1093/mnras/sty3402. 1811.01949 .\n- 317. Hummels, C. B. et al. The Impact of Enhanced Halo Resolution on the Simulated Circumgalactic Medium. ApJ 882 , 156 (2019). DOI 10.3847/1538-4357/ab378f. 1811.12410 .\n- 318. Rasia, E. et al. Cool Core Clusters from Cosmological Simulations. ApJ 813 , L17 (2015). DOI 10.1088/2041-8205/813/1/L17. 1509.04247 .\n- 319. Planelles, S. et al. On the role of AGN feedback on the thermal and chemodynamical properties of the hot intracluster medium. MNRAS 438 , 195-216 (2014). DOI 10.1093/mnras/stt2141. 1311.0818 .\n- 320. Biffi, V. et al. On the Nature of Hydrostatic Equilibrium in Galaxy Clusters. ApJ 827 , 112 (2016). DOI 10.3847/0004-637X/827/2/112. 1606.02293 .\n- 321. Vogelsberger, M. et al. The uniformity and time-invariance of the intra-cluster metal distribution in galaxy clusters from the IllustrisTNG simulations. MNRAS 474 , 2073-2093 (2018). DOI 10.1093/mnras/stx2955. 1707.05318 .\n- 322. Barnes, D. J. et al. A census of cool-core galaxy clusters in IllustrisTNG. MNRAS 481 , 18091831 (2018). DOI 10.1093/mnras/sty2078. 1710.08420 .\n- 323. Meneux, B. et al. The VIMOS-VLT Deep Survey (VVDS). The dependence of clustering on galaxy stellar mass at z ˜ 1. A&A 478 , 299-310 (2008). DOI 10.1051/0004-6361:20078182. 0706.4371 .\n- 324. Foucaud, S. et al. Clustering properties of galaxies selected in stellar mass: breaking down the link between luminous and dark matter in massive galaxies from z = 0 to z = 2. MNRAS 406 , 147-164 (2010). DOI 10.1111/j.1365-2966.2010.16682.x. 1003.2755 .\n- 325. Wake, D. A. et al. Galaxy Clustering in the NEWFIRM Medium Band Survey: The Relationship Between Stellar Mass and Dark Matter Halo Mass at 1 < z < 2. ApJ 728 , 46 (2011). DOI 10.1088/0004-637X/728/1/46. 1012.1317 .\n- 326. Artale, M. C. et al. Small-scale galaxy clustering in the eagle simulation. MNRAS 470 , 17711787 (2017). DOI 10.1093/mnras/stx1263. 1611.05064 .\n- 327. Marulli, F. et al. The VIMOS Public Extragalactic Redshift Survey (VIPERS) . Luminosity and stellar mass dependence of galaxy clustering at 0.5 < z < 1.1. A&A 557 , A17 (2013). DOI 10.1051/0004-6361/201321476. 1303.2633 .\n- 328. Shen, S. et al. The size distribution of galaxies in the Sloan Digital Sky Survey. MNRAS 343 , 978-994 (2003). DOI 10.1046/j.1365-8711.2003.06740.x. astro-ph/0301527 .\n- 329. Kormendy, J. & Ho, L. C. Coevolution (Or Not) of Supermassive Black Holes and Host Galaxies. ARA&A 51 , 511-653 (2013). DOI 10.1146/annurev-astro-082708-101811. 1304.7762 .\n- 330. Tremonti, C. A. et al. The Origin of the Mass-Metallicity Relation: Insights from 53,000 Star-forming Galaxies in the Sloan Digital Sky Survey. ApJ 613 , 898-913 (2004). DOI 10.1086/423264. astro-ph/0405537 .\n- 331. Dav'e, R., Rafieferantsoa, M. H., Thompson, R. J. & Hopkins, P. F. MUFASA: Galaxy star formation, gas, and metal properties across cosmic time. MNRAS 467 , 115-132 (2017). DOI 10.1093/mnras/stx108. 1610.01626 .\n- 332. De Rossi, M. E., Bower, R. G., Font, A. S., Schaye, J. & Theuns, T. Galaxy metallicity scaling relations in the EAGLE simulations. MNRAS 472 , 3354-3377 (2017). DOI 10.1093/mnras/stx2158. 1704.00006 .\n- 333. Torrey, P. et al. The evolution of the mass-metallicity relation and its scatter in IllustrisTNG. MNRAS 484 , 5587-5607 (2019). DOI 10.1093/mnras/stz243. 1711.05261 .\n- 334. Torrey, P. et al. Similar star formation rate and metallicity variability time-scales drive the fundamental metallicity relation. MNRAS 477 , L16-L20 (2018). DOI 10.1093/mnrasl/sly031. 1711.11039 .\n- 335. Trayford, J. W. et al. Colours and luminosities of z = 0.1 galaxies in the EAGLE simulation. MNRAS 452 , 2879-2896 (2015). DOI 10.1093/mnras/stv1461. 1504.04374 .\n- 336. Nelson, D. et al. First results from the IllustrisTNG simulations: the galaxy colour bimodality. MNRAS 475 , 624-647 (2018). DOI 10.1093/mnras/stx3040. 1707.03395 .\n- 337. Kauffmann, G., White, S. D. M. & Guiderdoni, B. The Formation and Evolution of Galaxies Within Merging Dark Matter Haloes. MNRAS 264 , 201 (1993). DOI 10.1093/mnras/264.1.201.\n- 338. Somerville, R. S. & Primack, J. R. Semi-analytic modelling of galaxy formation: the local Universe. MNRAS 310 , 1087-1110 (1999). DOI 10.1046/j.1365-8711.1999.03032.x. astro-ph/9802268 .\n- 339. Bower, R. G. et al. Breaking the hierarchy of galaxy formation. MNRAS 370 , 645-655 (2006). DOI 10.1111/j.1365-2966.2006.10519.x. astro-ph/0511338 .\n- 340. Croton, D. J. et al. The many lives of active galactic nuclei: cooling flows, black holes and the luminosities and colours of galaxies. MNRAS 365 , 11-28 (2006). DOI 10.1111/j.13652966.2005.09675.x. astro-ph/0508046 .\n- 341. Guo, Q. et al. From dwarf spheroidals to cD galaxies: simulating the galaxy population in a L CDM cosmology. MNRAS 413 , 101-131 (2011). DOI 10.1111/j.1365-2966.2010.18114.x. 1006.0106 .\n- 342. Fall, S. M. & Efstathiou, G. Formation and rotation of disc galaxies with haloes. MNRAS 193 , 189-206 (1980).\n- 343. Mo, H. J., Mao, S. & White, S. D. M. The formation of galactic discs. MNRAS 295 , 319-336 (1998). DOI 10.1046/j.1365-8711.1998.01227.x. astro-ph/9707093 .\n- 344. Scannapieco, C., White, S. D. M., Springel, V. & Tissera, P. B. The formation and survival of discs in a L CDM universe. MNRAS 396 , 696-708 (2009). DOI 10.1111/j.13652966.2009.14764.x. 0812.0976 .\n- 345. Balogh, M. L., Pearce, F. R., Bower, R. G. & Kay, S. T. Revisiting the cosmic cooling crisis. MNRAS 326 , 1228-1234 (2001). DOI 10.1111/j.1365-8711.2001.04667.x. arXiv: astro-ph/0104041 .\n- 346. Brook, C. B. et al. Hierarchical formation of bulgeless galaxies: why outflows have low angular momentum. MNRAS 415 , 1051-1060 (2011). DOI 10.1111/j.1365-2966.2011.18545.x. 1010.1004 .\n- 347. Okamoto, T., Eke, V. R., Frenk, C. S. & Jenkins, A. Effects of feedback on the morphology of galaxy discs. MNRAS 363 , 1299-1314 (2005). DOI 10.1111/j.1365-2966.2005.09525.x. astro-ph/0503676 .\n- 348. Brooks, A. M. et al. Interpreting the Evolution of the Size-Luminosity Relation for Disk Galaxies from Redshift 1 to the Present. ApJ 728 , 51 (2011). DOI 10.1088/0004-637X/728/1/51. 1011.0432 .\n- 349. Guedes, J., Callegari, S., Madau, P. & Mayer, L. Forming Realistic Late-type Spirals in a L CDM Universe: The Eris Simulation. ApJ 742 , 76 (2011). DOI 10.1088/0004637X/742/2/76. 1103.6030 .\n- 350. Aumer, M., White, S. D. M., Naab, T. & Scannapieco, C. Towards a more realistic population of bright spiral galaxies in cosmological simulations. MNRAS 434 , 3142-3164 (2013). DOI 10.1093/mnras/stt1230. 1304.1559 .\n- 351. Marinacci, F., Pakmor, R. & Springel, V. The formation of disc galaxies in high-resolution moving-mesh cosmological simulations. MNRAS 437 , 1750-1775 (2014). DOI 10.1093/mnras/stt2003. 1305.5360 .\n- 352. Wetzel, A. R. et al. Reconciling Dwarf Galaxies with L CDM Cosmology: Simulating a Realistic Population of Satellites around a Milky Way-mass Galaxy. ApJ 827 , L23 (2016). DOI 10.3847/2041-8205/827/2/L23. 1602.05957 .\n- 353. Pakmor, R., Marinacci, F. & Springel, V. Magnetic Fields in Cosmological Simulations of Disk Galaxies. ApJ 783 , L20 (2014). DOI 10.1088/2041-8205/783/1/L20. 1312.2620 .\n- 354. Beck, R. Magnetic fields in spiral galaxies. Astron. & Astropysics Rev. 24 , 4 (2015). DOI 10.1007/s00159-015-0084-4. 1509.04522 .\n- 355. Rieder, M. & Teyssier, R. A small-scale dynamo in feedback-dominated galaxies as the origin of cosmic magnetic fields - I. The kinematic phase. MNRAS 457 , 1722-1738 (2016). DOI 10.1093/mnras/stv2985. 1506.00849 .\n- 356. Rieder, M. & Teyssier, R. A small-scale dynamo in feedback-dominated galaxies - II. The saturation phase and the final magnetic configuration. MNRAS 471 , 2674-2686 (2017). DOI 10.1093/mnras/stx1670. 1704.05845 .  \n- 357. Rieder, M. & Teyssier, R. A small-scale dynamo in feedback-dominated galaxies - III. Cosmological simulations. MNRAS 472 , 4368-4373 (2017). DOI 10.1093/mnras/stx2276. 1708.01486 .\n- 358. Naab, T. et al. The ATLAS 3 D project - XXV. Two-dimensional kinematic analysis of simulated galaxies and the cosmological origin of fast and slow rotators. MNRAS 444 , 3357-3387 (2014). DOI 10.1093/mnras/stt1919. 1311.0284 .\n- 359. Kobayashi, C. GRAPE-SPH chemodynamical simulation of elliptical galaxies - II. Scaling relations and the fundamental plane. MNRAS 361 , 1216-1226 (2005). DOI 10.1111/j.13652966.2005.09248.x. astro-ph/0506094 .\n- 360. Feldmann, R., Carollo, C. M. & Mayer, L. The Hubble Sequence in Groups: The Birth of the Early-type Galaxies. ApJ 736 , 88 (2011). DOI 10.1088/0004-637X/736/2/88. 1008.3386 .\n- 361. Kobayashi, C. GRAPE-SPH chemodynamical simulation of elliptical galaxies - I. Evolution of metallicity gradients. MNRAS 347 , 740-758 (2004). DOI 10.1111/j.1365-2966.2004.07258.x. astro-ph/0310160 .\n- 362. Oser, L., Ostriker, J. P., Naab, T., Johansson, P. H. & Burkert, A. The Two Phases of Galaxy Formation. ApJ 725 , 2312-2323 (2010). DOI 10.1088/0004-637X/725/2/2312. 1010.1381 .\n- 363. Huang, S., Ho, L. C., Peng, C. Y., Li, Z.-Y. & Barth, A. J. Fossil Evidence for the Two-phase Formation of Elliptical Galaxies. ApJ 768 , L28 (2013). DOI 10.1088/2041-8205/768/2/L28. 1304.2299 .\n- 364. Rodriguez-Gomez, V. et al. The stellar mass assembly of galaxies in the Illustris simulation: growth by mergers and the spatial distribution of accreted stars. MNRAS 458 , 2371-2390 (2016). DOI 10.1093/mnras/stw456. 1511.08804 .\n- 365. Clauwens, B., Schaye, J., Franx, M. & Bower, R. G. The three phases of galaxy formation. MNRAS 478 , 3994-4009 (2018). DOI 10.1093/mnras/sty1229. 1711.00030 .\n- 366. Lagos, C. d. P. et al. The connection between mass, environment, and slow rotation in simulated galaxies. MNRAS 476 , 4327-4345 (2018). DOI 10.1093/mnras/sty489. 1712.01398 .\n- 367. Schulze, F. et al. Kinematics of simulated galaxies - I. Connecting dynamical and morphological properties of early-type galaxies at different redshifts. MNRAS 480 , 4636-4658 (2018). DOI 10.1093/mnras/sty2090. 1802.01583 .\n- 368. Springel, V. et al. Simulations of the formation, evolution and clustering of galaxies and quasars. Nat. 435 , 629-636 (2005b). DOI 10.1038/nature03597. arXiv:astro-ph/ 0504097 .\n- 369. Boylan-Kolchin, M., Springel, V., White, S. D. M., Jenkins, A. & Lemson, G. Resolving cosmic structure formation with the Millennium-II Simulation. MNRAS 398 , 1150-1164 (2009). DOI 10.1111/j.1365-2966.2009.15191.x. 0903.3041 .\n- 370. Teyssier, R. et al. Full-sky weak-lensing simulation with 70 billion particles. A&A 497 , 335341 (2009). DOI 10.1051/0004-6361/200810657. 0807.3651 .\n- 371. Klypin, A. A., Trujillo-Gomez, S. & Primack, J. Dark Matter Halos in the Standard Cosmological Model: Results from the Bolshoi Simulation. ApJ 740 , 102 (2011). DOI 10.1088/0004637X/740/2/102. 1002.3660 .\n- 372. Alimi, J.-M. et al. DEUS Full Observable LCDM Universe Simulation: the numerical challenge. arXiv e-prints arXiv:1206.2838 (2012). 1206.2838 .\n- 373. Prada, F., Klypin, A. A., Cuesta, A. J., Betancort-Rijo, J. E. & Primack, J. Halo concentrations in the standard L cold dark matter cosmology. MNRAS 423 , 3018-3030 (2012). DOI 10.1111/j.1365-2966.2012.21007.x. 1104.5130 .\n- 374. Ishiyama, T. et al. The n 2 GC simulations: Quantifying the dark side of the universe in the Planck cosmology. PASJ 67 , 61 (2015). DOI 10.1093/pasj/psv021. 1412.2860 .\n- 375. Diemand, J. et al. Clumps and streams in the local dark matter distribution. Nat. 454 , 735-738 (2008). DOI 10.1038/nature07153. 0805.1244 .\n- 376. Stadel, J. et al. Quantifying the heart of darkness with GHALO - a multibillion particle simulation of a galactic halo. MNRAS 398 , L21-L25 (2009). DOI 10.1111/j.1745-3933.2009.00699.x. 0808.2981 .\n- 377. Libeskind, N. I. et al. Constrained simulations of the Local Group: on the radial distribution of substructures. MNRAS 401 , 1889-1897 (2010). DOI 10.1111/j.1365-2966.2009.15766.x. 0909.4423 .\n- 378. Garrison-Kimmel, S., Boylan-Kolchin, M., Bullock, J. S. & Lee, K. ELVIS: Exploring the Local Volume in Simulations. MNRAS 438 , 2578-2596 (2014). DOI 10.1093/mnras/stt2377. 1310.6746 .\n- 379. Hellwing, W. A. et al. The Copernicus Complexio: a high-resolution view of the small-scale Universe. MNRAS 457 , 3492-3509 (2016). DOI 10.1093/mnras/stw214. 1505.06436 .\n- 380. Dubois, Y. et al. Dancing in the dark: galactic properties trace spin swings along the cosmic web. MNRAS 444 , 1453-1468 (2014). DOI 10.1093/mnras/stu1227. 1402.1165 .\n- 381. Khandai, N. et al. The MassiveBlack-II simulation: the evolution of haloes and galaxies to z ∼ 0. MNRAS 450 , 1349-1374 (2015). DOI 10.1093/mnras/stv627. 1402.0888 .\n- 382. Feng, Y. et al. The BlueTides simulation: first galaxies and reionization. MNRAS 455 , 27782791 (2016). DOI 10.1093/mnras/stv2484. 1504.06619 .\n- 383. Dav'e, R., Thompson, R. & Hopkins, P. F. MUFASA: galaxy formation simulations with meshless hydrodynamics. MNRAS 462 , 3265-3284 (2016). DOI 10.1093/mnras/stw1862. 1604.01418 .\n- 384. McCarthy, I. G., Schaye, J., Bird, S. & Le Brun, A. M. C. The BAHAMAS project: calibrated hydrodynamical simulations for large-scale structure cosmology. MNRAS 465 , 2936-2965 (2017). DOI 10.1093/mnras/stw2792. 1603.02702 .\n- 385. Tremmel, M. et al. The Romulus cosmological simulations: a physical approach to the formation, dynamics and accretion models of SMBHs. MNRAS 470 , 1121-1139 (2017). DOI 10.1093/mnras/stx1160. 1607.02151 .\n- 386. Ceverino, D. et al. Radiative feedback and the low efficiency of galaxy formation in lowmass haloes at high redshift. MNRAS 442 , 1545-1559 (2014). DOI 10.1093/mnras/stu956. 1307.0943 .\n- 387. Sawala, T. et al. The APOSTLE simulations: solutions to the Local Group's cosmic puzzles. MNRAS 457 , 1931-1943 (2016). DOI 10.1093/mnras/stw145. 1511.01098 .\n- 388. Barnes, D. J. et al. The redshift evolution of massive galaxy clusters in the MACSIS simulations. MNRAS 465 , 213-233 (2017). DOI 10.1093/mnras/stw2722. 1607.04569 .\n- 389. Cui, W. et al. The Three Hundred project: a large catalogue of theoretically modelled galaxy clusters for cosmological and astrophysical applications. MNRAS 480 , 2898-2915 (2018). DOI 10.1093/mnras/sty2111. 1809.04622 .\n- 390. Henden, N. A., Puchwein, E., Shen, S. & Sijacki, D. The FABLE simulations: a feedback model for galaxies, groups, and clusters. MNRAS 479 , 5385-5412 (2018). DOI 10.1093/mnras/sty1780. 1804.05064 .\n- 391. Tremmel, M. et al. Introducing ROMULUSC: a cosmological simulation of a galaxy cluster with an unprecedented resolution. MNRAS 483 , 3336-3362 (2019). DOI 10.1093/mnras/sty3336. 1806.01282 .\n- 392. Bullock, J. S. & Boylan-Kolchin, M. Small-Scale Challenges to the L CDMParadigm. ARA&A 55 , 343-387 (2017). DOI 10.1146/annurev-astro-091916-055313. 1707.04256 .\n- 393. Klypin, A., Kravtsov, A. V., Valenzuela, O. & Prada, F. Where Are the Missing Galactic Satellites? ApJ 522 , 82-92 (1999). DOI 10.1086/307643. astro-ph/9901240 .\n- 394. Zavala, J. et al. The Velocity Function in the Local Environment from L CDM and L WDM Constrained Simulations. ApJ 700 , 1779-1793 (2009). DOI 10.1088/0004-637X/700/2/1779. 0906.0585 .\n- 395. Papastergis, E., Martin, A. M., Giovanelli, R. & Haynes, M. P. The Velocity Width Function of Galaxies from the 40% ALFALFA Survey: Shedding Light on the Cold Dark Matter Overabundance Problem. ApJ 739 , 38 (2011). DOI 10.1088/0004-637X/739/1/38. 1106.0710 .\n- 396. Klypin, A., Karachentsev, I., Makarov, D. & Nasonova, O. Abundance of field galaxies. MNRAS 454 , 1798-1810 (2015). DOI 10.1093/mnras/stv2040. 1405.4523 .\n- 397. de Blok, W. J. G. & McGaugh, S. S. The dark and visible matter content of low surface brightness disc galaxies. MNRAS 290 , 533-552 (1997). DOI 10.1093/mnras/290.3.533. astro-ph/9704274 .\n- 398. Walker, M. G. & Pe˜narrubia, J. A Method for Measuring (Slopes of) the Mass Profiles of Dwarf Spheroidal Galaxies. ApJ 742 , 20 (2011). DOI 10.1088/0004-637X/742/1/20. 1108.2404 .\n- 399. Boylan-Kolchin, M., Bullock, J. S. & Kaplinghat, M. Too big to fail? The puzzling darkness of massive Milky Way subhaloes. MNRAS 415 , L40-L44 (2011). DOI 10.1111/j.17453933.2011.01074.x. 1103.0007 .\n- 400. Papastergis, E., Giovanelli, R., Haynes, M. P. & Shankar, F. Is there a 'too big to fail' problem in the field? A&A 574 , A113 (2015). DOI 10.1051/0004-6361/201424909. 1407.4665 .\n- 401. Oman, K. A. et al. The unexpected diversity of dwarf galaxy rotation curves. MNRAS 452 , 3650-3665 (2015). DOI 10.1093/mnras/stv1504. 1504.01437 .\n- 402. Navarro, J. F., Eke, V. R. & Frenk, C. S. The cores of dwarf galaxy haloes. MNRAS 283 , L72-L78 (1996). DOI 10.1093/mnras/283.3.L72. astro-ph/9610187 .\n- 403. Governato, F. et al. Cuspy no more: how outflows affect the central dark matter and baryon distribution in L cold dark matter galaxies. MNRAS 422 , 1231-1240 (2012). DOI 10.1111/j.13652966.2012.20696.x. 1202.0554 .\n- 404. O˜norbe, J. et al. Forged in FIRE: cusps, cores and baryons in low-mass dwarf galaxies. MNRAS 454 , 2092-2106 (2015). DOI 10.1093/mnras/stv2072. 1502.02036 .\n- 405. Chan, T. K. et al. The impact of baryonic physics on the structure of dark matter haloes: the view from the FIRE cosmological simulations. MNRAS 454 , 2981-3001 (2015). DOI 10.1093/mnras/stv2165. 1507.02282 .\n- 406. Read, J. I., Agertz, O. & Collins, M. L. M. Dark matter cores all the way down. MNRAS 459 , 2573-2590 (2016). DOI 10.1093/mnras/stw713. 1508.04143 .\n- 407. Green, A. M., Hofmann, S. & Schwarz, D. J. The power spectrum of SUSY-CDM on subgalactic scales. MNRAS 353 , L23-L27 (2004). DOI 10.1111/j.1365-2966.2004.08232.x. astro-ph/0309621 .\n- 408. Bode, P., Ostriker, J. P. & Turok, N. Halo Formation in Warm Dark Matter Models. ApJ 556 , 93-107 (2001). DOI 10.1086/321541. astro-ph/0010389 .\n- 409. Carlson, E. D., Machacek, M. E. & Hall, L. J. Self-interacting dark matter. ApJ 398 , 43-52 (1992). DOI 10.1086/171833.\n- 410. Spergel, D. N. & Steinhardt, P. J. Observational Evidence for Self-Interacting Cold Dark Matter. Phys. Rev. Lett. 84 , 3760-3763 (2000). DOI 10.1103/PhysRevLett.84.3760. astro-ph/ 9909386 .\n- 411. Vogelsberger, M., Sijacki, D., Kereˇs, D., Springel, V. & Hernquist, L. Moving mesh cosmology: numerical techniques and global statistics. MNRAS 425 , 3024-3057 (2012). DOI 10.1111/j.1365-2966.2012.21590.x. 1109.1281 .\n- 412. Cyr-Racine, F.-Y. et al. ETHOS - an effective theory of structure formation: From dark particle physics to the matter distribution of the Universe. Phys. Rev. D 93 , 123527 (2016). DOI 10.1103/PhysRevD.93.123527. 1512.05344 .\n- 413. Vogelsberger, M. et al. ETHOS-aneffective theory of structure formation: dark matter physics as a possible explanation of the small-scale CDM problems. MNRAS 460 , 1399-1416 (2016). DOI 10.1093/mnras/stw1076. 1512.05349 .\n- 414. Vogelsberger, M., Zavala, J. & Loeb, A. Subhaloes in self-interacting galactic dark matter haloes. MNRAS 423 , 3740-3752 (2012). DOI 10.1111/j.1365-2966.2012.21182.x. 1201. 5892 .  \n- 415. Peter, A. H. G., Rocha, M., Bullock, J. S. & Kaplinghat, M. Cosmological simulations with self-interacting dark matter - II. Halo shapes versus observations. MNRAS 430 , 105-120 (2013). DOI 10.1093/mnras/sts535. 1208.3026 .\n- 416. Fry, A. B. et al. All about baryons: revisiting SIDM predictions at small halo masses. MNRAS 452 , 1468-1479 (2015). DOI 10.1093/mnras/stv1330. 1501.00497 .\n- 417. Elbert, O. D. et al. Core formation in dwarf haloes with self-interacting dark matter: no finetuning necessary. MNRAS 453 , 29-37 (2015). DOI 10.1093/mnras/stv1470. 1412.1477 .\n- 418. Hui, L., Ostriker, J. P., Tremaine, S. & Witten, E. Ultralight scalars as cosmological dark matter. Phys. Rev. D 95 , 043541 (2017). DOI 10.1103/PhysRevD.95.043541. 1610.08297 .\n- 419. Lee, J.-W. & Koh, I.-G. Galactic halos as boson stars. Phys. Rev. D 53 , 2236-2239 (1996). DOI 10.1103/PhysRevD.53.2236. hep-ph/9507385 .\n- 420. Hu, W., Barkana, R. & Gruzinov, A. Fuzzy Cold Dark Matter: The Wave Properties of Ultralight Particles. Phys. Rev. Lett. 85 , 1158-1161 (2000). DOI 10.1103/PhysRevLett.85.1158. astro-ph/0003365 .\n- 421. Peebles, P. J. E. Fluid Dark Matter. ApJ 534 , L127-L129 (2000). DOI 10.1086/312677. astro-ph/0002495 .\n- 422. Chavanis, P.-H. Mass-radius relation of Newtonian self-gravitating Bose-Einstein condensates with short-range interactions. I. Analytical results. Phys. Rev. D 84 , 043531 (2011). DOI 10.1103/PhysRevD.84.043531. 1103.2050 .\n- 423. Su'arez, A., Robles, V. H. & Matos, T. A Review on the Scalar Field/Bose-Einstein Condensate Dark Matter Model. In Moreno Gonz'alez, C., Madriz Aguilar, J. E. & Reyes Barrera, L. M. (eds.) Accelerated Cosmic Expansion , vol. 38 of Astrophysics and Space Science Proceedings , 107 (2014). DOI 10.1007/978-3-319-02063-1 9. 1302.0903 .\n- 424. Matos, T., V'azquez-Gonz'alez, A. & Maga˜na, J. phi squared as dark matter. MNRAS 393 , 1359-1369 (2009). DOI 10.1111/j.1365-2966.2008.13957.x. 0806.0683 .\n- 425. Lundgren, A. P., Bondarescu, M., Bondarescu, R. & Balakrishna, J. Lukewarm Dark Matter: Bose Condensation of Ultralight Particles. ApJ 715 , L35-L39 (2010). DOI 10.1088/20418205/715/1/L35. 1001.0051 .\n- 426. Robles, V. H. & Matos, T. Exact Solution to Finite Temperature SFDM: Natural Cores without Feedback. ApJ 763 , 19 (2013). DOI 10.1088/0004-637X/763/1/19. 1207.5858 .\n- 427. Seidel, E. & Suen, W.-M. Dynamical evolution of boson stars: Perturbing the ground state. Phys. Rev. D 42 , 384-403 (1990). DOI 10.1103/PhysRevD.42.384.\n- 428. Sin, S.-J. Late-time phase transition and the galactic halo as a Bose liquid. Phys. Rev. D 50 , 3650-3654 (1994). DOI 10.1103/PhysRevD.50.3650. hep-ph/9205208 .\n- 429. Mocz, P. & Succi, S. Numerical solution of the nonlinear Schrodinger equation using smoothedparticle hydrodynamics. Phys. Rev. D 91 , 053304 (2015). DOI 10.1103/PhysRevE.91.053304. 1503.03869 .\n- 430. Nori, M. & Baldi, M. AX-GADGET: a new code for cosmological simulations of Fuzzy Dark Matter and Axion models. MNRAS 478 , 3935-3951 (2018). DOI 10.1093/mnras/sty1224. 1801.08144 .\n- 431. Copeland, E. J., Sami, M. & Tsujikawa, S. Dynamics of Dark Energy. Int. J. Mod. Phys. D 15 , 1753-1935 (2006). DOI 10.1142/S021827180600942X. hep-th/0603057 .\n- 432. Baldi, M. Dark Energy simulations. Phys. Dark Universe 1 , 162-193 (2012). DOI 10.1016/j.dark.2012.10.004. 1210.6650 .\n- 433. Linder, E. V. & Jenkins, A. Cosmic structure growth and dark energy. MNRAS 346 , 573-583 (2003). DOI 10.1046/j.1365-2966.2003.07112.x. astro-ph/0305286 .\n- 434. Grossi, M. & Springel, V. The impact of early dark energy on non-linear structure formation. MNRAS 394 , 1559-1574 (2009). DOI 10.1111/j.1365-2966.2009.14432.x. 0809.3404 .\n- 435. Penzo, C., Macci'o, A. V., Casarini, L., Stinson, G. S. & Wadsley, J. Dark MaGICC: the effect of dark energy on disc galaxy formation. Cosmology does matter. MNRAS 442 , 176-186 (2014). DOI 10.1093/mnras/stu857. 1401.3338 .\n- 436. Sefusatti, E. & Vernizzi, F. Cosmological structure formation with clustering quintessence. J. Cosmol. Astro-Particle Phys. 2011 , 047 (2011). DOI 10.1088/1475-7516/2011/03/047. 1101.1026 .\n- 437. Amendola, L. & Tocchini-Valentini, D. Stationary dark energy: The present universe as a global attractor. Phys. Rev. D 64 , 043509 (2001). DOI 10.1103/PhysRevD.64.043509. astro-ph/0011243 .\n- 438. Macci'o, A. V., Quercellini, C., Mainini, R., Amendola, L. & Bonometto, S. A. Coupled dark energy: Parameter constraints from N-body simulations. Phys. Rev. D 69 , 123516 (2004). DOI 10.1103/PhysRevD.69.123516. astro-ph/0309671 .\n- 439. Baldi, M., Pettorino, V., Robbers, G. & Springel, V. Hydrodynamical N-body simulations of coupled dark energy cosmologies. MNRAS 403 , 1684-1702 (2010). DOI 10.1111/j.13652966.2009.15987.x. 0812.3901 .\n- 440. Li, B. & Barrow, J. D. N-body simulations for coupled scalar-field cosmology. Phys. Rev. D 83 , 024007 (2011). DOI 10.1103/PhysRevD.83.024007. 1005.4231 .\n- 441. Li, B. & Barrow, J. D. On the effects of coupled scalar fields on structure formation. MNRAS 413 , 262-270 (2011). DOI 10.1111/j.1365-2966.2010.18130.x. 1010.3748 .\n- 442. Milgrom, M. A modification of the Newtonian dynamics as a possible alternative to the hidden mass hypothesis. ApJ 270 , 365-370 (1983). DOI 10.1086/161130.\n- 443. Famaey, B. & McGaugh, S. S. Modified Newtonian Dynamics (MOND): Observational Phenomenology and Relativistic Extensions. Living Rev. Relativ. 15 , 10 (2012). DOI 10.12942/lrr2012-10. 1112.3960 .\n- 444. Bekenstein, J. D. Relativistic gravitation theory for the modified Newtonian dynamics paradigm. Phys. Rev. D 70 , 083509 (2004). DOI 10.1103/PhysRevD.70.083509. astro-ph/ 0403694 .\n- 445. Skordis, C. Generalizing tensor-vector-scalar cosmology. Phys. Rev. D 77 , 123502 (2008). DOI 10.1103/PhysRevD.77.123502. 0801.1985 .\n- 446. McGaugh, S. S., Lelli, F. & Schombert, J. M. Radial Acceleration Relation in Rotationally Supported Galaxies. Phys. Rev. Lett. 117 , 201101 (2016). DOI 10.1103/PhysRevLett.117.201101. 1609.05917 .\n- 447. Lelli, F., McGaugh, S. S., Schombert, J. M. & Pawlowski, M. S. One Law to Rule Them All: The Radial Acceleration Relation of Galaxies. ApJ 836 , 152 (2017). DOI 10.3847/15384357/836/2/152. 1610.08981 .\n- 448. Keller, B. W. & Wadsley, J. W. L CDMisConsistent with SPARC Radial Acceleration Relation. ApJ 835 , L17 (2017). DOI 10.3847/2041-8213/835/1/L17. 1610.06183 .\n- 449. Ludlow, A. D. et al. Mass-Discrepancy Acceleration Relation: A Natural Outcome of Galaxy Formation in Cold Dark Matter Halos. Phys. Rev. Lett. 118 , 161103 (2017). DOI 10.1103/PhysRevLett.118.161103. 1610.07663 .\n- 450. Dutton, A. A., Macci'o, A. V., Obreja, A. & Buck, T. NIHAO - XVIII. Origin of the MOND phenomenology of galactic rotation curves in a L CDM universe. MNRAS 485 , 1886-1899 (2019). DOI 10.1093/mnras/stz531. 1902.06751 .\n- 451. Llinares, C. Simulation techniques for modified gravity. Int. J. Mod. Phys. D 27 , 1848003 (2018). DOI 10.1142/S0218271818480036.\n- 452. Li, B., Zhao, G.-B., Teyssier, R. & Koyama, K. ECOSMOG: an Efficient COde for Simulating MOdified Gravity. J. Cosmol. Astropart. Phys. 1 , 051 (2012). DOI 10.1088/14757516/2012/01/051. 1110.1379 .\n- 453. Puchwein, E., Baldi, M. & Springel, V. Modified-Gravity-GADGET: a new code for cosmological hydrodynamical simulations of modified gravity models. MNRAS 436 , 348-360 (2013). DOI 10.1093/mnras/stt1575. 1305.2418 .\n- 454. Llinares, C., Mota, D. F. & Winther, H. A. ISIS: a new N-body cosmological code with scalar fields based on RAMSES. Code presentation and application to the shapes of clusters. A&A 562 , A78 (2014). DOI 10.1051/0004-6361/201322412. 1307.6748 .\n- 455. Brax, P., Davis, A.-C., Li, B., Winther, H. A. & Zhao, G.-B. Systematic simulations of modified gravity: symmetron and dilaton models. J. Cosmol. Astro-Particle Phys. 2012 , 002 (2012). DOI 10.1088/1475-7516/2012/10/002. 1206.3568 .\n- 456. Barreira, A., Li, B., Hellwing, W. A., Baugh, C. M. & Pascoli, S. Nonlinear structure formation in the cubic Galileon gravity model. J. Cosmol. Astro-Particle Phys. 2013 , 027 (2013). DOI 10.1088/1475-7516/2013/10/027. 1306.3219 .\n- 457. Winther, H. A. et al. Modified gravity N-body code comparison project. MNRAS 454 , 42084234 (2015). DOI 10.1093/mnras/stv2253. 1506.06384 .\n- 458. Fontanot, F., Puchwein, E., Springel, V. & Bianchi, D. Semi-analytic galaxy formation in f(R)-gravity cosmologies. MNRAS 436 , 2672-2679 (2013). DOI 10.1093/mnras/stt1763. 1307.5065 .\n- 459. Jennings, E., Baugh, C. M., Li, B., Zhao, G.-B. & Koyama, K. Redshift-space distortions in f(R) gravity. MNRAS 425 , 2128-2143 (2012). DOI 10.1111/j.1365-2966.2012.21567.x. 1205.2698 .\n- 460. Naik, A. P., Puchwein, E., Davis, A.-C. & Arnold, C. Imprints of Chameleon f(R) gravity on Galaxy rotation curves. MNRAS 480 , 5211-5225 (2018). DOI 10.1093/mnras/sty2199. 1805.12221 .\n- 461. Lombriser, L. & Pe˜narrubia, J. How chameleons core dwarfs with cusps. Phys. Rev. D 91 , 084022 (2015). DOI 10.1103/PhysRevD.91.084022. 1407.7862 .\n- 462. Arnold, C., Leo, M. & Li, B. Realistic simulations of galaxy formation in f(R) modified gravity. Nat. Astron. 386 (2019). DOI 10.1038/s41550-019-0823-y. 1907.02977 .\n- 463. Davis, A.-C., Lim, E. A., Sakstein, J. & Shaw, D. J. Modified gravity makes galaxies brighter. Phys. Rev. D 85 , 123006 (2012). DOI 10.1103/PhysRevD.85.123006. 1102.5278 ."}