`lyapower.CLASS`¶

Linear power spectra and transfer functions (CLASS / cosmoprimo).

lyapower.CLASS ¶

Linear power spectra and transfer functions from CLASS / cosmoprimo.

Wraps the Boltzmann code CLASS (via classy) in :class:MyClass to compute and write linear matter power spectra P(k) and transfer functions T_i(k) in the CLASS or CAMB conventions (as needed by MP-Gadget / MP-Genic / 2LPT / CosmicIC initial-condition codes), with helpers for grid scans and massive neutrinos. :class:CosmoprimoInterface provides the same linear P(k) (and its no-wiggle version) through the cosmoprimo façade.

Convention: wavenumber k in h/Mpc, power spectrum in (h^-1 Mpc)^3.

@author: cravoux

MyClass ¶

Bases: object

A configured CLASS model producing power spectra and transfer functions.

Source code in lyapower/CLASS.py

class MyClass(object):
    """A configured CLASS model producing power spectra and transfer functions."""

    def __init__(self, pwd, settings, hierarchy="EQUI", class_version="public"):
        """Instantiate and configure the underlying CLASS model.

        Args:
            pwd (str): Working directory for outputs.
            settings (dict): CLASS input parameters.
            hierarchy (str, optional): Neutrino mass hierarchy label.
            class_version (str, optional): ``"public"`` (``classy``) or
                ``"pk1D"`` (``classy_pk1D``).
        """
        if class_version == "public":
            from classy import Class
        elif class_version == "pk1D":
            from classy_pk1D import Class

        self.pwd = pwd
        self.settings = settings
        self.model = Class()
        self.model.set(settings)

        self.define_output_format()

        self.hierarchy = hierarchy

    def define_output_format(self):
        """Set ``self.output_format`` from the settings (default ``"class"``)."""
        if "format" not in self.settings.keys():
            self.output_format = "class"
        else:
            self.output_format = self.settings["format"]

    def write_pk_tk(self, z, name, kmin=-4, kmax=3, nb_points=2000, output=True, verbose=True):
        """return Power Spectra with k in h.Mpc-1 and P in h3.Mpc-3 in the convention of CAMB and CLASS (the one needed for MP-Gadget).
        return Transfer function in the convention of CLASS (needed for MP-Genic) or CAMB (needed for 2LPT) !!!
        note : . CLASS format for MP-GENIC, CAMB format for 2LPT, CAMB format for COSMICIC
               . for 2LPT & COSMICIC: no header! k vector need to be equidistant in log-k space and length of pk must be the same than tk
        """
        self.model.compute()
        Power, sigma_8 = self.write_pk(
            name, z, kmin=kmin, kmax=kmax, nb_points=nb_points,output=output, verbose=verbose
        )
        Transfer = self.write_tk(name, z,output=output)
        return (Power, Transfer, sigma_8)

    def write_pk_tk_grid_variation(
        self,
        z,
        name,
        varying_dictionary,
        derived_params=None,
        kmin=-4,
        kmax=3,
        nb_points=2000,
        output=True,
    ):
        """return Power Spectra with k in h.Mpc-1 and P in h3.Mpc-3 in the convention of CAMB and CLASS (the one needed for MP-Gadget).
        return Transfer function in the convention of CLASS (needed for MP-Genic) or CAMB (needed for 2LPT) !!!
        note : . CLASS format for MP-GENIC, CAMB format for 2LPT, CAMB format for COSMICIC
               . for 2LPT & COSMICIC: no header! k vector need to be equidistant in log-k space and length of pk must be the same than tk
        """
        dict_deriv = {}
        if derived_params is not None:
            for i in range(len(derived_params)):
                dict_deriv[derived_params[i]] = []
        list_key = list(varying_dictionary.keys())
        for i in range(len(varying_dictionary[list_key[0]])):
            dict_to_set = {
                list_key[j]: varying_dictionary[list_key[j]][i]
                for j in range(len(list_key))
            }
            self.model.set(dict_to_set)
            self.model.compute()
            Power, sigma_8 = self.write_pk(
                name, z, kmin=kmin, kmax=kmax, nb_points=nb_points, output= output
            )
            Transfer = self.write_tk(name, z,output=output)
            if derived_params is not None:
                par = self.model.get_current_derived_parameters(derived_params)
                print(par.keys())
                for key in par.keys():
                    dict_deriv[key].append(par[key])
        return (Power, Transfer, sigma_8, dict_deriv)

    def write_pk_tk_neutrino_mass(
        self,
        z,
        name,
        kmin=-4,
        kmax=3,
        nb_points=2000,
        nb_neutrinos=3,
        nb_neutrinos_massif=1,
        neutrino_mass=None,
        method="cbnu",
        output=True,
    ):
        """return Power Spectra with k in h.Mpc-1 and P in h3.Mpc-3 in the convention of CAMB and CLASS (the one needed for MP-Gadget).
        return Transfer function in the convention of CLASS (needed for MP-Genic) or CAMB (needed for 2LPT) !!!
        """
        m_ncdm = ""
        Omega_ncdm = ""
        for j in range(nb_neutrinos):
            m_ncdm = m_ncdm + "0.0,"
            Omega_ncdm = Omega_ncdm + "0.0,"
        m_ncdm = m_ncdm[0:-1]
        Omega_ncdm = Omega_ncdm[0:-1]
        self.model.set(
            {"N_ncdm": nb_neutrinos, "m_ncdm": m_ncdm, "Omega_ncdm": Omega_ncdm}
        )
        self.model.compute()
        if neutrino_mass is not None:
            omega_cdm_ref = self.model.Omega0_cdm()
            omega_nu_ref = self.model.Omega_nu
            m_ncdm = ""
            for j in range(nb_neutrinos_massif):
                m_ncdm = m_ncdm + str(neutrino_mass / nb_neutrinos_massif) + ","
            for j in range(nb_neutrinos - nb_neutrinos_massif):
                m_ncdm = m_ncdm + "0.0,"
            print(m_ncdm)
            m_ncdm = m_ncdm[0:-1]
            if method == "cbnu":
                self.model.set(
                    {
                        "N_ncdm": nb_neutrinos,
                        "m_ncdm": m_ncdm,
                        "omega_cdm": str(
                            (
                                omega_cdm_ref
                                - self.get_omeganu_from_mass(neutrino_mass)
                                + omega_nu_ref
                            )
                            * self.model.h() ** 2
                        ),
                    }
                )
            if method == "cb":
                self.model.set({"N_ncdm": nb_neutrinos, "m_ncdm": m_ncdm})
            self.model.compute()

        Power, sigma_8 = self.write_pk(
            name, z, kmin=kmin, kmax=kmax, nb_points=nb_points, output=output
        )
        Transfer = self.write_tk(name, z,output=output)
        return (Power, Transfer, sigma_8)

    def write_pk(
        self, name, z, kmin=-4, kmax=3, nb_points=2000, output=True, header_output=True, verbose=True
    ):
        """Compute and (optionally) write the linear power spectrum P(k).

        Args:
            name (str): Output file base name.
            z (float): Redshift.
            kmin, kmax (float, optional): log10 wavenumber bounds (h/Mpc).
            nb_points (int, optional): Number of k samples.
            output (bool, optional): Write the ``pk_{name}_z{z}.dat`` file.
            header_output (bool, optional): Include a descriptive header.
            verbose (bool, optional): Print the cosmology summary.

        Returns:
            tuple: ``(Power, sigma_8)`` — an ``(N, 2)`` ``[k, P(k)]`` array and
            sigma_8.
        """
        if self.output_format == "class":
            (k_space, Pk, sigma_8) = self.compute_power_spectrum(
                z, kmin=kmin, kmax=kmax, nb_points=nb_points, verbose=verbose
            )  # k in h/Mpc
            header = header_pk_class.format(z, kmin, kmax, nb_points)
        elif self.output_format == "camb":
            tr = self.model.get_transfer(z, output_format="camb")
            k = np.array(tr["k (h/Mpc)"])
            kpower = np.logspace(
                *(np.log10(k[[0, -1]])) * (1 - 1e-5), len(k)
            )  # k in h/Mpc
            (k_space, Pk, sigma_8) = self.compute_power_spectrum(
                z, k_array=kpower, verbose=verbose
            )
            header = header_pk_camb.format(
                z, np.min(k_space), np.max(k_space), len(k_space)
            )
        Power = np.stack([k_space, Pk], axis=1)
        if output:
            if header_output:
                np.savetxt("pk_{}_z{}.dat".format(name, z), Power, header=header)
            else:
                np.savetxt("pk_{}_z{}.dat".format(name, z), Power)
        return (Power, sigma_8)

    def write_tk(self, name, z, output=True, header_output=True):
        """Compute and (optionally) write the transfer functions T_i(k).

        Args:
            name (str): Output file base name.
            z (float): Redshift.
            output (bool, optional): Write the ``tk_{name}_z{z}.dat`` file.
            header_output (bool, optional): Include a descriptive header.

        Returns:
            numpy.ndarray: The transfer-function table (CLASS or CAMB format).
        """
        if self.output_format == "class":
            Transfer_dict = self.model.get_transfer(z, output_format="class")
            Transfer = np.stack([Tk for key, Tk in Transfer_dict.items()], axis=1)
            header = header_tk_class.format(
                z, np.min(Transfer[:, 0]), np.max(Transfer[:, 0]), len(Transfer[:, 0])
            )
        elif self.output_format == "camb":
            Transfer_dict = self.model.get_transfer(z, output_format="camb")
            Transfer = np.stack([Tk for key, Tk in Transfer_dict.items()], axis=1)
            Transfer = self.interp_to_equidistant_log_space_camb_format(Transfer)
            header = header_tk_camb.format(
                z, np.min(Transfer[:, 0]), np.max(Transfer[:, 0]), len(Transfer[:, 0])
            )
        list_key = [key for key, Tk in Transfer_dict.items()]
        head = ""
        for i in range(len(list_key)):
            head = head + "            " + str(i + 1) + ":" + list_key[i]
        header = header + head
        if output:
            if header_output:
                np.savetxt(
                    "tk_{}_z{}.dat".format(name, z),
                    Transfer,
                    header=header,
                    delimiter="    ",
                )
            else:
                np.savetxt("tk_{}_z{}.dat".format(name, z), Transfer, delimiter="    ")
        return Transfer

    def compute_power_spectrum(
        self, z, kmin=-4, kmax=3, nb_points=2000, k_array=None, verbose=True
    ):
        """Evaluate the linear power spectrum on a k-grid (h-normalised).

        Args:
            z (float): Redshift.
            kmin, kmax (float, optional): log10 wavenumber bounds (h/Mpc).
            nb_points (int, optional): Number of k samples (if ``k_array`` None).
            k_array (numpy.ndarray, optional): Explicit k grid (h/Mpc).
            verbose (bool, optional): Print the cosmology summary.

        Returns:
            tuple: ``(k_space, Pk, sigma_8)`` with ``Pk`` in ``(h^-1 Mpc)^3``.
        """
        if k_array is None:
            k_space = np.logspace(kmin, kmax, num=nb_points)
        else:
            k_space = k_array
        sigma_8 = self.model.sigma(8 / self.model.h(), z)
        if verbose:
            print("Omega matter = " + str(self.model.Omega_m()))
            print("Omega lambda = " + str(self.model.Omega_Lambda()))
            print("Omega baryon = " + str(self.model.Omega_b()))
            print("Omega dm = " + str(self.model.Omega0_cdm()))
            print("Omega k = " + str(self.model.Omega0_k()))
            print(
                "Omega sum = "
                + str(
                    self.model.Omega0_cdm() + self.model.Omega_b() + self.model.Omega_nu
                )
            )
            print("Omega nu = " + str(self.model.Omega_nu))
            print("Omega rad = " + str(self.model.Omega_g()))
            print("sigma 8 at (z={}) = {}".format(z, sigma_8))
        Pk = []
        h = self.model.h()
        for k in k_space:
            Pk.append(self.model.pk(k * h, z) * h**3)
        return (k_space, Pk, sigma_8)

    def interp_to_equidistant_log_space_camb_format(self, Transfer):
        """Resample a transfer-function table onto a log-equidistant k grid.

        Required by the CAMB output format (2LPT / CosmicIC expect equidistant
        log-k).

        Args:
            Transfer (numpy.ndarray): Transfer table (first column is k).

        Returns:
            numpy.ndarray: The resampled table.
        """
        k_space = Transfer[:, 0]
        k_log_space = np.logspace(
            *(np.log10(k_space[[0, -1]])) * (1 - 1e-5), len(k_space)
        )  # k in h/Mpc
        for i in range(1, Transfer.shape[-1]):
            interp = interp1d(k_space, Transfer[:, i], bounds_error=True)
            Transfer[:, i] = interp(k_log_space)
        Transfer[:, 0] = k_log_space
        return Transfer

    def sigmaR_all_species_at(self, R, z):
        """Compute sigma(R, z) over all species.

        Args:
            R (float): Smoothing radius (Mpc/h).
            z (float): Redshift.

        Returns:
            float: The RMS density fluctuation sigma(R, z).
        """
        self.model.compute()
        return self.model.sigma(R / self.model.h(), z)

    def free_structure(self):
        """Free the CLASS internal structures (``struct_cleanup``)."""
        self.model.struct_cleanup()

    def close_class(self):
        """Empty the CLASS model (release its parameters)."""
        self.model.empty()

    def get_omeganu_from_mass(self, mass):
        """Convert a summed neutrino mass (eV) to ``omega_nu = m / (93.14 h^2)``.

        Args:
            mass (float): Summed neutrino mass (eV).

        Returns:
            float: ``omega_nu`` (physical density parameter).
        """
        return mass / (93.14 * (self.model.h()) ** 2)

    def get_mass_from_omeganu(self, omeganu):
        """Convert ``omega_nu`` back to a summed neutrino mass (eV).

        Args:
            omeganu (float): Physical neutrino density parameter.

        Returns:
            float: Summed neutrino mass (eV).
        """
        return omeganu * (93.14) * (self.model.h()) ** 2

