CLASS MANUAL
|
#include "primordial.h"
Functions | |
int | primordial_spectrum_at_k (struct primordial *ppm, int index_md, enum linear_or_logarithmic mode, double input, double *output) |
int | primordial_init (struct precision *ppr, struct perturbs *ppt, struct primordial *ppm) |
int | primordial_free (struct primordial *ppm) |
int | primordial_indices (struct perturbs *ppt, struct primordial *ppm) |
int | primordial_get_lnk_list (struct primordial *ppm, double kmin, double kmax, double k_per_decade) |
int | primordial_analytic_spectrum_init (struct perturbs *ppt, struct primordial *ppm) |
int | primordial_analytic_spectrum (struct primordial *ppm, int index_md, int index_ic1_ic2, double k, double *pk) |
int | primordial_inflation_potential (struct primordial *ppm, double phi, double *V, double *dV, double *ddV) |
int | primordial_inflation_hubble (struct primordial *ppm, double phi, double *H, double *dH, double *ddH, double *dddH) |
int | primordial_inflation_indices (struct primordial *ppm) |
int | primordial_inflation_solve_inflation (struct perturbs *ppt, struct primordial *ppm, struct precision *ppr) |
int | primordial_inflation_analytic_spectra (struct perturbs *ppt, struct primordial *ppm, struct precision *ppr, double *y_ini) |
int | primordial_inflation_spectra (struct perturbs *ppt, struct primordial *ppm, struct precision *ppr, double *y_ini) |
int | primordial_inflation_one_wavenumber (struct perturbs *ppt, struct primordial *ppm, struct precision *ppr, double *y_ini, int index_k) |
int | primordial_inflation_one_k (struct primordial *ppm, struct precision *ppr, double k, double *y, double *dy, double *curvature, double *tensor) |
int | primordial_inflation_find_attractor (struct primordial *ppm, struct precision *ppr, double phi_0, double precision, double *y, double *dy, double *H_0, double *dphidt_0) |
int | primordial_inflation_evolve_background (struct primordial *ppm, struct precision *ppr, double *y, double *dy, enum target_quantity target, double stop, short check_epsilon, enum integration_direction direction, enum time_definition time) |
int | primordial_inflation_check_potential (struct primordial *ppm, double phi, double *V, double *dV, double *ddV) |
int | primordial_inflation_check_hubble (struct primordial *ppm, double phi, double *H, double *dH, double *ddH, double *dddH) |
int | primordial_inflation_get_epsilon (struct primordial *ppm, double phi, double *epsilon) |
int | primordial_inflation_find_phi_pivot (struct primordial *ppm, struct precision *ppr, double *y, double *dy) |
int | primordial_inflation_derivs (double tau, double *y, double *dy, void *parameters_and_workspace, ErrorMsg error_message) |
int | primordial_external_spectrum_init (struct perturbs *ppt, struct primordial *ppm) |
Documented primordial module.
Julien Lesgourgues, 24.08.2010
This module computes the primordial spectra. It can be used in different modes: simple parametric form, evolving inflaton perturbations, etc. So far only the mode corresponding to a simple analytic form in terms of amplitudes, tilts and runnings has been developed.
The following functions can be called from other modules:
int primordial_spectrum_at_k | ( | struct primordial * | ppm, |
int | index_md, | ||
enum linear_or_logarithmic | mode, | ||
double | input, | ||
double * | output | ||
) |
Primordial spectra for arbitrary argument and for all initial conditions.
This routine evaluates the primordial spectrum at a given value of k by interpolating in the pre-computed table.
When k is not in the pre-computed range but the spectrum can be found analytically, it finds it. Otherwise returns an error.
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 primordial_init() has been called before, and primordial_free() has not been called yet.
ppm | Input: pointer to primordial structure containing tabulated primordial spectrum |
index_md | Input: index of mode (scalar, tensor, ...) |
mode | Input: linear or logarithmic |
input | Input: wavenumber in 1/Mpc (linear mode) or its logarithm (logarithmic mode) |
output | Output: for each pair of initial conditions, primordial spectra P(k) in ![]() |
Summary:
int primordial_init | ( | struct precision * | ppr, |
struct perturbs * | ppt, | ||
struct primordial * | ppm | ||
) |
This routine initializes the primordial structure (in particular, it computes table of primordial spectrum values)
ppr | Input: pointer to precision structure (defines method and precision for all computations) |
ppt | Input: pointer to perturbation structure (useful for knowing k_min, k_max, etc.) |
ppm | Output: pointer to initialized primordial structure |
Summary:
expression for alpha_s comes from:
ns_2 = (lnpk_plus-lnpk_pivot)/(dlnk)+1
ns_1 = (lnpk_pivot-lnpk_minus)/(dlnk)+1
alpha_s = dns/dlnk = (ns_2-ns_1)/dlnk = (lnpk_plus-lnpk_pivot-lnpk_pivot+lnpk_minus)/(dlnk)/(dlnk)
expression for beta_s:
ppm->beta_s = (alpha_plus-alpha_minus)/dlnk = (lnpk_plusplus-2.*lnpk_plus+lnpk_pivot - (lnpk_pivot-2.*lnpk_minus+lnpk_minusminus)/pow(dlnk,3)
int primordial_free | ( | struct primordial * | ppm | ) |
This routine frees all the memory space allocated by primordial_init().
To be called at the end of each run.
ppm | Input: pointer to primordial structure (which fields must be freed) |
int primordial_indices | ( | struct perturbs * | ppt, |
struct primordial * | ppm | ||
) |
This routine defines indices and allocates tables in the primordial structure
ppt | Input: pointer to perturbation structure |
ppm | Input/output: pointer to primordial structure |
int primordial_get_lnk_list | ( | struct primordial * | ppm, |
double | kmin, | ||
double | kmax, | ||
double | k_per_decade | ||
) |
This routine allocates and fills the list of wavenumbers k
ppm | Input/output: pointer to primordial structure |
kmin | Input: first value |
kmax | Input: last value that we should encompass |
k_per_decade | Input: number of k per decade |
int primordial_analytic_spectrum_init | ( | struct perturbs * | ppt, |
struct primordial * | ppm | ||
) |
This routine interprets and stores in a condensed form the input parameters in the case of a simple analytic spectra with amplitudes, tilts, runnings, in such way that later on, the spectrum can be obtained by a quick call to the routine primordial_analytic_spectrum(()
ppt | Input: pointer to perturbation structure |
ppm | Input/output: pointer to primordial structure |
int primordial_analytic_spectrum | ( | struct primordial * | ppm, |
int | index_md, | ||
int | index_ic1_ic2, | ||
double | k, | ||
double * | pk | ||
) |
This routine returns the primordial spectrum in the simple analytic case with amplitudes, tilts, runnings, for each mode (scalar/tensor...), pair of initial conditions, and wavenumber.
ppm | Input/output: pointer to primordial structure |
index_md | Input: index of mode (scalar, tensor, ...) |
index_ic1_ic2 | Input: pair of initial conditions (ic1, ic2) |
k | Input: wavenumber in same units as pivot scale, i.e. in 1/Mpc |
pk | Output: primordial power spectrum A (k/k_pivot)^(n+...) |
int primordial_inflation_potential | ( | struct primordial * | ppm, |
double | phi, | ||
double * | V, | ||
double * | dV, | ||
double * | ddV | ||
) |
This routine encodes the inflaton scalar potential
ppm | Input: pointer to primordial structure |
phi | Input: background inflaton field value in units of Mp |
V | Output: inflaton potential in units of ![]() |
dV | Output: first derivative of inflaton potential wrt the field |
ddV | Output: second derivative of inflaton potential wrt the field |
int primordial_inflation_hubble | ( | struct primordial * | ppm, |
double | phi, | ||
double * | H, | ||
double * | dH, | ||
double * | ddH, | ||
double * | dddH | ||
) |
This routine encodes the function
ppm | Input: pointer to primordial structure |
phi | Input: background inflaton field value in units of Mp |
H | Output: Hubble parameters in units of Mp |
dH | Output: ![]() |
ddH | Output: ![]() |
dddH | Output: ![]() |
int primordial_inflation_indices | ( | struct primordial * | ppm | ) |
This routine defines indices used by the inflation simulator
ppm | Input/output: pointer to primordial structure |
int primordial_inflation_solve_inflation | ( | struct perturbs * | ppt, |
struct primordial * | ppm, | ||
struct precision * | ppr | ||
) |
Main routine of inflation simulator. Its goal is to check the background evolution before and after the pivot value phi=phi_pivot, and then, if this evolution is suitable, to call the routine primordial_inflation_spectra().
ppt | Input: pointer to perturbation structure |
ppm | Input/output: pointer to primordial structure |
ppr | Input: pointer to precision structure |
Summary:
int primordial_inflation_analytic_spectra | ( | struct perturbs * | ppt, |
struct primordial * | ppm, | ||
struct precision * | ppr, | ||
double * | y_ini | ||
) |
Routine for the computation of an analytic apporoximation to the the primordial spectrum. In general, should be used only for comparing with exact numerical computation performed by primordial_inflation_spectra().
ppt | Input: pointer to perturbation structure |
ppm | Input/output: pointer to primordial structure |
ppr | Input: pointer to precision structure |
y_ini | Input: initial conditions for the vector of background/perturbations, already allocated and filled |
Summary
int primordial_inflation_spectra | ( | struct perturbs * | ppt, |
struct primordial * | ppm, | ||
struct precision * | ppr, | ||
double * | y_ini | ||
) |
Routine with a loop over wavenumbers for the computation of the primordial spectrum. For each wavenumber it calls primordial_inflation_one_wavenumber()
ppt | Input: pointer to perturbation structure |
ppm | Input/output: pointer to primordial structure |
ppr | Input: pointer to precision structure |
y_ini | Input: initial conditions for the vector of background/perturbations, already allocated and filled |
int primordial_inflation_one_wavenumber | ( | struct perturbs * | ppt, |
struct primordial * | ppm, | ||
struct precision * | ppr, | ||
double * | y_ini, | ||
int | index_k | ||
) |
Routine coordinating the computation of the primordial spectrum for one wavenumber. It calls primordial_inflation_one_k() to integrate the perturbation equations, and then it stores the result for the scalar/tensor spectra.
ppt | Input: pointer to perturbation structure |
ppm | Input/output: pointer to primordial structure |
ppr | Input: pointer to precision structure |
y_ini | Input: initial conditions for the vector of background/perturbations, already allocated and filled |
index_k | Input: index of wavenumber to be considered |
Summary
int primordial_inflation_one_k | ( | struct primordial * | ppm, |
struct precision * | ppr, | ||
double | k, | ||
double * | y, | ||
double * | dy, | ||
double * | curvature, | ||
double * | tensor | ||
) |
Routine integrating the background plus perturbation equations for each wavenumber, and returning the scalar and tensor spectrum.
ppm | Input: pointer to primordial structure |
ppr | Input: pointer to precision structure |
k | Input: Fourier wavenumber |
y | Input: running vector of background/perturbations, already allocated and initialized |
dy | Input: running vector of background/perturbation derivatives, already allocated |
curvature | Output: curvature perturbation |
tensor | Output: tensor perturbation |
Summary:
int primordial_inflation_find_attractor | ( | struct primordial * | ppm, |
struct precision * | ppr, | ||
double | phi_0, | ||
double | precision, | ||
double * | y, | ||
double * | dy, | ||
double * | H_0, | ||
double * | dphidt_0 | ||
) |
Routine searching for the inflationary attractor solution at a given phi_0, by iterations, with a given tolerance. If no solution found within tolerance, returns error message. The principle is the following. The code starts integrating the background equations from various values of phi, corresponding to earlier and earlier value before phi_0, and separated by a small arbitrary step size, corresponding roughly to 1 e-fold of inflation. Each time, the integration starts with the initial condition (slow-roll prediction). If the found value of
in phi_0 is stable (up to the parameter "precision"), the code considers that there is an attractor, and stops iterating. If this process does not converge, it returns an error message.
ppm | Input: pointer to primordial structure |
ppr | Input: pointer to precision structure |
phi_0 | Input: field value at which we wish to find the solution |
precision | Input: tolerance on output values (if too large, an attractor will always considered to be found) |
y | Input: running vector of background variables, already allocated and initialized |
dy | Input: running vector of background derivatives, already allocated |
H_0 | Output: Hubble value at phi_0 for attractor solution |
dphidt_0 | Output: field derivative value at phi_0 for attractor solution |
int primordial_inflation_evolve_background | ( | struct primordial * | ppm, |
struct precision * | ppr, | ||
double * | y, | ||
double * | dy, | ||
enum target_quantity | target, | ||
double | stop, | ||
short | check_epsilon, | ||
enum integration_direction | direction, | ||
enum time_definition | time | ||
) |
Routine integrating background equations only, from initial values stored in y, to a final value (if target = aH, until aH = aH_stop; if target = phi, till phi = phi_stop; if target = end_inflation, until (here t = proper time)). In output, y contains the final background values. In addition, if check_epsilon is true, the routine controls at each step that the expansion is accelerated and that inflation holds (wepsilon>1), otherwise it returns an error. Thanks to the last argument, it is also possible to specify whether the integration should be carried forward or backward in time. For the inflation_H case, only a 1st order differential equation is involved, so the forward and backward case can be done exactly without problems. For the inflation_V case, the equation of motion is 2nd order. What the module will do in the backward case is to search for an approximate solution, corresponding to the (first-order) attractor inflationary solution. This approximate backward solution is used in order to estimate some initial times, but the approximation made here will never impact the final result: the module is written in such a way that after using this approximation, the code always computes (and relies on) the exact forward solution.
ppm | Input: pointer to primordial structure |
ppr | Input: pointer to precision structure |
y | Input/output: running vector of background variables, already allocated and initialized |
dy | Input: running vector of background derivatives, already allocated |
target | Input: whether the goal is to reach a given aH or ![]() |
stop | Input: the target value of either aH or ![]() |
check_epsilon | Input: whether we should impose inflation (epsilon>1) at each step |
direction | Input: whether we should integrate forward or backward in time |
time | Input: definition of time (proper or conformal) |
int primordial_inflation_check_potential | ( | struct primordial * | ppm, |
double | phi, | ||
double * | V, | ||
double * | dV, | ||
double * | ddV | ||
) |
Routine checking positivity and negative slope of potential. The negative slope is an arbitrary choice. Currently the code can only deal with monotonic variations of the inflaton during inflation. So the slope had to be always negative or always positive... we took the first option.
ppm | Input: pointer to primordial structure |
phi | Input: field value where to perform the check |
V | Output: inflaton potential in units of ![]() |
dV | Output: first derivative of inflaton potential wrt the field |
ddV | Output: second derivative of inflaton potential wrt the field |
int primordial_inflation_check_hubble | ( | struct primordial * | ppm, |
double | phi, | ||
double * | H, | ||
double * | dH, | ||
double * | ddH, | ||
double * | dddH | ||
) |
Routine checking positivity and negative slope of . The negative slope is an arbitrary choice. Currently the code can only deal with monotonic variations of the inflaton during inflation. And H can only decrease with time. So the slope
has to be always negative or always positive... we took the first option: phi increases, H decreases.
ppm | Input: pointer to primordial structure |
phi | Input: field value where to perform the check |
H | Output: Hubble parameters in units of Mp |
dH | Output: ![]() |
ddH | Output: ![]() |
dddH | Output: ![]() |
int primordial_inflation_get_epsilon | ( | struct primordial * | ppm, |
double | phi, | ||
double * | epsilon | ||
) |
Routine computing the first slow-roll parameter epsilon
ppm | Input: pointer to primordial structure |
phi | Input: field value where to compute epsilon |
epsilon | Output: result |
int primordial_inflation_find_phi_pivot | ( | struct primordial * | ppm, |
struct precision * | ppr, | ||
double * | y, | ||
double * | dy | ||
) |
Routine searching phi_pivot when a given amount of inflation is requested.
ppm | Input/output: pointer to primordial structure |
ppr | Input: pointer to precision structure |
y | Input: running vector of background variables, already allocated and initialized |
dy | Input: running vector of background derivatives, already allocated |
Summary:
int primordial_inflation_derivs | ( | double | tau, |
double * | y, | ||
double * | dy, | ||
void * | parameters_and_workspace, | ||
ErrorMsg | error_message | ||
) |
Routine returning derivative of system of background/perturbation variables. Like other routines used by the generic integrator (background_derivs, thermodynamics_derivs, perturb_derivs), this routine has a generic list of arguments, and a slightly different error management, with the error message returned directly in an ErrMsg field.
tau | Input: time (not used explicitly inside the routine, but requested by the generic integrator) |
y | Input/output: running vector of background variables, already allocated and initialized |
dy | Input: running vector of background derivatives, already allocated |
parameters_and_workspace | Input: all necessary input variables apart from y |
error_message | Output: error message |
int primordial_external_spectrum_init | ( | struct perturbs * | ppt, |
struct primordial * | ppm | ||
) |
This routine reads the primordial spectrum from an external command, and stores the tabulated values. The sampling of the k's given by the external command is preserved.
Author: Jesus Torrado (torra) Date: 2013-12-20 doca cho@l oren tz.le iden univ. nl
ppt | Input/output: pointer to perturbation structure |
ppm | Input/output: pointer to primordial structure |
Summary: