CLASS MANUAL
|
#include "spectra.h"
Functions | |
int | spectra_cl_at_l (struct spectra *psp, double l, double *cl_tot, double **cl_md, double **cl_md_ic) |
int | spectra_pk_at_z (struct background *pba, struct spectra *psp, enum linear_or_logarithmic mode, double z, double *output_tot, double *output_ic) |
int | spectra_pk_at_k_and_z (struct background *pba, struct primordial *ppm, struct spectra *psp, double k, double z, double *pk_tot, double *pk_ic) |
int | spectra_pk_nl_at_z (struct background *pba, struct spectra *psp, enum linear_or_logarithmic mode, double z, double *output_tot) |
int | spectra_pk_nl_at_k_and_z (struct background *pba, struct primordial *ppm, struct spectra *psp, double k, double z, double *pk_tot) |
int | spectra_tk_at_z (struct background *pba, struct spectra *psp, double z, double *output) |
int | spectra_tk_at_k_and_z (struct background *pba, struct spectra *psp, double k, double z, double *output) |
int | spectra_init (struct precision *ppr, struct background *pba, struct perturbs *ppt, struct primordial *ppm, struct nonlinear *pnl, struct transfers *ptr, struct spectra *psp) |
int | spectra_free (struct spectra *psp) |
int | spectra_indices (struct background *pba, struct perturbs *ppt, struct transfers *ptr, struct primordial *ppm, struct spectra *psp) |
int | spectra_cls (struct background *pba, struct perturbs *ppt, struct transfers *ptr, struct primordial *ppm, struct spectra *psp) |
int | spectra_compute_cl (struct background *pba, struct perturbs *ppt, struct transfers *ptr, struct primordial *ppm, struct spectra *psp, int index_md, int index_ic1, int index_ic2, int index_l, int cl_integrand_num_columns, double *cl_integrand, double *primordial_pk, double *transfer_ic1, double *transfer_ic2) |
int | spectra_k_and_tau (struct background *pba, struct perturbs *ppt, struct spectra *psp) |
int | spectra_pk (struct background *pba, struct perturbs *ppt, struct primordial *ppm, struct nonlinear *pnl, struct spectra *psp) |
int | spectra_sigma (struct background *pba, struct primordial *ppm, struct spectra *psp, double R, double z, double *sigma) |
int | spectra_matter_transfers (struct background *pba, struct perturbs *ppt, struct spectra *psp) |
int | spectra_output_tk_data (struct background *pba, struct perturbs *ppt, struct spectra *psp, enum file_format output_format, double z, int number_of_titles, double *data) |
Documented spectra module
Julien Lesgourgues, 25.08.2010
This module computes the anisotropy and Fourier power spectra 's given the transfer and Bessel functions (for anisotropy spectra), the source functions (for Fourier spectra) and the primordial spectra.
The following functions can be called from other modules:
int spectra_cl_at_l | ( | struct spectra * | psp, |
double | l, | ||
double * | cl_tot, | ||
double ** | cl_md, | ||
double ** | cl_md_ic | ||
) |
Anisotropy power spectra 's for all types, modes and initial conditions.
This routine evaluates all the 's at a given value of l by interpolating in the pre-computed table. When relevant, it also sums over all initial conditions for each mode, and over all modes.
This function can be called from whatever module at whatever time, provided that spectra_init() has been called before, and spectra_free() has not been called yet.
psp | Input: pointer to spectra structure (containing pre-computed table) |
l | Input: multipole number |
cl_tot | Output: total ![]() |
cl_md | Output: ![]() |
cl_md_ic | Output: ![]() |
Summary:
int spectra_pk_at_z | ( | struct background * | pba, |
struct spectra * | psp, | ||
enum linear_or_logarithmic | mode, | ||
double | z, | ||
double * | output_tot, | ||
double * | output_ic | ||
) |
Matter power spectrum for arbitrary redshift and for all initial conditions.
This routine evaluates the matter power spectrum at a given value of z by interpolating in the pre-computed table (if several values of z have been stored) or by directly reading it (if it only contains values at z=0 and we want P(k,z=0))
Can be called in two modes: linear or logarithmic.
One little subtlety: in case of several correlated initial conditions, the cross-correlation spectrum can be negative. Then, in logarithmic mode, the non-diagonal elements contain the cross-correlation angle (from -1 to 1) instead of
This function can be called from whatever module at whatever time, provided that spectra_init() has been called before, and spectra_free() has not been called yet.
pba | Input: pointer to background structure (used for converting z into tau) |
psp | Input: pointer to spectra structure (containing pre-computed table) |
mode | Input: linear or logarithmic |
z | Input: redshift |
output_tot | Output: total matter power spectrum P(k) in ![]() |
output_ic | Output: for each pair of initial conditions, matter power spectra P(k) in ![]() |
Summary:
int spectra_pk_at_k_and_z | ( | struct background * | pba, |
struct primordial * | ppm, | ||
struct spectra * | psp, | ||
double | k, | ||
double | z, | ||
double * | pk_tot, | ||
double * | pk_ic | ||
) |
Matter power spectrum for arbitrary wavenumber, redshift and initial condition.
This routine evaluates the matter power spectrum at a given value of k and z by interpolating in a table of all P(k)'s computed at this z by spectra_pk_at_z() (when kmin <= k <= kmax), or eventually by using directly the primordial spectrum (when 0 <= k < kmin): the latter case is an approximation, valid when kmin << comoving Hubble scale today. Returns zero when k=0. Returns an error when k<0 or k > kmax.
This function can be called from whatever module at whatever time, provided that spectra_init() has been called before, and spectra_free() has not been called yet.
pba | Input: pointer to background structure (used for converting z into tau) |
ppm | Input: pointer to primordial structure (used only in the case 0 < k < kmin) |
psp | Input: pointer to spectra structure (containing pre-computed table) |
k | Input: wavenumber in 1/Mpc |
z | Input: redshift |
pk_tot | Output: total matter power spectrum P(k) in ![]() |
pk_ic | Output: for each pair of initial conditions, matter power spectra P(k) in ![]() |
Summary:
int spectra_pk_nl_at_z | ( | struct background * | pba, |
struct spectra * | psp, | ||
enum linear_or_logarithmic | mode, | ||
double | z, | ||
double * | output_tot | ||
) |
Non-linear total matter power spectrum for arbitrary redshift.
This routine evaluates the non-linear matter power spectrum at a given value of z by interpolating in the pre-computed table (if several values of z have been stored) or by directly reading it (if it only contains values at z=0 and we want P(k,z=0))
Can be called in two modes: linear or logarithmic.
This function can be called from whatever module at whatever time, provided that spectra_init() has been called before, and spectra_free() has not been called yet.
pba | Input: pointer to background structure (used for converting z into tau) |
psp | Input: pointer to spectra structure (containing pre-computed table) |
mode | Input: linear or logarithmic |
z | Input: redshift |
output_tot | Output: total matter power spectrum P(k) in ![]() |
Summary:
int spectra_pk_nl_at_k_and_z | ( | struct background * | pba, |
struct primordial * | ppm, | ||
struct spectra * | psp, | ||
double | k, | ||
double | z, | ||
double * | pk_tot | ||
) |
Non-linear total matter power spectrum for arbitrary wavenumber and redshift.
This routine evaluates the matter power spectrum at a given value of k and z by interpolating in a table of all P(k)'s computed at this z by spectra_pk_nl_at_z() (when kmin <= k <= kmax), or eventually by using directly the primordial spectrum (when 0 <= k < kmin): the latter case is an approximation, valid when kmin << comoving Hubble scale today. Returns zero when k=0. Returns an error when k<0 or k > kmax.
This function can be called from whatever module at whatever time, provided that spectra_init() has been called before, and spectra_free() has not been called yet.
pba | Input: pointer to background structure (used for converting z into tau) |
ppm | Input: pointer to primordial structure (used only in the case 0 < k < kmin) |
psp | Input: pointer to spectra structure (containing pre-computed table) |
k | Input: wavenumber in 1/Mpc |
z | Input: redshift |
pk_tot | Output: total matter power spectrum P(k) in ![]() |
Summary:
int spectra_tk_at_z | ( | struct background * | pba, |
struct spectra * | psp, | ||
double | z, | ||
double * | output | ||
) |
Matter transfer functions for arbitrary redshift and for all initial conditions.
This routine evaluates the matter transfer functions at a given value of z by interpolating in the pre-computed table (if several values of z have been stored) or by directly reading it (if it only contains values at z=0 and we want )
This function can be called from whatever module at whatever time, provided that spectra_init() has been called before, and spectra_free() has not been called yet.
pba | Input: pointer to background structure (used for converting z into tau) |
psp | Input: pointer to spectra structure (containing pre-computed table) |
z | Input: redshift |
output | Output: matter transfer functions |
Summary:
int spectra_tk_at_k_and_z | ( | struct background * | pba, |
struct spectra * | psp, | ||
double | k, | ||
double | z, | ||
double * | output | ||
) |
Matter transfer functions for arbitrary wavenumber, redshift and initial condition.
This routine evaluates the matter transfer functions at a given value of k and z by interpolating in a table of all 's computed at this z by spectra_tk_at_z() (when kmin <= k <= kmax). Returns an error when k<kmin or k > kmax.
This function can be called from whatever module at whatever time, provided that spectra_init() has been called before, and spectra_free() has not been called yet.
pba | Input: pointer to background structure (used for converting z into tau) |
psp | Input: pointer to spectra structure (containing pre-computed table) |
k | Input: wavenumber in 1/Mpc |
z | Input: redshift |
output | Output: matter transfer functions |
Summary:
int spectra_init | ( | struct precision * | ppr, |
struct background * | pba, | ||
struct perturbs * | ppt, | ||
struct primordial * | ppm, | ||
struct nonlinear * | pnl, | ||
struct transfers * | ptr, | ||
struct spectra * | psp | ||
) |
This routine initializes the spectra structure (in particular, computes table of anisotropy and Fourier spectra )
ppr | Input: pointer to precision structure |
pba | Input: pointer to background structure (will provide H, Omega_m at redshift of interest) |
ppt | Input: pointer to perturbation structure |
ptr | Input: pointer to transfer structure |
ppm | Input: pointer to primordial structure |
pnl | Input: pointer to nonlinear structure |
psp | Output: pointer to initialized spectra structure |
Summary:
int spectra_free | ( | struct spectra * | psp | ) |
This routine frees all the memory space allocated by spectra_init().
To be called at the end of each run, only when no further calls to spectra_cls_at_l(), spectra_pk_at_z(), spectra_pk_at_k_and_z() are needed.
psp | Input: pointer to spectra structure (which fields must be freed) |
int spectra_indices | ( | struct background * | pba, |
struct perturbs * | ppt, | ||
struct transfers * | ptr, | ||
struct primordial * | ppm, | ||
struct spectra * | psp | ||
) |
This routine defines indices and allocates tables in the spectra structure
pba | Input: pointer to background structure |
ppt | Input: pointer to perturbation structure |
ptr | Input: pointer to transfers structure |
ppm | Input: pointer to primordial structure |
psp | Input/output: pointer to spectra structure |
int spectra_cls | ( | struct background * | pba, |
struct perturbs * | ppt, | ||
struct transfers * | ptr, | ||
struct primordial * | ppm, | ||
struct spectra * | psp | ||
) |
This routine computes a table of values for all harmonic spectra 's, given the transfer functions and primordial spectra.
pba | Input: pointer to background structure |
ppt | Input: pointer to perturbation structure |
ptr | Input: pointer to transfers structure |
ppm | Input: pointer to primordial structure |
psp | Input/Output: pointer to spectra structure |
Summary:
int spectra_compute_cl | ( | struct background * | pba, |
struct perturbs * | ppt, | ||
struct transfers * | ptr, | ||
struct primordial * | ppm, | ||
struct spectra * | psp, | ||
int | index_md, | ||
int | index_ic1, | ||
int | index_ic2, | ||
int | index_l, | ||
int | cl_integrand_num_columns, | ||
double * | cl_integrand, | ||
double * | primordial_pk, | ||
double * | transfer_ic1, | ||
double * | transfer_ic2 | ||
) |
This routine computes the 's for a given mode, pair of initial conditions and multipole, but for all types (TT, TE...), by convolving the transfer functions with the primordial spectra.
pba | Input: pointer to background structure |
ppt | Input: pointer to perturbation structure |
ptr | Input: pointer to transfers structure |
ppm | Input: pointer to primordial structure |
psp | Input/Output: pointer to spectra structure (result stored here) |
index_md | Input: index of mode under consideration |
index_ic1 | Input: index of first initial condition in the correlator |
index_ic2 | Input: index of second initial condition in the correlator |
index_l | Input: index of multipole under consideration |
cl_integrand_num_columns | Input: number of columns in cl_integrand |
cl_integrand | Input: an allocated workspace |
primordial_pk | Input: table of primordial spectrum values |
transfer_ic1 | Input: table of transfer function values for first initial condition |
transfer_ic2 | Input: table of transfer function values for second initial condition |
int spectra_k_and_tau | ( | struct background * | pba, |
struct perturbs * | ppt, | ||
struct spectra * | psp | ||
) |
This routine computes the values of k and tau at which the matter power spectra and the matter transfer functions
will be stored.
pba | Input: pointer to background structure (for z to tau conversion) |
ppt | Input: pointer to perturbation structure (contain source functions) |
psp | Input/Output: pointer to spectra structure |
Summary:
int spectra_pk | ( | struct background * | pba, |
struct perturbs * | ppt, | ||
struct primordial * | ppm, | ||
struct nonlinear * | pnl, | ||
struct spectra * | psp | ||
) |
This routine computes a table of values for all matter power spectra P(k), given the source functions and primordial spectra.
pba | Input: pointer to background structure (will provide H, Omega_m at redshift of interest) |
ppt | Input: pointer to perturbation structure (contain source functions) |
ppm | Input: pointer to primordial structure |
pnl | Input: pointer to nonlinear structure |
psp | Input/Output: pointer to spectra structure |
Summary:
int spectra_sigma | ( | struct background * | pba, |
struct primordial * | ppm, | ||
struct spectra * | psp, | ||
double | R, | ||
double | z, | ||
double * | sigma | ||
) |
This routine computes sigma(R) given P(k) (does not check that k_max is large enough)
pba | Input: pointer to background structure |
ppm | Input: pointer to primordial structure |
psp | Input: pointer to spectra structure |
z | Input: redshift |
R | Input: radius in Mpc |
sigma | Output: variance in a sphere of radius R (dimensionless) |
int spectra_matter_transfers | ( | struct background * | pba, |
struct perturbs * | ppt, | ||
struct spectra * | psp | ||
) |
This routine computes a table of values for all matter power spectra P(k), given the source functions and primordial spectra.
pba | Input: pointer to background structure (will provide density of each species) |
ppt | Input: pointer to perturbation structure (contain source functions) |
psp | Input/Output: pointer to spectra structure |
Summary:
int spectra_output_tk_data | ( | struct background * | pba, |
struct perturbs * | ppt, | ||
struct spectra * | psp, | ||
enum file_format | output_format, | ||
double | z, | ||
int | number_of_titles, | ||
double * | data | ||
) |