wale.CosmologyModel module

class wale.CosmologyModel.Cosmology_function(h, Oc, Ob, w=-1.0, wa=0.0, **kwargs)[source]

Bases: object

A class for initializing and handling cosmological computations using PyCCL.

This class encapsulates key cosmological parameters and provides methods to compute quantities such as the Hubble parameter, comoving distances, redshift inversions, non-linear matter power spectrum, and lensing weights. Internally, it constructs a pyccl.Cosmology object for use with the Core Cosmology Library (CCL).

Methods

get_H(z)

Compute the Hubble parameter H(z) in km/s/Mpc.

get_chi(z)

Compute the comoving radial distance χ(z) in Mpc.

get_lensing_weight(chis, chisource, **kwargs)

Dispatch method to compute lensing weights depending on source distribution.

get_lensing_weight_array(chis, chi_source)

Compute lensing weight array W(chi) for a single source plane at chi_source.

get_lensing_weight_array_nz(chis, z_nz, n_z)

Compute lensing weight W(chi) using a redshift distribution n(z), without relying on pyccl.

get_nonlinear_pk(z[, ks])

Compute the non-linear matter power spectrum P(k, z) using HALOFIT.

get_z_from_chi(chi_target[, z_max, tol])

Invert χ(z) to find z such that χ(z)=chi_target.

get_H(z)[source]

Compute the Hubble parameter H(z) in km/s/Mpc.

Parameters:

z (float or array_like) – Redshift(s) at which to evaluate H(z).

Returns:

H – Hubble parameter at redshift z, in units of km/s/Mpc.

Return type:

float or ndarray

get_chi(z)[source]

Compute the comoving radial distance χ(z) in Mpc.

This function numerically integrates the inverse of the dimensionless Hubble parameter E(z) = H(z)/H0 to compute:

χ(z) = ∫₀ᶻ (c / H(z’)) dz’

Parameters:

z (float or array_like) – Redshift(s) at which to evaluate the comoving distance.

Returns:

chi – Comoving radial distance(s) to redshift z, in units of Mpc.

Return type:

float or ndarray

Notes

  • Uses scipy.integrate.quad for numerical integration.

  • Assumes a flat universe with constant c = 299792.458 km/s.

  • The integration is done individually for each redshift if an array is given.

get_lensing_weight(chis, chisource, **kwargs)[source]

Dispatch method to compute lensing weights depending on source distribution.

Parameters:
  • chis (array_like) – Comoving radial distances where lensing weights are evaluated.

  • chisource (float) – Source-plane comoving distance (used only if no n(z) is provided).

  • **kwargs

    Optional arguments: - nz_file : str, optional

    If provided, uses redshift distribution from file (not yet implemented).

Returns:

  • z_values (ndarray) – Redshifts corresponding to chis.

  • lensing_weight (ndarray) – Lensing weights W(chi) using either a single plane or n(z).

get_lensing_weight_array(chis, chi_source)[source]

Compute lensing weight array W(chi) for a single source plane at chi_source.

Parameters:
  • chis (array_like) – Comoving radial distances (chi) at which to compute lensing weights.

  • chi_source (float) – Comoving distance to the source plane (single redshift).

Returns:

  • z_values (ndarray) – Redshift values corresponding to chis, computed via inversion.

  • lensing_weight (ndarray) – Lensing weight W(chi) evaluated at each input chi.

get_lensing_weight_array_nz(chis, z_nz, n_z)[source]

Compute lensing weight W(chi) using a redshift distribution n(z), without relying on pyccl.

Parameters:
  • chis (array_like) – Comoving radial distances (chi) at which to compute lensing weights.

  • z_nz (array_like) – Redshift values at which the distribution n(z) is defined.

  • n_z (array_like) – Source redshift distribution n(z), not necessarily normalized.

Returns:

  • z_values (ndarray) – Redshift values corresponding to the input chis, estimated via interpolation.

  • lensing_weight (ndarray) – Lensing weight values W(chi) computed using the n(z) distribution.

get_nonlinear_pk(z, ks=None)[source]

Compute the non-linear matter power spectrum P(k, z) using HALOFIT.

Parameters:
  • z (float) – Redshift at which to evaluate the power spectrum.

  • ks (array_like, optional) – Wavenumber values in h/Mpc. If None, uses the default self.k array.

Returns:

Pnl – Non-linear matter power spectrum at redshift z for the given ks, in (Mpc/h)^3.

Return type:

ndarray

get_z_from_chi(chi_target, z_max=10.0, tol=1e-08)[source]

Invert χ(z) to find z such that χ(z)=chi_target.

Parameters:
  • chi_target (float or array_like) – Comoving distance(s) in Mpc.

  • z_max (float) – Maximum redshift to search within (must be large enough that χ(z_max) > chi_target).

  • tol (float) – Desired precision in z.

Returns:

z – Redshift(s) corresponding to chi_target.

Return type:

float or ndarray