CLASS MANUAL
|
#include "transfer.h"
Functions | |
int | transfer_functions_at_q (struct transfers *ptr, int index_md, int index_ic, int index_tt, int index_l, double q, double *transfer_function) |
int | transfer_init (struct precision *ppr, struct background *pba, struct thermo *pth, struct perturbs *ppt, struct nonlinear *pnl, struct transfers *ptr) |
int | transfer_free (struct transfers *ptr) |
int | transfer_indices_of_transfers (struct precision *ppr, struct perturbs *ppt, struct transfers *ptr, double q_period, double K, int sgnK) |
int | transfer_get_l_list (struct precision *ppr, struct perturbs *ppt, struct transfers *ptr) |
int | transfer_get_q_list (struct precision *ppr, struct perturbs *ppt, struct transfers *ptr, double q_period, double K, int sgnK) |
int | transfer_get_k_list (struct perturbs *ppt, struct transfers *ptr, double K) |
int | transfer_get_source_correspondence (struct perturbs *ppt, struct transfers *ptr, int **tp_of_tt) |
int | transfer_source_tau_size (struct precision *ppr, struct background *pba, struct perturbs *ppt, struct transfers *ptr, double tau_rec, double tau0, int index_md, int index_tt, int *tau_size) |
int | transfer_compute_for_each_q (struct precision *ppr, struct background *pba, struct perturbs *ppt, struct transfers *ptr, int **tp_of_tt, int index_q, int tau_size_max, double tau_rec, double ***pert_sources, double ***pert_sources_spline, struct transfer_workspace *ptw) |
int | transfer_interpolate_sources (struct perturbs *ppt, struct transfers *ptr, int index_q, int index_md, int index_ic, int index_type, double *pert_source, double *pert_source_spline, double *interpolated_sources) |
int | transfer_sources (struct precision *ppr, struct background *pba, struct perturbs *ppt, struct transfers *ptr, double *interpolated_sources, double tau_rec, int index_q, int index_md, int index_tt, double *sources, double *tau0_minus_tau, double *w_trapz, int *tau_size_out) |
int | transfer_selection_function (struct precision *ppr, struct perturbs *ppt, struct transfers *ptr, int bin, double z, double *selection) |
int | transfer_dNdz_analytic (struct transfers *ptr, double z, double *dNdz, double *dln_dNdz_dz) |
int | transfer_selection_sampling (struct precision *ppr, struct background *pba, struct perturbs *ppt, struct transfers *ptr, int bin, double *tau0_minus_tau, int tau_size) |
int | transfer_lensing_sampling (struct precision *ppr, struct background *pba, struct perturbs *ppt, struct transfers *ptr, int bin, double tau0, double *tau0_minus_tau, int tau_size) |
int | transfer_source_resample (struct precision *ppr, struct background *pba, struct perturbs *ppt, struct transfers *ptr, int bin, double *tau0_minus_tau, int tau_size, int index_md, double tau0, double *interpolated_sources, double *sources) |
int | transfer_selection_times (struct precision *ppr, struct background *pba, struct perturbs *ppt, struct transfers *ptr, int bin, double *tau_min, double *tau_mean, double *tau_max) |
int | transfer_selection_compute (struct precision *ppr, struct background *pba, struct perturbs *ppt, struct transfers *ptr, double *selection, double *tau0_minus_tau, double *w_trapz, int tau_size, double *pvecback, double tau0, int bin) |
int | transfer_compute_for_each_l (struct transfer_workspace *ptw, struct precision *ppr, struct perturbs *ppt, struct transfers *ptr, int index_q, int index_md, int index_ic, int index_tt, int index_l, double l, double q_max_bessel, radial_function_type radial_type) |
int | transfer_integrate (struct perturbs *ppt, struct transfers *ptr, struct transfer_workspace *ptw, int index_q, int index_md, int index_tt, double l, int index_l, double k, radial_function_type radial_type, double *trsf) |
int | transfer_limber (struct transfers *ptr, struct transfer_workspace *ptw, int index_md, int index_q, double l, double q, radial_function_type radial_type, double *trsf) |
int | transfer_limber_interpolate (struct transfers *ptr, double *tau0_minus_tau, double *sources, int tau_size, double tau0_minus_tau_limber, double *S) |
int | transfer_limber2 (int tau_size, struct transfers *ptr, int index_md, int index_k, double l, double k, double *tau0_minus_tau, double *sources, radial_function_type radial_type, double *trsf) |
Documented transfer module.
Julien Lesgourgues, 28.07.2013
This module has two purposes:
Hence the following functions can be called from other modules:
Note that in the standard implementation of CLASS, only the pre-computed values of the transfer functions are used, no interpolation is necessary; hence the routine transfer_functions_at_q() is actually never called.
int transfer_functions_at_q | ( | struct transfers * | ptr, |
int | index_md, | ||
int | index_ic, | ||
int | index_tt, | ||
int | index_l, | ||
double | q, | ||
double * | transfer_function | ||
) |
Transfer function at a given wavenumber q.
For a given mode (scalar, vector, tensor), initial condition, type (temperature, polarization, lensing, etc) and multipole, computes the transfer function for an arbitrary value of q by interpolating between pre-computed values of q. This function can be called from whatever module at whatever time, provided that transfer_init() has been called before, and transfer_free() has not been called yet.
Wavenumbers are called q in this module and k in the perturbation module. In flat universes k=q. In non-flat universes q and k differ through , where m=0,1,2 for scalar, vector, tensor. q should be used throughout the transfer module, excepted when interpolating or manipulating the source functions S(k,tau) calculated in the perturbation module: for a given value of q, this should be done at the corresponding k(q).
ptr | Input: pointer to transfer structure |
index_md | Input: index of requested mode |
index_ic | Input: index of requested initial condition |
index_tt | Input: index of requested type |
index_l | Input: index of requested multipole |
q | Input: any wavenumber |
transfer_function | Output: transfer function |
Summary:
int transfer_init | ( | struct precision * | ppr, |
struct background * | pba, | ||
struct thermo * | pth, | ||
struct perturbs * | ppt, | ||
struct nonlinear * | pnl, | ||
struct transfers * | ptr | ||
) |
This routine initializes the transfers structure, (in particular, computes table of transfer functions )
Main steps:
ppr | Input: pointer to precision structure |
pba | Input: pointer to background structure |
pth | Input: pointer to thermodynamics structure |
ppt | Input: pointer to perturbation structure |
pnl | Input: pointer to nonlinear structure |
ptr | Output: pointer to initialized transfers structure |
Summary:
int transfer_free | ( | struct transfers * | ptr | ) |
This routine frees all the memory space allocated by transfer_init().
To be called at the end of each run, only when no further calls to transfer_functions_at_k() are needed.
ptr | Input: pointer to transfers structure (which fields must be freed) |
int transfer_indices_of_transfers | ( | struct precision * | ppr, |
struct perturbs * | ppt, | ||
struct transfers * | ptr, | ||
double | q_period, | ||
double | K, | ||
int | sgnK | ||
) |
This routine defines all indices and allocates all tables in the transfers structure
Compute list of (k, l) values, allocate and fill corresponding arrays in the transfers structure. Allocate the array of transfer function tables.
ppr | Input: pointer to precision structure |
ppt | Input: pointer to perturbation structure |
ptr | Input/Output: pointer to transfer structure |
q_period | Input: order of magnitude of the oscillation period of transfer functions |
K | Input: spatial curvature (in absolute value) |
sgnK | Input: spatial curvature sign (open/closed/flat) |
Summary:
This routine defines the number and values of multipoles l for all modes.
ppr | Input: pointer to precision structure |
ppt | Input: pointer to perturbation structure |
ptr | Input/Output: pointer to transfers structure containing l's |
Summary:
int transfer_get_q_list | ( | struct precision * | ppr, |
struct perturbs * | ppt, | ||
struct transfers * | ptr, | ||
double | q_period, | ||
double | K, | ||
int | sgnK | ||
) |
This routine defines the number and values of wavenumbers q for each mode (goes smoothly from logarithmic step for small q's to linear step for large q's).
ppr | Input: pointer to precision structure |
ppt | Input: pointer to perturbation structure |
ptr | Input/Output: pointer to transfers structure containing q's |
q_period | Input: order of magnitude of the oscillation period of transfer functions |
K | Input: spatial curvature (in absolute value) |
sgnK | Input: spatial curvature sign (open/closed/flat) |
This routine infers from the q values a list of corresponding k values for each mode.
ppt | Input: pointer to perturbation structure |
ptr | Input/Output: pointer to transfers structure containing q's |
K | Input: spatial curvature |
int transfer_get_source_correspondence | ( | struct perturbs * | ppt, |
struct transfers * | ptr, | ||
int ** | tp_of_tt | ||
) |
This routine defines the correspondence between the sources in the perturbation and transfer module.
ppt | Input: pointer to perturbation structure |
ptr | Input: pointer to transfers structure containing l's |
tp_of_tt | Input/Output: array with the correspondence (allocated before, filled here) |
Summary:
int transfer_source_tau_size | ( | struct precision * | ppr, |
struct background * | pba, | ||
struct perturbs * | ppt, | ||
struct transfers * | ptr, | ||
double | tau_rec, | ||
double | tau0, | ||
int | index_md, | ||
int | index_tt, | ||
int * | tau_size | ||
) |
the code makes a distinction between "perturbation sources" (e.g. gravitational potential) and "transfer sources" (e.g. total density fluctuations, obtained through the Poisson equation, and observed with a given selection function).
This routine computes the number of sampled time values for each type of transfer sources.
ppr | Input: pointer to precision structure |
pba | Input: pointer to background structure |
ppt | Input: pointer to perturbation structure |
ptr | Input: pointer to transfers structure |
tau_rec | Input: recombination time |
tau0 | Input: time today |
index_md | Input: index of the mode (scalar, tensor) |
index_tt | Input: index of transfer type |
tau_size | Output: pointer to number of sampled times |
int transfer_compute_for_each_q | ( | struct precision * | ppr, |
struct background * | pba, | ||
struct perturbs * | ppt, | ||
struct transfers * | ptr, | ||
int ** | tp_of_tt, | ||
int | index_q, | ||
int | tau_size_max, | ||
double | tau_rec, | ||
double *** | pert_sources, | ||
double *** | pert_sources_spline, | ||
struct transfer_workspace * | ptw | ||
) |
Summary:
int transfer_interpolate_sources | ( | struct perturbs * | ppt, |
struct transfers * | ptr, | ||
int | index_q, | ||
int | index_md, | ||
int | index_ic, | ||
int | index_type, | ||
double * | pert_source, | ||
double * | pert_source_spline, | ||
double * | interpolated_sources | ||
) |
This routine interpolates sources for each mode, initial condition and type (of perturbation module), to get them at the right values of k, using the spline interpolation method.
ppt | Input: pointer to perturbation structure |
ptr | Input: pointer to transfers structure |
index_q | Input: index of wavenumber |
index_md | Input: index of mode |
index_ic | Input: index of initial condition |
index_type | Input: index of type of source (in perturbation module) |
pert_source | Input: array of sources |
pert_source_spline | Input: array of second derivative of sources |
interpolated_sources | Output: array of interpolated sources (filled here but allocated in transfer_init() to avoid numerous reallocation) |
Summary:
int transfer_sources | ( | struct precision * | ppr, |
struct background * | pba, | ||
struct perturbs * | ppt, | ||
struct transfers * | ptr, | ||
double * | interpolated_sources, | ||
double | tau_rec, | ||
int | index_q, | ||
int | index_md, | ||
int | index_tt, | ||
double * | sources, | ||
double * | tau0_minus_tau, | ||
double * | w_trapz, | ||
int * | tau_size_out | ||
) |
The code makes a distinction between "perturbation sources" (e.g. gravitational potential) and "transfer sources" (e.g. total density fluctuations, obtained through the Poisson equation, and observed with a given selection function).
This routine computes the transfer source given the interpolated perturbation source, and copies it in the workspace.
ppr | Input: pointer to precision structure |
pba | Input: pointer to background structure |
ppt | Input: pointer to perturbation structure |
ptr | Input: pointer to transfers structure |
interpolated_sources | Input: interpolated perturbation source |
tau_rec | Input: recombination time |
index_q | Input: index of wavenumber |
index_md | Input: index of mode |
index_tt | Input: index of type of (transfer) source |
sources | Output: transfer source |
tau0_minus_tau | Output: values of (tau0-tau) at which source are sample |
w_trapz | Output: trapezoidal weights for integration over tau |
tau_size_out | Output: pointer to size of previous two arrays, converted to double |
Summary:
int transfer_selection_function | ( | struct precision * | ppr, |
struct perturbs * | ppt, | ||
struct transfers * | ptr, | ||
int | bin, | ||
double | z, | ||
double * | selection | ||
) |
Arbitrarily normalized selection function dN/dz(z,bin)
ppr | Input: pointer to precision structure |
ppt | Input: pointer to perturbation structure |
ptr | Input: pointer to transfers structure |
bin | Input: redshift bin number |
z | Input: one value of redshift |
selection | Output: pointer to selection function |
int transfer_dNdz_analytic | ( | struct transfers * | ptr, |
double | z, | ||
double * | dNdz, | ||
double * | dln_dNdz_dz | ||
) |
Analytic form for dNdz distribution, from arXiv:1004.4640
ptr | Input: pointer to transfer structure |
z | Input: redshift |
dNdz | Output: density per redshift, dN/dZ |
dln_dNdz_dz | Output: dln(dN/dz)/dz, used optionally for the source evolution |
int transfer_selection_sampling | ( | struct precision * | ppr, |
struct background * | pba, | ||
struct perturbs * | ppt, | ||
struct transfers * | ptr, | ||
int | bin, | ||
double * | tau0_minus_tau, | ||
int | tau_size | ||
) |
For sources that need to be multiplied by a selection function, redefine a finer time sampling in a small range
ppr | Input: pointer to precision structure |
pba | Input: pointer to background structure |
ppt | Input: pointer to perturbation structure |
ptr | Input: pointer to transfers structure |
bin | Input: redshift bin number |
tau0_minus_tau | Output: values of (tau0-tau) at which source are sample |
tau_size | Output: pointer to size of previous array |
int transfer_lensing_sampling | ( | struct precision * | ppr, |
struct background * | pba, | ||
struct perturbs * | ppt, | ||
struct transfers * | ptr, | ||
int | bin, | ||
double | tau0, | ||
double * | tau0_minus_tau, | ||
int | tau_size | ||
) |
For lensing sources that need to be convolved with a selection function, redefine the sampling within the range extending from the tau_min of the selection function up to tau0
ppr | Input: pointer to precision structure |
pba | Input: pointer to background structure |
ppt | Input: pointer to perturbation structure |
ptr | Input: pointer to transfers structure |
bin | Input: redshift bin number |
tau0 | Input: time today |
tau0_minus_tau | Output: values of (tau0-tau) at which source are sample |
tau_size | Output: pointer to size of previous array |
int transfer_source_resample | ( | struct precision * | ppr, |
struct background * | pba, | ||
struct perturbs * | ppt, | ||
struct transfers * | ptr, | ||
int | bin, | ||
double * | tau0_minus_tau, | ||
int | tau_size, | ||
int | index_md, | ||
double | tau0, | ||
double * | interpolated_sources, | ||
double * | sources | ||
) |
For sources that need to be multiplied by a selection function, redefine a finer time sampling in a small range, and resample the perturbation sources at the new value by linear interpolation
ppr | Input: pointer to precision structure |
pba | Input: pointer to background structure |
ppt | Input: pointer to perturbation structure |
ptr | Input: pointer to transfers structure |
bin | Input: redshift bin number |
tau0_minus_tau | Output: values of (tau0-tau) at which source are sample |
tau_size | Output: pointer to size of previous array |
index_md | Input: index of mode |
tau0 | Input: time today |
interpolated_sources | Input: interpolated perturbation source |
sources | Output: resampled transfer source |
int transfer_selection_times | ( | struct precision * | ppr, |
struct background * | pba, | ||
struct perturbs * | ppt, | ||
struct transfers * | ptr, | ||
int | bin, | ||
double * | tau_min, | ||
double * | tau_mean, | ||
double * | tau_max | ||
) |
For each selection function, compute the min, mean and max values of conformal time (associated to the min, mean and max values of redshift specified by the user)
ppr | Input: pointer to precision structure |
pba | Input: pointer to background structure |
ppt | Input: pointer to perturbation structure |
ptr | Input: pointer to transfers structure |
bin | Input: redshift bin number |
tau_min | Output: smallest time in the selection interval |
tau_mean | Output: time corresponding to z_mean |
tau_max | Output: largest time in the selection interval |
int transfer_selection_compute | ( | struct precision * | ppr, |
struct background * | pba, | ||
struct perturbs * | ppt, | ||
struct transfers * | ptr, | ||
double * | selection, | ||
double * | tau0_minus_tau, | ||
double * | w_trapz, | ||
int | tau_size, | ||
double * | pvecback, | ||
double | tau0, | ||
int | bin | ||
) |
Compute and normalize selection function for a set of time values
ppr | Input: pointer to precision structure |
pba | Input: pointer to background structure |
ppt | Input: pointer to perturbation structure |
ptr | Input: pointer to transfers structure |
selection | Output: normalized selection function |
tau0_minus_tau | Input: values of (tau0-tau) at which source are sample |
w_trapz | Input: trapezoidal weights for integration over tau |
tau_size | Input: size of previous two arrays |
pvecback | Input: allocated array of background values |
tau0 | Input: time today |
bin | Input: redshift bin number |
int transfer_compute_for_each_l | ( | struct transfer_workspace * | ptw, |
struct precision * | ppr, | ||
struct perturbs * | ppt, | ||
struct transfers * | ptr, | ||
int | index_q, | ||
int | index_md, | ||
int | index_ic, | ||
int | index_tt, | ||
int | index_l, | ||
double | l, | ||
double | q_max_bessel, | ||
radial_function_type | radial_type | ||
) |
This routine computes the transfer functions ) as a function of wavenumber k for a given mode, initial condition, type and multipole l passed in input.
For a given value of k, the transfer function is inferred from the source function (passed in input in the array interpolated_sources) and from Bessel functions (passed in input in the bessels structure), either by convolving them along tau, or by a Limber approximation. This elementary task is distributed either to transfer_integrate() or to transfer_limber(). The task of this routine is mainly to loop over k values, and to decide at which k_max the calculation can be stopped, according to some approximation scheme designed to find a compromise between execution time and precision. The approximation scheme is defined by parameters in the precision structure.
ptw | Input: pointer to transfer_workspace structure (allocated in transfer_init() to avoid numerous reallocation) |
ppr | Input: pointer to precision structure |
ppt | Input: pointer to perturbation structure |
ptr | Input/output: pointer to transfers structure (result stored there) |
index_q | Input: index of wavenumber |
index_md | Input: index of mode |
index_ic | Input: index of initial condition |
index_tt | Input: index of type of transfer |
index_l | Input: index of multipole |
l | Input: multipole |
q_max_bessel | Input: maximum value of argument q at which Bessel functions are computed |
radial_type | Input: type of radial (Bessel) functions to convolve with |
Summary:
int transfer_integrate | ( | struct perturbs * | ppt, |
struct transfers * | ptr, | ||
struct transfer_workspace * | ptw, | ||
int | index_q, | ||
int | index_md, | ||
int | index_tt, | ||
double | l, | ||
int | index_l, | ||
double | k, | ||
radial_function_type | radial_type, | ||
double * | trsf | ||
) |
This routine computes the transfer functions ) for each mode, initial condition, type, multipole l and wavenumber k, by convolving the source function (passed in input in the array interpolated_sources) with Bessel functions (passed in input in the bessels structure).
ppt | Input: pointer to perturbation structure |
ptr | Input: pointer to transfers structure |
ptw | Input: pointer to transfer_workspace structure (allocated in transfer_init() to avoid numerous reallocation) |
index_q | Input: index of wavenumber |
index_md | Input: index of mode |
index_tt | Input: index of type |
l | Input: multipole |
index_l | Input: index of multipole |
k | Input: wavenumber |
radial_type | Input: type of radial (Bessel) functions to convolve with |
trsf | Output: transfer function ![]() |
Summary:
int transfer_limber | ( | struct transfers * | ptr, |
struct transfer_workspace * | ptw, | ||
int | index_md, | ||
int | index_q, | ||
double | l, | ||
double | q, | ||
radial_function_type | radial_type, | ||
double * | trsf | ||
) |
This routine computes the transfer functions ) for each mode, initial condition, type, multipole l and wavenumber k, by using the Limber approximation, i.e by evaluating the source function (passed in input in the array interpolated_sources) at a single value of tau (the Bessel function being approximated as a Dirac distribution).
ptr | Input: pointer to transfers structure |
ptw | Input: pointer to transfer workspace structure |
index_md | Input: index of mode |
index_q | Input: index of wavenumber |
l | Input: multipole |
q | Input: wavenumber |
radial_type | Input: type of radial (Bessel) functions to convolve with |
trsf | Output: transfer function ![]() |
Summary:
int transfer_limber_interpolate | ( | struct transfers * | ptr, |
double * | tau0_minus_tau, | ||
double * | sources, | ||
int | tau_size, | ||
double | tau0_minus_tau_limber, | ||
double * | S | ||
) |
int transfer_limber2 | ( | int | tau_size, |
struct transfers * | ptr, | ||
int | index_md, | ||
int | index_k, | ||
double | l, | ||
double | k, | ||
double * | tau0_minus_tau, | ||
double * | sources, | ||
radial_function_type | radial_type, | ||
double * | trsf | ||
) |
This routine computes the transfer functions ) for each mode, initial condition, type, multipole l and wavenumber k, by using the Limber approximation at order two, i.e as a function of the source function and its first two derivatives at a single value of tau
tau_size | Input: size of conformal time array |
ptr | Input: pointer to transfers structure |
index_md | Input: index of mode |
index_k | Input: index of wavenumber |
l | Input: multipole |
k | Input: wavenumber |
tau0_minus_tau | Input: array of values of (tau_today - tau) |
sources | Input: source functions |
radial_type | Input: type of radial (Bessel) functions to convolve with |
trsf | Output: transfer function ![]() |
Summary: