-
PDF
- Split View
-
Views
-
Cite
Cite
Daniele Sorini, John A Peacock, Lucas Lombriser, The impact of the cosmological constant on past and future star formation, Monthly Notices of the Royal Astronomical Society, Volume 535, Issue 2, December 2024, Pages 1449–1474, https://doi.org/10.1093/mnras/stae2236
- Share Icon Share
ABSTRACT
We present an extended analytical model for cosmic star formation, with the aim of investigating the impact of cosmological parameters on the star formation history within the |$\Lambda$|CDM paradigm. Constructing an ensemble of flat |$\Lambda$|CDM models where the cosmological constant varies between |$\Lambda = 0$| and |$10^5$| times the observed value, |$\Lambda _{\rm obs}$|, we find that the fraction of cosmic baryons that are converted into stars over the entire history of the universe peaks at |$\sim$| 27 per cent for |$0.01 \lesssim \Lambda /\Lambda _{\rm obs} \lesssim 1$|. We explain, from first principles, that the decline of this asymptotic star formation efficiency for lower and higher values of |$\Lambda$| is driven, respectively, by the astrophysics of star formation, and by the suppression of cosmic structure formation. However, the asymptotic efficiency declines slowly as |$\Lambda$| increases, falling below 5 per cent only for |$\Lambda \gt 100 \, \Lambda _{\rm obs}$|. Making the minimal assumption that the probability of generating observers is proportional to this efficiency, and following Weinberg in adopting a flat prior on |$\Lambda$|, the median posterior value of |$\Lambda$| is |$539 \, \Lambda _{\rm obs}$|. Furthermore, the probability of observing |$\Lambda \le \Lambda _{\rm obs}$| is only 0.5 per cent. Although this work has not considered recollapsing models with |$\Lambda \lt 0$|, the indication is thus that |$\Lambda _{\rm obs}$| appears to be unreasonably small compared to the predictions of the simplest multiverse ensemble. This poses a challenge for anthropic reasoning as a viable explanation for cosmic coincidences and the apparent fine-tuning of the Universe: either the approach is invalid or more parameters than |$\Lambda$| alone must vary within the ensemble.
1 INTRODUCTION
Galaxies are the most obvious tracers of large-scale structure in the Universe, and understanding how they emerge and evolve throughout cosmic history has been a longstanding goal of cosmology and astrophysics. The |$\Lambda$|CDM paradigm approaches this via the gravitational collapse of dark matter haloes from an initial quasi-homogeneous Universe. This process is well understood, thanks to early analytical models of halo assembly (Lacey & Cole 1993), which were subsequently validated numerically with full N-body simulations (Springel et al. 2005; Klypin, Trujillo-Gomez & Primack 2011; Angulo et al. 2012; Fosalba et al. 2015). However, many of the detailed baryonic processes that govern the build-up of galaxies, such as gas accretion, star formation, and outflows, are not yet fully understood. A successful theory of galaxy formation needs to account for these processes in order to predict the main observed properties of galaxies. In particular, reproducing the observed efficiency of star formation, both locally within individual galaxies and globally over a cosmologically representative volume, constitutes a crucial test for any model of galaxy formation (see the review by Madau & Dickinson 2014). In this paper, we consider predictions for this efficiency and how it depends on cosmological parameters, specifically the cosmological constant.
Early fully analytical models of star formation were based on simple prescriptions for the cooling time for the gas within galaxies and the typical time-scale for converting it into stars (e.g. Hernquist & Springel 2003), while following the growth of structures with well-established analytical forms for the halo mass function (Press & Schechter 1974; Sheth & Tormen 1999, 2002). However, these insightful models neglected more sophisticated mechanisms such as feedback processes driven by stars or active galactic nuclei (AGNs). With a slightly different approach, White & Frenk (1991) implemented the effect of subgalactic physics through a series of approximate formulae, while keeping the treatment of structure formation nearly analytical. This seminal work set the stage for a sequence of refinements of the modelling (e.g. Kauffmann, White & Guiderdoni 1993; Cole et al. 1994, 2000; Guiderdoni et al. 1998; Kauffmann et al. 1999 ). Further extensions included the assembly of the central black hole, which enabled a description of the co-evolution of galaxies and quasars (e.g. Kauffmann & Haehnelt 2000; Somerville et al. 2008; Henriques et al. 2015; Lacey et al. 2016). Other semi-analytical techniques followed the formation of dark matter haloes in full N-body simulations, coupled with analytical recipes for the baryonic physics (e.g. Croton et al. 2006; Henriques et al. 2020).
Full hydrodynamical simulations incorporate several physical processes into the modelling and follow the evolution of baryons and dark matter from first principles. Several cosmological simulations (e.g. Schaye et al. 2010; Dubois et al. 2014; Hopkins et al. 2014; Vogelsberger et al. 2014; Lukić et al. 2015; Schaye et al. 2015; Davé, Thompson & Hopkins 2016; McCarthy et al. 2017; Pillepich et al. 2018; Davé et al. 2019), often based on different numerical approaches (Springel, Yoshida & White 2001a; Springel et al. 2005, 2021; Almgren et al. 2013; Bryan et al. 2014; Hopkins 2015 ), generally managed to reproduce a plethora of observations to a satisfactory level of accuracy. However, feedback processes are implemented via various different numerical ‘subgrid’ prescriptions, whose parameters are tuned to reproduce certain observations. It is therefore important to ask if the predictions of these calculations are robust, or whether they are fine-tuned to our Universe and have a high sensitivity to model parameters. While a substantial body of literature has focused on the effect on star formation history of varying the subgrid parameters (see the review by Somerville & Davé 2015), the impact of changing the cosmological parameters has historically received less attention. This may partially reflect the tremendous progress in constraining the parameters of the |$\Lambda$|CDM paradigm (e.g. Planck Collaboration VI 2020; but see also Verde, Treu & Riess 2019 for a discussion on recent tensions). However, studying the predictions of hydrodynamical simulations regarding the star formation history in counter-factual cosmological models would constitute an interesting ‘stress-test’ of our understanding of galaxy formation.
There is also a more fundamental reason to explore this question. Although the |$\Lambda$|CDM model is highly successful in many ways, it has theoretical issues that are hard to ignore (see e.g. the review by Bull et al. 2016). While the cosmological constant |$\Lambda$| is required to explain the accelerating expansion of the universe, there is no consensus on its physical meaning. If |$\Lambda$| is a manifestation of the energy of the quantum vacuum, then one can estimate its value by integrating over the zero-point energy of all possible modes. Adopting a cut-off energy |$E_{\rm v}$| yields a vacuum density of |$E_{\rm v}^4$| in natural units; but the observed |$\Lambda$| requires a cut-off at rather low energies, |$E_{\rm v}\sim 1\, {\rm meV}$|, in gross conflict with the scale of new physics at perhaps 10 TeV or above (Martin 2012). This common calculation is in fact deeply flawed as it is non-relativistic and yields the wrong vacuum equation of state (Koksma & Prokopec 2011): a more sophisticated calculation yields a vacuum density of order |$M^4$|, in terms of the particle mass. However, since particles exist with M near to the TeV scale (0.17 TeV for the top quark), the discrepancy in the estimated |$\Lambda$| is much the same as in the naive approach. For a more detailed discussion of this ‘cosmological constant problem’ see, for example, Weinberg (2000a) or Abel, Bryan & Norman (2002).
Many attempts have been made to move beyond the idea of |$\Lambda$| representing a simple vacuum density, and to hypothesize some more general ‘dark energy’ contributing to an effective |$\Lambda$| that may vary with time. Ratra & Peebles (1988) proposed that a dynamical scalar field could cause the accelerating expansion of the universe. While the scalar field proposed by Ratra & Peebles (1988) was not motivated by fundamental physics, later models tried to connect it to extensions of the standard model of particle physics. In practice, the kinetic or potential terms of the Lagrangian of the scalar field depend on some fundamental mass scale (e.g. Zlatev, Wang & Steinhardt 1999; Armendariz-Picon, Mukhanov & Steinhardt 2001). An alternative view is that the cosmic acceleration in fact shows the need for some modified theory of gravity (see e.g. Joyce, Lombriser & Schmidt 2016). Other approaches include mechanisms that prevent vacuum energy from gravitating (Kaloper & Padilla 2014). However, all these models effectively move the problem of the value of the |$\Lambda$| to the fine-tuning of some other parameter of the theory, such as the mass scale associated with the scalar field (see e.g. Amendola & Tsujikawa 2010). There is also a more radical position asserting that the debate on the physical nature of the cosmological constant is moot, and that |$\Lambda$| is simply a fundamental constant emerging within the theory of General Relativity (Bianchi & Rovelli 2010).
Regardless of the physical interpretation of the cosmological constant, the oddly small non-null value of |$\Lambda$| gives rise to a number of coincidences that cry out for an explanation. Perhaps the most well-known of these is that the Universe became |$\Lambda$|-dominated at |$z\approx 0.39$|, near to the formation time of the Sun. We therefore appear to live near the unique era when |$\Lambda$| transitions from being negligible to dominating the Universe: this is the ‘why-now problem’ (Velten, vom Marttens & Zimdahl 2014). A further time coincidence is that recombination occurred at the same time of baryon-radiation equality. In addition, Lombriser & Smer-Barreto (2017) pointed out that the equality between |$\Lambda$| and radiation occurs around the mid-point of cosmic reionization. All these time-scale puzzles are in principle distinct from the cosmological constant problem described earlier (but see also the discussion in Lombriser 2023).
A possible explanation for the ‘why now’ problem was suggested by Weinberg (1987). He noted that if the cosmological constant had been much larger than observed, the accelerating expansion of the Universe would have set in at earlier times, freezing out the growth of structure before galaxy-scale haloes had been able to form. Thus the star formation in galaxies that is necessary for the creation of observers would not occur if |$\Lambda$| was substantially larger than the observed value. Such an argument is an example of anthropic reasoning (Carter 1974), which effectively considers the existence of observers (such as ourselves) as a ‘data point’ and explores the implications for the cosmological parameters conditional on this information.
To make this Bayesian argument (probability of |$\Lambda$| given that it is observed), we need there to be some physical mechanism that allows |$\Lambda$| to vary. It is also common to invoke a multiverse: an ensemble of different universes. The probability calculus is the same whether or not the members of the ensemble actually exist, or merely have the potential to do so. However, the idea of a concrete multiverse underlying Weinberg’s argument was given stronger motivation by the theory of inflation. Here, a multiverse of causally disconnected ‘bubble universes’ arises in models of stochastic inflation where inflation proceeds eternally (Vilenkin 1983; Linde 1986; Guth 2007; Freivogel 2011). All causally disconnected bubbles evolve as independent universes, each characterized by a different set of constants, including |$\Lambda$|. The advantage of this picture is that for a sufficiently large ensemble there are guaranteed to exist universes suitable for the formation of structure. Thus, that would explain the existence of our Universe, no matter how atypical it is within the ensemble.
Anthropic arguments often encounter significant resistance, with many physicists arguing that efforts should be focused on finding solutions to cosmological puzzles from first principles (e.g. Kane, Perry & Zytkow 2002). Of course, such efforts should always be pursued. However, we can note that anthropic approaches pervade other fields of astronomy, without generating controversy. A prime example is the concept of circumstellar ‘habitable zone’, which is defined based on the conditions that can sustain life on a planet (see e.g. Kasting, Whitmire & Reynolds 1993). Of course, unlike with exoplanets, there is only one Universe that can actually be observed (although see Aguirre & Kozaczuk 2013; Wainwright et al. 2014; Johnson et al. 2016). As such, anthropic arguments cannot be tested in the Galilean sense ingrained in the scientific method. However, they can still be tested, since we can in principle predict the probability distribution of observed values of cosmological parameters over the ensemble. This prediction can be compared with the single ‘data point’ of our Universe, and a sufficiently large deviation from the average serves to rule out the model on which the prediction was based.
After the initial formulation given by Weinberg (1987), investigations of the anthropic approach became progressively more refined. Weinberg (1989) extended his argument by noting that observer selection does not require |$\Lambda$| to be exactly zero, and thus a small non-zero (in principle even negative) value of |$\Lambda$| would be anthropically predicted. A critical element of this argument is the idea that the prior distribution of |$\Lambda$| should be flat because |$\Lambda =0$| is not a spatial point, and therefore any continuous prior must be treatable as a constant in a narrow range near zero. Efstathiou (1995) revisited the argument by making more detailed estimates of the abundance of galaxies in universes with |$\Lambda \ge 0$|, subject to the constraint that the observed temperature of the cosmic microwave background (CMB) is equal to |$2.73 \, \rm K$|. He concluded that anthropic arguments may explain a value of |$\Lambda$| close to the one that is observed. Later works allowed for more than one cosmological parameter to vary. Garriga, Livio & Vilenkin (1999) considered simultaneous variations of |$\Lambda$| and the density contrast at the time of recombination. Peacock (2007) explored anthropic arguments treating both the CMB temperature and |$\Lambda$| as free parameters, and further allowing |$\Lambda$| to assume negative values. Bousso & Leichenauer (2010) also considered the variation of multiple cosmological parameters, as well as negative values of |$\Lambda$|. In their work, the generation of observers is tied to the global efficiency of star formation in the different members of the multiverse ensemble. The star formation history is predicted with a semi-analytical model, under the assumption that the astrophysics of star formation does not vary throughout the multiverse (Bousso & Leichenauer 2009). Other semi-analytical work by Sudoh et al. (2017) showed that incorporating the astrophysics of galaxy formation in anthropic reasoning can affect the range of the anthropically favoured values of |$\Lambda$| by almost an order of magnitude.
More recently, progress in numerical power has enabled the testing of anthropic reasoning with full hydrodynamical simulations, exploring the simulated past and future star formation history for a wide range of |$\Lambda$| (Barnes et al. 2018; Salcido et al. 2018; Oh et al. 2022). While such simulations include several astrophysical processes (albeit in an approximate parametrized form), it is hard to probe a wide parameter space, or to explore the far future of the universe (beyond |$\sim$| 100 Gyr cosmic time), without a massive commitment of computational resources.
Analytic models of star formation therefore represent an attractive complementary approach as they are not subject to the same limitations as hydrodynamical simulations. While inevitably simplified in terms of astrophysics, they are an efficient technique that can offer a more intuitive picture of the evolution of star formation (e.g. Rasera & Teyssier 2006; Davé, Finlator & Oppenheimer 2012; Sharma & Theuns 2019; Salcido, Bower & Theuns 2020; Fukugita & Kawasaki 2022). Recently, Sorini & Peacock (2021) generalized the Hernquist & Springel (2003) model of cosmic star formation such that it can be applied to arbitrarily large times. In this work, we further adapt the Sorini & Peacock (2021) formalism to explore the impact of different values of |$\Lambda$| on past and future star formation. The main result of our study will be that values of |$\Lambda$| in the range |$0.01\lesssim \Lambda / \Lambda _{\rm obs} \lesssim \Lambda _{\rm obs}$| maximize the global efficiency of star formation. However, the interesting question is how effectively star formation is truncated by large values of |$\Lambda$|. If this suppression is not effective until |$\Lambda$| is vastly beyond the observed value, then Weinberg’s flat prior will mean that the observed Universe risks being an implausibly rare outlier. The main aim of this paper is to quantify just how unusual our Universe is in this respect.
This manuscript is organized as follows. In Section 2, we give an overview of the formalism in Sorini & Peacock (2021) and improve it by introducing extra features. In Section 3, we explain how we adapt it to |$\Lambda$|CDM models with arbitrary non-negative values of |$\Lambda$| and discuss the impact of changing the cosmological constant on the halo mass function and on the efficiency of star formation within haloes. In Section 4, we present our predictions for the long-term efficiency of star formation for different values of |$\Lambda$|. In Section 5, we discuss the implications of our results for anthropic reasoning. We mainly focus on the cosmological constant problem but also marginally consider the why-now problem. We also compare our results with previous literature on the subject and discuss the limitations of our model. Finally, in Section 6 we summarize the main conclusions of our work and discuss the future developments of our line of research.
Throughout this work, unless otherwise indicated, units of distance are understood to be proper units. Comoving units are designated with a ‘c’ prefix (e.g. ckpc, cMpc).
2 FORMALISM
In this work, we aim to exploit the analytical model for cosmic star formation developed by Sorini & Peacock (2021) – hereafter SP21. We therefore summarize the main aspects of the formalism in Section 2.1. In Section 2.2, we will show how we extended the SP21 model to obtain greater accuracy in the predictions of the late-time behaviour of star formation. As we will explain, this generalization will be crucial in answering the main scientific questions addressed in this work.
2.1 Summary of the SP21 model
The SP21 model predicts the cosmic star formation rate density (CSFRD) from first principles, generalizing the seminal work by Hernquist & Springel (2003) – hereafter HS03. The basic idea of the formalism is that the CSFRD is obtained by integrating the star formation rate (SFR) in all haloes within a given comoving volume, weighted by the halo multiplicity function:
where the |$\bar{\rho }_0$| is the comoving mean matter density of the universe and |$\mathrm{ d}F(M, \, z)/\mathrm{ d}\ln M$| is the halo multiplicity function, with |$F(M, \, z)$| being the collapsed mass fraction in haloes with total mass |$\gt M$|. These quantities encode the impact of background cosmology and of the growth of large-scale structure on the CSFRD. The specifically astrophysical component of the CSFRD is encapsulated in the term |$s(M, \, z)={\rm SFR}/M$|, which represents the average SFR over a population of haloes of mass M at redshift z, normalized by the total halo mass M. We will refer to |$s(M, \, z)$| as ‘normalized SFR’ (nSFR).
In this work, we follow SP21’s choice of modelling |$F(M, \, z)$| via the Sheth–Tormen formalism (Sheth & Tormen 1999, 2002). We also adopt the same definitions of the virial radius, mass, and temperature as in SP21. The virial radius R of a halo of virial mass M at redshift z is defined as the radius of the sphere that contains a matter density equal to |$\Delta \, \rho _{\rm c}(z)$|, where |$\rho _{\rm c}(z)$| is the critical density of the universe and |$\Delta$| a suitable multiplying factor:
SP21 adopted |$\Delta =200$|. One can then define a characteristic virial velocity
where G is the universal gravitational constant. We additionally define the virial temperature T such that
where |$k_{\rm B}$| is the Boltzmann constant and |$\mu$| is the mean molecular weight. We assume that |$\mu \approx 0.6 m_{\rm p}$|, which is a valid approximation for a fully ionized plasma of primordial composition. From equations (2)–(4), it follows that
Although the CSFRD in equation (1) is expressed as an integral over the halo mass, it may be at times more convenient to switch the integration variable to the virial temperature. In particular, within the HS03 formalism, the nSFR is more naturally expressed for a population of haloes at given T, as we shall now explain.
At any given z, the nSFR is set by a characteristic time-scale. At high redshift, haloes are denser and hence gas cooling is more efficient. The bottleneck of star formation is therefore represented by the typical time-scale that governs the conversion of cool gas into stars. This ‘average gas consumption time-scale’, |$\langle t_* \rangle$|, is assumed to be a constant, with no dependence on the virial temperature of the halo or on the cosmological epoch. This assumption was originally introduced by HS03 following the results of full hydrodynamical cosmological simulations (Springel & Hernquist 2003a, b), but really asserts that |$\langle t_* \rangle$| is a time-scale set by the local microphysics of molecular clouds. The same prescription was adopted by SP21, who determined a physically reasonable range of |$\langle t_* \rangle$| by matching their predictions of the Kennicutt–Schmidt law (Kennicutt 1998) with measurements by Genzel et al. (2010). The value of |$\langle t_* \rangle$| that best reproduces the observations was found to be |$2.39\, \rm Gyr$|, but values in the range |$(1.63-3.87)\, \rm Gyr$| are consistent with the data within |$3\sigma$|. SP21 showed that values of |$\langle t_* \rangle$| within this range can predict the CSFRD at high redshift within a factor of two, and we adopt the same prescription for the average high-redshift nSFR of a halo population at given virial temperature:
where x is the fraction of cold gas clouds, |$\beta$| is the mass fraction of massive (|$\gt 8 \, \rm {\rm M}_{\odot }$|) short-lived stars, and |$f_{\rm gas}(T,\, z)$| is the mass fraction of gas within the haloes considered. Following the same choice as in HS03 and SP21, we set |$x=0.95$|. This value is again motivated by the subgrid model for star formation by Springel & Hernquist (2003b), which was utilized in hydrodynamical simulations by Springel & Hernquist (2003a). As in SP21, we set |$\beta = 0.21$|. This is determined by assuming a Chabrier (2003) initial mass function (IMF) over a stellar mass range between |$0.1 \,$| and |$80 \, \rm {\rm M}_{\odot }$|.
At low redshift, haloes are less dense and gas cooling thus proceeds more slowly. The nSFR is then no longer dictated by an internal time-scale, but is limited by the supply of new cold gas, which is set by the cooling time-scale |$t_{\rm cool}$|. This quantity depends on the local gas density, which is assumed to follow a spherically symmetric power-law profile:
where |$M_{\rm gas}$| is the total gas mass enclosed in the halo and the slope |$\eta$| is a free parameter of the model. Clearly, one must have |$\eta \lt 3$| to prevent the halo gas mass from diverging. We also assume that |$\eta \gt 0$|, so that density falls with radius. The assumption of a power-law gas density profile is generally supported by full cosmological simulations (e.g. Sorini et al. 2024b), but more complex models may provide a more realistic physical picture (see e.g. Mathews & Prochaska 2017; Prochaska & Zheng 2019; Khrykin et al. 2024; Sorini et al. 2024a).
The cooling time at distance r from the centre of the halo is given by
where |$n_{\rm H}$| is the number density of hydrogen and |$\mathcal {C}(T)$| is the cooling function. As in HS03 and SP21, we assume a primordial cooling function (Sutherland & Dopita 1993), effectively neglecting metal cooling. However, this limitation is expected to alter the CSFRD at low redshift in a manner that does not significantly affect conclusions regarding the long-term star formation history (see the discussion in SP21 and HS03), so that metallicity evolution is unimportant for the scope of this work. The cooling function depends purely on temperature, so there is some advantage in considering the evolution of the nSFR for haloes with given T rather than given M.
The cooling rate of gas within haloes of virial temperature T at a given time is then estimated by following the expansion of a cooling front from the centre of the halo outwards. At any time t, the cooling front reaches the cooling radius |$r_{\rm cool}(t)$|, defined by the criterion |$t_{\rm cool}(r_{\rm cool}(t))=t$|. The gas mass |$M_{\rm cool}$| within |$r_{\rm cool}$| cools down and remains cool thereafter. The cooling rate is then determined by solving the equation
Here, the meaning of the time t needs to be defined quite carefully. In principle, this would be the time since the formation of the halo, i.e. since the last major merger. In the matter-dominated era, low-mass haloes have a life span comparable to the age of the universe, whereas high-mass haloes beyond the exponential cut-off survive for shorter times. However, massive haloes are rare, so one can effectively make the reasonable assumption that t (and |$t_{\rm cool}$|) are comparable to the cosmic time. However, it has been suggested that |$t_{\rm cool}$| should be of the order of the dynamical time of the halo |$t_{\rm dyn} = R/V$|, as it is on this time-scale that the gas profile reacts to pressure losses, and hence should set the extent of the cooling radius (Springel et al. 2001b). Because this prescription provides good agreement with simulations (Yoshida et al. 2002), HS03 solved equation (9) imposing |$t_{\rm cool} = t_{\rm dyn}$|. However, SP21 noted that this assumption breaks down in the far future, when the universe becomes |$\Lambda$|-dominated. In this regime, merging eventually ceases and haloes are isolated. Therefore, haloes can in principle cool down for arbitrarily large times. Thus, SP21 adopted the prescription
where
By construction, equation (10) yields |$t_{\rm cool} \approx f_{\rm dyn} t_{\rm dyn}$| for |$t \ll t_{\rm dyn}$|, where |$f_{\rm dyn}$| is a constant of order unity, and |$t_{\rm cool} \approx t$| for |$t \gg t_{\rm dyn}$|. SP21 showed that the values of the parameters |$f_{\rm dyn}$| and m have minimal impact on the predicted CSFRD. We will then adopt their fiducial values, |$f_{\rm dyn}=1$| and |$m=2$|.
In a flat |$\Lambda$|CDM universe, the correspondence between redshift and cosmic time can be obtained analytically outside the radiation-dominated era:
where |$H_0$|, |$\Omega _{\rm m}$|, and |$\Omega _{\Lambda }$| are the Hubble parameter and density parameters of matter and |$\Lambda$| at |$z=0$|, respectively. A principal interest of this paper is the impact of changing the value of |$\Lambda$|, which then alters the expansion history of the universe and changes all cosmological parameters. We discuss below how to handle these changes in a consistent fashion. Combining equation (12) with equations (9) and (10), one can then obtain the cooling rate as a function of the virial temperature and redshift. Using the cooling rate as a proxy for the low-redshift nSFR, SP21 found
where
and |$\tilde{S}(T)$| is a temperature-dependent proportionality factor with the units of an SFR per unit mass.
Clearly, to determine the nSFR both at high and low redshift, it is crucial to have an expression for |$f_{\rm gas}(T, \, z)$|. SP21 adopt the approximation |$f_{\rm gas}(T,\, z)\approx f_{\rm b,\, halo} (T,\, z)$|, where |$f_{\rm b,\, halo} (T,\, z)$| is the baryon mass fraction in haloes of temperature T at redshift z. This is in turn determined by studying the balance between the gravitational force and outward pressure exerted by supernova-driven winds on a gas parcel within the virial radius, as a function of its distance from the centre. One can then define a critical radius within which baryons are bound to the halo and escape otherwise. The question is whether such a critical radius falls within the cooling radius, which is the region of the halo that can produce stars in the simplified picture considered by SP21 and HS03. SP21 show that above a certain critical virial temperature |$T_{\rm crit}(z)$| the critical radius falls outside the cooling radius, hence all haloes retain their cosmic share of baryons. Below |$T_{\rm crit}(z)$|, this approach yields a baryon mass fraction enclosed in the halo that scales as a power of its virial temperature:
The transition between the two-temperature regimes is smoothed with a suitable analytic function. Equation (15) effectively provides a prediction for the baryonic Tully–Fisher relationship (bTFR), an empirical correlation between the baryonic mass and total mass of haloes. Comparing the predicted bTFR with observations of the bTFR by Lelli, McGaugh & Schombert (2016), SP21 conclude that values of |$\eta$| in the range |$1.9-2.4$| are compatible with the data within |$3 \sigma$|. The same values are simultaneously compatible, within a similar level of precision, with the Genzel et al. (2010) observations of the Kennicutt–Schmidt relationship.
The final expression of the nSFR as a function of time is then obtained by simply interpolating the high-redshift (equation 6) and low-redshift solutions (equation 13) with a smooth function, such that
SP21 set |$m=2$|, although they showed that the impact of this parameter on the predicted CSFRD is minimal. Equation (16), combined with the Sheth–Tormen formalism for the halo multiplicity function, allows one to obtain an analytic expression for the CSFRD given by equation (1).
2.2 Improving the accuracy of the SFR at late times
As explained in Section 2.1, one simplifying assumption adopted in SP21 is that the baryon mass fraction in haloes is roughly equal to the gas mass fraction. This approximation can be easily justified at high redshift, shortly after the onset of star formation. However, this assumption can be too strong at later times, when the stellar mass fraction of haloes is not negligible with respect to their gas mass fraction (e.g. McGaugh et al. 2010). This could become an even more important issue when considering the future of the universe. Since the major focus of this work is understanding the impact of |$\Lambda$| on both the past and future star formation history, we need to first extend the SP21 formalism by keeping |$f_{\rm gas}(T,\, z)$| and |$f_{\rm b,\, halo} (T,\, z)$| distinct.
One would be tempted to compute |$f_{\rm gas}$| from equation (15), replacing the cosmic baryon mass fraction with |$f_{\rm b}-f_*(T, \, z)$|, where |$f_*(T, \, z)$| is the stellar mass fraction for a halo with virial temperature T at redshift z. The stellar mass fraction is related to the time integral of the nSFR provided by equation (16). The trouble is that the nSFR depends on |$f_{\rm gas}$|, while |$f_{\rm gas}$| in turn depends on the nSFR.
To break this circularity, we adopt the following iterative method to calculate the nSFR:
we assume that, at any fixed virial temperature, |$f_{\rm gas}(T,\, z)\approx f_{\rm b,\, halo} (T,\, z)$| at a sufficiently high redshift |$z_{\rm in}$|, and compute the corresponding nSFR following the original SP21 model;
- we obtain the stellar mass fraction as a function of redshift by integrating the nSFR over time, or, in terms of redshift:(17)$$\begin{eqnarray} f_*(T,\, z) = \frac{1}{M(T, \, z)} \int ^{z_{\rm in}} _{z} \frac{M(T,\, z^{\prime }) \, s(T, \, z^{\prime })}{(1+z^{\prime }) H(z^{\prime })} \mathrm{ d}z^{\prime } \, ; \end{eqnarray}$$
we compute the cooling radius as a function of redshift from the criterion |$t_{\rm cool}(r_{\rm cool}(t))=t$| explained in Section 2.1; the cooling time is given by equation (8), with the gas density given by equation (7), where we further impose that |$M_{ \rm gas}(T, \, z) = (f_{\rm b}-f_*(T,\, z)) M$|. We then compare the evolution of the cooling radius with that of the critical radius described in the previous section to determine the critical temperature |$T_{\rm crit} (z)$| (see the appendix in SP21 for more details);
we update |$f_{\rm gas}(T,\, z)$| via equation (15), where we now replace |$f_{\rm b}$| with |$f_{\rm b}-f_*(T, \, z)$|;
we use the new gas mass fraction to re-calculate the nSFR via equations (6), (13), and (16);
we repeat the above protocol restarting from point (ii).
To ensure that the iterative procedure described above produces converging results, we examined the nSFR, stellar mass fraction, and critical temperature computed at each iteration. We plot the evolution of these quantities in Fig. 1. The nSFR and stellar mass fraction refer to a halo with virial temperature |$T=10^5 \, \rm K$| in our |$\Lambda$|CDM Universe. For all quantities plotted, convergence is achieved very quickly – within five iterations. As expected, a non-null stellar mass fraction reduces the gas available for star formation, so that the nSFR is reduced. This in turns makes the cooling radius smaller: haloes that were at the critical temperature in the SP21 work are therefore now above the critical temperature so that the critical temperature is lowered with respect to the original formalism.

Redshift evolution of the nSFR (left panel) and stellar mass fraction (middle panel) for a halo with virial temperature |$T=10^5 \, \rm K$|, as predicted at every iteration of our updated version of Sorini & Peacock (2021) described in Section 2.2. The right panel shows the prediction for the redshift evolution of the critical temperature. The black lines labelled ‘iteration 0’ refer to the predictions of the original SP21 model with |$f_*(T,\, z)=0$|; thus no black line appears in the middle panel, which adopts a logarithmic scale on the y-axis. For all quantities plotted in this figure, full convergence is achieved after five iterations.
As we can see in the middle panel of Fig. 1, the stellar mass fraction asymptotes to a constant in the far future. That is a consequence of the decay in the nSFR, which tends to zero in the limit of infinitely large cosmic time. We will discuss this in detail in Section 4.1. The stellar mass fraction is a monotonic function of time, because the nSFR is always positive. There is no explicit mechanism in the model that removes stars that have gone supernova or became stellar remnants. The contribution of short-lived and massive stars to the nSFR is removed through the |$\beta$| factor in equations (6) and (13), but all other stars are ‘eternal’ once formed. This is a reasonably good approximation, since the time-scale for white dwarfs to turn into black dwarfs has been estimated to be |$\sim$| |$10^5\, \rm Gyr$| (Dyson 1979; Adams & Laughlin 1997; Caplan 2020); we will discuss the subsequent implications for anthropic reasoning in Section 5.1.
We verified that convergence is achieved in the range of interest for this study, that is |$10^4 \, \mathrm{K} \lt T \lt 10^{10} \, \rm K$| (as explained in SP21, the contribution of haloes with |$T\gt 10^{10} \, \rm K$| to the CSFRD is negligible; we verified that this holds also for the cosmological models considered in this work). At higher temperatures, convergence is faster: at |$T=10^7 \, \rm K$|, two iterations are sufficient. For all quantities shown in Fig. 1, we chose |$z_{\rm in}=20$| as the initial redshift, but we verified that the exact value is unimportant as far as our analysis is concerned (see Section 4.1 for a more detailed discussion).
To assess whether the stellar mass fractions predicted by our revised SP21 method are physically sensible, let us consider the mass of the Milky Way. Recent estimates point to a stellar mass of |$(6.08 \pm 1.14) \times 10^{10} \, \rm {\rm M}_{\odot }$| (Licquia & Newman 2015) and a virial mass of |$1.54_{-0.44}^{+0.75} \times 10^{12} \, \rm {\rm M}_{\odot }$| (Watkins et al. 2019). Therefore, the stellar-to-total mass ratio is expected to be in the range |$2.2-6.6~{{\ \rm per\ cent}}$|. Applying equation (5), we can derive the virial temperatures corresponding to the virial mass found by Watkins et al. (2019), and then compute the corresponding stellar mass fraction at |$z=0$| via equation (17), obtaining 6.7 per cent. Thus, the predictions of our model are compatible with the upper bound of the stellar-to-total halo mass range provided by the observations. Considering the simplifications made in the formalism, this is a reassuring result.
The key quantity under study in this paper is the CSFRD, so it is important to see how this is affected by the changes in normalized SFR shown in Fig. 1, which result from our generalization of the SP21 model. In the upper-left panel of Fig. 2 we show the CSFRD in the redshift range |$0\lt z\lt 10$|, computed both with the original and revised formalism. It is clear that our procedure yields a slightly lower CSFRD at low redshift, with little change in the slope at late times. In the upper-right panel, we show the behaviour of the two models in the future of the Universe. As future times correspond to negative redshifts, and the CSFRD approaches zero as time tends to infinity, the scale of both axes switches from logarithmic to linear when moving from the upper-left to the upper-right panel. We quantify the relative difference between the original SP21 model and our updated formalism in the lower panels of Fig. 2 (dashed black line). The differences are at the per cent level for |$z\gt 5$|, and then grow progressively, reaching |$\sim$| 10 per cent at the peak of star formation and |$\sim$| 20 per cent at |$z=0$|. The discrepancy grows slowly but steadily in the far future of the Universe, stopping short of 30 per cent in the limit |$t\rightarrow \infty$|. Such differences are subdominant with respect to the uncertainties of the parameters of the model. The original SP21 model is thus still suitable for rapidly predicting the past star formation history of our Universe, but our generalization becomes indispensable for the future star formation history. This is even more important when one considers larger values of the cosmological constant: we verified that the relative difference with respect to the SP21 model grows as |$\Lambda$| increases.

Upper panels: CSFRD, as computed with our extension of the SP21 model (dashed black line) and as given by the original SP21 model (dot–dashed red line). In the left panel, we further show with the empirical fit (dotted black line) to a compilation of observational data provided by Madau & Dickinson (2014) (grey data points). The right panel shows the behaviour of the SP21 model and our extended formalism in the future of the Universe. As future times correspond to negative redshifts, and the CSFRD becomes asymptotically null at arbitrarily large times, both axes switch from logarithmic to linear scale when moving from the left to the right panel. The grey hatched shaded area in the right panel excludes the non-physical region corresponding to |$\dot{\rho }_*\lt 0$|. Lower panels: Relative difference between the two models showed above (dashed black line), and between our formalism and the fit from Madau & Dickinson (2014) (dotted black line). The latter comparison is available only for positive redshifts (i.e. past cosmic times), as the fit represents a purely empirical fit to the data, and not a predictive theoretical model. The SP21 model agrees well with the improved method introduced in this work at high redshift, but overestimates the CSFRD at lower redshift and in the future. In the redshift range |$0\lt z\lt 10$|, both the original SP21 model and our updated formalism agree with the empirical fit within a factor of |$\sim \, 2.$|
In the upper-left panel of Fig. 2, we also show the compilation of observational data provided by Madau & Dickinson 2014, alongside their empirical fit to the data. Both data and fitting function have been properly rescaled from a Salpeter (1955) to a Chabrier (2003) IMF, for consistency with the SP21 formalism and the model presented in this work (the data points were not rescaled in this way in the SP21 paper). The relative difference between the SP21 model and the Madau & Dickinson (2014) fit are quantified in the lower-left panel, confirming that the SP21 predictions recover the fitting function within a factor of |$\sim 2$| for |$0\lt z\lt 10$|. Since our updated formalism deviates from SP21 by at most 20 per cent in the same redshift range, we can conclude that our model also matches the Madau & Dickinson (2014) fit within a factor of |$\sim$| 2. This level of agreement is remarkable, considering the simplicity of the SP21 model. As a reference, sophisticated hydrodynamic cosmological simulations, while providing an overall better agreement with the data, exhibit discrepancies of a factor of |$\sim$| 2 or larger around cosmic noon (McCarthy et al. 2017; Davé et al. 2019) and within a factor of |$\sim$| 1.6 at |$z\lesssim 3$| (Salcido et al. 2018). However, other state-of-the-art cosmological simulations, such as IllustrisTNG (Pillepich et al. 2018), achieve a significantly more accurate match with observations of the CSFRD (Weinberger et al. 2017). Thus, whereas the SP21 model manages to broadly reproduce both observed and simulated trends of the CSFRD (Scharré, Sorini & Davé 2024), there is certainly room for further improvement. A detailed discussion on the limitations of the SP21 model and the resulting impact on the predicted CSFRD can be found in Sorini & Peacock (2021).
It is also worth noting that the data in the Madau & Dickinson (2014) compilation might be affected by biases introduced by the parametric models underlying the estimation of the CSFRD from observations of galaxy spectral energy distributions (Carnall et al. 2019). Additionally, subsequent estimates of the CSFRD deviate from the Madau & Dickinson (2014) fit, both at high (Gruppioni et al. 2015; Rowan-Robinson et al. 2016) and low (Gruppioni et al. 2015) redshifts. Other measurements exhibited a lower normalization for the CSFRD at the peak of star formation (Liu et al. 2018), or a plateau at cosmic noon (Traina et al. 2024). The spread of some of the data points from these more recent measurements is generally comparable to, or even larger than, the discrepancy between the SP21 model and the Madau & Dickinson (2014) fit. Given this range of predictions, the precision of the CSFRD given by the SP21 model seems satisfactory for the purpose of the present investigation.
As for the performance of the model in predicting the future star formation history, we note that the shape of the future CSFRD is essentially dictated by the evolution of the cooling radius and of the critical temperature (see Appendix A for a detailed derivation of the asymptotic scaling of the CSFRD). The normalization of the CSFRD for |$z\lt 0$| will follow the overall normalization of the CSFRD at earlier times. Thus, given a factor of |$\sim$| 2 uncertainty in the past star formation history at |$z\approx 0$|, we can reasonably expect a similar systematic imprecision in the future CSFRD. As it will become apparent in Section 4.1, the future behaviour of the CSFRD may contribute significantly to the overall production of stars in the entire history of the universe. However, as we will show, the differences in the total stellar mass produced in the universe for different values of the cosmological constant can span several orders of magnitude, thus our main conclusions will be largely unaffected by an uncertainty of a factor 2 in the CSFRD.
We also highlight that our analysis is relative, i.e. based on the ratio of the star formation histories produced in different cosmologies, and will not depend on the absolute value of the star formation efficiency in a single cosmology. Of course, our posterior for |$\Lambda$| will still be subject to accepting the underlying simplified astrophysical model for cosmic star formation. We discuss the impact of these limitations of our study, and prospects for their future amelioration, in Section 5.4.
3 IMPACT OF CHANGING |$\Lambda$|
Having introduced the basic formalism in the previous section, we will now explain how varying |$\Lambda$| would affect the star formation history. Changing the cosmological constant obviously affects the evolution of large-scale structure; this is the subject of Section 3.1. However, a different value of |$\Lambda$| can also directly affect the cooling rate within haloes. The astrophysical impact of changing |$\Lambda$| will be discussed in Section 3.2.
3.1 Cosmological impact
To begin with, we need to understand how the empirical |$z=0$| cosmological parameters such as |$H_0$| are altered when we vary |$\Lambda$|. We will restrict our discussion to arbitrary flat|$\Lambda$|CDM universes only – implicitly assuming an inflationary multiverse. In addition, all the dimensionless parameters of the |$\Lambda$|CDM model are assumed to be unchanged: the horizon-scale amplitude |$A_s$|, the ratio of baryonic and CDM densities, and the ratio of photon and baryon number densities. The consequence is that all model universes under study are indistinguishable copies of our own Universe at very high redshifts where |$\Lambda$| is dynamically unimportant. However, in endowing them with different values of |$\Lambda$|, the histories of these copies become diverse at later times. The reason for this restriction is the measure problem: for |$\Lambda$| we can appeal to Weinberg’s argument for a flat prior, but if other parameters were to vary then we have little idea what the appropriate prior would be. Therefore, we concentrate on this simplest ensemble.
It is convenient to quantify these alternative models using the standard set of parameters for the |$\Lambda$|CDM model: |$\Omega _{\rm m}$|, |$\Omega _{\Lambda }$|, etc. We will add the subscript ‘ref’ when we indicate the corresponding parameters in our own Universe. This requires a definition of the ‘present’, which we take to be the point at which the CMB temperature equals the standard value (|$T_0=2.7255 \, \rm K$|; see Planck Collaboration VI 2020). We normalize the scale factor such that |$a=1$| at this point, so that the radiation density in all models scales |$\propto a^{-4}$|, with the identical constant of proportionality in all cases. However, changing the value of |$\Lambda$| means that the time corresponding to |$T=T_0$| and the Hubble parameter at that point take values that are different from those in the observed Universe, as we now explain.
We start by noting that the Friedmann equation for a flat universe applies in all cases, so at times relevant for star formation we have
We preserve the matter density and re-scale the vacuum density by a factor |$\alpha _{\Lambda } = \Lambda /\Lambda _{\rm ref}$|, where |$\Lambda _{\rm ref}$| is the cosmological constant in the reference cosmology, so that |$\rho _{\rm m} = \rho _{\rm m, \, ref}$| and |$\rho _{\Lambda } = \alpha _{\Lambda } \rho _{\Lambda , \, \rm ref}$| in all universes. The Hubble parameter in the generic universe of the ensemble can be then be re-written as
Writing the ratio of the two terms on the right hand side of equations (18) and (19), and remembering that flatness always ensures |$\Omega _\Lambda =1-\Omega _m$|, we deduce
However, the unchanged matter density tells us that |$H_0^2=\Omega _{\rm m,\ ref} H_{0,\rm ref}^2/ \Omega _{\rm m}$|, so that
Hence, we readily obtain the new Hubble parameter and density parameter at |$a=1$|, for a given scaling of |$\Lambda$|. The cosmic baryon mass fraction is taken as unchanged. Similar reasoning was given by Oh et al. (2022).
This argument applies for any sign of |$\alpha _\Lambda$|, but rescaling to negative values complicates things. The expansion ceases at some maximum value of a and the universe subsequently undergoes recollapse. In principle, the SP21 code should be able to handle such changes, although in practice it will require amendment in order to cope with a non-monotonic |$a(t)$| relation. However, in any case, the physical effect on the CSFRD of a negative |$\Lambda$| is qualitatively different to that of a positive value. In the latter case, a large value of |$\Lambda$| is expected to suppress star formation because structure growth ceases before galaxy-scale haloes can be generated. The main aim of this paper is to make a quantitative estimate of this suppression. In recollapsing models, however, growth does not freeze out and the relevant question is how much star formation can occur in the limited time before the big crunch. We aim to consider this question elsewhere.
In Fig. 3, we show the resulting redshift evolution of the Hubble parameter for different values of |$\alpha _{\Lambda }$|. As expected, for larger values of |$\Lambda$| the Hubble parameter reaches its asymptotic value |$H\rightarrow H_0 \sqrt{\Omega _{\Lambda }}$| at earlier redshifts, and this value increases with |$\Lambda$|. This reflects the fact that the early phase of all models is Einstein–de Sitter (EdS), with |$H\propto a^{-3/2}$|, while H freezes out at the value of the EdS model at the point where |$\Lambda$| comes to dominate.

Evolution of the Hubble parameter |$H(z)$| for different re-scalings of the cosmological constant. For larger values of |$\Lambda$|, the asymptotic value of the Hubble constant in the far future of the universe (negative redshift) increases. In an EdS universe, the Hubble constant becomes indefinitely small as the cosmic time tends to infinity.
Calculating the new value of |$\sigma _8$| and its scaling with |$\alpha _{\Lambda }$| is more complicated: |$\sigma _8$| is affected both by the new |$\Omega _{\rm m}$|, which changes the power spectrum shape and evolution, and by the new Hubble constant, which changes the scale |$8 \, h^{-1} \, \rm cMpc$|. These effects can be allowed for by recalling that conditions at very high z are the same in all models, so that |$\sigma (8\, h_{\rm ref}^{-1}{\rm cMpc})$| is known at |$z\gg 1$|. We then compute the scale-dependent correction numerically, using the code camb (Lewis, Challinor & Lasenby 2000; Lewis & Challinor 2011). This dependence of |$\sigma _8$| on |$\Lambda$| was studied by Oh et al. (2022), who showed that as the universe approaches an EdS solution, |$\sigma _8$| exhibits an asymptotic behaviour (see their fig. 1). For a Planck-2018 cosmology, the asymptotic value of |$\sigma _8$| for low values of |$\Lambda$| is 0.6786 (Oh et al. 2022). If one increases the value of |$\Lambda$| with respect to the fiducial one, |$\sigma _8$| reaches a maximum of |$\sigma _8 \sim 0.9$| at around |$\alpha _{\Lambda }\sim 8$|. For larger values of |$\Lambda$|, it keeps declining. As a reference, for |$\alpha _{\Lambda }=1000$| one obtains |$\sigma _8 \approx 0.5$|.
We are now fully equipped to study the impact of |$\Lambda$| on the CSFRD. The most obvious effect of altering |$\Lambda$| is a modification of the cosmological term in the integrand of equation (1), i.e. the halo mass function. A larger cosmological constant causes the |$\Lambda$|-dominated accelerating phase to begin earlier. Consequently, the freeze-out of structure formation would also occur earlier, thus shifting the truncation of the halo mass function to smaller halo masses, and decreasing the merger rate. This is exactly what we observe in Fig. 4, which plots the halo multiplicity function for different values of |$\alpha _{\Lambda }$| and redshift. For a cosmological constant 100 times larger than the observed value (dashed lines), the cut-off in the HMF drops by two orders of magnitude at |$z=0$|. For reduced values of |$\Lambda$|, the trend is opposite. At any given redshift, models with |$\alpha _{\Lambda }\lt 1$| overproduce massive haloes relative to the reference cosmology. In the limit of an EdS universe (not plotted in the figure), there is no freeze-out and the halo mass cut-off becomes arbitrarily large at sufficiently great times.

Sheth–Tormen halo mass function at different redshift, represented by different colours, and for different values of |$\Lambda$| (different line styles). For a given redshift, larger values of |$\Lambda$| result in a mass cut at smaller halo masses. In an EdS universe, the mass cut grows indefinitely in the future (|$z\lt 0$|). This is not the case for |$\Lambda \gt 0$|, where the mass cut approaches an asymptotic value as |$z\rightarrow -1$| (i.e. |$t\rightarrow \infty$|). This is a consequence of the freeze out of structure formation in the |$\Lambda$|CDM model.
To summarize, lower positive values of |$\Lambda$| promote structure formation. However, it would be premature to assume that the probability of generating observers is solely dependent on the collapsed mass fraction in the universe, and thus to conclude that the EdS model is anthropically favoured. Rather, we must consider the effect of |$\Lambda$| on the astrophysical processes that govern cosmological star formation. This will be the focus of the next section.
3.2 Astrophysical impact
Altering |$\Lambda$| affects the astrophysics of star formation because both the gas mass fraction in haloes and the gas cooling rate depend on the virial quantities of haloes, which in turn are sensitive to the background cosmology via the Hubble parameter (see Section 2). In this section, we will illustrate these points in detail.
In this work, we will fix the astrophysical parameters to the ‘best-fitting’ choice in SP21, i.e.: |$\eta =1.9$|, |$\langle t_* \rangle = 3.87\, \rm Gyr$| and |$T_{\rm min} = 10^{4.5} \, \rm K$|.
3.2.1 Cooling rate
As we recalled in Section 2, within the SP21 model the low-redshift SFR is set by the gas cooling rate. This is in turn determined by the extent of the cooling radius |$r_{\rm cool}$|. We will now analyse how the virial temperature of the halo affects the cooling radius, and how this evolves depending on the cosmological model. A particular emphasis will be given to the asymptotic behaviour of the cooling radius for arbitrarily large cosmic times.
We can express the cooling radius in terms of the cooling time by obtaining the gas density at |$r_{\rm cool}$| via equation (7), and then inserting it into equation (8):
To obtain an explicit dependence of |$r_{\rm cool}$| on redshift, we can simply replace |$t_{\rm cool}$| with the prescription for the ‘effective cooling time’ given by equation (10), making use of the cosmic time-redshift relationship defined by equation (12). We can further apply the definition |$M_{\rm gas} = f_{\rm gas}(T,\, z) M$|, and then express all virial quantities in terms of T through equations (2)–(5).
We can therefore follow the redshift evolution of |$r_{\rm cool}$|; this is shown in Fig. 5 for different virial temperatures and different values of |$\Lambda$|.

Redshift evolution of the cooling radius |$r_{\rm cool}$|, in units of the virial radius, at fixed virial temperature. Different colours refer to different values of |$\Lambda$|, while different line styles represent different virial temperatures. In an EdS universe, the cooling radius shrinks indefinitely in the far future. For |$\Lambda \gt 0$|, it reaches a minimum value, and then slowly resumes increasing. This is a consequence of the freeze-out of structure formation: see Section 3.2.1 for further details.
For a fixed value of |$\Lambda$|, |$r_{\rm cool}/R$| decreases with increasing virial temperature. In haloes above the critical temperature |$T_{\rm crit}$| (see Section 2), this behaviour can be readily understood from equation (23). If |$T\gt T_{\rm crit}$|, then the gas mass fraction in haloes reaches its maximum value, and no longer depends on the virial temperature. The cooling radius is then simply proportional to |$(\mathcal {C}(T)/T)^{1/\eta }$|. For |$T\gtrsim 10^6 \, \rm K$|, a Sutherland & Dopita (1993) cooling function for an H/He plasma with primordial abundances scales approximately as |$\propto T^{1/2}$|. Therefore, |$r_{\rm cool}/R \propto T^{-1/2 \eta }$|, so the cooling radius decreases with temperature in the regime considered. For |$T\lt T_{\rm crit}$|, there is an extra temperature-dependence carried by the |$f_{\rm gas}(T, \, z)$| factor, but we verified that this does not qualitatively impact the scaling of |$r_{\rm cool}$| with temperature.
We can now examine the redshift evolution of |$r_{\rm cool}/R$| for a fixed virial temperature, for different cosmological models. To begin with, we restrict the discussion to universes with |$\Lambda \gt 0$|. At high redshift, the effective cooling time is approximately equal to the dynamical time (see equation 10). Thus, the factor containing the cosmic time in equation (23) is very close to unity. Also, the critical temperature drops at sufficiently high redshift (see SP21 and Fig. 6). Therefore, at early times, we have that |$f_{\rm gas} \approx f_{\rm b, \, halo} \approx f_{\rm b}$|, so that the evolution of the cooling radius is dictated primarily by the Hubble parameter in equation (23). In this regime, the universe is dominated by matter, thus the cooling radius scales as |$H(z)^{1/\eta } \propto (1+z)^{3/2\eta }$|. Recalling that we set |$\eta =1.9$|, we would then expect that the redshift-evolution of |$r_{\rm cool}/R$| to be a power law with index close to |$3/4$| at high redshift, where all models are matter dominated. This is exactly what we observe in Fig. 5.

Redshift evolution of the critical virial temperature |$T_{\rm crit}$| above which the baryon mass fraction in haloes saturates to the cosmic value |$f_{\rm b}=\Omega _{\rm b}/\Omega _{\rm m}$|, according to the SP21 model, for different values of the cosmological constant (see colour bar). The scaling of |$T_{\rm crit}$| with redshift mirrors the relationship between cosmic time and redshift. In an EdS universe, the long-term behaviour resembles a power law. As |$\Lambda$| increases, |$T_{\rm crit}(z)$| deviates from a power law at progressively earlier redshifts (see Section 3.2.2 for further details).
Conversely, as the universe becomes dominated by |$\Lambda$|, the effective cooling time given by equation (10) approaches the cosmic time, so that the cooling radius scales as |$\propto (f_{\rm gas} \,\, t)^{1/\eta }$|. To determine the evolution of the cooling radius, we then need to understand how |$f_{\rm gas}$| scales with redshift or, equivalently, cosmic time. Because for |$\Lambda \gt 0$| the critical temperature is monotonically increasing with cosmic time after the point of |$\Lambda$|-domination (see Fig. 6), there is always a sufficiently late cosmic time after which the virial temperature of a halo drops below |$T_{\rm crit}(z)$|. Thereafter, the gas mass fraction within haloes of a given virial temperature scales with |$T/T_{\rm crit}(z)$| as in equation (15).
The redshift evolution of |$T_{\rm crit}(z)$| needs to be computed numerically, but we can find an analytical approximation by simplifying the Sutherland & Dopita (1993) cooling function with a piece-wise power law. Under this assumption, the critical temperature itself scales as a power of cosmic time, hence so does the cooling radius (see Appendix A2 for details). We verified that the index of the power law is positive and below unity, explaining why the cooling radius increases after the onset of |$\Lambda$|-domination, while its time derivative flattens out to zero. This is indeed the trend that we observe in Fig. 5 at |$z\lt 0$| for all models with |$\Lambda \gt 0$|. For a fixed redshift and virial temperature and after the point of |$\Lambda$|-domination, |$r_{\rm cool}/R$| becomes larger as |$\Lambda$| increases. To prevent the cooling radius from overflowing beyond the virial radius, we cap the ratio |$r_{\rm cool}/R$| to unity with a suitably smooth function, similar to the one used for the nSFR in equation (16).
The behaviour that we just described descends from the definition of the boundaries of haloes: for a given virial temperature, it follows from equations (2)–(5) that both the virial radius and virial mass of a halo scale as |$1/H(z)$|. Thus, in a universe with a larger value of |$\Lambda$|, haloes of a given temperature become more concentrated (see Fig. 3). Therefore, the gas density becomes overall larger, and consequently cooling becomes more efficient. However, in the far future of the universe structure formation freezes out, and haloes become isolated. Eventually, the SP21 framework assumes that haloes will expel their baryon content through feedback processes, progressively becoming devoid of gas and unable to form stars (see Fig. 1).
To summarize, we saw that in a universe with |$\Lambda \gt 0$| the cooling radius decreases as a power of time at high redshift, but that it eventually switches to evolve in the opposite direction. The point of turnaround between these two regimes occurs at the transition to |$\Lambda$|-domination, when the effective cooling time in equation (10) switches from being approximated by the dynamical time to tracking the cosmic time. As can be seen in Fig. 5, increasing |$\Lambda$| moves the turnaround of the cooling radius to earlier redshifts. In contrast, there is no such turnaround in an EdS universe, and the cooling radius (in units of the virial radius) shrinks indefinitely as time goes by. Because there is no lower limit to the Hubble parameter, the size of haloes with a given virial temperature can become arbitrarily large. This makes the diffuseness of haloes and consequent inefficiency of gas cooling more dramatic than in any universe with positive |$\Lambda$|.
3.2.2 Gas mass fraction in haloes
As explained in Section 2, the SP21 model assumes that the baryon mass fraction in haloes scales with the virial temperature below a certain critical threshold |$T_{\rm crit}$|, and saturates to |$f_{\rm b}$| above this value. In our improved variation laid out in Section 2.2, the gas mass in haloes above the critical temperature saturates to |$f_{\rm b}-f_*(T, \, z)$|, where |$f_*(T, \, z)$| is the stellar mass fraction of a halo of virial temperature T at redshift z. In either case, the physical significance of the critical temperature is the same.
In haloes with |$T\lt T_{\rm crit}$|, the gaseous component is gravitationally bound only within a certain critical radius |$r_{\rm crit}$|. Beyond this scale, the momentum injected by supernova-driven winds into the surrounding gas is sufficient to overcome gravity and eject the gas. In practice then, |$T_{\rm crit}$| is the temperature at which |$r_{\rm crit}$| coincides with the boundary of the star-forming region within the halo. At high redshift, this is the entirety of the halo, so that |$T_{\rm crit}$| is determined by the condition |$r_{\rm crit}\lt R$|. At low redshift, star formation is cooling-driven, and the relevant condition is |$r_{\rm crit}\lt r_{\rm cool}$|.
Since the cosmological model affects both the virial radius and the cooling radius (see Section 3.2.1), the critical temperature is also sensitive to the cosmological parameters. In particular, in Fig. 6 we show the redshift evolution of |$T_{\rm crit}$| for different values of |$\alpha _{\Lambda }$|. The late-time behaviour of the critical temperature can be understood starting from the condition |$r_{\rm crit}=r_{\rm cool}$|, which defines |$T_{\rm crit}$|. As discussed in the previous section, the cooling radius explicitly depends on the cooling time; in the far future of the universe, |$r_{\rm cool}$| corresponds to the point at which the cooling time is approximately equal to the cosmic time.
In a universe with |$\Lambda \gt 0$|, under the reasonable approximation of the cooling function as a piece-wise power law, the critical temperature scales as a power of cosmic time in the limit |$t \rightarrow \infty$| (see Appendix A2 for details). We verified that the index of this power law is positive in all physically relevant cases. Thus, for |$\Lambda \gt 0$|, the critical temperature increases at arbitrarily large times. The mapping between the scale factor and cosmic time becomes exponential for |$t\rightarrow \infty$| in models with |$\Lambda \gt 0$|. Therefore, in the far future, the critical temperature scales with redshift as a power of |$-\ln (1+z)$|, and that is why the increase of |$T_{\rm crit}$| with redshift is rather slow (see Fig. 6).
In an EdS universe, we observe a qualitatively different behaviour for the critical temperature, which decreases rather than increasing at late times. This can be understood by considering that for a fixed virial temperature, the virial radius scales |$\propto H(z)^{-1} = (1+z)^{-3/2}$|, hence diverging for |$z\rightarrow -1$|. Therefore, all haloes of a given T become arbitrarily large in the far future, so that their baryonic mass fraction must eventually match the cosmic value. Within the SP21 formalism, this means that all haloes must asymptotically exceed the critical temperature. This condition is achieved if |$T_{\rm crit} \rightarrow 0$| for |$t \rightarrow \infty$|, in accord with our numerical results.
At the high-redshift end of Fig. 6, all models start out with the same |$T_{\rm crit}(z)$| relation, because all universes are matter dominated at early times. As we saw in Section 3.2.1, |$r_{\rm cool}$| in the far future increases for larger values of |$\Lambda$| (see Fig. 5). The critical temperature increases accordingly for |$\alpha _{\Lambda } \le 100$|. For |$\alpha _{\Lambda }\gt 100$|, the critical temperature appears to invert this trend, and although it still increases at later times, its value is overall lower than in universes with a smaller cosmological constant. The reason is that in universes with |$\alpha _{\Lambda } \gg 1$|, haloes reach overall higher densities: for a fixed halo mass M (or, equivalently, for a fixed virial temperature at future times), the virial radius scales as |$R \propto M^{1/3} H(z)^{-2/3} \propto M^{1/3} \alpha _{\Lambda }^{-1/3}$|. A higher gas density makes the cooling process more efficient (see equation 8), so that the cooling radius extends to the virial radius. This means that in practice the condition for determining the critical temperature, |$r_{\rm crit} \lt r_{\rm cool}$|, translates to |$r_{\rm crit}\lt R$|. It easily follows that |$T_{\rm crit}(z) \propto H(z)^{-2}$|, regardless of the value of |$\eta$|. At late times, |$H(z)^2 \propto \alpha _{\Lambda }$|, therefore the critical temperature diminishes as |$\Lambda$| increases (see Fig. 3).
The behaviour of |$T_{\rm crit}$| shown in Fig. 6 directly affects the nSFR both at high and low redshift, as it affects the evolution of the baryon mass fraction in haloes. Piecing together the discussion in Section 3.2.1 and in this section, we gain two important physical insights into the astrophysical impact of |$\Lambda$| on the SFR. For small values of |$\Lambda$|, there is a larger gas mass fraction available in haloes at earlier times, but a progressively smaller fraction will cool down at later times. For larger values of |$\Lambda$|, there is an overall smaller mass fraction of gas that resides within haloes at later times, but a larger fraction will cool down as time goes by. These features will play a prominent role in shaping the dependence of the cosmic star formation efficiency as a function of |$\Lambda$|, as we will discuss in Sections 4.1 and 4.2.
4 RESULTS
4.1 Long-term efficiency of star formation
In this section, we will show our results for the star formation history for different values of |$\Lambda$|, focusing especially on the long-term behaviour of star formation. Determining the asymptotic behaviour of the nSFR will tell us whether there is a well-defined final cumulative efficiency of cosmic star formation – both in individual haloes and in the universe as a whole.
Let us then start with the evolution of the nSFR. In the original SP21 formalism, the nSFR is most naturally expressed for a fixed virial temperature, as that is the independent variable in the cooling function. Nevertheless, as we explained earlier in Section 2, the Sheth–Tormen formalism for the clustering of haloes allows us to convert straightforwardly between the fixed-T and the fixed-M views, making use of the definitions expressed in equations (2)–(5). For the remainder of this paper, the fixed-M view will be more convenient.
In the upper panels in Fig. 7, we show the time evolution of the nSFR for populations of haloes at given virial mass. From left to right, we show the results for a typical dwarf galaxy (|$10^{10} \, \rm {\rm M}_{\odot }$|), a Milky Way-like galaxy (|$10^{12} \, \rm {\rm M}_{\odot }$|), and a galaxy cluster (|$10^{14} \, \rm {\rm M}_{\odot }$|), respectively. The main independent variable in the plot is cosmic time; the upper x-axis shows the redshift that would be measured by an observer in the reference Universe for a given cosmic time. Therefore, it applies only to the case |$\alpha _{\Lambda }=1$|.

Upper panels: Time evolution of the SFR, normalized by the total halo mass. Each panel refers to a population of haloes with fixed virial mass M, as annotated in the figure. Different colours refer to different values of |$\Lambda$| in units of the observed value |$\Lambda _{\rm obs}$|, as reported in the colour bar. The lower x-axis reports the cosmic time, while the upper x-axis the corresponding redshift that would be observed in our Universe (|$\Lambda = \Lambda _{\rm obs}$|). The break in the upper x-axis indicates that the ticks are omitted for |$z\gt -1+10^{-10}$|, since the exponential growth of the scale factor in the future makes it increasingly hard to represent the redshift on a meaningful scale. Regardless of the halo mass, the normalized SFR tends to zero in the limit of |$t\rightarrow \infty$|. Lower panels: Fraction of the stellar mass produced up to a certain cosmic time, with respect to the total halo mass, for haloes with fixed virial mass. This is the time integral of the quantity shown in the corresponding upper panels. The colour coding and the upper and lower x-axes are the same as in the upper panels. In all haloes and cosmological models considered, the cumulative stellar mass fraction reaches an asymptotic value for |$t \rightarrow \infty$|.
We note that for all halo masses and all cosmologies, the nSFR tends to a universal value at sufficiently early times. This happens because, for a given virial mass, the virial temperature increases at earlier times, while the critical temperature is lower (see Fig. 6). There will thus always be a sufficiently early time when the baryon mass fraction within haloes saturates to |$f_{\rm b}$|. Equation (6) then tells us that the nSFR reaches a fixed constant value in this regime.
Since all universes in our ensemble are identical and matter-dominated at a sufficiently high redshift, they are naturally indistinguishable at early times. But for all halo masses considered, the nSFR drops steeply in the future, and the departure from the EdS solution occurs earlier for universes with a larger value of |$\Lambda$|. The asymptotic behaviour can be deduced by studying each factor in the r.h.s. of equation (13). As discussed in Section 3.2, the main challenge is understanding the late-time behaviour of the critical temperature, which needs to be determined numerically. However, under the reasonable assumption of a piece-wise cooling function, this asymptotic behaviour can be approximated analytically with a power of time. The index of the power law depends on the virial temperature of the halo and on the slope of the gas density profile. It is then possible to prove that the time integral of the nSFR at fixed mass (i.e. the stellar mass fraction) is convergent (although the argument is relatively detailed, and is presented in Appendix A). This convergence holds both in an EdS universe and for |$\Lambda \gt 0$|, and is seen in Fig. 7, where the cumulative stellar mass fraction reaches a plateau at late times. We plot this asymptotic stellar mass fraction, |$f_{*\, \infty }$|, as a function of halo mass in Fig. 8.

Asymptotic value (i.e. in the limit |$t\rightarrow \infty$|) of the stellar mass fraction for a halo population with given total mass M. Different colours refer to different values of |$\Lambda$|, as indicated in the colour bar. Even for large values of |$\Lambda$|, massive haloes can still be very efficient (see Section 4.1 for the detailed explanation).
For all values of |$\Lambda$| considered, we notice that |$f_{*\, \infty }$| increases with mass for very low masses. For |$\alpha _{\Lambda } \lt 1$|, the efficiency then peaks at |$M\approx 10^{11} \, \rm {\rm M}_{\odot }$|, and then decreases for larger halo masses, with the drop being faster in an EdS universe. However, even for still relatively small values of |$\Lambda$|, such as |$\alpha _{\Lambda }=0.1$|, the asymptotic efficiency decreases more slowly after reaching its maximum, and for |$M\approx 10^{18} \, \rm {\rm M}_{\odot }$| it is almost one order of magnitude larger than in an EdS universe. For |$\Lambda$| at the observed level or above, the initial decline in efficiency above |$M=10^{11} \, \rm {\rm M}_{\odot }$| is rather slow, with a plateau of nearly constant efficiency, followed by an eventual sharper drop. In the case of |$\alpha _{\Lambda }=10^4$|, |$f_{* \, \infty }$| exhibits little variation over the nine decades of halo mass shown in Fig. 8. In such an extreme cosmology, the freeze out of structure formation occurs very early and the critical temperature is so low (see Fig. 6) that most haloes retain their cosmic share of baryons; thus they are always efficient at forming stars, regardless of their mass.
The absence of a sharp peak in the stellar mass fraction and subsequent steep decline at the higher mass end might seem puzzling, at least for the case |$\alpha =1$|, as it appears to contradict the current paradigm of galaxy formation (see e.g. Behroozi et al. 2019). However, we stress that Fig. 8 shows the asymptotic stellar mass fraction (i.e. in the limit |$t\rightarrow \infty$|), whereas the consensus on the mass dependence of the stellar mass fraction is based on observations and models of the past star formation history. We verified that our formalism would also predict a sharper peak and a steeper decline of the stellar mass fraction for |$M\gtrsim 10^{12} \, {\rm M}_{\odot }$| at |$z=0$|, in qualitative agreement with the accepted framework for galaxy formation. Therefore, Fig. 8 tells us that, according to our model, high-mass haloes produce a larger fraction of their asymptotic stellar mass in the future of the universe compared to low-mass haloes. Higher mass haloes appear to be especially efficient at forming stars in cosmologies with larger values of |$\Lambda$|, where haloes of a given virial mass are overall denser, hence more efficient at cooling (see Section 3.2.2). While this argument provides us with a basic understanding of the trend observed in Fig. 5, one should bear in mind that gas cooling is not the only relevant process for star formation, even at late times. While we incorporate the SP21 model of supernova-driven winds (see their Appendix A) in our calculations, we do not include an explicit self-regulation mechanism due to AGNs. This may have a noteworthy impact on the future SFR. We will discuss the subsequent implications for our conclusions in Section 5.3.
It is important to bear in mind that Fig. 8 ignores the cosmology dependence of halo multiplicities. We then need to weight the efficiencies shown in Fig. 8 by the HMF. If we then further integrate over halo mass, that gives us the stellar mass density (SMD) of all stars formed up to time t from the time corresponding to the onset of star formation |$t_{\rm in}$|. This quantity is by definition the time integral of the CSFRD or, in terms of redshift:
where |$\dot{\rho }_*(z)$| is given by equation (1). To study the asymptotic behaviour of |$\rho _*(z)$|, we should thus study the convergence of the integrand, which depends both on the nSFR and the HMF. We have already discussed the convergence of the nSFR for |$\Lambda \gt 0$| above, and we now focus on the behaviour of the HMF for |$t\rightarrow \infty$|. In a universe with |$\Lambda \gt 0$|, the asymptotic behaviour of |$\sigma (M, \, z)$| is a non-null constant, for any M. Therefore, the integral in equation (24) converges as long as the stellar mass fraction converges – and we have already seen that this quantity does indeed asymptote to a finite constant in the far future. In an EdS model, the HMF continues to evolve indefinitely, but in fact the total stellar density still converges (see Appendix A1 for the detailed argument).
The asymptotic study of the CSFRD described above also allows us to define physically motivated analytic approximations for how the CSFRD approaches its asymptotic value in the far future (see Appendix A for details). The results of this exercise are shown in the upper panel of Fig. 9, plotting the stellar density as a function of time for various values of |$\alpha _{\Lambda }$|. In an EdS universe, the asymptotic SMD is reached at |$t\sim 100\, \rm Gyr$|; but as |$\Lambda$| is increased from zero the asymptotic SMD becomes larger, because now the cooling radius will not shrink indefinitely at late times, hence increasing the star formation efficiency of haloes. However, for |$\alpha _{\Lambda } \gt 0.1$|, the trend is reversed. This might seem somewhat at odds with Fig. 8, which shows that the stellar mass fraction within haloes keeps increasing for |$\alpha _{\Lambda } \gt 0.1$|. However, Fig. 8 focuses on isolated haloes, without accounting for their multiplicity. Instead, the stellar production needs to be weighted by the HMF when computing the SMD. Thus, Fig. 9 tells us that the lower asymptotic SMD for large values of |$\Lambda$| is a consequence of the suppression of structure due to the cut-off of the HMF occurring at lower masses.

Upper panel: Cumulative stellar mass produced in a unit comoving volume up to a certain cosmic time. Different lines correspond to different values of |$\Lambda$|, as indicated in the colour bar. The upper x-axis shows the redshift that corresponds to the cosmic time, for the reference cosmology. Therefore, the correspondence between the upper and lower x-axes applies only to the case |$\Lambda = \Lambda _{\rm obs}$|. The break in the upper x-axis indicates that the ticks are omitted for |$z\gt -1+10^{-10}$|, since the exponential growth of the scale factor in the future makes it increasingly hard to represent the redshift on a meaningful scale. Lower panel: Same as in the upper panel, but normalized to the total stellar mass density produced over the entire history of the universe. The stellar mass fraction produced up to present time in our Universe is |$\sim$| 32 per cent, hence we are typical observers. However, in general the mid-point of star formation history does not coincide with |$\Lambda$|-domination in other universes, meaning that the ‘why-now’ problem is specific to |$\Lambda = \Lambda _{\rm obs}$| (see discussion in Section 5.1.1).
In the lower panel of Fig. 9, we normalize the SMD by the asymptotic value. For the observed value of |$\Lambda$| (green line), about 32 per cent of the stellar mass that will ever exist in our Universe has already been formed. In the original SP21 paper, this quantity was quoted to be close to 50 per cent, but in that case the SMD was computed differently. Namely, SP21 did not go through the detailed convergence study undertaken here, but simply integrated the CSFRD up to |$z=-0.9999$|, which corresponds to |$t\sim 200\, \rm Gyr$| in our Universe; Fig. 9 teaches us that this is not always sufficiently late to reach the true asymptotic stellar mass density.
For small positive |$\Lambda$|, the fractional SMD in the lower panel of Fig. 9 departs from the EdS solution much earlier than for the SMD evolution shown in the upper panel, reflecting the larger asymptotic values reached for |$0 \lt \Lambda \lt 0.1$|. Thus, similar values of the SMD at early times can still differ considerably in fractional terms. Instead, for large values of |$\Lambda$|, haloes are more efficient at producing stars, therefore the asymptotic SMD is reached earlier. Fig. 8 shows that for |$\Lambda \lesssim 10^3$|, star formation is driven by higher mass haloes. For extreme cosmologies such as |$\alpha _{\Lambda }=10^4$|, haloes are at nearly full star formation efficiency irrespective of their mass, but the cut-off of the HMF at late times occurs at |$M\sim 10^9 \, \rm {\rm M}_{\odot }$|. Recalling that we consider only haloes with |$T\gt 10^{4.5} \, \rm K$| as eligible for star formation, this means that stars are produced only by haloes in the narrow mass range |$10^8 \, \mathrm{M}_{\odot } \lesssim M \lesssim 10^9 \, \rm {\rm M}_{\odot }$|. Therefore, it takes longer to reach the asymptotic SMD.
In conclusion, if |$\Lambda$| is higher, halo formation ceases at early times: haloes are unable to attain masses as large as in our Universe (see Section 3.1), and are less efficient at forming stars. The dominant fraction of stars is formed in the future (i.e. at negative redshift), when haloes are isolated and efficient at forming stars due to their larger cooling rate with respect to their counterparts in our own Universe (see Section 3.2.1). This long-term future of star formation is expected to play an important role in determining the asymptotic efficiency of star formation in universes with a large cosmological constant.
4.2 Cosmic star formation efficiency
In the previous section, we showed that the comoving stellar mass density produced in a universe with |$\Lambda \ge 0$| asymptotes to a constant in the limit |$t \rightarrow \infty$|. We now want to understand how this asymptotic density, |$\rho _{* \, \infty }$|, varies as a function of |$\Lambda$|. In practice, we consider the cosmic star formation efficiency, defined as |$\varepsilon = \rho _{* \, \infty }/\bar{\rho }_{\rm b}$|, where |$\bar{\rho }_{\rm b}$| is the mean comoving baryon density in the universe.
We plot this cosmic stellar efficiency in Fig. 10, where we see the striking result that the peak of cosmic star formation efficiency occurs within one order of magnitude of |$\Lambda _{\rm obs}$|. At the peak efficiency, |$\Lambda \approx 0.1 \Lambda _{\rm obs}$|, the model predicts that about 27 per cent of the baryonic mass in the universe will eventually be converted into stars. In our Universe, this figure is 23 per cent. For larger values of |$\Lambda$|, the cosmic efficiency declines, becoming negligibly small as |$\alpha _{\Lambda }$| approaches |$10^5$|.

Right panel: Fraction of baryonic mass (per unit comoving volume) that is converted into stars over the entire history of the universe, as a function of |$\Lambda$|. This global stellar efficiency peaks around 1/10 of the observed value of |$\Lambda$| (thin black vertical line), and becomes negligibly small for large values of |$\Lambda$|. This is due to the suppression of cosmological structure formation caused by the earlier and larger acceleration of the expansion of the universe. The horizontal thin black line represents the global stellar efficiency in an EdS universe, which is smaller than in our Universe. The decrease in the cosmic star formation efficiency at small positive |$\Lambda$| is driven by astrophysical rather than cosmological factors: haloes of a given virial mass become larger as |$\Lambda$| decreases, and the consequent reduction in internal gas density diminishes the gas cooling rate, which in turn determines the SFR at low redshift (see discussion in Section 4.2).
It might seem surprising that the cosmic stellar efficiency does not fall monotonically as |$\Lambda$| increases. A larger cosmological constant would raise the acceleration of the cosmic expansion at earlier times, hence suppressing the gravitational collapse of galaxy-scale haloes; higher values of |$\Lambda$| thus correspond to smaller fractions of matter in collapsed structures (e.g. Weinberg 1987; Martel, Shapiro & Weinberg 1998). However, the picture is more complicated when considering the impact of |$\Lambda$| on the fraction of baryonic rather than total matter within haloes (see e.g. the discussion in Page 2011). Furthermore, in Fig. 10 we are concerned with the stellar component only, which is subject to the complex interplay of astrophysical processes and the evolution of the large-scale structure of the universe.
The suppression of the star formation efficiency at high |$\Lambda$| cannot be attributed to astrophysical processes, since Fig. 8 shows that high-mass haloes remain highly efficient at producing stars, even for large values of |$\Lambda$|. Therefore, the small |$\varepsilon$| found for large |$\Lambda$| must be caused by the cut-off of the HMF at small halo masses (see Section 3.1). Ultimately, this is a consequence of the earlier freeze-out that results from a large cosmological constant. In short, our model predicts that haloes are generally efficient at forming stars when |$\Lambda$| is high, but their number density is then so low that the overall effect is a reduction in the cosmic star formation efficiency.
For small values of |$\Lambda$|, the opposite is true. There is now stronger clustering and haloes can attain larger masses. However, for a given virial mass, lower values of |$\Lambda$| correspond to a larger virial radius, given that |$R \propto M^{1/3} H_0^{-2/3}$| and |$H_0$| is smaller for lower values of |$\alpha _{\Lambda }$| (see equation 21 and discussion in Section 3.2.2). Thus, a halo of a given virial mass has a lower internal density in a universe with a smaller cosmological constant. This reflects the continued progress of mergers, which dilute the internal density of haloes to a multiple of the cosmic density, up to the point where |$\Lambda$| dominates. The reduced gas density results in a smaller cooling rate, making star formation less efficient at late times. The cosmic stellar efficiency at low values of |$\Lambda$| in Fig. 10 is therefore suppressed for astrophysical reasons.
These astrophysics-driven and cosmology-driven suppression effects on star formation balance out around |$\alpha _\Lambda \sim 0.1$|. Fig. 10 suggests that the observed value of |$\Lambda$| is somewhat peculiar, in the sense that it is close to the optimal value for maximum cosmic star formation efficiency. However, this is not the same as saying that it is the most likely value of |$\Lambda$| for observers to experience. In the next section, we will address the implications of Fig. 10 for such anthropic issues.
5 DISCUSSION
5.1 Implications for anthropic arguments
5.1.1 Typicality and the ‘why now’ problem
The first aspect of anthropic arguments that we will address concerns the so-called why-now problem, i.e. the fact that we live at an usual epoch when the density of matter and dark energy are of the same order of magnitude. The ‘weak’ anthropic approach to this puzzle asserts that typical observers should live when the age of the universe corresponds roughly to the lifetime of Sun-like stars, assuming these to be necessary for forming and sustaining habitable planets (Carter 1974). We can test this principle to some extent, by asking how typical we are as observers. Consider the fraction of stellar mass density produced up to present time in our own Universe only (green line in Fig. 9): this predicts that roughly 32 per cent of all the stars that will ever exist are already in place, so it appears that the present time is indeed typical, lying neither at the beginning nor at the end of the cosmic star formation history.
However, this argument neglects the question of how long intelligent life may to take to emerge following the formation of a star. This is an issue on which cosmologists are unlikely to have a reliable opinion. However, we can note that if the emergence of life is rare, then indeed an offset of the order of the typical stellar lifetime might be expected, and this is the basis of the standard anthropic argument for the current age of the Universe. Alternatively, we can follow Peacock (2007) and eliminate the biological element by noting that the Sun formed at |$z\approx 0.45$|. From our results, we find that |$\sim \,$|25 per cent of the asymptotic stellar density was in place at that point, so the Sun is a slight ‘early adopter’ but by no means atypical.
The Sun thus formed almost exactly at the point of matter-|$\Lambda$| equality, and the interesting question is whether this coincidence would still have held for a typical star if |$\Lambda$| had taken a very different value. It is apparent from the predictions in Fig. 9 that the observed coincidence is not to be expected in all cases. For |$\alpha _\Lambda =100$|, the 50 per cent point of stellar production is reached at |$t\approx 20\, \rm Gyr$|. That corresponds to |$1+z\approx 10^{-6}$| in that universe, at which point the density from |$\Lambda$| exceeds that of matter by a factor of approximately |$10^{14}$|. We also note that the predicted CSFRD is not so very different from the observed one when we carry out the prediction for the EdS universe. In that case, there is still a well-defined typical era of star formation, which is clearly not dictated by the cosmological constant.
Our conclusion is then that the why-now puzzle is specific to the actual value of |$\Lambda$|. Weak anthropic selection that biases observers towards living at a special time in the universe would not have produced such a coincidence in general. We are therefore driven to ask if there is some more general selection effect that favours the observation of this particular value of |$\Lambda$|.
5.1.2 Ensembles and the observed value of |$\Lambda$|
The multiverse approach to |$\Lambda$| considers an ensemble of different universes, each with its own cosmological constant, and argues that only the universes that have a small enough value of |$\Lambda$| would be able to host observers. As was emphasized by Efstathiou (1995), this situation can be analysed in a Bayesian framework. If we consider the existence of observers as the ‘data’, then the posterior on |$\Lambda$| is
In the equation above, |$p(\rm O \vert \Lambda)$| is the likelihood, i.e. the probability of observers to emerge for a given value of |$\Lambda$|, and |$p(\Lambda$|) is the prior on |$\Lambda$|.
Both the prior on |$\Lambda$| and the likelihood are non-trivial to determine. The main conceptual difficulty is that one would need to choose a well-defined measure for weighting different members of the ensemble. However, for an inflationary multiverse, the bubble universes are formally of infinite volume, and it is therefore not clear how to define the measure (see, e.g., Wenmackers 2024). However, if |$\Lambda$| is the sole parameter that varies across the ensemble, then the problem can be evaded following Weinberg (1987, 1989, 2000a, b), who pointed out that there is no known physical mechanism that preferences |$\Lambda =0$|. Therefore, the prior on |$\Lambda$| should be flat at least in a neighbourhood of |$\Lambda =0$|. We will adopt Weinberg’s convention of a flat prior, and extend it to the full range of values of |$\Lambda$| considered in our work, i.e. |$0 \le \alpha _\Lambda \le 10^5$|.
Regarding the likelihood |$p(\rm O \vert \Lambda)$|, we will adhere to the view expressed in the previous section that stars are a precondition for observers. Thus, we will simply assume that the likelihood is proportional to the cosmic star formation efficiency |$\varepsilon$|. We stress that we are agnostic concerning the time at which observers appear, unlike in the previous section: |$\varepsilon$| is an integrated quantity over the entire history of the universe.
We can now write the posterior probability density in the interval |$(\Lambda ,\, \Lambda + \mathrm{d}\Lambda)$| as:
or, more conveniently, as the distribution of |$\ln \Lambda$|:
This probability distribution is shown in the left panel of Fig. 11 (see also Sorini 2022). We notice immediately that the peak corresponds to |$\alpha _{\Lambda }\approx 800$|: the most likely order of magnitude of |$\Lambda$|, weighted by observer number across the multiverse, is vastly greater than the observed value. Furthermore, |$\alpha _{\Lambda }=1$| (indicated with a vertical black line in the plot) appears to be an outlier of the distribution. To quantify this more precisely, we also show the corresponding cumulative probability distribution in the right panel of Fig. 11 – emphasizing once again that this neglects the entirely distinct case of recollapsing models. We then find that the probability of observing a cosmological constant equal to or below the value of |$\Lambda$| observed in our Universe is around 0.5 per cent. Thus, with our stated assumptions of (1) the extension to the SP21 model of cosmic star formation presented in this work; (2) a flat prior on |$\Lambda$|; and (3) assuming a direct proportionality between the cosmic star-forming efficiency and the likelihood of generating observers at a given |$\Lambda$|, the simplest multiverse framework is not able to account naturally for the small observed value of |$\Lambda$|.

Left panel: probability of generating observers in a given universe, at any point in time, as a function of |$\Lambda$|, assuming a direct proportionality to the stellar efficiency shown in Fig. 10 (right panel) and imposing a flat prior on |$\Lambda$|. The peak occurs around |$\Lambda \approx 800 \Lambda _{\rm obs}$|. The observed value of |$\Lambda$| (thin black vertical line) is an outlier in the distribution. Right panel: cumulative probability function of the distribution shown in the left panel. The probability that |$\Lambda \le \Lambda _{\rm obs}$| is |$\sim$| 0.5 per cent, meaning that this value is not anthropically favoured.
5.2 Alternative anthropic assumptions
The above failure to account for the magnitude of |$\Lambda$| can mean one of three things: (i) there is no multiverse; (ii) members of the ensemble vary in a more complicated way than simply altering |$\Lambda$|; and (iii) some of our basic assumptions are inappropriate.
Undoubtedly, it is interesting to consider a more complex multiverse ensemble, where more fundamental parameters are varied at the same time (up to 31: Tegmark et al. 2006). Indeed, other works have studied the simultaneous impact of varying two or three cosmological parameters on the star formation history (Cline, Frey & Holder 2008; Bousso & Hall 2013). For instance, a universe with both a larger value of |$\Lambda$| and a larger amount of matter may well be evolving in a similar fashion as our own Universe, as the increased acceleration of the expansion induced by a higher cosmological constant would be compensated by a stronger clustering of matter. More radically, we may consider altering other fundamental parameters, such as the gravitational constant. Increasing G would obviously favour the collapse of structures, and may counteract the effect of a larger |$\Lambda$|. At the same time, such a change would also affect the structure of haloes and subsequently the cooling rate. If one changes other coupling constants, such as the fine structure constant, then the picture becomes increasingly more complex. That would affect the very structure of matter, and modify cooling rates in unusual ways. A comprehensive discussion on how star formation would proceed in such exotic scenarios is provided by Adams (2019). We will explore these generalizations and other matters in future work.
The problem with all these more elaborate ensembles is the difficulty in deciding on the prior for the additional parameters. We have adopted Weinberg’s prior on |$\Lambda$|, arguing that it is a special case. Nevertheless, it is interesting to consider what would happen if Weinberg’s argument was in error, since a prior that gives more weight to lower values of |$\Lambda$| could radically change the conclusions. A prior that is flat in |$\ln (\alpha _{\Lambda })$| would yield by construction a posterior for |$\Lambda$| that is identical to the cosmic efficiency |$\varepsilon$|, properly renormalized. However, a logarithmic prior cannot be applied if we allow for |$\Lambda =0$|, let alone if we accept that |$\Lambda \lt 0$| is possible. A power-law prior on |$\Lambda$| would be allowed if we insist that |$\Lambda \ge 0$|, and we show posteriors on |$\Lambda$| for various power-law slopes in Fig. 12. As the power law becomes steeper, the probability distribution becomes broader, and the peak moves towards smaller values of |$\alpha _{\Lambda }$|.
![Same as in the left panel of Fig. 11, but for different priors on $\Lambda$. The red line corresponds to the flat prior, and is therefore the same quantity plotted in the left panel of Fig. 11. The other lines refer to a power-law prior, as reported in the legend. A power-law prior with an index in the range $(-1, -0.9]$ would strongly favour the observed value of $\Lambda$, but of course such a prior would need to be motivated by a theory and not chosen a posteriori in order to validate anthropic reasoning (see discussion in the main text).](https://oup.silverchair-cdn.com/oup/backfile/Content_public/Journal/mnras/535/2/10.1093_mnras_stae2236/1/m_stae2236fig12.jpeg?Expires=1750605377&Signature=T19d-R8I6irknft54-rg24v5WK445pM7Dp4g7k2Dc5NCbmTohEHDVMOJCuIR3AMpcNWgq1dTH8C3iG6Q2RNfeRRz~aOjUXtgkszuLqCj3z5nRK7uTgkCYyewWXWjlmERFxPuovWeUYD9ynrfpyWIHUp1UMAZH7tjoKmOmM2nw9YZZYvrk4idVjeTv4EXhAAbEu4-2KYMnsMLGFgeeFVsnBUNHAlbC9Rz6Jjb8c2cmd9z-LLSBXd2C~XwSinkgduhJ3vJiNEEpGAq~3MbFGO6NW-s5l5IeNCik2PXFW-vNfqS6JDJIlvV8oaQDmtE7MsK6z~PQ5ScvdZ315JPJ-QIuQ__&Key-Pair-Id=APKAIE5G5CRDK6RD3PGA)
Same as in the left panel of Fig. 11, but for different priors on |$\Lambda$|. The red line corresponds to the flat prior, and is therefore the same quantity plotted in the left panel of Fig. 11. The other lines refer to a power-law prior, as reported in the legend. A power-law prior with an index in the range |$(-1, -0.9]$| would strongly favour the observed value of |$\Lambda$|, but of course such a prior would need to be motivated by a theory and not chosen a posteriori in order to validate anthropic reasoning (see discussion in the main text).
Of particular interest is the case of a prior in the form |$\alpha _{\Lambda }^{-0.75}$|. If one assumes that the cosmological constant is set by the energy of the vacuum |$E_{\rm v}$|, then |$\Lambda \propto E_{\rm v}^4$|. One might then opt for a flat prior on |$E_{\rm v}$|, corresponding to a prior on |$\Lambda$| in the form |$\alpha _{\Lambda }^{-0.75}$|. In this case our model would yield a probability of observing |$\Lambda \le \Lambda _{\rm obs}$| of 39 per cent, very well consistent with observation. However, the prior should be motivated by a priori theoretical considerations, rather than being chosen in order to achieve by construction some preferred conclusion. At present, we do not have a compelling argument that would justify a uniform prior in |$E_{\rm v}$|, nor for the associated assumption that |$\Lambda \gt 0$|.
There have been several attempts at solving the cosmological constant problem that predict a non-negative value. Supersymmetry generates a strictly positive vacuum energy density if it is spontaneously broken (see e.g. Martin 2012). Kaloper & Padilla (2014) proposed a reformulation of General Relativity whereby all vacuum energy is sequestered from the matter sector. This prevents the vacuum energy from sourcing spatial curvature and predicts a positive cosmological constant, which is expressed as the space-time average of the trace of the energy–momentum tensor. However, the model requires the universe to end in order to produce a positive |$\Lambda$|. This is not necessary in a later generalization of the sequestering model by Lombriser (2019b), who revisited the work of Kaloper & Padilla (2014) by explicitly modelling the effect of collapsed structures. A positive |$\Lambda$| is also predicted by models attempting to solve the cosmological constant problem by varying the Planck mass (Lombriser 2019a; Sobral-Blanco & Lombriser 2020, 2021), at least in their simplest incarnation. The basic idea is to vary the Einstein–Hilbert action with respect to both the metric and the Planck mass over a manifold. The resulting cosmological constant can be expressed in terms of the average energy–momentum tensor, and can be interpreted as a ‘backreaction’ of structures on the expanding universe. Assuming the typical size of collapsed structures, the model predicts |$\Omega _{\Lambda }=0.704$|, very close to the observed value. Recently, Gaztañaga (2021) proposed that |$\Lambda$| represents a zero-action boundary term in the Einstein–Hilbert action. It can be derived that |$\Lambda$| is proportional to the average density of the universe and is thus a non-negative quantity.
There are, however, other strong theoretical grounds for considering a negative cosmological constant. In supergravity, the vacuum energy can be negative (Martin 2012). In the landscape of string theory, negative values of |$\Lambda$| emerge naturally, and predicting a small positive value remains a challenge (see e.g. Lerh Feng et al. 2021, Demirtas et al. 2022). These scenarios would then not be compatible with a power-law or logarithmic prior on |$\Lambda$|.
Lastly, our conclusions would certainly be subject to change if we were to alter the model of cosmic star formation. We next discuss how the limitations of our present modelling could affect our conclusions in Section 5.3. In Section 5.4, we will then make a detailed comparison of our findings with other related works in the literature.
5.3 Limitations of our study
This paper is based on the SP21 model of cosmological star formation. This aimed to retain the physical transparency of the original Hernquist & Springel (2003) model, which comes at the price of some simplifications. Chief among these is the fact that the model formalism does not explicitly include AGN feedback, but only supernova-driven winds. Nevertheless, the parameters underlying the stellar feedback model in the SP21 formalism are tuned to reproduce observations of the baryonic Tully–Fisher relationship (McGaugh et al. 2000) and of the Kennicutt–Schmidt relationship (Kennicutt 1998) in our Universe, which does harbour AGN feedback. Thus, the effect of AGN feedback is implicitly allowed for, as in the case of cosmological simulations that contain stellar feedback only, but which calibrate the parameters of their subgrid prescriptions to reproduce observed data. This approach was followed, for example, by Oh et al. (2022), who analysed the impact of |$\Lambda$| in a suite of hydrodynamical simulations that extended up to |$t=97 \, \rm Gyr$|, including only stellar feedback. They argued that any energy fed back into the environment by black holes at late times will happen on a time-scale correlated to the feedback from star formation. Then, even without including an AGN feedback model, their feedback parameters would capture the global energy input that would be expected from both star formation and black holes. The empirical correlation of stellar mass and black hole mass in galaxies is commonly taken as evidence for co-evolution of these two populations (see e.g. Kormendy & Ho 2013; Heckman & Best 2014), and the neglect of explicit AGN feedback here is justifiable as long as co-evolution holds exactly.
However, several hydrodynamical simulations do indicate that AGN-driven jets are an important mechanism in the quenching of star formation after |$z=2$| (see e.g. Vogelsberger et al. 2013; Weinberger et al. 2017; Pillepich et al. 2018; Sorini et al. 2022; Scharré et al. 2024). While other groups found a more limited impact of AGN feedback on the low-redshift CSFRD (McCarthy et al. 2017; Salcido et al. 2018; see also the discussion in Salcido et al. 2020), the general consensus is that AGN feedback should be explicitly included in realistic models of the star formation history. Indeed, the total energy input from feedback processes is not the sole important factor in determining the quenching of star formation. Cosmological simulations show that the way such energy propagates, for example through diffuse winds or collimated AGN-driven jets, has a potentially significant impact on the star formation history and the diverse gaseous phases within haloes (e.g. Weinberger et al. 2017; Sorini et al. 2022; Scharré et al. 2024; Yang et al. 2024). However, such considerations are beyond the scope of this work, and their incorporation would require a major expansion of the SP21 model.
As SP21 pointed out, the introduction of AGN feedback in their formalism could plausibly diminish the low-z CSFRD and hence also the long-term efficiency of star formation. Because most of the stellar mass in universes with large values of |$\Lambda$| is formed in the far future (see Fig. 9), AGN feedback may suppress the cosmic star formation efficiency at the high-|$\Lambda$| end more rapidly. This would in turn increase the probability of observing a value of |$\Lambda$| that is |$\le \Lambda _{\rm obs}$|, reducing the discrepancy with observation. It remains to be seen whether the explicit inclusion of AGN feedback would significantly alter the conclusions in our work regarding the anthropically favoured value of |$\Lambda$|, or whether its impact would be sub-dominant. We leave this important analysis, which requires incorporating a first-principles model of AGN feedback within our formalism, for future work.
Other potentially important physical ingredients omitted from the SP21 model are chemical evolution and metal enrichment of the interstellar and intergalactic media. Following the evolution of the metallicity of stars and gas within and outside haloes would alter the cooling function over time. At later times, and reasonably also in the future of the universe, the metallicity would increase. For virial temperatures above |$T=10^4 \, \rm K$|, which are those of interest for the SP21 model, the cooling function becomes larger as the metallicity increases (see Sutherland & Dopita 1993). Hence, in the future of the universe the cooling rate would become higher, increasing the star formation efficiency.
We also acknowledge that more exotic effects might become significant in the far future of the universe. For instance, gas cooling via inverse Compton scattering of CMB photons might become an important additional cooling mechanism in the far future (Adams & Laughlin 1997; Bousso & Leichenauer 2010; Bousso & Hall 2013). As a further example, Adams & Laughlin (1997) argued that gravitational instabilities within galaxies could deviate the trajectory of stars from their orbits, so that stars either collapse towards the centre of the galaxy, or drift away from it, becoming isolated and unbound objects. In the former scenario, there can be a secondary spark of star formation, triggered by the clash between brown dwarfs around the galactic centre. Such exotic possible channels of star formation occur in the extreme far future, |$t\approx 10^{15}\, \rm Gyr$|. The ejection of stars from the galaxy due to the accumulation of stellar encounters would occur on a smaller, but still colossal, time-scale, |$t\approx 10^{10}-10^{11}\, \rm Gyr$| (Adams & Laughlin 1997). These processes occur well after the point where our calculations of star formation have converged, and there is no reason to believe that they will dominate over previous activity in terms of total stellar production.
One can also ask whether it is reasonable to believe that ‘conventional star formation’ proceeds all the way up to the |$10^{10}\, \rm Gyr$| time-scale represented in Fig. 9. Adams & Laughlin (1997) argued that, based on the life-time of typical red dwarfs, ‘conventional’ star formation should halt on a time-scale of |$t\sim 10^4\, \rm Gyr$|. If we imposed such a cut-off in the star formation history, that would lead to a lower star formation efficiency in universes with lower |$\Lambda$|. On the contrary, it would not affect as significantly universes with higher |$\Lambda$|, where the bulk of star formation occurs at earlier times (see Fig. 9). Therefore, the posterior on |$\Lambda$| would be more skewed towards higher values, hence making the anthropic prediction even less consistent with observation. However, arguments on the time-scales of different channels of star formation in the far future of the universe are themselves fairly uncertain. We therefore opted for not introducing an additional ill-defined parameter in our model, and simply integrated the star formation history up to |$t \rightarrow \infty$|.
Throughout this work, we have assumed that generation of observers is simply proportional to the cosmic star formation efficiency. More realistically, we might consider what is known about statistics of exoplanets as a function of stellar mass. This variation would not matter provided the IMF was invariant, but given that the typical density and temperature of star-forming gas is likely to change with epoch, there is ample scope for the IMF to alter. We have neglected this complication, as have most authors in this area.
Another option could be to weight the CSFRD by the metallicity of the stars formed. Indeed, Dayal et al. (2015) argued that the stellar mass, total metal mass, and SFR all conspire in enhancing the habitability of galaxies such as metal-rich giant elliptical galaxies. Furthermore, accounting for the effect of radiation on the emergence of life on potentially habitable planets (e.g. Totani et al. 2019; Gobat et al. 2021) could bring in the structure of the galaxy as a relevant factor. Indeed, recent works defined a ‘galactic habitable zone’ (Gobat & Hong 2016), in analogy to the circumstellar habitable zone (Huang 1959; Kasting et al. 1993). However, these variations in the generation of observers within and between galaxies are not a concern for the main focus of this paper, which is on the global production of observers summed over the universe and over all cosmic time. Therefore, it seems reasonable to neglect such issues for the present.
In conclusion, our model can certainly be improved in different ways, and we hope to address some of them in the future. However, despite its simplicity the model is capable of capturing the main physical processes that are relevant for a sound description of the star formation history in the universe, both in the past and the future, and in different cosmologies.
5.4 Comparison with previous work
5.4.1 Posterior on |$\Lambda$| and anthropic considerations
In this section, we will discuss our results on the posterior distribution on |$\Lambda$| and the implications for anthropic arguments in the context of previous relevant work in the literature.
Oh et al. (2022) ran a suite of hydrodynamical simulations based on the enzo code (Bryan et al. 2014), for an EdS universe and for |$\Lambda$|CDM universes with different values of |$\Lambda$|, up to 100 times the observed value. The simulations included a subgrid model for stellar feedback (Oh et al. 2020), but no AGN feedback. This model was shown to yield physically sensible results even when run up to a cosmic time of |$100\, \rm Gyr$| (Oh et al. 2021). Using the same assumptions adopted in this work concerning the probability of generating observers in a given universe, Oh et al. (2022) found that the median of the posterior distribution on |$\Lambda$| corresponds to |$8.0 \, \Lambda _{\rm obs}$|, with 95 per cent of the distribution spanning the range |$\alpha _{\Lambda } = 0.32-105$|, and thus consistent with the observed |$\Lambda$|. They find that the probability that |$\Lambda \le \Lambda _{\rm obs}$| is 12.9 per cent. Oh et al. (2022) also consider the traditional assumption for the anthropic weighting, based on the single-scale approximation for the collapsed fraction of matter in galaxy-scale haloes (see Peacock 2007). Under this assumption, the median of the posterior distribution on |$\Lambda$| rises to |$12.6 \, \Lambda _{\rm obs}$|, and the corresponding probability of an observer measuring |$\Lambda \le \Lambda _{\rm obs}$| is 8.2 per cent. The overall conclusions from the detailed modelling in Oh et al. (2022) are thus not hugely different from the simple single-scale approach. These conclusions are qualitatively in accord with our own, to the extent that the observed value of |$\Lambda$| is a low outlier in the ensemble. However, we find an unacceptably small probability of |$\Lambda \le \Lambda _{\rm obs}$| – whereas Oh et al. (2022) find a probability that, while small, does not strongly rule out the multiverse hypothesis.
An earlier example of a large ‘multiverse simulation’ is given by Barnes et al. (2018). They ran a suite of variants of the EAGLE hydrodynamical simulation (Schaye et al. 2015), for values of |$\alpha _{\Lambda }$| in the range |$0-300$|, until cosmic time |$t=20.7\, \rm Gyr$|. These simulations include recipes for both stellar feedback and AGN feedback. Barnes et al. (2018) consider three possible models for the generation of observers:
Observers are generated following the formation of stars, after a fixed delay of |$5\, \rm Gyr$|, which would represent the time-scale for the evolution of intelligent life;
The rate of generation of observers at a given time around a certain stellar population is proportional to the fraction of stars that are still on the main sequence, to account for the fact that life appears to be a rare phenomenon in the universe, and therefore should occur towards the end of the main-sequence lifetime of a star (Carter 1983; Barrow & Tipler 1986);
The rate of generation of observers at a given time is proportional to the metallicity of the fraction of stars that are still on the main sequence for the stellar population considered, to reflect the argument that life would more easily emerge on rocky planets, which would be more likely to form around stars with higher metallicity (Gonzalez 1997; Fischer & Valenti 2005; but see also Buchhave & Latham 2015; Wang & Fischer 2015).
The main impact of these different choices is to change when the observers form, rather than the total eventual numbers, which is our main concern. The results for all these models are reported in Table 1, in correspondence of the cells in the second column, labelled ‘star formation + delay’, ‘star formation + lifetime’, and ‘star formation + metals’, respectively. For each of the aforementioned models, Barnes et al. (2018) further considered three possible ways of constructing a measure in the multiverse:
a mass-weighted measure, whereby a given mass element in the universe can inhabit a region with any value of |$\Lambda$| with equal probability;
a causal-patch measure, in which regions of the universe that have equal comoving volume of their causal patch at a certain cosmic time host the same number of observers, for any value of |$\Lambda$| considered;
a causal-diamond measure, analogous to the previous measure, except that the comoving volume of the region (at given cosmic time t) that is enclosed by a photon that departs from a world line at reheating and returns at t is considered instead of the volume of the causal patch.
It is important to note that in all cases described above Barnes et al. (2018) consider a flat prior on |$\Lambda$|. Thus their posterior on |$\Lambda$| is the product of the measure times the weighting based on star formation. The mass-weighted measure is essentially what we adopt in this work.
The results that Barnes et al. (2018) obtain for all observer generation models and all multiverse measures are listed in Table 1. The mass-weighted measures yields a median of the posterior distribution on |$\Lambda$| of 45 or 59, depending on the observer generation prescription. The 84|$^{\rm th}$| percentile of the distribution can reach values as large as |$194\, \Lambda _{\rm obs}$|. On the other hand, the causal patch and causal diamond measures move the median |$\Lambda$| to much smaller values, and always below |$\Lambda _{\rm obs}$|. Consequently, the mass-weighted measure would predict a significantly smaller probability of observing |$\Lambda \le \Lambda _{\rm obs}$| (1.9 per cent–2.5 per cent) than the other two measures (73 per cent–93 per cent). Thus, the choice of the measure appears to be crucial.
Median, |$16^{\rm th}$| and |$84^{\rm th}$| percentiles of the posterior distribution on |$\alpha _{\Lambda }$| for different choices of the prior on |$\Lambda$| or of the measure of observers in the multiverse, and for different models. The table also shows the probability of observers measuring a value of |$\Lambda$| equal to or smaller than |$\Lambda _{\rm obs}$|. Unless otherwise explicitly stated in the main text, results from other works are quoted directly from the respective manuscripts, or derived from fitting formulae therein.
Star formation model . | Observers model . | Prior/measure . | |$16^{\rm th}$| percentile . | Median . | |$84^{\rm th}$| percentile . | |$p(\Lambda \le \Lambda _{\rm obs})$| . |
---|---|---|---|---|---|---|
This work | Star formation | Flat in |$\alpha _{\Lambda }$| | 73 | 539 | 2306 | 0.52 per cent |
|$\alpha _{\Lambda }^{-0.1}$| | 44 | 389 | 1868 | 1.0 per cent | ||
|$\alpha _{\Lambda }^{-0.75}$| | 0.041 | 2.8 | 84 | 39 per cent | ||
|$\alpha _{\Lambda }^{-0.9}$| | |$6.9 \times 10^{-4}$| | 0.18 | 12 | 65 per cent | ||
Oh et al. (2022) | Star formation | Single-scale approximation | 2.3 | 12.6 | 43 | 8.2 per cent |
Simulation | 1.4 | 8.0 | 34 | 12.9 per cent | ||
Barnes et al. (2018) | Star formation + delay | Mass weighted | 10 | 59 | 194 | 1.9 per cent |
Causal patch | 0.04 | 0.34 | 0.96 | 86 per cent | ||
Causal diamond | 0.13 | 0.65 | 3.15 | 73 per cent | ||
Star formation + lifetime | Mass weighted | 10 | 59 | 194 | 1.9 per cent | |
Causal patch | 0.009 | 0.089 | 0.849 | 90 per cent | ||
Causal diamond | 0.01 | 0.25 | 0.96 | 86 per cent | ||
Star formation + metals | Mass weighted | 8 | 45 | 163 | 2.5 per cent | |
Causal patch | 0.004 | 0.07 | 0.71 | 93 per cent | ||
Causal diamond | 0.01 | 0.17 | 0.87 | 90 per cent |
Star formation model . | Observers model . | Prior/measure . | |$16^{\rm th}$| percentile . | Median . | |$84^{\rm th}$| percentile . | |$p(\Lambda \le \Lambda _{\rm obs})$| . |
---|---|---|---|---|---|---|
This work | Star formation | Flat in |$\alpha _{\Lambda }$| | 73 | 539 | 2306 | 0.52 per cent |
|$\alpha _{\Lambda }^{-0.1}$| | 44 | 389 | 1868 | 1.0 per cent | ||
|$\alpha _{\Lambda }^{-0.75}$| | 0.041 | 2.8 | 84 | 39 per cent | ||
|$\alpha _{\Lambda }^{-0.9}$| | |$6.9 \times 10^{-4}$| | 0.18 | 12 | 65 per cent | ||
Oh et al. (2022) | Star formation | Single-scale approximation | 2.3 | 12.6 | 43 | 8.2 per cent |
Simulation | 1.4 | 8.0 | 34 | 12.9 per cent | ||
Barnes et al. (2018) | Star formation + delay | Mass weighted | 10 | 59 | 194 | 1.9 per cent |
Causal patch | 0.04 | 0.34 | 0.96 | 86 per cent | ||
Causal diamond | 0.13 | 0.65 | 3.15 | 73 per cent | ||
Star formation + lifetime | Mass weighted | 10 | 59 | 194 | 1.9 per cent | |
Causal patch | 0.009 | 0.089 | 0.849 | 90 per cent | ||
Causal diamond | 0.01 | 0.25 | 0.96 | 86 per cent | ||
Star formation + metals | Mass weighted | 8 | 45 | 163 | 2.5 per cent | |
Causal patch | 0.004 | 0.07 | 0.71 | 93 per cent | ||
Causal diamond | 0.01 | 0.17 | 0.87 | 90 per cent |
Median, |$16^{\rm th}$| and |$84^{\rm th}$| percentiles of the posterior distribution on |$\alpha _{\Lambda }$| for different choices of the prior on |$\Lambda$| or of the measure of observers in the multiverse, and for different models. The table also shows the probability of observers measuring a value of |$\Lambda$| equal to or smaller than |$\Lambda _{\rm obs}$|. Unless otherwise explicitly stated in the main text, results from other works are quoted directly from the respective manuscripts, or derived from fitting formulae therein.
Star formation model . | Observers model . | Prior/measure . | |$16^{\rm th}$| percentile . | Median . | |$84^{\rm th}$| percentile . | |$p(\Lambda \le \Lambda _{\rm obs})$| . |
---|---|---|---|---|---|---|
This work | Star formation | Flat in |$\alpha _{\Lambda }$| | 73 | 539 | 2306 | 0.52 per cent |
|$\alpha _{\Lambda }^{-0.1}$| | 44 | 389 | 1868 | 1.0 per cent | ||
|$\alpha _{\Lambda }^{-0.75}$| | 0.041 | 2.8 | 84 | 39 per cent | ||
|$\alpha _{\Lambda }^{-0.9}$| | |$6.9 \times 10^{-4}$| | 0.18 | 12 | 65 per cent | ||
Oh et al. (2022) | Star formation | Single-scale approximation | 2.3 | 12.6 | 43 | 8.2 per cent |
Simulation | 1.4 | 8.0 | 34 | 12.9 per cent | ||
Barnes et al. (2018) | Star formation + delay | Mass weighted | 10 | 59 | 194 | 1.9 per cent |
Causal patch | 0.04 | 0.34 | 0.96 | 86 per cent | ||
Causal diamond | 0.13 | 0.65 | 3.15 | 73 per cent | ||
Star formation + lifetime | Mass weighted | 10 | 59 | 194 | 1.9 per cent | |
Causal patch | 0.009 | 0.089 | 0.849 | 90 per cent | ||
Causal diamond | 0.01 | 0.25 | 0.96 | 86 per cent | ||
Star formation + metals | Mass weighted | 8 | 45 | 163 | 2.5 per cent | |
Causal patch | 0.004 | 0.07 | 0.71 | 93 per cent | ||
Causal diamond | 0.01 | 0.17 | 0.87 | 90 per cent |
Star formation model . | Observers model . | Prior/measure . | |$16^{\rm th}$| percentile . | Median . | |$84^{\rm th}$| percentile . | |$p(\Lambda \le \Lambda _{\rm obs})$| . |
---|---|---|---|---|---|---|
This work | Star formation | Flat in |$\alpha _{\Lambda }$| | 73 | 539 | 2306 | 0.52 per cent |
|$\alpha _{\Lambda }^{-0.1}$| | 44 | 389 | 1868 | 1.0 per cent | ||
|$\alpha _{\Lambda }^{-0.75}$| | 0.041 | 2.8 | 84 | 39 per cent | ||
|$\alpha _{\Lambda }^{-0.9}$| | |$6.9 \times 10^{-4}$| | 0.18 | 12 | 65 per cent | ||
Oh et al. (2022) | Star formation | Single-scale approximation | 2.3 | 12.6 | 43 | 8.2 per cent |
Simulation | 1.4 | 8.0 | 34 | 12.9 per cent | ||
Barnes et al. (2018) | Star formation + delay | Mass weighted | 10 | 59 | 194 | 1.9 per cent |
Causal patch | 0.04 | 0.34 | 0.96 | 86 per cent | ||
Causal diamond | 0.13 | 0.65 | 3.15 | 73 per cent | ||
Star formation + lifetime | Mass weighted | 10 | 59 | 194 | 1.9 per cent | |
Causal patch | 0.009 | 0.089 | 0.849 | 90 per cent | ||
Causal diamond | 0.01 | 0.25 | 0.96 | 86 per cent | ||
Star formation + metals | Mass weighted | 8 | 45 | 163 | 2.5 per cent | |
Causal patch | 0.004 | 0.07 | 0.71 | 93 per cent | ||
Causal diamond | 0.01 | 0.17 | 0.87 | 90 per cent |
However, the measure is not a degree of freedom that can be adjusted at will (see the discussion in Barnes et al. 2018). We confess that we find the motivation for the causal measures obscure: in our own Universe there are definitely pairs of galaxies that are causally disconnected (those in opposite directions on the sky at high redshift), and yet each can potentially host observers. Barnes et al. (2018) mention theoretical motivations for choosing the causal patch and diamond measures in the context of the quantum xeroxing paradox (Susskind, Thorlacius & Uglum 1993; Bousso 2006; Bousso, Freivogel & Yang 2006) and holographic probability in eternal inflation (Bousso 2006). However, it seems to us that the results are biased, and systematically exclude potential observers at early times when the horizon is small. We will therefore compare our predictions with the results in Barnes et al. (2018) arising from their mass-weighted measure only.
We can see that our results predict a median value of |$\Lambda$| that is |$9-12$| times larger than that of Barnes et al. (2018), and our probability for |$\Lambda \le \Lambda _{\rm obs}$| is about |$\sim$| |$4-5$| times smaller. However, we integrate the CSFRD all the way to |$t\rightarrow \infty$|, whereas Barnes et al. (2018) halt at |$20.7\, \rm Gyr$| in all cases. Nevertheless, we verified that even if we stopped our integration at this same early time, the probability of measuring |$\Lambda \le \Lambda _{\rm obs}$| increases only slightly to 0.7 per cent. This increase does not qualitatively affect our conclusions.
Anthropic calculations have also been pursued using semi-analytical methods. Sudoh et al. (2017) adopted the semi-analytical model of galaxy formation by Nagashima & Yoshii (2004), coupled with a Monte Carlo calculation for the merger history of dark matter haloes. Adopting a flat prior on |$\Lambda$| and using stellar mass density as a proxy for observer weighting, they obtained a median value of |$\Lambda$| in the range |$\alpha _{\Lambda } = 9.1-12$|, depending on the exact details of their model. Overall, the distribution is broad, and their probability of observing |$\Lambda \le \Lambda _{\rm obs}$| is between 6.7 per cent and 9.7 per cent. Such values are in agreement with the results by Oh et al. (2022), but are in tension with our conclusions and with Barnes et al. (2018). It is not surprising that Sudoh et al. (2017) obtain a lower median |$\alpha _{\Lambda }$| compared with us. Rather than integrating the mass of cool gas that is converted into stars over the entire history of the universe, they estimate the number of observers by calculating the total stellar mass density at the fixed cosmic time |$t=15\, \rm Gyr$|. The authors argue that this is a reasonable approximation, but it is not an assumption that is supported by our model: rather, we predict that in universes with a large value of |$\Lambda$| the majority of stars form at |$t \gg 15\, \rm Gyr$| (see e.g. Fig. 9).
To summarize, a common feature of all work on the simplest anthropic ensemble is that the posterior probability distribution on |$\Lambda$| is broad, with a median that is substantially larger than |$\Lambda _{\rm obs}$|. The exact value is heavily dependent on the underlying star formation model and the probability measure, but in all works considered the observed value of the cosmological constant seems to be located in the tail of the posterior for |$\Lambda$| (except for the causal patch and diamond measures in Barnes et al. 2018). Nevertheless, in the previous literature discussed here the probability of measuring |$\Lambda \le \Lambda _{\rm obs}$| is still of the order of a few per cent or even |$\sim$| 10 per cent, thus failing to reject the multiverse model. In contrast, our results pose a challenge in this respect, given that this probability is only 0.5 per cent. The viability of Weinberg’s simple anthropic explanation for |$\Lambda$| thus depends on which of these calculations is closest to a correct physical prediction of star formation in counter-factual universes. All existing studies have limitations, and ours is definitely no exception. However, in yielding a result that is so far from the observed level of |$\Lambda$|, it raises the motivation for continuing to seek a more realistic treatment of this problem.
5.4.2 Dependence of the cosmic stellar efficiency on |$\Lambda$|
The conclusions on the posterior distribution of |$\Lambda$| reached by the different works discussed in the previous section are predicated on the dependence of the comic star formation efficiency on |$\Lambda$|. Analysing this key quantity is informative regardless of any anthropic consideration, as it explores the predictions of different cosmic star formation models outside the canonical values of the underlying cosmological parameters. In this section, we will discuss similarities and discrepancies of our results on the cosmic star formation efficiency (Fig. 10) with respect to the results from direct simulations.
Fig. 13 shows the cosmic star formation efficiency predicted by our model, i.e. the same as in Fig. 10, in comparison with the efficiencies obtained with the Eagle simulations by Barnes et al. (2018) and the Enzo simulations by Oh et al. (2022). For |$\Lambda \gt \Lambda _{\rm obs}$|, there is general qualitative agreement, with all works considered showing a suppression of the cosmic star formation efficiency for high values of |$\Lambda$|. Our model predicts a peak at |$\Lambda \approx 0.1 \, \Lambda _{\rm obs}$|, and lower efficiencies as we approach an EdS universe. By contrast, the fitting functions to the numerical results appear to indicate a monotonic increase in the asymptotic cosmic star formation efficiency as |$\Lambda \rightarrow 0$|. While this trend is qualitatively different from our result, no simulations actually cover the range |$0\lt \Lambda \lt \Lambda _{\rm obs}$|. It would be interesting to run further simulations in this range in order to establish whether the simulated dependence of the cosmic star formation efficiency on |$\Lambda$| is genuinely monotonic.

Same as Fig. 10, for the model presented in this work (red line), and the cosmological simulations by Barnes et al. (2018) and Oh et al. (2022) (blue circles and purple squares, respectively). The corresponding fitting formulae to the numerical results, proposed by Oh et al. (2022), are represented with blue and purple solid lines, respectively. Integrating the CSFRD given by our model up to the same finite cosmic time as in the simulations (|$t_{\rm max}=20.7 \, \rm Gyr$| and |$t_{\rm max}=97 \, \rm Gyr$| for Barnes et al. 2018 and Oh et al. 2022, respectively) improves the agreement with the numerical results (dashed and solid salmon lines, respectively). The match with Barnes et al. (2018) is further improved by properly changing the slope of the gas density profiles within haloes, |$\eta$|, which is a free parameter in our model (dotted salmon line; see discussion in Section 5.4.2 for details).
Nevertheless, assuming that the fitting functions proposed by Oh et al. (2022) do represent the actual trend of the efficiency in the range |$0\lt \Lambda \lt \Lambda _{\rm obs}$|, it is worth speculating on the possible origins of the different result that we obtain. We should, however, first make sure that we make a fair comparison across all models considered. The simulations are inevitably stopped at some finite future time; for our cosmology and an EdS universe, this maximum time is |$t_{\rm max}=97\, \rm Gyr$| and |$t_{\rm max} = 20.7 \, \rm Gyr$| in Barnes et al. (2018) and Oh et al. (2022), respectively. To account for the star formation history at arbitrarily large times, Barnes et al. (2018) extrapolate the time evolution of the CSFRD beyond |$t_{\rm max}$| with a decaying exponential, and conclude that the contribution from times |$t\gt t_{\rm max}$| to the total cosmic star formation efficiency is negligible. Therefore, for a fairer comparison between our model and the Oh et al. (2022) fitting formula to the Barnes et al. (2018) simulations (blue line in Fig. 13), we should re-compute our cosmic star formation efficiencies by integrating the CSFRD up to |$t_{\rm max}=20.7\, \rm Gyr$|, and not in the limit |$t\rightarrow \infty$|.
The result is represented by the dashed salmon line in Fig. 13, which considerably improves the agreement with the fitting formula to the Barnes et al. (2018) results. Our model recovers the cosmic efficiency in the Barnes et al. (2018) simulations within 50 per cent for an EdS universe, and yields an excellent match at the high-|$\Lambda$| end. The largest discrepancy occurs for |$\Lambda \approx 0.3 \Lambda _{\rm obs}$|, and amounts to less than a factor of |$\sim$| 3. This deviation is smaller than the discrepancy between Barnes et al. (2018) and Oh et al. (2022). However, the non-monotonic dependence of the cosmic star formation efficiency on |$\Lambda$| for |$\Lambda \lt \Lambda _{\rm obs}$| does not disappear by simply imposing a cut-off in the time integral of the CSFRD. However, the peak becomes much less prominent, and the discrepancy with the numerical results is reduced.
If we apply the same reasoning for a fair comparison between our model and the Oh et al. (2022) simulations, we reach overall similar conclusions. To estimate the integral of the CSFRD in the limit |$t\rightarrow \infty$|, Oh et al. (2022) extrapolated their numerical results beyond |$t_{\rm max}$| by fitting their CSFRD with the same functional shape used by Madau & Dickinson (2014). This analytical formula allowed us to subtract the contribution from cosmic times |$t\gt t_{\rm max}$| from the cosmic star formation efficiencies computed by Oh et al. (2022), and we verified that the change would be in any case negligible. We can then simply cut the time integral of the CSFRD at |$t_{\rm max}=97\, \rm Gyr$|. The result is given by the salmon solid line in Fig. 13. The peak is less conspicuous, and the agreement with the fit to the numerical results (purple line) improves. However, the decay at large values of |$\Lambda$| is less steep compared to Oh et al. (2022) work, and the efficiency in the EdS universe is recovered within a factor of |$\sim$| 4.5.
Without adding any new physics to our model, we can check whether a change in the underlying parameters would yield a better match with cosmological simulations. This can be done by considering the late-time behaviour of the CSFRD. As mentioned earlier, Oh et al. (2022) showed that their CSFRDs can be well described with the same fitting formula adopted by Madau & Dickinson (2014). We verified that this is the case for the Barnes et al. (2018) CSFRDs too. At late times, the Madau & Dickinson (2014) curve can be in turn approximated with a power-law scaling of the CSFRD with |$(1+z)$|. Our model also predicts a power-law scaling in an EdS universe at late times, with an index that depends on the slope of the gas density profile within haloes, |$\eta$| (see Appendix A1). We can therefore find the value of |$\eta$| that yields the same late-time power-law scaling in an EdS universe as in Barnes et al. (2018) and Oh et al. (2022).
In the first case, we find that the desired value is |$\eta =2.75$|. We then run our model for such value of |$\eta$|, while leaving all other parameters the same. We integrate the CSFRD for any value of |$\Lambda$| up to |$t_{\rm max} = 20.7 \, \rm Gyr$|, and obtain the cosmic star formation efficiency. Moving from |$\eta =1.9$| to |$\eta =2.75$| further reduces the prominence of the peak in the cosmic star formation efficiency. However, the overall agreement with Barnes et al. (2018) does not significantly improve. Also, the value |$\eta =2.75$|, while mathematically acceptable, would fall outside the range that provides a consistency within |$3\sigma$| with observational constraints from the Kennicutt–Schmidt relationship and the baryonic Tully–Fisher relationship (see SP21 for details). If we apply the same logic to the Oh et al. (2022) simulations, we would find |$\eta =3.21$|. This is not even mathematically acceptable in our formalism, as values of |$\eta \gt 3$| would cause the gas mass of haloes to diverge (see equation 7). Therefore, we cannot fine tune |$\eta$| to match the late-time behaviour of the CSFRD found by Oh et al. (2022) for an EdS universe.
To summarize, the qualitative discrepancy between the cosmic star formation efficiency predicted by our model and by analytical fits to cosmological simulations is considerably mitigated once a cut-off time is taken into consideration, for a fair comparison with the numerical results. Where possible, the significance of the peak in the cosmic star formation efficiency predicted by our model can be further reduced by regulating the value of |$\eta$| to match the late-time behaviour of the simulated CSFRD in an EdS universe. However, the peak does not disappear completely, therefore some change in the physics included in our model may be necessary in order to obtain a cosmic star formation efficiency that monotonically decreases with |$\Lambda$|.
The residual qualitative discrepancies between our model and numerical simulation occur in the range |$0\lt \Lambda \lt \Lambda _{\rm obs}$|. Within our formalism, the asymptotic stellar mass produced in a halo of a fixed virial mass |$M\lesssim 10^{11}\, \rm {\rm M}_{\odot }$| does not vary appreciably with the cosmological constant for |$\Lambda \lt \Lambda _{\rm obs}$|. However, for the same cosmologies, the stellar mass produced in high-mass haloes (|$M\gtrsim 10^{11}\, \rm {\rm M}_{\odot }$|) can vary by up to one order of magnitude (see Fig. 8). Therefore, the treatment of star formation within high-mass haloes is probably the main reason behind the discrepancies between our model and simulations. In this respect, an obvious missing piece of physics in our model is AGN feedback, as discussed in Section 5.3. An explicit AGN feedback mechanism could reduce the contribution of high-mass haloes to the star formation history, and may bring the cosmic efficiency at |$0\lt \Lambda \lt \Lambda _{\rm obs}$| closer to that of an EdS universe.
It is worth noting that Oh et al. (2022) predict a decreasing efficiency even in the absence of AGN feedback. There may therefore be other ways to achieve the same qualitative trend without (or in addition to) AGN feedback. For example, one may include chemical evolution models, or an explicit shutdown of star formation at very late times, as detailed in Section 5.3. We intend to address these issues in future work. At the same time, we point out that it would be very informative to run cosmological simulations in the range |$0\lt \Lambda \lt \Lambda _{\rm obs}$|, which is the range where our predictions differ from the fits to the numerical results by Barnes et al. (2018) and Oh et al. (2022).
6 CONCLUSIONS AND PERSPECTIVES
We have improved the SP21 analytical model for cosmic star formation (Sorini & Peacock 2021), which was itself based on the work by Hernquist & Springel (2003). We have then used this model to predict the cosmic star formation history in flat |$\Lambda$|CDM cosmological models where the cosmological constant is varied between |$\Lambda =0$| and |$\Lambda =10^5 \, \Lambda _{\rm obs}$|. This ensemble of universes was designed so that all members are indistinguishable at high redshift, but later part company as their different cosmological constants come to dominate.
The CSFRD is obtained by integrating the SFR within haloes over the halo mass function. We assumed that the same astrophysical processes operate in each member of the ensemble as in our own universe. At high redshift, the SFR is determined by the typical time taken to convert cool gas into stars, whereas at low redshift it is set by the cooling time. We were then able to compute the total stellar mass generated by a given time, as a fraction of the total baryonic mass, per comoving volume. We showed that this cumulative star formation efficiency tends to a well-defined limit as |$t\rightarrow \infty$|, and we determined the dependence of this asymptotic efficiency on the value of |$\Lambda$|. We then investigated the implications of our results for anthropic attempts to explain the observed value of |$\Lambda$|. Our main conclusions are as follows:
The SP21 formalism predicts a converging total comoving stellar mass density as the cosmic time tends to infinity. In our universe, about 32 per cent of the eventual stellar mass has been formed by |$z=0$| (Fig. 9). Thus, we are typical observers in terms of the cosmological epoch in which we live, which is marginally |$\Lambda$|-dominated. However, we find that this coincidence between typical eras of star formation and |$\Lambda$|-domination does not hold for other values of |$\Lambda$|; hence single-universe anthropic arguments do not account for the ‘why-now’ problem.
Our model predicts that the cosmic star formation efficiency peaks at around 27 per cent for |$\Lambda \approx 0.1 \Lambda _{\rm obs}$| (Fig. 10). In an Einstein–de Sitter universe with |$\Lambda =0$|, the efficiency is approximately a factor of 3 lower than this; but the efficiency is much more strongly suppressed at high values of |$\Lambda$|. The decline is, however, relatively slow, and the efficiency only becomes negligibly small as |$\Lambda$| approaches |$10^5 \, \Lambda _{\rm obs}$|.
For a given |$\Lambda$|, there is a characteristic halo mass above which the stellar efficiency within haloes is maximized (Fig. 8). A halo with a given virial mass is smaller in universes with larger |$\Lambda$|, so that its typical density is higher. This in turn makes gas cooling, and subsequently star formation, more efficient. Thus, the suppression of the cosmic star formation efficiency at high |$\Lambda$| must be cosmological in nature, and not astrophysical: it is induced by the cut-off in the halo mass function occurring at lower halo masses. Conversely, the reason for the reduced cosmic star formation efficiency at low |$\Lambda$| is astrophysical: haloes are more diffuse as the cosmological model approximates an Einstein–de Sitter universe, so that the cooling rate is reduced. Thus, according to our model, astrophysics favours star formation in those haloes that exist when |$\Lambda$| is high, whereas cosmology increases the abundance of haloes when |$\Lambda$| is low. The observed value of the cosmological constant finds itself near the sweet spot of these two opposite trends (Fig. 10). However, that does not imply that the generation of observers is maximized around the same value of |$\Lambda$|.
Assuming a flat prior on |$\Lambda$| and that the observer weighting of a universe with a given |$\Lambda$| is proportional to the cosmic star formation efficiency in that universe, the peak of the posterior on |$\Lambda$| occurs at |$\Lambda \approx 800 \, \Lambda _{\rm obs}$| (Fig. 11). The median of the distribution is |$\Lambda =539$|, and the |$16^{\rm th}-84^{\rm th}$| percentile range corresponds to |$\Lambda = 73-2306$| (Table 1). The probability that an observer measures a value of the cosmological constant |$\Lambda \le \Lambda _{\rm obs}$| is |$\sim$| 0.5 per cent. Thus, Weinberg’s anthropic reasoning would be disfavoured by the SP21 model, under the given assumptions.
It is striking that the observed value of |$\Lambda$| is within one order of magnitude of the optimal value that maximizes the cosmic star formation efficiency. We also succeed in providing an explanation for why this is the case, based on calculations from first principles. However, we also found that this asymptotic efficiency of star formation declines rather slowly as |$\Lambda$| increases, falling below 1/10 of the peak value only for |$\Lambda \gt 500\, \Lambda _{\rm obs}$| (Fig. 10). Hence, with Weinberg’s assumption of a uniform prior on |$\Lambda$|, it is inevitable that typical observers would expect to experience values of |$\Lambda$| much greater than the one we observe, and this is why we find that values of |$\Lambda$| as small as observed are highly unlikely. In principle this rules out the simple anthropic ensemble considered here, subject to the assumptions of our model. But we have reviewed alternative calculations and found that some give a less severe disagreement – although there is unanimity that the observed |$\Lambda$| is unusually small at some level.
Future progress on this question requires more detailed and realistic modelling: we will endeavour to generalize our model of cosmic star formation by including features such as chemical evolution of the interstellar/intergalactic medium and AGN feedback. In particular, the absence of an explicit AGN feedback mechanism may lead to an overproduction of stars in massive haloes, which, according to our model, are especially efficient in universes with a high value of |$\Lambda$|. Whether or not this is important depends on the applicability of exact co-evolution between black holes and stellar content in galaxies. If this concept holds exactly, the modelled effective energy feedback from star formation may be able to capture the effect of AGN feedback. However, if co-evolution does not hold exactly, then the inclusion of separate AGN feedback in our formalism might alter the resulting conclusions on the anthropically favoured range of |$\Lambda$|. We plan to address this question in future work.
One would also like to see further direct calculations via simulations. However, hydrodynamical simulations are expensive to run, especially beyond |$z=0$|, thus it is important to invest the computational resources to study models that promise to bring the most valuable physical insights. Our analytical model enables a fast exploration of the dependence of star formation history on the cosmological parameters, so it can aid in designing a suite of simulations exploring cosmological models that differ from the fiducial |$\Lambda$|CDM cosmology in especially interesting ways.
Eventually, we can hope for a convergence of these studies and a consensus on the viability of Weinberg’s anthropic explanation for |$\Lambda$|. However, even if the verdict is negative, this would be far from a rejection of anthropic reasoning in general. This manuscript focuses on changing only one fundamental parameter, i.e. the cosmological constant – and even there we have restricted ourselves to the case of |$\Lambda \ge 0$|, leaving the whole topic of recollapsing models for a further investigation. Beyond this, one should in principle vary multiple parameters simultaneously, and explore their joint impact on star formation history and the generation of observers. We are already working on this kind of study (Lombriser et al., in preparation).
ACKNOWLEDGEMENTS
We thank the anonymous reviewers for helpful comments. The authors thank Sownak Bose and Boon Kiat Oh for useful discussions, and Ed Wheeler for comments on a draft of this manuscript. DS and LL acknowledge support from the Swiss National Science Foundation Professorship grant no. 202671. DS and JAP were supported by the European Research Council, under grant no. 670193 and by the Science and Technology Facilities Council consolidated grant no. RA5496. This work made extensive use of the NASA Astrophysics Data System and of the astro-ph preprint archive at arXiv.org. For the purpose of open access, the author has applied a Creative Commons Attribution (CC BY) licence to any Author Accepted Manuscript version arising from this submission.
DS dedicates this work to the memory of his grandmother Linda.
DATA AVAILABILITY
The code that we developed to produce the results presented in this manuscript is undergoing further extension, and we plan to release it publicly as part of follow-up work. Until then, the code and the output data published in this manuscript may be released upon reasonable request to the corresponding author.
REFERENCES
APPENDIX A: DERIVATION OF THE ASYMPTOTIC BEHAVIOUR OF STAR FORMATION HISTORY
In Section 4.1, we showed that the stellar mass fraction in haloes with fixed mass M asymptotes to a constant in the far future of the universe. This means that the corresponding nSFR, and consequently the CSFRD, tend to zero in the limit |$t\rightarrow \infty$|. These conclusions result from the numerical implementation of the iterative procedure that we introduced to extend the SP21 analytic model (see Section 2.2). In this Appendix, we will show that, under well-motivated approximations, it is possible to derive the asymptotic behaviour of the nSFR and CSFRD analytically. Crucially, this will enable us to prove that the SMD, which is the time integral of the CSFRD, converges in the limit |$t\rightarrow \infty$|.
Since the CSFRD depends on the HMF and the nSFR at fixed M (see equation 1), we will first need to study the limit of these two quantities for arbitrarily large cosmic times. Once we obtain the asymptotic behaviour of the CSFRD, it will be straightforward to study the convergence of the SMD. We will need to distinguish between the two qualitatively different cases of |$\Lambda =0$| and |$\Lambda \gt 0$|. All computations in this section have been verified with the software mathematica (Wolfram Research, Inc. 2023).
A1 Einstein–de Sitter universe
A1.1 Normalized SFR
The evolution of the nSFR of a halo with fixed virial temperature T in the low-redshift regime is given by equation (13). In an EdS universe, the effective cooling time coincides with the dynamical time at any cosmic epoch (see equation (10)). Therefore, equation (13) can be simplified to
switching from fixed virial temperature to fixed virial mass as the independent variable. Thus, the first term in the r.h.s of the equation above acquires a redshift dependence, which can be determined by replacing the virial temperature in SP21 equation (25) with the expression given by equation (5). Recalling the definition of |$\tilde{S}(T)$| from SP21,
it becomes obvious that this task is complicated by the presence of the cooling function |$\mathcal {C}(T)$|, which has a non-trivial dependence on temperature. However, we can reasonably approximate the Sutherland & Dopita (1993) cooling function for a plasma of primordial composition, in the temperature range of our interest (|$T\gt 10^{4.5}\, \rm K$|), with a piece-wise power law:
Within this approximation, and considering that in an EdS universe |$H(z)\propto (1+z)^{3/2}$|, it follows that
Let us now turn to the second factor in the r.h.s. of equation (A1). As discussed in Section 3.2.2, the critical temperature decreases with cosmic time in an EdS universe. Recalling that for the computation of the CSFRD we only consider haloes above a certain virial temperature threshold (in this work, |$10^{4.5} \, \rm K$|), it is apparent that all haloes that contribute to the CSFRD will exceed the critical temperature at a sufficiently late time. Thus, |$f_{\rm gas}$| will asymptotically reach a non-null constant, so that the only remaining redshift dependence in the r.h.s. of equation (A1) is |$H(z)^{3/\eta } \propto (1+z)^{9/2 \eta }$|. Combining it with equation (A4), it follows that |$s_{\rm low}(M,\, z) \propto (1+z)^{\zeta }$|, with
The nSFR then converges as |$z \rightarrow -1$| (or, equivalently, |$t \rightarrow \infty$|), as long as |$\zeta \gt -1$|. In the physically relevant range |$3/2\lt \eta \lt 3$|, this condition is satisfied for all values of |$\theta$| quoted in equation (A3).
A1.2 Halo mass function
Within the Sheth–Tormen formalism (Sheth & Tormen 1999, 2002), the halo multiplicity function is
where |$\sigma _0(M)$| is the linear-theory fractional variance of matter density fluctuations averaged over spheres containing a mass M at |$z=0$|, and N is a normalization constant chosen such that all cosmic mass is contained in haloes. We also defined
where |$D(z)$| is the growth function, which, in an EdS universe, scales as |$D(z) \propto (1+z)^{-1}$|. Therefore, the scaling of y with redshift is |$y\propto (1+z)$|.
It follows that the redshift dependence of the halo multiplicity function is given by
In the limit |$z \rightarrow -1$|, the scaling above reduces to |$\mathrm{d}F/\mathrm{d}\ln M \propto (1+z)$|.
A1.3 Stellar mass density
The convergence of the CSFRD depends on the asymptotic behaviour of the integrand in equation (1). As explained in Section 3.2.2, the critical temperature becomes arbitrarily low at large times, therefore the relevant regime in our piece-wise approximation of the cooling function is |$T_{\rm crit} \lt 10^{4.5} \, \rm K$|, which corresponds to |$\theta =10.75$|. If we further set |$\eta =1.9$|, we obtain that |$\zeta \approx 6.1$|. We saw in Section A1.2 that the asymptotic redshift dependence of the halo multiplicity function in equation (1) is |$\propto (1+z)$|. Therefore, the integrand scales approximately as |$(1+z)^{\zeta +1} \approx (1+z)^{7.1}$|. This is then the redshift-dependence of the CSFRD in the far future.
The SMD is simply the time integral of the CSFRD. Changing the integration variable from cosmic time to redshift, the SMD is defined by equation (24). For an EdS model, the integrand in that equation is approximately proportional to |$\dot{\rho }_*(z) (1+z)^{-5/2} \propto (1+z)^{4.6}$|, which is convergent for |$z\rightarrow -1$|. This proves that the mass density of stars formed over the entire history of an EdS universe is finite. We therefore used the scaling |$\rho _*(z) \propto (1+z)^{\zeta -1/2} \approx (1+z)^{5.6}$| to extrapolate the SMD in the far future and compute the asymptotic SMD. We verified that the analytic approximation matches the numerical results at |$z\lesssim -0.99$|.
A2 Positive cosmological constant
To study the convergence of the SMD in a universe with |$\Lambda \gt 0$|, we will still follow a similar logic to the one outlined in the previous section. We will highlight the differences in the reasoning that stem from having a positive rather than null cosmological constant.
A2.1 Normalized SFR
In the case |$\Lambda \gt 0$|, it is more convenient to consider the asymptotic behaviour of the nSFR in terms of cosmic time rather than redshift. This can be straightforwardly obtained by adapting equation (15) in SP21 such that haloes are taken at a fixed virial mass instead of a fixed virial temperature:
In a universe with |$\Lambda \gt 0$| the Hubble parameter asymptotes to a positive constant, so we can ignore any |$H(z)$| factor in the equation above for our convergence study. From equation (4), it also follows that haloes with a fixed mass will have an approximately constant virial temperature in the far future, so we can disregard any redshift dependence originating from the relationship between virial temperature and virial mass. Finally, we notice that our prescription for the cooling time (equation 10) yields |$t_{\rm cool} \sim t$| in the far future. Therefore, we can immediately remove the last factor on the r.h.s. of equation (A9). Its asymptotic behaviour in the limit |$t \rightarrow \infty$| is therefore
We then need to study the scaling of |$f_{\rm gas}$| with cosmic time in this limit.
Our numerical results in Fig. 6 suggest that the critical temperature increases with time once |$\Lambda$| dominates. Therefore, there will eventually be a time when the virial temperature falls below the critical temperature. At this point, the scaling for |$f_{\rm gas}$| is set by the critical temperature – see the case |$T\lt T_{\rm crit}$| in equation (15). We then need to determine the asymptotic time dependence of the critical temperature.
As mentioned in Section 3.2.1, the critical temperature is the virial temperature above which haloes retain a baryon mass fraction equal to the cosmic value |$f_{\rm b}$|. At low redshift, this temperature is the one that satisfies the following condition (see SP21):
where we make the implicit assumption that the stellar mass fraction asymptotes to a constant in the far future. This assumption will be validated a posteriori, when we will show that the formalism does predict a converging stellar mass fraction.
Taking |$t_{\rm cool} \approx t$| and approximating again the cooling function as a piece-wise power law (see equation A3), we obtain:
It can be verified that, if |$1.5\lt \eta \lt 3$|, the exponent in the equation above is positive for |$\theta \lesssim 1.4$|. Our numerical results show that for |$\eta =1.9$|, the critical temperature is above |$10^5 \, \rm K$| at late times in universes with |$\Lambda \gt 0$| (see Fig. 6). In this temperature range, equation (A3) returns |$\theta = -1$| or |$\theta =0.5$|. In either case, this satisfies the conditions for which the exponent in equation (A12) is positive. Therefore, the critical temperature keeps increasing with time.
One might wonder whether this result is specific to the value of |$\eta =1.9$| adopted here. However, SP21 showed that, for |$\Lambda = \Lambda _{\rm obs}$|, a different value of |$\eta$| would change the late-time critical temperature only within a factor of |$\sim$| |$2-3$|. Thus, our considerations on the asymptotic behaviour of |$T_{\rm crit}$| in a universe with |$\Lambda \gt 0$| should not be altered if we changed the value of |$\eta$|.
Inserting the time scaling of the critical temperature above in equation (23), it can be easily seen that the cooling radius also scales with t as a power law in the far future, with a positive index smaller than unity. This proves that the cooling radius increases with cosmic time, but does so progressively more slowly, as its time derivative flattens to zero as |$t\rightarrow \infty$| (see the discussion in Section 3.2.1).
If we now insert the scaling given by equation (A12) in the |$T\lt T_{\rm crit}$| case of equation (15), and use the resulting time dependence for |$f_{\rm gas}$| in equation (A10), we deduce that |$s_{\rm low}(T(M, \, z), z) \propto t^{\xi }$|, with
We find that, for the value |$\eta =1.9$| that we are adopting, |$\xi \lt -1$| for |$-0.3 \lesssim \theta \lesssim 1.5$|. This includes |$\theta =0.5$|, which corresponds to |$T_{\rm crit}\gt 10^6 \, \rm K$|, a condition that is certainly satisfied asymptotically. Hence, our analysis would predict that the nSFR converges in the limit |$t \rightarrow \infty$|.
Having thoroughly tested our results by a comparison of numerics with the approximate analytical approach presented in this Appendix, we are therefore confident of the robustness of our findings regarding the convergence of the total stellar mass within haloes.
A2.2 Stellar mass density
In a universe with |$\Lambda \gt 0$|, the growth function asymptotes to a constant, so that the halo multiplicity function follows the same asymptotic behaviour for a fixed halo mass (see also Section A1.2). The time dependence of the CSFRD is then simply obtained by integrating the nSFR over all halo masses, weighted by the asymptotic multiplicity function. However, there is a complication: not all haloes lie below the critical temperature threshold at any fixed cosmic time. Some will therefore follow the time-dependence for the nSFR found in the previous section, which applies to the case |$T\lt T_{\rm crit}$|, while others will yet have to enter that regime. Therefore, the asymptotic time dependence of the nSFR cannot be simply extracted from the integral in equation (1), and the late-time behaviour of the CSFRD will deviate from that.
We empirically verified that the a good analytical approximation of the late-time CSFRD is given by
where |$t_{\rm l}$| is a sufficiently late cosmic time and |$\phi$| is the exponent that minimizes the residuals with the CSFRD that we computed numerically. In the above, we defined |$\rm U = 1 \rm {\rm M}_{\odot } \, yr^{-1} \, cMpc^{-3}$|, so that the argument of the logarithms is dimensionless. We verified that our fitting formula yields an accuracy of a few per cent at |$1+z \lesssim 10^{-5}$|, for all cosmological models considered. The values of the fitting parameters depended very weakly on the choice of |$t_{\rm l}$|.
The fitting function allows us to calculate the asymptotic SMD analytically, as we shall now show. The SMD at any time |$t\gt t_{\rm l}$| is given by:
where we defined |$P = \ln (\dot{\rho }_*(t_{\rm l})/ \rm U)$|. With the substitution |$t^{\prime } \rightarrow \tilde{t} = -P (t^{\prime }/t_{\rm l})^{\phi }$|, we obtain
where |$\gamma _{\rm l}$| is the lower incomplete gamma function. We note that because the CSFRD decays with time, one can always choose a large enough |$t_{\rm l}$| such that |$P\lt 0$|. For |$t\gt t_{\rm l}$|, the logarithm of the CSFRD is negative and increases in absolute value, so that |$\ln \dot{\rho }_*(t)/\ln \dot{\rho }_*(t_{\rm l})$| increases with time. Thus, |$\phi \gt 0$| and the lower incomplete gamma function in equation (A16) is always well defined. In the limit |$t\rightarrow \infty$|, the above expression converges to
where |$\gamma _{\rm u}$| is the upper incomplete gamma function. This is the expression that we used to compute the asymptotic SMD in universes with |$\Lambda \gt 0$|.