Cosmology¶
This module is an implementation of the standard FLRW cosmology with a number of dark energy models including \(\Lambda CDM\), wCDM, and varying dark energy equations of state. The cosmology object models the contributions from dark matter, baryons, curvature, photons, neutrinos, and dark energy.
Basics¶
In Colossus, the cosmology is set globally, and all functions respect that global cosmology. Colossus does not set a default cosmology, meaning that the user must set a cosmology before using any cosmological functions or any other functions that rely on the Cosmology module. This documentation contains coding examples of the most common operations. Much more extensive code samples can be found in the Tutorials.
Setting and getting cosmologies¶
First, we import the cosmology module:
from colossus.cosmology import cosmology
Setting a cosmology is almost always achieved with the setCosmology()
function, which can be
used in multiple ways:
Set one of the predefined cosmologies:
cosmology.setCosmology('planck18')
Set one of the predefined cosmologies, but overwrite certain parameters:
cosmology.setCosmology('planck18', {'print_warnings': False})
Add a new cosmology to the global list of available cosmologies. This has the advantage that the new cosmology can be set from anywhere in the code. Only the main cosmological parameters are mandatory, all other parameters can be left to their default values:
params = {'flat': True, 'H0': 67.2, 'Om0': 0.31, 'Ob0': 0.049, 'sigma8': 0.81, 'ns': 0.95} cosmology.addCosmology('myCosmo', params) cosmo = cosmology.setCosmology('myCosmo')
Set a new cosmology without adding it to the global list of available cosmologies:
params = {'flat': True, 'H0': 67.2, 'Om0': 0.31, 'Ob0': 0.049, 'sigma8': 0.81, 'ns': 0.95} cosmo = cosmology.setCosmology('myCosmo', params)
Set a selfsimilar cosmology with a powerlaw power spectrum of a certain slope, and the default settings set in the
powerlaw
cosmology:cosmo = cosmology.setCosmology('powerlaw_2.60')
Whichever way a cosmology is set, the current cosmology is stored in a global variable and can be obtained at any time:
cosmo = cosmology.getCurrent()
For more extensive examples, please see the Tutorials.
Changing and switching cosmologies¶
The current cosmology can also be set to an already existing cosmology object, for example when switching between cosmologies:
cosmo1 = cosmology.setCosmology('WMAP9')
cosmo2 = cosmology.setCosmology('planck18')
cosmology.setCurrent(cosmo1)
The user can change the cosmological parameters of an existing cosmology object at runtime, but
MUST call the update function checkForChangedCosmology()
directly after the changes. This function ensures that the parameters are consistent (e.g.,
flatness) and that no outdated cached quantities are used:
cosmo = cosmology.setCosmology('WMAP9')
cosmo.Om0 = 0.31
cosmo.checkForChangedCosmology()
Only userdefined cosmological parameters, that is, parameters that can be passed to the constructor
of the cosmology object, can be changed in this way. Changing other internal variables of the class
can have unintended consequences! If derived parameters (such as Ok0 or Onu0) are changed, those
changes will simply be overwritten when
checkForChangedCosmology()
is called.
Converting to and from Astropy cosmologies¶
Colossus can easily interface with the cosmology module of the popular
Astropy code. Astropy cosmology objects can be converted to Colossus
cosmologies with the fromAstropy()
function:
import astropy.cosmology
params = dict(H0 = 70, Om0 = 0.27, Ob0 = 0.0457, Tcmb0 = 2.7255, Neff = 3.046)
sigma8 = 0.82
ns = 0.96
astropy_cosmo = astropy.cosmology.FlatLambdaCDM(**params)
colossus_cosmo = cosmology.fromAstropy(astropy_cosmo, sigma8, ns, name = 'my_cosmo')
The name
parameter is not necessary if a name is set in the Astropy cosmology. The sigma8
and ns
parameters must be set by the user because the Astropy cosmology does not contain them
(because it does not compute power spectrumrelated quantities). The conversion supports the
LambdaCDM
, FlatLambdaCDM
, wCDM
, FlatwCDM
, w0waCDM
, and Flatw0waCDM
Astropy cosmology classes. Conversely, to convert a Colossus cosmology to Astropy, we simply use
the toAstropy()
function:
colossus_cosmo = cosmology.setCosmology('WMAP9')
astropy_cosmo = colossus_cosmo.toAstropy()
Naturally, both conversion functions will fail if astropy.cosmology cannot be imported.
Summary of getter and setter functions¶

Set a cosmology. 

Add a set of cosmological parameters to the global list. 

Set the current global cosmology to a cosmology object. 
Get the current global cosmology. 


Convert an Astropy cosmology object to Colossus and set it as current cosmology. 
Standard cosmologies¶
The following sets of cosmological parameters can be chosen using the
setCosmology()
function:
ID 
Paper 
Location 
Explanation 

planck18only 
Table 2 
Bestfit, Planck only (column 5) 

planck18 
Table 2 
Bestfit with BAO (column 6) 

planck15only 
Table 4 
Bestfit, Planck only (column 2) 

planck15 
Table 4 
Bestfit with ext (column 6) 

planck13only 
Table 2 
Bestfit, Planck only 

planck13 
Table 5 
Bestfit with BAO etc. 

WMAP9only 
Table 2 
Max. likelihood, WMAP only 

WMAP9ML 
Table 2 
Max. likelihood, with eCMB, BAO and H0 

WMAP9 
Table 4 
Bestfit, with eCMB, BAO and H0 

WMAP7only 
Table 1 
Max. likelihood, WMAP only 

WMAP7ML 
Table 1 
Max. likelihood, with BAO and H0 

WMAP7 
Table 1 
Bestfit, with BAO and H0 

WMAP5only 
Table 1 
Max. likelihood, WMAP only 

WMAP5ML 
Table 1 
Max. likelihood, with BAO and SN 

WMAP5 
Table 1 
Bestfit, with BAO and SN 

WMAP3ML 
Table 2 
Max.likelihood, WMAP only 

WMAP3 
Table 5 
Best fit, WMAP only 

WMAP1ML 
Table 1/4 
Max.likelihood, WMAP only 

WMAP1 
Table 7/4 
Best fit, WMAP only 

illustris 
– 
Cosmology of the Illustris simulation 

bolshoi 
– 
Cosmology of the Bolshoi simulation 

multidarkplanck 
Table 1 
Cosmology of the MultidarkPlanck simulations 