define_output_format ¶

define_output_format()

Set self.output_format from the settings (default "class").

Source code in lyapower/CLASS.py

def define_output_format(self):
    """Set ``self.output_format`` from the settings (default ``"class"``)."""
    if "format" not in self.settings.keys():
        self.output_format = "class"
    else:
        self.output_format = self.settings["format"]

write_pk_tk ¶

write_pk_tk(z, name, kmin=-4, kmax=3, nb_points=2000, output=True, verbose=True)

return Power Spectra with k in h.Mpc-1 and P in h3.Mpc-3 in the convention of CAMB and CLASS (the one needed for MP-Gadget). return Transfer function in the convention of CLASS (needed for MP-Genic) or CAMB (needed for 2LPT) !!! note : . CLASS format for MP-GENIC, CAMB format for 2LPT, CAMB format for COSMICIC . for 2LPT & COSMICIC: no header! k vector need to be equidistant in log-k space and length of pk must be the same than tk

Source code in lyapower/CLASS.py

def write_pk_tk(self, z, name, kmin=-4, kmax=3, nb_points=2000, output=True, verbose=True):
    """return Power Spectra with k in h.Mpc-1 and P in h3.Mpc-3 in the convention of CAMB and CLASS (the one needed for MP-Gadget).
    return Transfer function in the convention of CLASS (needed for MP-Genic) or CAMB (needed for 2LPT) !!!
    note : . CLASS format for MP-GENIC, CAMB format for 2LPT, CAMB format for COSMICIC
           . for 2LPT & COSMICIC: no header! k vector need to be equidistant in log-k space and length of pk must be the same than tk
    """
    self.model.compute()
    Power, sigma_8 = self.write_pk(
        name, z, kmin=kmin, kmax=kmax, nb_points=nb_points,output=output, verbose=verbose
    )
    Transfer = self.write_tk(name, z,output=output)
    return (Power, Transfer, sigma_8)

write_pk_tk_grid_variation ¶

write_pk_tk_grid_variation(z, name, varying_dictionary, derived_params=None, kmin=-4, kmax=3, nb_points=2000, output=True)

return Power Spectra with k in h.Mpc-1 and P in h3.Mpc-3 in the convention of CAMB and CLASS (the one needed for MP-Gadget). return Transfer function in the convention of CLASS (needed for MP-Genic) or CAMB (needed for 2LPT) !!! note : . CLASS format for MP-GENIC, CAMB format for 2LPT, CAMB format for COSMICIC . for 2LPT & COSMICIC: no header! k vector need to be equidistant in log-k space and length of pk must be the same than tk

Source code in lyapower/CLASS.py

def write_pk_tk_grid_variation(
    self,
    z,
    name,
    varying_dictionary,
    derived_params=None,
    kmin=-4,
    kmax=3,
    nb_points=2000,
    output=True,
):
    """return Power Spectra with k in h.Mpc-1 and P in h3.Mpc-3 in the convention of CAMB and CLASS (the one needed for MP-Gadget).
    return Transfer function in the convention of CLASS (needed for MP-Genic) or CAMB (needed for 2LPT) !!!
    note : . CLASS format for MP-GENIC, CAMB format for 2LPT, CAMB format for COSMICIC
           . for 2LPT & COSMICIC: no header! k vector need to be equidistant in log-k space and length of pk must be the same than tk
    """
    dict_deriv = {}
    if derived_params is not None:
        for i in range(len(derived_params)):
            dict_deriv[derived_params[i]] = []
    list_key = list(varying_dictionary.keys())
    for i in range(len(varying_dictionary[list_key[0]])):
        dict_to_set = {
            list_key[j]: varying_dictionary[list_key[j]][i]
            for j in range(len(list_key))
        }
        self.model.set(dict_to_set)
        self.model.compute()
        Power, sigma_8 = self.write_pk(
            name, z, kmin=kmin, kmax=kmax, nb_points=nb_points, output= output
        )
        Transfer = self.write_tk(name, z,output=output)
        if derived_params is not None:
            par = self.model.get_current_derived_parameters(derived_params)
            print(par.keys())
            for key in par.keys():
                dict_deriv[key].append(par[key])
    return (Power, Transfer, sigma_8, dict_deriv)

write_pk_tk_neutrino_mass ¶

write_pk_tk_neutrino_mass(z, name, kmin=-4, kmax=3, nb_points=2000, nb_neutrinos=3, nb_neutrinos_massif=1, neutrino_mass=None, method='cbnu', output=True)

return Power Spectra with k in h.Mpc-1 and P in h3.Mpc-3 in the convention of CAMB and CLASS (the one needed for MP-Gadget). return Transfer function in the convention of CLASS (needed for MP-Genic) or CAMB (needed for 2LPT) !!!

Source code in lyapower/CLASS.py

def write_pk_tk_neutrino_mass(
    self,
    z,
    name,
    kmin=-4,
    kmax=3,
    nb_points=2000,
    nb_neutrinos=3,
    nb_neutrinos_massif=1,
    neutrino_mass=None,
    method="cbnu",
    output=True,
):
    """return Power Spectra with k in h.Mpc-1 and P in h3.Mpc-3 in the convention of CAMB and CLASS (the one needed for MP-Gadget).
    return Transfer function in the convention of CLASS (needed for MP-Genic) or CAMB (needed for 2LPT) !!!
    """
    m_ncdm = ""
    Omega_ncdm = ""
    for j in range(nb_neutrinos):
        m_ncdm = m_ncdm + "0.0,"
        Omega_ncdm = Omega_ncdm + "0.0,"
    m_ncdm = m_ncdm[0:-1]
    Omega_ncdm = Omega_ncdm[0:-1]
    self.model.set(
        {"N_ncdm": nb_neutrinos, "m_ncdm": m_ncdm, "Omega_ncdm": Omega_ncdm}
    )
    self.model.compute()
    if neutrino_mass is not None:
        omega_cdm_ref = self.model.Omega0_cdm()
        omega_nu_ref = self.model.Omega_nu
        m_ncdm = ""
        for j in range(nb_neutrinos_massif):
            m_ncdm = m_ncdm + str(neutrino_mass / nb_neutrinos_massif) + ","
        for j in range(nb_neutrinos - nb_neutrinos_massif):
            m_ncdm = m_ncdm + "0.0,"
        print(m_ncdm)
        m_ncdm = m_ncdm[0:-1]
        if method == "cbnu":
            self.model.set(
                {
                    "N_ncdm": nb_neutrinos,
                    "m_ncdm": m_ncdm,
                    "omega_cdm": str(
                        (
                            omega_cdm_ref
                            - self.get_omeganu_from_mass(neutrino_mass)
                            + omega_nu_ref
                        )
                        * self.model.h() ** 2
                    ),
                }
            )
        if method == "cb":
            self.model.set({"N_ncdm": nb_neutrinos, "m_ncdm": m_ncdm})
        self.model.compute()

    Power, sigma_8 = self.write_pk(
        name, z, kmin=kmin, kmax=kmax, nb_points=nb_points, output=output
    )
    Transfer = self.write_tk(name, z,output=output)
    return (Power, Transfer, sigma_8)

write_pk ¶

write_pk(name, z, kmin=-4, kmax=3, nb_points=2000, output=True, header_output=True, verbose=True)

Compute and (optionally) write the linear power spectrum P(k).

Parameters:

Name	Type	Description	Default
`name`	`str`	Output file base name.	required
`z`	`float`	Redshift.	required
`kmin, kmax`	`float`	log10 wavenumber bounds (h/Mpc).	required
`nb_points`	`int`	Number of k samples.	`2000`
`output`	`bool`	Write the `pk_{name}_z{z}.dat` file.	`True`
`header_output`	`bool`	Include a descriptive header.	`True`
`verbose`	`bool`	Print the cosmology summary.	`True`

Returns:

Name	Type	Description
`tuple`		`(Power, sigma_8)` — an `(N, 2)` `[k, P(k)]` array and
		sigma_8.

Source code in lyapower/CLASS.py

def write_pk(
    self, name, z, kmin=-4, kmax=3, nb_points=2000, output=True, header_output=True, verbose=True
):
    """Compute and (optionally) write the linear power spectrum P(k).

    Args:
        name (str): Output file base name.
        z (float): Redshift.
        kmin, kmax (float, optional): log10 wavenumber bounds (h/Mpc).
        nb_points (int, optional): Number of k samples.
        output (bool, optional): Write the ``pk_{name}_z{z}.dat`` file.
        header_output (bool, optional): Include a descriptive header.
        verbose (bool, optional): Print the cosmology summary.

    Returns:
        tuple: ``(Power, sigma_8)`` — an ``(N, 2)`` ``[k, P(k)]`` array and
        sigma_8.
    """
    if self.output_format == "class":
        (k_space, Pk, sigma_8) = self.compute_power_spectrum(
            z, kmin=kmin, kmax=kmax, nb_points=nb_points, verbose=verbose
        )  # k in h/Mpc
        header = header_pk_class.format(z, kmin, kmax, nb_points)
    elif self.output_format == "camb":
        tr = self.model.get_transfer(z, output_format="camb")
        k = np.array(tr["k (h/Mpc)"])
        kpower = np.logspace(
            *(np.log10(k[[0, -1]])) * (1 - 1e-5), len(k)
        )  # k in h/Mpc
        (k_space, Pk, sigma_8) = self.compute_power_spectrum(
            z, k_array=kpower, verbose=verbose
        )
        header = header_pk_camb.format(
            z, np.min(k_space), np.max(k_space), len(k_space)
        )
    Power = np.stack([k_space, Pk], axis=1)
    if output:
        if header_output:
            np.savetxt("pk_{}_z{}.dat".format(name, z), Power, header=header)
        else:
            np.savetxt("pk_{}_z{}.dat".format(name, z), Power)
    return (Power, sigma_8)

write_tk ¶

write_tk(name, z, output=True, header_output=True)

Compute and (optionally) write the transfer functions T_i(k).

Parameters:

Name	Type	Description	Default
`name`	`str`	Output file base name.	required
`z`	`float`	Redshift.	required
`output`	`bool`	Write the `tk_{name}_z{z}.dat` file.	`True`
`header_output`	`bool`	Include a descriptive header.	`True`

Returns:

Type	Description
	numpy.ndarray: The transfer-function table (CLASS or CAMB format).

Source code in lyapower/CLASS.py

def write_tk(self, name, z, output=True, header_output=True):
    """Compute and (optionally) write the transfer functions T_i(k).

    Args:
        name (str): Output file base name.
        z (float): Redshift.
        output (bool, optional): Write the ``tk_{name}_z{z}.dat`` file.
        header_output (bool, optional): Include a descriptive header.

    Returns:
        numpy.ndarray: The transfer-function table (CLASS or CAMB format).
    """
    if self.output_format == "class":
        Transfer_dict = self.model.get_transfer(z, output_format="class")
        Transfer = np.stack([Tk for key, Tk in Transfer_dict.items()], axis=1)
        header = header_tk_class.format(
            z, np.min(Transfer[:, 0]), np.max(Transfer[:, 0]), len(Transfer[:, 0])
        )
    elif self.output_format == "camb":
        Transfer_dict = self.model.get_transfer(z, output_format="camb")
        Transfer = np.stack([Tk for key, Tk in Transfer_dict.items()], axis=1)
        Transfer = self.interp_to_equidistant_log_space_camb_format(Transfer)
        header = header_tk_camb.format(
            z, np.min(Transfer[:, 0]), np.max(Transfer[:, 0]), len(Transfer[:, 0])
        )
    list_key = [key for key, Tk in Transfer_dict.items()]
    head = ""
    for i in range(len(list_key)):
        head = head + "            " + str(i + 1) + ":" + list_key[i]
    header = header + head
    if output:
        if header_output:
            np.savetxt(
                "tk_{}_z{}.dat".format(name, z),
                Transfer,
                header=header,
                delimiter="    ",
            )
        else:
            np.savetxt("tk_{}_z{}.dat".format(name, z), Transfer, delimiter="    ")
    return Transfer

compute_power_spectrum ¶

compute_power_spectrum(z, kmin=-4, kmax=3, nb_points=2000, k_array=None, verbose=True)

Evaluate the linear power spectrum on a k-grid (h-normalised).

Parameters:

Name	Type	Description	Default
`z`	`float`	Redshift.	required
`kmin, kmax`	`float`	log10 wavenumber bounds (h/Mpc).	required
`nb_points`	`int`	Number of k samples (if `k_array` None).	`2000`
`k_array`	`ndarray`	Explicit k grid (h/Mpc).	`None`
`verbose`	`bool`	Print the cosmology summary.	`True`

Returns:

Name	Type	Description
`tuple`		`(k_space, Pk, sigma_8)` with `Pk` in `(h^-1 Mpc)^3`.

Source code in lyapower/CLASS.py

def compute_power_spectrum(
    self, z, kmin=-4, kmax=3, nb_points=2000, k_array=None, verbose=True
):
    """Evaluate the linear power spectrum on a k-grid (h-normalised).

    Args:
        z (float): Redshift.
        kmin, kmax (float, optional): log10 wavenumber bounds (h/Mpc).
        nb_points (int, optional): Number of k samples (if ``k_array`` None).
        k_array (numpy.ndarray, optional): Explicit k grid (h/Mpc).
        verbose (bool, optional): Print the cosmology summary.

    Returns:
        tuple: ``(k_space, Pk, sigma_8)`` with ``Pk`` in ``(h^-1 Mpc)^3``.
    """
    if k_array is None:
        k_space = np.logspace(kmin, kmax, num=nb_points)
    else:
        k_space = k_array
    sigma_8 = self.model.sigma(8 / self.model.h(), z)
    if verbose:
        print("Omega matter = " + str(self.model.Omega_m()))
        print("Omega lambda = " + str(self.model.Omega_Lambda()))
        print("Omega baryon = " + str(self.model.Omega_b()))
        print("Omega dm = " + str(self.model.Omega0_cdm()))
        print("Omega k = " + str(self.model.Omega0_k()))
        print(
            "Omega sum = "
            + str(
                self.model.Omega0_cdm() + self.model.Omega_b() + self.model.Omega_nu
            )
        )
        print("Omega nu = " + str(self.model.Omega_nu))
        print("Omega rad = " + str(self.model.Omega_g()))
        print("sigma 8 at (z={}) = {}".format(z, sigma_8))
    Pk = []
    h = self.model.h()
    for k in k_space:
        Pk.append(self.model.pk(k * h, z) * h**3)
    return (k_space, Pk, sigma_8)

interp_to_equidistant_log_space_camb_format ¶

interp_to_equidistant_log_space_camb_format(Transfer)

Resample a transfer-function table onto a log-equidistant k grid.

Required by the CAMB output format (2LPT / CosmicIC expect equidistant log-k).

Parameters:

Name	Type	Description	Default
`Transfer`	`ndarray`	Transfer table (first column is k).	required

Returns:

Type	Description
	numpy.ndarray: The resampled table.

Source code in lyapower/CLASS.py

def interp_to_equidistant_log_space_camb_format(self, Transfer):
    """Resample a transfer-function table onto a log-equidistant k grid.

    Required by the CAMB output format (2LPT / CosmicIC expect equidistant
    log-k).

    Args:
        Transfer (numpy.ndarray): Transfer table (first column is k).

    Returns:
        numpy.ndarray: The resampled table.
    """
    k_space = Transfer[:, 0]
    k_log_space = np.logspace(
        *(np.log10(k_space[[0, -1]])) * (1 - 1e-5), len(k_space)
    )  # k in h/Mpc
    for i in range(1, Transfer.shape[-1]):
        interp = interp1d(k_space, Transfer[:, i], bounds_error=True)
        Transfer[:, i] = interp(k_log_space)
    Transfer[:, 0] = k_log_space
    return Transfer

sigmaR_all_species_at ¶

sigmaR_all_species_at(R, z)

Compute sigma(R, z) over all species.

Parameters:

Name	Type	Description	Default
`R`	`float`	Smoothing radius (Mpc/h).	required
`z`	`float`	Redshift.	required

Returns:

Name	Type	Description
`float`		The RMS density fluctuation sigma(R, z).

Source code in lyapower/CLASS.py

def sigmaR_all_species_at(self, R, z):
    """Compute sigma(R, z) over all species.

    Args:
        R (float): Smoothing radius (Mpc/h).
        z (float): Redshift.

    Returns:
        float: The RMS density fluctuation sigma(R, z).
    """
    self.model.compute()
    return self.model.sigma(R / self.model.h(), z)

free_structure ¶

free_structure()

Free the CLASS internal structures (struct_cleanup).

Source code in lyapower/CLASS.py

def free_structure(self):
    """Free the CLASS internal structures (``struct_cleanup``)."""
    self.model.struct_cleanup()

close_class ¶

close_class()

Empty the CLASS model (release its parameters).

Source code in lyapower/CLASS.py

def close_class(self):
    """Empty the CLASS model (release its parameters)."""
    self.model.empty()

get_omeganu_from_mass ¶

get_omeganu_from_mass(mass)

Convert a summed neutrino mass (eV) to omega_nu = m / (93.14 h^2).

Parameters:

Name	Type	Description	Default
`mass`	`float`	Summed neutrino mass (eV).	required

Returns:

Name	Type	Description
`float`		`omega_nu` (physical density parameter).

Source code in lyapower/CLASS.py

def get_omeganu_from_mass(self, mass):
    """Convert a summed neutrino mass (eV) to ``omega_nu = m / (93.14 h^2)``.

    Args:
        mass (float): Summed neutrino mass (eV).

    Returns:
        float: ``omega_nu`` (physical density parameter).
    """
    return mass / (93.14 * (self.model.h()) ** 2)

get_mass_from_omeganu ¶

get_mass_from_omeganu(omeganu)

Convert omega_nu back to a summed neutrino mass (eV).

Parameters:

Name	Type	Description	Default
`omeganu`	`float`	Physical neutrino density parameter.	required

Returns:

Name	Type	Description
`float`		Summed neutrino mass (eV).

Source code in lyapower/CLASS.py

def get_mass_from_omeganu(self, omeganu):
    """Convert ``omega_nu`` back to a summed neutrino mass (eV).

    Args:
        omeganu (float): Physical neutrino density parameter.

    Returns:
        float: Summed neutrino mass (eV).
    """
    return omeganu * (93.14) * (self.model.h()) ** 2

CosmoprimoInterface ¶

Linear power spectra via the cosmoprimo façade (CLASS engine).

Source code in lyapower/CLASS.py

class CosmoprimoInterface:
    """Linear power spectra via the cosmoprimo façade (CLASS engine)."""

    def __init__(self, pwd, settings):
        """Store the working directory and cosmology settings.

        Args:
            pwd (str): Working directory.
            settings (dict): Cosmology parameters for cosmoprimo.
        """
        self.pwd = pwd
        self.settings = settings

    def Pl_class_cosmoprimo(self, k_array, z):
        """Linear matter power spectrum from cosmoprimo (CLASS engine).

        Args:
            k_array (numpy.ndarray): Wavenumbers (h/Mpc).
            z (float): Redshift.

        Returns:
            numpy.ndarray: ``P(k, z)`` in ``(h^-1 Mpc)^3``.
        """
        cosmo = Cosmology(engine="class", **self.settings)
        fo = Fourier(cosmo, engine="class")
        pk = fo.pk_interpolator()
        return pk(k_array, z=z)

    def Pl_class_cosmoprimo_no_wiggle(self, k_array, z):
        """No-wiggle (BAO-smoothed) linear power spectrum from cosmoprimo.

        Args:
            k_array (numpy.ndarray): Wavenumbers (h/Mpc).
            z (float): Redshift.

        Returns:
            numpy.ndarray: The BAO-filtered ``P(k, z)`` in ``(h^-1 Mpc)^3``.
        """
        cosmo = Cosmology(engine="class", **self.settings)
        fo = Fourier(cosmo, engine="class")
        pk = fo.pk_interpolator()
        pknow = PowerSpectrumBAOFilter(
            pk, engine="wallish2018"
        ).smooth_pk_interpolator()
        return pknow(k_array, z=z)

Pl_class_cosmoprimo ¶

Pl_class_cosmoprimo(k_array, z)

Linear matter power spectrum from cosmoprimo (CLASS engine).

Parameters:

Name	Type	Description	Default
`k_array`	`ndarray`	Wavenumbers (h/Mpc).	required
`z`	`float`	Redshift.	required

Returns:

Type	Description
	numpy.ndarray: `P(k, z)` in `(h^-1 Mpc)^3`.

Source code in lyapower/CLASS.py

def Pl_class_cosmoprimo(self, k_array, z):
    """Linear matter power spectrum from cosmoprimo (CLASS engine).

    Args:
        k_array (numpy.ndarray): Wavenumbers (h/Mpc).
        z (float): Redshift.

    Returns:
        numpy.ndarray: ``P(k, z)`` in ``(h^-1 Mpc)^3``.
    """
    cosmo = Cosmology(engine="class", **self.settings)
    fo = Fourier(cosmo, engine="class")
    pk = fo.pk_interpolator()
    return pk(k_array, z=z)

Pl_class_cosmoprimo_no_wiggle ¶

Pl_class_cosmoprimo_no_wiggle(k_array, z)

No-wiggle (BAO-smoothed) linear power spectrum from cosmoprimo.

Parameters:

Name	Type	Description	Default
`k_array`	`ndarray`	Wavenumbers (h/Mpc).	required
`z`	`float`	Redshift.	required

Returns:

Type	Description
	numpy.ndarray: The BAO-filtered `P(k, z)` in `(h^-1 Mpc)^3`.

Source code in lyapower/CLASS.py

def Pl_class_cosmoprimo_no_wiggle(self, k_array, z):
    """No-wiggle (BAO-smoothed) linear power spectrum from cosmoprimo.

    Args:
        k_array (numpy.ndarray): Wavenumbers (h/Mpc).
        z (float): Redshift.

    Returns:
        numpy.ndarray: The BAO-filtered ``P(k, z)`` in ``(h^-1 Mpc)^3``.
    """
    cosmo = Cosmology(engine="class", **self.settings)
    fo = Fourier(cosmo, engine="class")
    pk = fo.pk_interpolator()
    pknow = PowerSpectrumBAOFilter(
        pk, engine="wallish2018"
    ).smooth_pk_interpolator()
    return pknow(k_array, z=z)

lyapower.CLASS¶

lyapower.CLASS ¶

MyClass ¶

define_output_format ¶

write_pk_tk ¶

write_pk_tk_grid_variation ¶

write_pk_tk_neutrino_mass ¶

write_pk ¶

write_tk ¶

compute_power_spectrum ¶

interp_to_equidistant_log_space_camb_format ¶

sigmaR_all_species_at ¶

free_structure ¶

close_class ¶

get_omeganu_from_mass ¶

get_mass_from_omeganu ¶

CosmoprimoInterface ¶

Pl_class_cosmoprimo ¶

Pl_class_cosmoprimo_no_wiggle ¶

`lyapower.CLASS`¶