millennium 
– 
Cosmology of the Millennium simulation 

EdS 
– 
– 
Einsteinde Sitter cosmology 
powerlaw 
– 
– 
Default settings for powerlaw cosms. 
Those cosmologies that refer to particular simulations (such as bolshoi
and millennium
) are
generally set to ignore relativistic species, i.e. photons and neutrinos, because they are not
modeled in the simulations. The EdS
cosmology refers to an Einsteinde Sitter model, i.e. a flat
cosmology with only dark matter and \(\Omega_{\rm m} = 1\).
Dark energy and curvature¶
All the default parameter sets above represent flat \(\Lambda CDM\) cosmologies, i.e. model dark energy as a cosmological constant and contain no curvature. To add curvature, the default for flatness must be overwritten, and the dark energy content of the universe must be set (which is otherwise computed from the matter and relativistic contributions):
params = cosmology.cosmologies['planck18']
params['flat'] = False
params['Ode0'] = 0.75
cosmo = cosmology.setCosmology('planck_curvature', params)
Multiple models for the dark energy equation of state parameter \(w(z)\) are implemented,
namely a cosmological constant (\(w=1\)), a constant \(w\), a linearly varying
\(w(z) = w_0 + w_a (1  a)\), and arbitrary usersupplied functions for \(w(z)\). To set,
for example, a linearly varying EOS, we change the de_model
parameter:
params = cosmology.cosmologies['planck18']
params['de_model'] = 'w0wa'
params['w0'] = 0.8
params['wa'] = 0.1
cosmo = cosmology.setCosmology('planck_w0wa', params)
We can implement more exotic models by supplying an arbitrary function:
def wz_func(z):
return 1.0 + 0.1 * z
params = cosmology.cosmologies['planck18']
params['de_model'] = 'user'
params['wz_function'] = wz_func
cosmo = cosmology.setCosmology('planck_wz', params)
Please note that the redshift range into the future is reduced from \(z = 0.995\)
(\(a = 200\)) to \(z = 0.9\) (\(a = 10\)) for the w0wa
and user
dark energy
models to avoid numerical issues.
Power spectrum models¶
By default, Colossus relies on fitting functions for the matter power spectrum which, in turn,
is the basis for the variance and correlation function. These models are implemented in the
power_spectrum
module, documented at the bottom of this file.
Derivatives and inverses¶
Almost all cosmology functions that are interpolated (e.g.,
age()
,
luminosityDistance()
or
sigma()
) can be evaluated as an nth derivative. Please note
that some functions are interpolated in log space, resulting in a logarithmic derivative, while
others are interpolated and differentiated in linear space. Please see the function documentations
below for details.
The derivative functions were not systematically tested for accuracy. Their accuracy will depend on how well the function in question is represented by the interpolating spline approximation. In general, the accuracy of the derivatives will be worse that the error quoted on the function itself, and get worse with the order of the derivative.
Furthermore, the inverse of interpolated functions can be evaluated by passing inverse = True
.
In this case, for a function y(x), x(y) is returned instead. Those functions raise an Exception if
the requested value lies outside the range of the interpolating spline.
The inverse and derivative flags can be combined to give the derivative of the inverse, i.e. dx/dy. Once again, please check the function documentation whether that derivative is in linear or logarithmic units.
Performance optimization and accuracy¶
This module is optimized for fast performance, particularly in computationally intensive
functions such as the correlation function. Almost all quantities are, by
default, tabulated, stored in files, and reloaded when the same cosmology is set again (see the
storage
module for details). For some rare applications (for example, MCMC chains
where functions are evaluated few times, but for a large number of cosmologies), the user can turn
this behavior off:
cosmo = cosmology.setCosmology('planck18', {'interpolation': False, 'persistence': ''})
For more details, please see the documentation of the interpolation
and persistence
parameters. In order to turn off the interpolation temporarily, the user can simply switch the
interpolation
parameter off:
cosmo.interpolation = False
Pk = cosmo.matterPowerSpectrum(k)
cosmo.interpolation = True
In this example, the power spectrum is evaluated directly without interpolation. The interpolation is fairly accurate (see specific notes in the function documentation), meaning that it is very rarely necessary to use the exact routines.
Module reference¶

class
cosmology.cosmology.
Cosmology
(name=None, flat=True, Om0=None, Ode0=None, Ob0=None, H0=None, sigma8=None, ns=None, de_model='lambda', w0=None, wa=None, wz_function=None, relspecies=True, Tcmb0=2.7255, Neff=3.046, power_law=False, power_law_n=0.0, interpolation=True, persistence='rw', print_info=False, print_warnings=True)¶ A cosmology is set via the parameters passed to the constructor. Any parameter whose default value is
None
must be set by the user. The easiest way to set these parameters is to use thesetCosmology()
function with one of the predefined sets of cosmological parameters listed above.In addition, the user can choose between different equations of state for dark energy, including an arbitrary \(w(z)\) function.
A few cosmological parameters are free in principle, but are well constrained and have a subdominant impact on the computations. For such parameters, default values are preset so that the user does not have to choose them manually. This includes the CMB temperature today (
Tcmb0
= 2.7255 K, Fixsen 2009) and the effective number of neutrino species (Neff
= 3.046, Planck Collaboration 2018). These values are compatible with the most recent observational measurements and can be changed by the user if necessary. Parameters
 name: str
A name for the cosmology, e.g.
WMAP9
or a userdefined name. If a userdefined set of cosmological parameters is used, it is advisable to use a name that does not represent any of the preset cosmologies. flat: bool
If flat, there is no curvature, \(\Omega_{\rm k} = 0\), and the dark energy content of the universe is computed as \(\Omega_{\rm de} = 1  \Omega_{\rm m}  \Omega_{\gamma}  \Omega_{\nu}\) where \(\Omega_{\rm m}\) is the density of matter (dark matter and baryons) in units of the critical density, \(\Omega_{\gamma}\) is the density of photons, and \(\Omega_{\nu}\) the density of neutrinos. If
flat == False
, theOde0
parameter must be passed. Om0: float
\(\Omega_{\rm m}\), the matter density in units of the critical density at z = 0 (includes all nonrelativistic matter, i.e., dark matter and baryons but not neutrinos).
 Ode0: float
\(\Omega_{\rm de}\), the dark energy density in units of the critical density at z = 0. This parameter is ignored if
flat == True
. Ob0: float
\(\Omega_{\rm b}\), the baryon density in units of the critical density at z = 0.
 H0: float
The Hubble constant in km/s/Mpc.
 sigma8: float
The normalization of the power spectrum, i.e. the variance when the field is filtered with a top hat filter of radius 8 Mpc/h. See the
sigma()
function for details on the variance. ns: float
The tilt of the primordial power spectrum.
 de_model: str
An identifier indicating which dark energy equation of state is to be used. The DE equation of state can either be a cosmological constant (
de_model = lambda
), a constant w (de_model = w0
, thew0
parameter must be set), a linear function of the scale factor according to the parameterization of Linder 2003 where \(w(z) = w_0 + w_a (1  a)\) (de_model = w0wa
, thew0
andwa
parameters must be set), or a function supplied by the user (de_model = user
). In the latter case, the w(z) function must be passed using thewz_function
parameter. w0: float
If
de_model == w0
, this variable gives the constant dark energy equation of state parameter w. Ifde_model == w0wa
, this variable gives the constant component w (seede_model
parameter). wa: float
If
de_model == w0wa
, this variable gives the varying component of w, otherwise it is ignored (seede_model
parameter). wz_function: function
If
de_model == user
, this field must give a function that represents the dark energy equation of state. This function must take z as its only input variable and return w(z). relspecies: bool
If
relspecies == False
, all relativistic contributions to the energy density of the universe (such as photons and neutrinos) are ignored. Ifrelspecies == True
, their energy densities are computed based on theTcmb0
andNeff
parameters. Tcmb0: float
The temperature of the CMB at z = 0 in Kelvin.
 Neff: float
The effective number of neutrino species.
 power_law: bool
Create a selfsimilar cosmology with a powerlaw matter power spectrum, \(P(k) = k^{\rm power\_law\_n}\).
 power_law_n: float
See
power_law
. interpolation: bool
By default, lookup tables are created for certain computationally intensive quantities, cutting down the computation times for future calculations. If
interpolation == False
, all interpolation is switched off. This can be useful when evaluating quantities for many different cosmologies (where computing the tables takes a prohibitively long time). However, many functions will be much slower if this setting isFalse
, and the derivatives and inverses will not work. Thus, please useinterpolation == False
only if absolutely necessary. persistence: str
By default, interpolation tables and other data are stored in a permanent file for each cosmology. This avoids recomputing the tables when the same cosmology is set again. However, if either read or write file access is to be avoided (for example in MCMC chains), the user can set this parameter to any combination of read (
'r'
) and write ('w'
), such as'rw'
(read and write, the default),'r'
(read only),'w'
(write only), or''
(no persistence). print_info: bool
Output information to the console.
 print_warnings: bool
Output warnings to the console.
Methods
Ez
(self, z)The Hubble parameter as a function of redshift, in units of \(H_0\).
Hz
(self, z)The Hubble parameter as a function of redshift.
Ob
(self, z)The baryon density of the universe, in units of the critical density.
Ode
(self, z)The dark energy density of the universe, in units of the critical density.
Ogamma
(self, z)The density of photons in the universe, in units of the critical density.
Ok
(self, z)The curvature density of the universe in units of the critical density.
Om
(self, z)The matter density of the universe, in units of the critical density.
Onu
(self, z)The density of neutrinos in the universe, in units of the critical density.
Or
(self, z)The density of relativistic species, in units of the critical density.
age
(self, z[, derivative, inverse])The age of the universe at redshift z.
angularDiameterDistance
(self, z[, derivative])The angular diameter distance to redshift z.
checkForChangedCosmology
(self)Check whether the cosmological parameters have been changed by the user.
comovingDistance
(self[, z_min, z_max, …])The transverse or lineofsight comoving distance.
correlationFunction
(self, R[, z, …])The linear mattermatter correlation function at radius R.
distanceModulus
(self, z)The distance modulus to redshift z in magnitudes.
filterFunction
(self, filt, k, R)The Fourier transform of certain filter functions.
getName
(self)Return the name of this cosmology.
growthFactor
(self, z[, derivative, inverse])The linear growth factor normalized to z = 0, \(D_+(z) / D_+(0)\).
growthFactorUnnormalized
(self, z)The linear growth factor, \(D_+(z)\).
hubbleTime
(self, z)The Hubble time, \(1/H(z)\).
lookbackTime
(self, z[, derivative, inverse])The lookback time since redshift z.
luminosityDistance
(self, z[, derivative, …])The luminosity distance to redshift z.
matterPowerSpectrum
(self, k[, z, model, …])The matter power spectrum at a scale k.
rho_b
(self, z)The baryon density of the universe at redshift z.
rho_c
(self, z)The critical density of the universe at redshift z.
rho_de
(self, z)The dark energy density of the universe at redshift z.
rho_gamma
(self, z)The photon density of the universe at redshift z.
rho_m
(self, z)The matter density of the universe at redshift z.
rho_nu
(self, z)The neutrino density of the universe at redshift z.
rho_r
(self, z)The density of relativistic species in the universe at redshift z.
sigma
(self, R, z[, j, filt, inverse, …])The rms variance of the linear density field on a scale R, \(\sigma(R)\).
soundHorizon
(self)The sound horizon at recombination.
toAstropy
(self)Create an equivalent Astropy cosmology object.
wz
(self, z)The dark energy equation of state parameter.

getName
(self)¶ Return the name of this cosmology.

checkForChangedCosmology
(self)¶ Check whether the cosmological parameters have been changed by the user. If there are changes, all precomputed quantities (e.g., interpolation tables) are discarded and recomputed if necessary.

toAstropy
(self)¶ Create an equivalent Astropy cosmology object.
This function throws an error if Astropy cannot be imported. Most standard Colossus cosmologies can be converted, exceptions are selfsimilar cosmologies with powerlaw spectra and userdefined dark energy equations of state.
 Returns
 cosmo_astropy: FLRW
An astropy.cosmology.FLRW class object.

Ez
(self, z)¶ The Hubble parameter as a function of redshift, in units of \(H_0\).
 Parameters
 z: array_like
Redshift; can be a number or a numpy array.
 Returns
 E: array_like
\(H(z) / H_0\); has the same dimensions as
z
.
See also
Hz
The Hubble parameter as a function of redshift.

Hz
(self, z)¶ The Hubble parameter as a function of redshift.
 Parameters
 z: array_like
Redshift; can be a number or a numpy array.
 Returns
 H: array_like
\(H(z)\) in units of km/s/Mpc; has the same dimensions as
z
.
See also
Ez
The Hubble parameter as a function of redshift, in units of \(H_0\).

wz
(self, z)¶ The dark energy equation of state parameter.
The EOS parameter is defined as \(w(z) = P(z) / \rho(z)\). Depending on its chosen functional form (see the
de_model
parameter toCosmology()
), w(z) can be 1, another constant, a linear function of a, or an arbitrary function chosen by the user. Parameters
 z: array_like
Redshift; can be a number or a numpy array.
 Returns
 w: array_like
\(w(z)\), has the same dimensions as
z
.

hubbleTime
(self, z)¶ The Hubble time, \(1/H(z)\).
 Parameters
 z: array_like
Redshift; can be a number or a numpy array.
 Returns
 tH: float
\(1/H\) in units of Gyr; has the same dimensions as
z
.
See also
lookbackTime
The lookback time since z.
age
The age of the universe at redshift z.

lookbackTime
(self, z, derivative=0, inverse=False)¶ The lookback time since redshift z.
The lookback time corresponds to the difference between the age of the universe at redshift z and today.
 Parameters
 z: array_like
Redshift, where \(0.995 < z < 200\); can be a number or a numpy array.
 derivative: int
If greater than 0, evaluate the nth derivative, \(d^nt/dz^n\).
 inverse: bool
If True, evaluate \(z(t)\) instead of \(t(z)\). In this case, the
z
field must contain the time(s) in Gyr.
 Returns
 t: array_like
The lookback time (or its derivative) since z in units of Gigayears; has the same dimensions as z.
See also
hubbleTime
The Hubble time, \(1/H_0\).
age
The age of the universe at redshift z.

age
(self, z, derivative=0, inverse=False)¶ The age of the universe at redshift z.
 Parameters
 z: array_like
Redshift, where \(0.995 < z < 200\); can be a number or a numpy array.
 derivative: int
If greater than 0, evaluate the nth derivative, \(d^nt/dz^n\).
 inverse: bool
If True, evaluate \(z(t)\) instead of \(t(z)\). In this case, the
z
field must contain the time(s) in Gyr.
 Returns
 t: array_like
The age of the universe (or its derivative) at redshift z in Gigayears; has the same dimensions as z.
See also
hubbleTime
The Hubble time, \(1/H_0\).
lookbackTime
The lookback time since z.

comovingDistance
(self, z_min=0.0, z_max=0.0, transverse=True)¶ The transverse or lineofsight comoving distance.
This function returns the comoving distance between two points. Depending on the chosen geometry, the output can have two different meanings. If
transverse = False
, the lineofsight distance is returned,\[d_{\rm com,los}(z) = \frac{c}{H_0} \int_{0}^{z} \frac{1}{E(z)} .\]However, if
transverse = False
, the function returns the comoving distance between two points separated by an angle of one radian atz_max
(ifz_min
is zero). This quantity depends on the spatial curvature of the universe,\[\begin{split}d_{\rm com,trans}(z) = \left\{ \begin{array}{ll} \frac{c/H_0}{\sqrt{\Omega_{\rm k,0}}} \sinh \left(\frac{\sqrt{\Omega_{\rm k,0}}}{c/H_0} d_{\rm com,los} \right) & \forall \, \Omega_{\rm k,0} > 0 \\ d_{\rm com,los} & \forall \, \Omega_{\rm k,0} = 0 \\ \frac{c/H_0}{\sqrt{\Omega_{\rm k,0}}} \sin \left(\frac{\sqrt{\Omega_{\rm k,0}}}{c/H_0} d_{\rm com,los} \right) & \forall \, \Omega_{\rm k,0} < 0 \\ \end{array} \right. \end{split}\]In Colossus, this distance is referred to as the “transverse comoving distance” (e.g., Hogg 1999), but a number of other terms are used in the literature, e.g., “comoving angular diameter distance” (Dodelson 2003), “comoving coordinate distance” (Mo et al. 2010),or “angular size distance” (Peebles 1993). The latter is not to be confused with the angular diameter distance.
Either
z_min
orz_max
can be a numpy array; in those cases, the samez_min
/z_max
is applied to all values of the other. If both are numpy arrays, they need to have the same dimensions, and the comoving distance returned corresponds to a series of differentz_min
andz_max
values.This function does not use interpolation (unlike the other distance functions) because it accepts both
z_min
andz_max
parameters which would necessitate a 2D interpolation. Thus, for fast evaluation, theluminosityDistance()
andangularDiameterDistance()
functions should be used. Parameters
 zmin: array_like
Redshift; can be a number or a numpy array.
 zmax: array_like
Redshift; can be a number or a numpy array.
 transverse: bool
Whether to return the transverse of lineofsight comoving distance. The two are the same in flat cosmologies.
 Returns
 d: array_like
The comoving distance in Mpc/h; has the same dimensions as
zmin
and/orzmax
.
See also
luminosityDistance
The luminosity distance to redshift z.
angularDiameterDistance
The angular diameter distance to redshift z.

luminosityDistance
(self, z, derivative=0, inverse=False)¶ The luminosity distance to redshift z.
 Parameters
 z: array_like
Redshift, where \(0.995 < z < 200\); can be a number or a numpy array.
 derivative: int
If greater than 0, evaluate the nth derivative, \(d^nD/dz^n\).
 inverse: bool
If True, evaluate \(z(D)\) instead of \(D(z)\). In this case, the
z
field must contain the luminosity distance in Mpc/h.
 Returns
 d: array_like
The luminosity distance (or its derivative) in Mpc/h; has the same dimensions as
z
.
See also
comovingDistance
The comoving distance between redshift \(z_{\rm min}\) and \(z_{\rm max}\).
angularDiameterDistance
The angular diameter distance to redshift z.

angularDiameterDistance
(self, z, derivative=0)¶ The angular diameter distance to redshift z.
The angular diameter distance is the transverse distance that, at redshift z, corresponds to an angle of one radian. Note that the inverse is not available for this function because it is not strictly increasing or decreasing with redshift, making its inverse multivalued.
 Parameters
 z: array_like
Redshift, where \(0.995 < z < 200\); can be a number or a numpy array.
 derivative: int
If greater than 0, evaluate the nth derivative, \(d^nD/dz^n\).
 Returns
 d: array_like
The angular diameter distance (or its derivative) in Mpc/h; has the same dimensions as
z
.
See also
comovingDistance
The comoving distance between redshift \(z_{\rm min}\) and \(z_{\rm max}\).
luminosityDistance
The luminosity distance to redshift z.

distanceModulus
(self, z)¶ The distance modulus to redshift z in magnitudes.
 Parameters
 z: array_like
Redshift; can be a number or a numpy array.
 Returns
 mu: array_like
The distance modulus in magnitudes; has the same dimensions as
z
.

soundHorizon
(self)¶ The sound horizon at recombination.
This function returns the sound horizon in Mpc/h, according to Eisenstein & Hu 1998, Equation 26. This fitting function is accurate to 2% where \(\Omega_{\rm b} h^2 > 0.0125\) and \(0.025 < \Omega_{\rm m} h^2 < 0.5\).
 Returns
 s: float
The sound horizon at recombination in Mpc/h.

rho_c
(self, z)¶ The critical density of the universe at redshift z.
 Parameters
 z: array_like
Redshift; can be a number or a numpy array.
 Returns
 rho_critical: array_like
The critical density in units of physical \(M_{\odot} h^2 / {\rm kpc}^3\); has the same dimensions as
z
.

rho_m
(self, z)¶ The matter density of the universe at redshift z.
 Parameters
 z: array_like
Redshift; can be a number or a numpy array.
 Returns
 rho_matter: array_like
The matter density in units of physical \(M_{\odot} h^2 / {\rm kpc}^3\); has the same dimensions as
z
.
See also
Om
The matter density of the universe, in units of the critical density.

rho_b
(self, z)¶ The baryon density of the universe at redshift z.
 Parameters
 z: array_like
Redshift; can be a number or a numpy array.
 Returns
 rho_baryon: array_like
The baryon density in units of physical \(M_{\odot} h^2 / {\rm kpc}^3\); has the same dimensions as
z
.
See also
Ob
The baryon density of the universe, in units of the critical density.

rho_de
(self, z)¶ The dark energy density of the universe at redshift z.
 Parameters
 z: array_like
Redshift; can be a number or a numpy array.
 Returns
 rho_de: float
The dark energy density in units of physical \(M_{\odot} h^2 / {\rm kpc}^3\); has the same dimensions as
z
.
See also
Ode
The dark energy density of the universe, in units of the critical density.

rho_gamma
(self, z)¶ The photon density of the universe at redshift z.
If
relspecies == False
, this function returns 0. Parameters
 z: array_like
Redshift; can be a number or a numpy array.
 Returns
 rho_gamma: array_like
The photon density in units of physical \(M_{\odot} h^2 / {\rm kpc}^3\); has the same dimensions as
z
.
See also
Ogamma
The density of photons in the universe, in units of the critical density.

rho_nu
(self, z)¶ The neutrino density of the universe at redshift z.
If
relspecies == False
, this function returns 0. Parameters
 z: array_like
Redshift; can be a number or a numpy array.
 Returns
 rho_nu: array_like
The neutrino density in units of physical \(M_{\odot} h^2 / {\rm kpc}^3\); has the same dimensions as
z
.
See also
Onu
The density of neutrinos in the universe, in units of the critical density.

rho_r
(self, z)¶ The density of relativistic species in the universe at redshift z.
This density is the sum of the photon and neutrino densities. If
relspecies == False
, this function returns 0. Parameters
 z: array_like
Redshift; can be a number or a numpy array.
 Returns
 rho_relativistic: array_like
The density of relativistic species in units of physical \(M_{\odot} h^2 / {\rm kpc}^3\); has the same dimensions as
z
.
See also
Or
The density of relativistic species in the universe, in units of the critical density.

Om
(self, z)¶ The matter density of the universe, in units of the critical density.
 Parameters
 z: array_like
Redshift; can be a number or a numpy array.
 Returns
 Omega_matter: array_like
Has the same dimensions as
z
.
See also
rho_m
The matter density of the universe at redshift z.

Ob
(self, z)¶ The baryon density of the universe, in units of the critical density.
 Parameters
 z: array_like
Redshift; can be a number or a numpy array.
 Returns
 Omega_baryon: array_like
Has the same dimensions as
z
.
See also
rho_b
The baryon density of the universe at redshift z.

Ode
(self, z)¶ The dark energy density of the universe, in units of the critical density.
 Parameters
 z: array_like
Redshift; can be a number or a numpy array.
 Returns
 Omega_de: array_like
Has the same dimensions as
z
.
See also
rho_de
The dark energy density of the universe at redshift z.

Ogamma
(self, z)¶ The density of photons in the universe, in units of the critical density.
 Parameters
 z: array_like
Redshift; can be a number or a numpy array.
 Returns
 Omega_gamma: array_like
Has the same dimensions as
z
.
See also
rho_gamma
The photon density of the universe at redshift z.

Onu
(self, z)¶ The density of neutrinos in the universe, in units of the critical density.
 Parameters
 z: array_like
Redshift; can be a number or a numpy array.
 Returns
 Omega_nu: array_like
Has the same dimensions as
z
.
See also
rho_nu
The neutrino density of the universe at redshift z.

Or
(self, z)¶ The density of relativistic species, in units of the critical density.
This function returns the sum of the densities of photons and neutrinos.
 Parameters
 z: array_like
Redshift; can be a number or a numpy array.
 Returns
 Omega_relativistic: array_like
Has the same dimensions as
z
.
See also
rho_r
The density of relativistic species in the universe at redshift z.

Ok
(self, z)¶ The curvature density of the universe in units of the critical density.
In a flat universe, \(\Omega_{\rm k} = 0\).
 Parameters
 z: array_like
Redshift; can be a number or a numpy array.
 Returns
 Omega_curvature: array_like
Has the same dimensions as
z
.

growthFactorUnnormalized
(self, z)¶ The linear growth factor, \(D_+(z)\).
The growth factor describes the linear evolution of over and underdensities in the dark matter density field. There are three regimes:
In the matterradiation regime, we use an approximate analytical formula (Equation 5 in Gnedin et al. 2011. If relativistic species are ignored, \(D_+(z) \propto a\).
In the matterdominated regime, \(D_+(z) \propto a\).
In the matterdark energy regime, we evaluate \(D_+(z)\) through integration as defined in Eisenstein & Hu 1999, Equation 8 (see also Heath 1977) for LCDM cosmologies. For cosmologies where \(w(z) \neq 1\), this expression is not valid and we instead solve the ordinary differential equation for the evolution of the growth factor (Equation 11 in Linder & Jenkins 2003).
At the transition between the integral and analytic approximation regimes, the two expressions do not quite match up, with differences of the order <1E3. in order to avoid a discontinuity, we introduce a transition regime where the two quantities are linearly interpolated.
The normalization is such that the growth factor approaches \(D_+(a) = a\) in the matterdominated regime. There are other normalizations of the growth factor (e.g., Percival 2005, Equation 15), but since we almost always care about the growth factor normalized to z = 0, the normalization does not matter too much (see the
growthFactor()
function). Parameters
 z: array_like
Redshift, where \(0.995 < z\); the high end of z is only limited by the validity of the analytical approximation mentioned above. Can be a number or a numpy array.
 Returns
 D: array_like
The linear growth factor; has the same dimensions as
z
.
Warning
This function directly evaluates the growth factor by integration or analytical approximation. In most cases, the
growthFactor()
function should be used since it interpolates and is thus much faster.See also
growthFactor
The linear growth factor normalized to z = 0, \(D_+(z) / D_+(0)\).

growthFactor
(self, z, derivative=0, inverse=False)¶ The linear growth factor normalized to z = 0, \(D_+(z) / D_+(0)\).
The growth factor describes the linear evolution of over and underdensities in the dark matter density field. This function is sped up through interpolation which barely degrades its accuracy, but if you wish to evaluate the exact integral or compute the growth factor for very high redshifts (z > 200), please use the
growthFactorUnnormalized()
function. Parameters
 z: array_like
Redshift, where \(0.995 < z < 200\); can be a number or a numpy array.
 derivative: int
If greater than 0, evaluate the nth derivative, \(d^nD_+/dz^n\).
 inverse: bool
If True, evaluate \(z(D_+)\) instead of \(D_+(z)\). In this case, the
z
field must contain the normalized growth factor.
 Returns
 D: array_like
The linear growth factor (or its derivative); has the same dimensions as
z
.
See also
growthFactorUnnormalized
The linear growth factor, \(D_+(z)\).

matterPowerSpectrum
(self, k, z=0.0, model='eisenstein98', path=None, derivative=False)¶ The matter power spectrum at a scale k.
By default, the power spectrum is computed using a model for the transfer function (see the
transferFunction()
function). The default Eisenstein & Hu 1998 approximation is accurate to about 5%, and the interpolation introduces errors significantly smaller than that.Alternatively, the user can supply the path to a file with a tabulated power spectrum using the
path
parameter. The file must contain two columns, \(\log_{10}(k)\) and \(\log_{10}(P)\) where k and P(k) are in the same units as in this function. This table is interpolated with a thirdorder spline. Note that the tabulated spectrum is normalized to the value if \(\sigma_8\) set in the cosmology.Also note that if a power spectrum is to be read from a file, the
persistence
parameter must allow for reading (though not necessarily writing) of files. Parameters
 k: array_like
The wavenumber k (in comoving h/Mpc), where \(10^{20} < k < 10^{20}\); can be a number or a numpy array. If a usersupplied table is used, the limits of that table apply.
 z: float
The redshift at which the power spectrum is evaluated, zero by default. If nonzero, the power spectrum is scaled with the linear growth factor, \(P(k, z) = P(k, 0) D_{+}^2(z)\).
 model: str
A model for the power spectrum (see the
cosmology.power_spectrum
module). If a tabulated power spectrum is used (seepath
parameter), this name must still be passed. Internally, the power spectrum is saved using this name, so the name must not overlap with any other models. path: str
A path to a file containing the power spectrum as a table, where the two columns are \(\log_{10}(k)\) (in comoving h/Mpc) and \(\log_{10}(P)\) (in \(({\rm Mpc}/h)^3\)).
 derivative: bool
If False, return P(k). If True, return \(d \log(P) / d \log(k)\).
 Returns
 Pk: array_like
The matter power spectrum; has the same dimensions as k and units of \(({\rm Mpc}/h)^3\). Alternatively, the dimensionless logarithmic derivative if
derivative == True
.
Warning
If a usersupplied power spectrum table is used, integrals over the power spectrum such as the variance and correlation function are integrated only within the limits of the given power spectrum. By default, Boltzmann codes return a relatively small range in wavenumber. Please increase this range if necessary, and check that the computed quantities are converged.

filterFunction
(self, filt, k, R)¶ The Fourier transform of certain filter functions.
The main use of the filter function is in computing the variance, please see the documentation of the
sigma()
function for details. This function is dimensionless, the input units are k in comoving h/Mpc and R in comoving Mpc/h. Available filters aretophat
,\[\tilde{W}_{\rm tophat} = \frac{3}{(kR)^3} \left[ \sin(kR)  kR \times \cos(kR) \right] \,,\]a
gaussian
filter,\[\tilde{W}_{\rm gaussian} = \exp \left[ \frac{(kR)^2}{2} \right] \,,\]and a
sharpk
filter,\[\tilde{W}_{\rm sharpk} = \Theta(1  kR) \,, \]where \(\Theta\) is the Heaviside step function.
 Parameters
 filt: str
Either
tophat
(a tophat filter in real space),sharpk
(a tophat filter in Fourier space), orgaussian
(a Gaussian in both real and Fourier space). k: float
A wavenumber k (in comoving h/Mpc).
 R: float
A radius R (in comoving Mpc/h).
 Returns
 filter: float
The value of the filter function.

sigma
(self, R, z, j=0, filt='tophat', inverse=False, derivative=False, kmin=None, kmax=None, ps_args={'model': 'eisenstein98', 'path': None})¶ The rms variance of the linear density field on a scale R, \(\sigma(R)\).
The variance and its higher moments are defined as the integral
\[\sigma^2(R,z) = \frac{1}{2 \pi^2} \int_0^{\infty} k^2 k^{2j} P(k,z) \tilde{W}(kR)^2 dk\]where \(\tilde{W}(kR)\) is the Fourier transform of the
filterFunction()
, and \(P(k,z) = D_+^2(z)P(k,0)\) is thematterPowerSpectrum()
. See the documentation offilterFunction()
for possible filters.By default, the power spectrum is computed using the transfer function approximation of Eisenstein & Hu 1998 (see the
cosmology.power_spectrum
module). With this approximation, the variance is accurate to about 2% or better (see the Colossus code paper for details). Using a tabulated power spectrum can make this computation more accurate, but please note that the limits of the corresponding table are used for the integration.Higher moments of the variance (such as \(\sigma_1\), \(\sigma_2\) etc) can be computed by setting j > 0 (see Bardeen et al. 1986). Furthermore, the logarithmic derivative \(d \log(\sigma) / d \log(R)\) can be evaluated by setting
derivative == True
. Parameters
 R: array_like
The radius of the filter in comoving Mpc/h, where \(10^{12} < R < 10^3\); can be a number or a numpy array.
 z: float
Redshift; for z > 0, \(\sigma(R)\) is multiplied by the linear growth factor.
 j: integer
The order of the integral. j = 0 corresponds to the variance, j = 1 to the same integral with an extra \(k^2\) term etc; see Bardeen et al. 1986 for mathematical details.
 filt: str
Either
tophat
,sharpk
orgaussian
(seefilterFunction()
). Higher moments (j > 0) can only be computed for the gaussian filter. inverse: bool
If True, compute \(R(\sigma)\) rather than \(\sigma(R)\).
 derivative: bool
If True, return the logarithmic derivative, \(d \log(\sigma) / d \log(R)\), or its inverse, \(d \log(R) / d \log(\sigma)\) if
inverse == True
. kmin: float
The lower limit of the variance integral in kspace. If
None
, the limit is determined automatically (it should be zero in principle, but will be set to a very small number depending on the type of power spectrum). Settingkmin
can be useful when considering finite simulation volumes where the largest scales (smallest kmodes) are not taken into account. Ifkmin
is notNone
and a tophat filter is used, the variance will oscillate at certain radii. Thus, the number of points in the interpolation table is automatically increased, but the solution is nevertheless unreliable close to the cutoff scale. kmax: float
The upper limit of the variance integral in kspace. If
None
, the limit is determined automatically (it should be infinity in principle, but will be set to a large number where the power spectrum has fallen off sufficiently that the integral is converged). Ifkmax
is notNone
and a tophat filter is used, the variance will oscillate at certain radii. Thus, the number of points in the interpolation table is automatically increased, but the solution is nevertheless unreliable close to the cutoff scale. ps_args: dict
Arguments passed to the
matterPowerSpectrum()
function.
 Returns
 sigma: array_like
The rms variance; has the same dimensions as
R
. If inverse and/or derivative are True, the inverse, derivative, or derivative of the inverse are returned. If j > 0, those refer to higher moments.
See also
matterPowerSpectrum
The matter power spectrum at a scale k.

correlationFunction
(self, R, z=0.0, derivative=False, ps_args={'model': 'eisenstein98', 'path': None})¶ The linear mattermatter correlation function at radius R.
The linear correlation function is defined as
\[\xi(R,z) = \frac{1}{2 \pi^2} \int_0^\infty k^2 P(k,z) \frac{\sin(kR)}{kR} dk\]where P(k) is the
matterPowerSpectrum()
. By default, the power spectrum is computed using the transfer function approximation of Eisenstein & Hu 1998 (see thecosmology.power_spectrum
module). With this approximation, the correlation function is accurate to ~5% over the range \(10^{2} < R < 200\) (see the Colossus code paper for details). Using a tabulated power spectrum can make this computation more accurate, but please note that the limits of the corresponding table are used for the integration. Parameters
 R: array_like
The radius in comoving Mpc/h; can be a number or a numpy array.
 z: float
Redshift; if nonzero, the correlation function is scaled with the linear growth factor, \(\xi(R, z) = \xi(R, 0) D_{+}^2(z)\).
 derivative: bool
If
derivative == True
, the linear derivative \(d \xi / d R\) is returned. ps_args: dict
Arguments passed to the
matterPowerSpectrum()
function.
 Returns
 xi: array_like
The correlation function, or its derivative; has the same dimensions as
R
.
See also
matterPowerSpectrum
The matter power spectrum at a scale k.

cosmology.cosmology.
setCosmology
(cosmo_name, params=None)¶ Set a cosmology.
This function provides a convenient way to create a cosmology object without setting the parameters of the
Cosmology
class manually. See the Basic Usage section for examples. Whichever way the cosmology is set, the global variable is updated so that thegetCurrent()
function returns the set cosmology. Parameters
 cosmo_name: str
The name of the cosmology. Can be the name of a preset cosmology, or another name in which case the
params
dictionary needs to be provided. params: dictionary
The parameters of the constructor of the
Cosmology
class. Not necessary ifcosmo_name
is the name of a preset cosmology.
 Returns
 cosmo: Cosmology
The created cosmology object.

cosmology.cosmology.
addCosmology
(cosmo_name, params)¶ Add a set of cosmological parameters to the global list.
After this function is executed, the new cosmology can be set using
setCosmology()
from anywhere in the code. Parameters
 cosmo_name: str
The name of the cosmology.
 params: dictionary
A set of parameters for the constructor of the Cosmology class.

cosmology.cosmology.
setCurrent
(cosmo)¶ Set the current global cosmology to a cosmology object.
Unlike
setCosmology()
, this function does not create a new cosmology object, but allows the user to set a cosmology object to be the current cosmology. This can be useful when switching between cosmologies, since many routines use thegetCurrent()
routine to obtain the current cosmology. Parameters
 cosmo: Cosmology
The cosmology object to be set as the global current cosmology.

cosmology.cosmology.
getCurrent
()¶ Get the current global cosmology.
This function should be used whenever access to the cosmology is needed. By using the globally set cosmology, there is no need to pass cosmology objects around the code. If no cosmology is set, this function raises an Exception that reminds the user to set a cosmology.
 Returns
 cosmo: Cosmology
The current globally set cosmology.

cosmology.cosmology.
fromAstropy
(astropy_cosmo, sigma8, ns, cosmo_name=None, **kwargs)¶ Convert an Astropy cosmology object to Colossus and set it as current cosmology.
This function throws an error if Astropy cannot be imported. The user needs to supply sigma8 and ns because Astropy does not include power spectrum calculations in its cosmology module. There are a few dark energy model implemented in Astropy that are not implemented in Colossus, namely “wpwaCDM” and “w0wzCDM”. The corresponding classes cannot be converted with this function.
 Parameters
 astropy_cosmo: FLRW
Astropy cosmology object of type FLRW or its derivative classes.
 sigma8: float
The normalization of the power spectrum, i.e. the variance when the field is filtered with a top hat filter of radius 8 Mpc/h. This parameter is not part of an Astropy cosmology but must be set in Colossus.
 ns: float
The tilt of the primordial power spectrum. This parameter is not part of an Astropy cosmology but must be set in Colossus.
 cosmo_name: str
The name of the cosmology. Can be
None
if a name is supplied in the Astropy cosmology, which is optional in Astropy. kwargs: dictionary
Additional parameters that are passed to the Colossus cosmology.
 Returns
 cosmo: Cosmology
The cosmology object derived from Astropy, which is also set globally.
Power spectrum¶
This module implements models for the matter power spectrum, or more exactly, for the transfer
function. Generally speaking, the transfer function should be evaluated using the
matterPowerSpectrum()
function. This module is automatically
imported with the cosmology module.
Power spectrum models¶
The following models are supported, and are listed in the models
dictionary. Their ID can
be passed as the model
parameter to the transferFunction()
function:
ID 
Reference 
Comment 

sugiyama95 
A semianalytical fitting function 

eisenstein98 
A semianalytical fitting function 

eisenstein98_zb 
The zerobaryon version, i.e., no BAO 
Module contents¶
Characteristics of power spectrum models. 

Dictionary containing a list of models. 


The transfer function. 

The transfer function according to Sugiyama 1995. 

The transfer function according to Eisenstein & Hu 1998. 

The zerobaryon transfer function according to Eisenstein & Hu 1998. 
Module reference¶

class
cosmology.power_spectrum.
PowerSpectrumModel
¶ Characteristics of power spectrum models.
This object is currently empty. The
models
dictionary contains one item of this class for each available model.

cosmology.power_spectrum.
models
= {'eisenstein98': <cosmology.power_spectrum.PowerSpectrumModel object>, 'eisenstein98_zb': <cosmology.power_spectrum.PowerSpectrumModel object>, 'sugiyama95': <cosmology.power_spectrum.PowerSpectrumModel object>}¶ Dictionary containing a list of models.
An ordered dictionary containing one
PowerSpectrumModel
entry for each model.

cosmology.power_spectrum.
transferFunction
(k, h, Om0, Ob0, Tcmb0, model='eisenstein98')¶ The transfer function.
The transfer function transforms the spectrum of primordial fluctuations into the linear power spectrum of the matter density fluctuations. The primordial power spectrum is usually described as a power law, leading to a power spectrum
\[P(k) = T(k)^2 k^{n_s}\]where P(k) is the matter power spectrum, T(k) is the transfer function, and \(n_s\) is the tilt of the primordial power spectrum. See the
Cosmology
class for further details on the cosmological parameters. Parameters
 k: array_like
The wavenumber k (in comoving h/Mpc); can be a number or a numpy array.
 h: float
The Hubble constant in units of 100 km/s/Mpc.
 Om0: float
\(\Omega_{\rm m}\), the matter density in units of the critical density at z = 0.
 Ob0: float
\(\Omega_{\rm b}\), the baryon density in units of the critical density at z = 0.
 Tcmb0: float
The temperature of the CMB at z = 0 in Kelvin.
 Returns
 Tk: array_like
The transfer function; has the same dimensions as
k
.

cosmology.power_spectrum.
modelSugiyama95
(k, h, Om0, Ob0, Tcmb0)¶ The transfer function according to Sugiyama 1995.
This function computes the Sugiyama 1995 approximation to the transfer function at a scale k, which is based on the Bardeen et al. 1986 formulation. Note that this approximation is not as accurate as the
eisenstein98
model, with deviations of about 1020% in the power spectrum, variance, and correlation function. Parameters
 k: array_like
The wavenumber k (in comoving h/Mpc); can be a number or a numpy array.
 h: float
The Hubble constant in units of 100 km/s/Mpc.
 Om0: float
\(\Omega_{\rm m}\), the matter density in units of the critical density at z = 0.
 Ob0: float
\(\Omega_{\rm b}\), the baryon density in units of the critical density at z = 0.
 Tcmb0: float
The temperature of the CMB at z = 0 in Kelvin.
 Returns
 Tk: array_like
The transfer function; has the same dimensions as
k
.

cosmology.power_spectrum.
modelEisenstein98
(k, h, Om0, Ob0, Tcmb0)¶ The transfer function according to Eisenstein & Hu 1998.
This function computes the Eisenstein & Hu 1998 approximation to the transfer function at a scale k. The code was adapted from Matt Becker’s cosmocalc code.
This function was tested against numerical calculations based on the CAMB code (Lewis et al. 2000) and found to be accurate to 5% or better up to k of about 100 h/Mpc (see the Colossus code paper for details).
 Parameters
 k: array_like
The wavenumber k (in comoving h/Mpc); can be a number or a numpy array.
 h: float
The Hubble constant in units of 100 km/s/Mpc.
 Om0: float
\(\Omega_{\rm m}\), the matter density in units of the critical density at z = 0.
 Ob0: float
\(\Omega_{\rm b}\), the baryon density in units of the critical density at z = 0.
 Tcmb0: float
The temperature of the CMB at z = 0 in Kelvin.
 Returns
 Tk: array_like
The transfer function; has the same dimensions as
k
.
See also
modelEisenstein98ZeroBaryon
The zerobaryon transfer function according to Eisenstein & Hu 1998.

cosmology.power_spectrum.
modelEisenstein98ZeroBaryon
(k, h, Om0, Ob0, Tcmb0)¶ The zerobaryon transfer function according to Eisenstein & Hu 1998.
This fitting function is significantly simpler than the full
modelEisenstein98()
version, and still approximates numerical calculations from a Boltzmann code to better than 10%, and almost as accurate when computing the variance or correlation function (see the Colossus code paper for details). Parameters
 k: array_like
The wavenumber k (in comoving h/Mpc); can be a number or a numpy array.
 h: float
The Hubble constant in units of 100 km/s/Mpc.
 Om0: float
\(\Omega_{\rm m}\), the matter density in units of the critical density at z = 0.
 Ob0: float
\(\Omega_{\rm b}\), the baryon density in units of the critical density at z = 0.
 Tcmb0: float
The temperature of the CMB at z = 0 in Kelvin.
 Returns
 Tk: array_like
The transfer function; has the same dimensions as
k
.
See also
modelEisenstein98
The transfer function according to Eisenstein & Hu 1998.