Skip to content

Commit

Permalink
Merge pull request #437 from LSSTDESC/relabel_spline_params
Browse files Browse the repository at this point in the history
Rename K_MIN_DEFAULT to K_MIN
  • Loading branch information
elisachisari authored Aug 17, 2018
2 parents e578d7f + 8d28c44 commit 3b5055c
Show file tree
Hide file tree
Showing 11 changed files with 41 additions and 34 deletions.
Binary file modified doc/0000-ccl_note/0000-ccl_note.pdf
Binary file not shown.
Binary file modified doc/0000-ccl_note/main.pdf
Binary file not shown.
20 changes: 11 additions & 9 deletions doc/0000-ccl_note/main.tex
Original file line number Diff line number Diff line change
Expand Up @@ -476,7 +476,7 @@ \subsubsection{Spline parameters \& the INI file}
\item {\tt A$\_$SPLINE$\_$MIN$\_$PK}: maximum value for the logarithmic part and minimum value for the linear part of the scale factor spline. The default is {\tt A$\_$SPLINE$\_$MIN$\_$PK}$=0.1$. This only applies for the 2D power spectrum splines.
\item {\tt K$\_$MAX$\_$SPLINE}: The maximum value of wavenumber considered in the spline. This is explained in more detail in the coming subsections. The default is {\tt K$\_$MAX$\_$SPLINE}$=50/$Mpc. In the case of the emulator's nonlinear power spectrum, the maximum wavenumber is {\tt K$\_$MAX$\_$SPLINE}$=5/$Mpc.
\item {\tt K$\_$MAX}: The maximum value of wavenumber when integrating over $k$. The default is {\tt K$\_$MAX}$=1000$/Mpc.
\item {\tt K$\_$MIN$\_$DEFAULT}: The minimum value of wavenumber when integrating over $k$. The default is {\tt K$\_$MIN$\_$DEFAULT}$=5 \times 10^{-5}$/Mpc.
\item {\tt K$\_$MIN}: The minimum value of wavenumber when integrating over $k$. The default is {\tt K$\_$MIN}$=5 \times 10^{-5}$/Mpc.
\item {\tt N$\_$K}: Number of bins per decade for the wavenumber. The default is {\tt N$\_$K}$=167$.
\item {\tt N$\_$K$\_$3DCOR}: Number of bins per decade for the wavenumber in the 3d correlation function calculation.
The default is {\tt N$\_$K$\_$3DCOR}~=~100,000.
Expand Down Expand Up @@ -511,11 +511,11 @@ \subsubsection{Spline parameters \& the INI file}
\subsubsection{Extrapolation for the nonlinear power spectrum}
\label{sec:NLextrapol}

The computation of the power spectrum from {\tt CLASS} can be significantly sped up by extrapolating in the range $k>$~{\tt K$\_$MAX$\_$SPLINE} and $k<$~{\tt K$\_$MIN}. In this section, we describe the implementation of the extrapolation and the accuracy attained. These tests are performed in a flat $\Lambda$CDM cosmology with $\Omega_c=0.25$, $\Omega_b=0.05$, $A_s=2.1\times10^{-9}$, $h=0.7$ and $n_s=0.96$.
The computation of the power spectrum from {\tt CLASS} can be significantly sped up by extrapolating in the range $k>$~{\tt K$\_$MAX$\_$SPLINE} and $k<$~{\tt K$\_$MIN$\_$SPLINE}. In this section, we describe the implementation of the extrapolation and the accuracy attained. These tests are performed in a flat $\Lambda$CDM cosmology with $\Omega_c=0.25$, $\Omega_b=0.05$, $A_s=2.1\times10^{-9}$, $h=0.7$ and $n_s=0.96$.

We first describe the extrapolation at high wavenumbers. The introduction of the parameter {\tt K$\_$MAX$\_$SPLINE} allows us to spline the matter power spectrum within the {\tt cosmo} structure up to that value of $k$ (in units of $1/$Mpc). A separate {\tt K$\_$MAX} parameter sets the limit for evaluation of the matter power spectrum. The range between {\tt K$\_$MAX$\_$SPLINE}$<k<${\tt K$\_$MAX} is evaluated by performing a second order Taylor expansion in $\ln k$ within the static routine {\tt ccl$\_$power$\_$extrapol$\_$highk}.

First, we compute the first and second derivative of the $\ln P(k,z)$ at $k_0={\rm \tt K\_MAX}-2\Delta\ln k$ by computing the numerical derivatives by finite differences using GSL. The fiducial choice for $\Delta\ln k$ is $10^{-2}$. We then apply a second order Taylor expansion to extrapolate the matter power spectrum to $k>${\tt K$\_$MAX$\_$SPLINE}. The Taylor expansion gives
First, we compute the first and second derivative of the $\ln P(k,z)$ at $k_0={\rm \tt K\_MAX\_SPLINE}-2\Delta\ln k$ by computing the numerical derivatives by finite differences using GSL. The fiducial choice for $\Delta\ln k$ is $10^{-2}$. We then apply a second order Taylor expansion to extrapolate the matter power spectrum to $k>${\tt K$\_$MAX$\_$SPLINE}. The Taylor expansion gives
%
\begin{equation}
\ln P(k,z) \simeq \ln P(k_0,z) + \frac{d\ln P}{d\ln k}(\ln k_0,z) (\ln k-\ln k_0) + \frac{1}{2} \frac{d^2\ln P}{d\ln k^2}(\ln k_0,z) (\ln k-\ln k_0)^2.
Expand All @@ -524,11 +524,13 @@ \subsubsection{Extrapolation for the nonlinear power spectrum}

The accuracy of this approximation is shown in Figure \ref{fig:NLextrapol} for redshifts $z=0$, $z=3$ and $z=20$. We compare the nonlinear matter power spectrum at these redshifts computed with the previously described approximation, to the matter power spectrum obtained by setting the power spectrum splines to high values. We find that for typical values of $\Delta \ln k=10^{-2}$ and {\tt K$\_$MAX$\_$SPLINE}$=50$/Mpc, $\ln P$ has converged to an accuracy that surpasses the expected impact of baryonic effects on the matter power spectrum at $k>10/$Mpc. (For an estimate of the impact of baryons on the total matter power spectrum, see \citealt{Schneider15}.) The lower {\tt K$\_$MAX$\_$SPLINE} is, the faster \ccl will run. The optimum choice of {\tt K$\_$MAX$\_$SPLINE} is left to the user for their particular application.

We also extrapolate the power spectrum at small wavenumbers within the static routine {\tt ccl\_power\_extrapol\_lowk}. In this case, the power spectrum below {\tt K$\_$MIN} is obtained by a power-law extrapolation with index $n_s$:
We also extrapolate the power spectrum at small wavenumbers within the static routine {\tt ccl\_power\_extrapol\_lowk}. In this case, the power spectrum below {\tt K$\_$MIN$\_$SPLINE} is obtained by a power-law extrapolation with index $n_s$:
\begin{equation}
\log P(k<{\tt K\_MIN},z) = \log P({\tt K\_MIN},z) + n_s (\log k-\log{\tt K\_MIN})
\log P(k<{\tt K\_MIN\_SPLINE},z) = \log P({\tt K\_MIN\_SPLINE},z) + n_s (\log k-\log{\tt K\_MIN\_SPLINE})
\end{equation}
The value adopted for {\tt K\_MIN} depends on the choice of power spectrum method. For CLASS and the nonlinear power spectrum, we adopt {\tt K\_MIN} that coincides with the smallest wavenumber output by CLASS, {\tt K\_MIN}$=7\times 10^{-6}$/Mpc.\footnote{For BBKS, the power spectrum is computed analytically at all $k$, there is no extrapolation. For the Eisenstein \& Hu implementation, the splines of the power spectrum span {\tt K\_MIN\_DEFAULT}$<k<${\tt K$\_$MAX$\_$SPLINE}, so there is only extrapolation at high $k$.} Note that this parameter is different from {\tt K\_MIN\_DEFAULT}, which sets the minimum $k$ for integrations and which is set by default to {\tt K\_MIN\_DEFAULT}$=5\times 10^{-5}$/Mpc. Hence, in practice, no extrapolation is occurring in this case, unless the user specifically asks for an output power spectra below {\tt K\_MIN\_DEFAULT} for their own purposes.
The value adopted for {\tt K\_MIN\_SPLINE} depends on the choice of power spectrum method and is not currently accessible to the user via the .ini file. For {\tt CLASS} and the nonlinear power spectrum, we adopt {\tt K\_MIN\_SPLINE} that coincides with the smallest wavenumber output by {\tt CLASS}, {\tt K\_MIN\_SPLINE}$=7\times 10^{-6}$/Mpc. Note that this parameter is different from {\tt K\_MIN}, which sets the minimum $k$ for integrations and which is set by default to {\tt K\_MIN}$=5\times 10^{-5}$/Mpc. Hence, in practice, no extrapolation is occurring in this case, unless the user specifically asks for an output power spectra below {\tt K\_MIN} for their own purposes.

For BBKS, the power spectrum is computed analytically at all $k$, there is no extrapolation. For the Eisenstein \& Hu implementation, the splines of the power spectrum span {\tt K\_MIN}$<k<${\tt K$\_$MAX$\_$SPLINE}, so there is only extrapolation at high $k$. For the nonlinear matter power spectrum from the emulator, {\tt K\_MIN\_SPLINE} and {\tt K\_MAX\_SPLINE} are set to fixed values that are determined from the range of validity of the emulator: {\tt K\_MIN\_SPLINE}$=10^{-3}$ Mpc$^{-1}$ and {\tt K\_MAX\_SPLINE}$=5$ Mpc$^{-1}$.


%------------------------
Expand All @@ -545,7 +547,7 @@ \subsubsection{Extrapolation for the linear power spectrum}

With the implementation described in the previous section, the power spectrum splines are initialized up to {\tt K$\_$MAX$\_$SPLINE}. This is also true for the linear matter power spectrum, which is used within \ccl in particular to obtain $\sigma_8$ (see Eq. \ref{eq:sigR}). We have tested here how the procedure described in the previous section affects the convergence of the linear matter power spectrum. The result is shown in Figure \ref{fig:NLextrapol}. For some applications that use the linear power spectrum, the user might need to increase the value of {\tt K$\_$MAX$\_$SPLINE}.

As in the previous section, the power spectrum at small wavenumber is extrapolated using a power-law. This extrapolation is performed below a fiducial value of {\tt K\_MIN\_DEFAULT}$=5\times 10^{-5}$.
As in the previous section, the power spectrum at small wavenumber is extrapolated using a power-law. This extrapolation is performed below a fiducial value of {\tt K\_MIN\_SPLINE} that coincides with the smallest wavenumber output by {\tt CLASS}, as in the case of the nonlinear power spectrum described above.

We have found that changing {\tt N$\_$A} to $200$, or changing the sampling of the wavenumber to $5000$ points, does not change the results presented in Figures \ref{fig:NLextrapol} in this section.

Expand All @@ -555,7 +557,7 @@ \subsubsection{Wishlist for the future}
\begin{itemize}
\item CAMB,
\item other emulators,
\item Halo model/HOD.
\item HOD.
\end{itemize}


Expand All @@ -566,7 +568,7 @@ \subsubsection{Normalization of the power spectrum}

In practice, there is only one argument that encodes the normalization. This is the argument ${\tt norm\_pk}$, which can be passed the power spectrum normalization parameterized by $\sigma_8$ or $A_\mathrm{s}$. As noted above, {\tt ccl$\_$parameters$\_$create} switches to $\sigma_8$ normalization if ${\tt norm\_pk} > 10^{-5}$, and to $A_{\mathrm s}$ normalization otherwise.

In the {\tt python} implementation, {\tt CCL} allows for either $\sigma 8$ or {\tt A\_s} to be passed as parameters.
In the {\tt python} implementation, {\tt CCL} allows for either $\sigma_8$ or {\tt A\_s} to be passed as parameters.

\subsection{Angular power spectra}
\label{sec:cl}
Expand Down
4 changes: 2 additions & 2 deletions examples/ccl_sample_baryons.c
Original file line number Diff line number Diff line change
Expand Up @@ -32,7 +32,7 @@ int main(int argc, char * argv[])
//The linear power spectrum is not changed when baryons are passed
/*printf("Linear matter PS\n");
printf("# k [1/Mpc],P(k,z=0),P(k,z=1),P(k,z=2),P(k,z=3)\n");
for (k = ccl_splines->K_MIN_DEFAULT; k<ccl_splines->K_MAX; k*=1.05) {
for (k = ccl_splines->K_MIN; k<ccl_splines->K_MAX; k*=1.05) {
p = ccl_linear_matter_power(cosmo, k,1.0, &status);
p1 = ccl_linear_matter_power(cosmo,k, a_at_z1,&status);
p2 = ccl_linear_matter_power(cosmo,k, a_at_z2,&status);
Expand All @@ -41,7 +41,7 @@ int main(int argc, char * argv[])
}*/
printf("# Total matter power spectrum\n");
printf("# k [1/Mpc],P(k,z=0),P(k,z=1),P(k,z=2),P(k,z=3)\n");
for (k = ccl_splines->K_MIN_DEFAULT; k<ccl_splines->K_MAX; k*=1.05) {
for (k = ccl_splines->K_MIN; k<ccl_splines->K_MAX; k*=1.05) {
p = ccl_nonlin_matter_power(cosmo, k,1.0,&status);
p1 = ccl_nonlin_matter_power(cosmo,k, a_at_z1,&status);
p2 = ccl_nonlin_matter_power(cosmo,k, a_at_z2,&status);
Expand Down
4 changes: 2 additions & 2 deletions examples/ccl_sample_power.c
Original file line number Diff line number Diff line change
Expand Up @@ -25,7 +25,7 @@ int main(int argc, char * argv[])
double a_at_z2=1./3.;
double a_at_z3=0.25;
if(cosmo->config.matter_power_spectrum_method==ccl_linear) {
for (k = ccl_splines->K_MIN_DEFAULT; k<ccl_splines->K_MAX; k*=1.05) {
for (k = ccl_splines->K_MIN; k<ccl_splines->K_MAX; k*=1.05) {
p = ccl_linear_matter_power(cosmo, k,1.0, &status);
p1 = ccl_linear_matter_power(cosmo,k, a_at_z1,&status);
p2 = ccl_linear_matter_power(cosmo,k, a_at_z2,&status);
Expand All @@ -35,7 +35,7 @@ int main(int argc, char * argv[])
}
else {
if(cosmo->config.matter_power_spectrum_method==ccl_halofit) {
for (k = ccl_splines->K_MIN_DEFAULT; k<ccl_splines->K_MAX; k*=1.05) {
for (k = ccl_splines->K_MIN; k<ccl_splines->K_MAX; k*=1.05) {
p = ccl_nonlin_matter_power(cosmo, k,1.0,&status);
p1 = ccl_nonlin_matter_power(cosmo,k, a_at_z1,&status);
p2 = ccl_nonlin_matter_power(cosmo,k, a_at_z2,&status);
Expand Down
2 changes: 1 addition & 1 deletion include/ccl_params.h
Original file line number Diff line number Diff line change
Expand Up @@ -32,7 +32,7 @@ typedef struct ccl_spline_params {
//k-splines and integrals
double K_MAX_SPLINE;
double K_MAX;
double K_MIN_DEFAULT;
double K_MIN;
int N_K;
int N_K_3DCOR;

Expand Down
2 changes: 1 addition & 1 deletion include/ccl_params.ini
Original file line number Diff line number Diff line change
Expand Up @@ -35,7 +35,7 @@ A_SPLINE_NLOG_PK=11
;These are in units of Mpc (no factor of h)
K_MAX_SPLINE=50
K_MAX=1e3
K_MIN_DEFAULT=5e-5
K_MIN=5e-5
; This is now the number of k per decade
N_K=167
; Number of k per decade for 3D-correlation
Expand Down
8 changes: 4 additions & 4 deletions src/ccl_cls.c
Original file line number Diff line number Diff line change
Expand Up @@ -896,13 +896,13 @@ static double *get_lkarr(ccl_cosmology *cosmo,CCL_ClWorkspace *w,
//First compute relevant k-range for this ell
double kmin,kmax,lkmin,lkmax;
if(l>w->l_limber) {
kmin=CCL_MAX(ccl_splines->K_MIN_DEFAULT,0.8*(l+0.5)/chimax);
kmin=CCL_MAX(ccl_splines->K_MIN,0.8*(l+0.5)/chimax);
kmax=CCL_MIN(ccl_splines->K_MAX,1.2*(l+0.5)/chimin);
}
else {
double xmin,xmax;
limits_bessel(l,CCL_FRAC_RELEVANT,&xmin,&xmax);
kmin=CCL_MAX(ccl_splines->K_MIN_DEFAULT,xmin/chimax);
kmin=CCL_MAX(ccl_splines->K_MIN,xmin/chimax);
kmax=CCL_MIN(ccl_splines->K_MAX,xmax/chimin);
//Cap by maximum meaningful argument of the Bessel function
kmax=CCL_MIN(kmax,2*(w->l_arr[w->n_ls-1]+0.5)/chimin); //Cap by 2 x inverse scale corresponding to l_max
Expand Down Expand Up @@ -1079,7 +1079,7 @@ static void get_k_interval(ccl_cosmology *cosmo,CCL_ClWorkspace *w,

if(l<w->l_limber) {
chimin=2*(l+0.5)/ccl_splines->K_MAX;
chimax=0.5*(l+0.5)/ccl_splines->K_MIN_DEFAULT;
chimax=0.5*(l+0.5)/ccl_splines->K_MIN;
}
else {
if(cut_low_1) {
Expand All @@ -1098,7 +1098,7 @@ static void get_k_interval(ccl_cosmology *cosmo,CCL_ClWorkspace *w,
}
else {
chimin=0.5*(l+0.5)/ccl_splines->K_MAX;
chimax=2*(l+0.5)/ccl_splines->K_MIN_DEFAULT;
chimax=2*(l+0.5)/ccl_splines->K_MIN;
}
}

Expand Down
3 changes: 1 addition & 2 deletions src/ccl_core.c
Original file line number Diff line number Diff line change
Expand Up @@ -64,7 +64,6 @@ void ccl_cosmology_read_config(void)
// Use default ini file
param_file = EXPAND_STR(__CCL_DATA_DIR__) "/ccl_params.ini";
}

if ((fconfig=fopen(param_file, "r")) == NULL) {
char msg[256];
snprintf(msg, 256, "ccl_core.c: Failed to open config file: %s", param_file);
Expand Down Expand Up @@ -103,7 +102,7 @@ void ccl_cosmology_read_config(void)
if(strcmp(var_name,"A_SPLINE_NLOG_PK")==0) ccl_splines->A_SPLINE_NLOG_PK=(int) var_dbl;
if(strcmp(var_name,"K_MAX_SPLINE")==0) ccl_splines->K_MAX_SPLINE=var_dbl;
if(strcmp(var_name,"K_MAX")==0) ccl_splines->K_MAX=var_dbl;
if(strcmp(var_name,"K_MIN_DEFAULT")==0) ccl_splines->K_MIN_DEFAULT=var_dbl;
if(strcmp(var_name,"K_MIN")==0) ccl_splines->K_MIN=var_dbl;
if(strcmp(var_name,"N_K")==0) ccl_splines->N_K=(int) var_dbl;
// 3dcorr parameters
if(strcmp(var_name,"N_K_3DCOR")==0) ccl_splines->N_K_3DCOR=(int) var_dbl;
Expand Down
4 changes: 2 additions & 2 deletions src/ccl_correlation.c
Original file line number Diff line number Diff line change
Expand Up @@ -452,9 +452,9 @@ void ccl_correlation_3d(ccl_cosmology *cosmo, double a,
double *k_arr,*pk_arr,*r_arr,*xi_arr;

//number of data points for k and pk array
N_ARR=(int)(ccl_splines->N_K_3DCOR*log10(ccl_splines->K_MAX/ccl_splines->K_MIN_DEFAULT));
N_ARR=(int)(ccl_splines->N_K_3DCOR*log10(ccl_splines->K_MAX/ccl_splines->K_MIN));

k_arr=ccl_log_spacing(ccl_splines->K_MIN_DEFAULT,ccl_splines->K_MAX,N_ARR);
k_arr=ccl_log_spacing(ccl_splines->K_MIN,ccl_splines->K_MAX,N_ARR);
if(k_arr==NULL) {
*status=CCL_ERROR_MEMORY;
strcpy(cosmo->status_message,"ccl_correlation.c: ccl_correlation_3d ran out of memory\n");
Expand Down
28 changes: 17 additions & 11 deletions src/ccl_power.c
Original file line number Diff line number Diff line change
Expand Up @@ -421,6 +421,7 @@ static void ccl_cosmology_compute_power_class(ccl_cosmology * cosmo, int * statu
return;
}

//These are the limits of the splining range
cosmo->data.k_min_lin=2*exp(sp.ln_k[0]);
cosmo->data.k_max_lin=ccl_splines->K_MAX_SPLINE;

Expand Down Expand Up @@ -493,6 +494,7 @@ static void ccl_cosmology_compute_power_class(ccl_cosmology * cosmo, int * statu
strcpy(cosmo->status_message ,"ccl_power.c: ccl_cosmology_compute_power_class(): K_MIN is less than CLASS's kmin. Not yet supported for nonlinear P(k).\n");
}

//These are the limits of the splining range
cosmo->data.k_min_nl=2*exp(sp.ln_k[0]);
cosmo->data.k_max_nl=ccl_splines->K_MAX_SPLINE;

Expand Down Expand Up @@ -802,8 +804,9 @@ static double eh_power(ccl_parameters *params,eh_struct *eh,double k,int wiggled

static void ccl_cosmology_compute_power_eh(ccl_cosmology * cosmo, int * status)
{
cosmo->data.k_min_lin = ccl_splines->K_MIN_DEFAULT;
cosmo->data.k_min_nl = ccl_splines->K_MIN_DEFAULT;
//These are the limits of the splining range
cosmo->data.k_min_lin = ccl_splines->K_MIN;
cosmo->data.k_min_nl = ccl_splines->K_MIN;
cosmo->data.k_max_lin = ccl_splines->K_MAX;
cosmo->data.k_max_nl = ccl_splines->K_MAX;
double kmin = cosmo->data.k_min_lin;
Expand Down Expand Up @@ -979,8 +982,9 @@ TASK: provide spline for the BBKS power spectrum with baryonic correction

static void ccl_cosmology_compute_power_bbks(ccl_cosmology * cosmo, int * status)
{
cosmo->data.k_min_lin=ccl_splines->K_MIN_DEFAULT;
cosmo->data.k_min_nl=ccl_splines->K_MIN_DEFAULT;
//These are the limits of the splining range
cosmo->data.k_min_lin=ccl_splines->K_MIN;
cosmo->data.k_min_nl=ccl_splines->K_MIN;
cosmo->data.k_max_lin=ccl_splines->K_MAX;
cosmo->data.k_max_nl=ccl_splines->K_MAX;
double kmin = cosmo->data.k_min_lin;
Expand Down Expand Up @@ -1222,9 +1226,10 @@ static void ccl_cosmology_compute_power_emu(ccl_cosmology * cosmo, int * status)
return;
}

//These are the limits of the splining range
cosmo->data.k_min_lin=2*exp(sp.ln_k[0]);
cosmo->data.k_max_lin=ccl_splines->K_MAX_SPLINE;
//CLASS calculations done - now allocate CCL splines
//CLASS calculations done - now allocate CCL splines
double kmin = cosmo->data.k_min_lin;
double kmax = ccl_splines->K_MAX_SPLINE;
//Compute nk from number of decades and N_K = # k per decade
Expand Down Expand Up @@ -1286,6 +1291,7 @@ static void ccl_cosmology_compute_power_emu(ccl_cosmology * cosmo, int * status)
}

//Now start the NL computation with the emulator
//These are the limits of the splining range
cosmo->data.k_min_nl=K_MIN_EMU;
cosmo->data.k_max_nl=K_MAX_EMU;
amin = A_MIN_EMU; //limit of the emulator
Expand Down Expand Up @@ -1394,15 +1400,15 @@ INPUT: ccl_cosmology * cosmo, a, k [1/Mpc]
TASK: extrapolate power spectrum at high k
*/
static double ccl_power_extrapol_highk(ccl_cosmology * cosmo, double k, double a,
gsl_spline2d * powerspl, double kmax, int * status)
gsl_spline2d * powerspl, double kmax_spline, int * status)
{
double log_p_1;
double deltak=1e-2; //step for numerical derivative;
double deriv_pk_kmid,deriv2_pk_kmid;
double lkmid;
double lpk_kmid;

lkmid = log(kmax)-2*deltak;
lkmid = log(kmax_spline)-2*deltak;

int gslstatus = gsl_spline2d_eval_e(powerspl, lkmid,a,NULL ,NULL ,&lpk_kmid);
if(gslstatus != GSL_SUCCESS) {
Expand Down Expand Up @@ -1432,16 +1438,16 @@ static double ccl_power_extrapol_highk(ccl_cosmology * cosmo, double k, double a

}

/*------ ROUTINE: ccl_power_extrapol_hxighk -----
/*------ ROUTINE: ccl_power_extrapol_lowk -----
INPUT: ccl_cosmology * cosmo, a, k [1/Mpc]
TASK: extrapolate power spectrum at low k
*/
static double ccl_power_extrapol_lowk(ccl_cosmology * cosmo, double k, double a,
gsl_spline2d * powerspl, double kmin, int * status)
gsl_spline2d * powerspl, double kmin_spline, int * status)
{
double log_p_1;
double deltak=1e-2; //safety step
double lkmin=log(kmin)+deltak;
double lkmin=log(kmin_spline)+deltak;
double lpk_kmin;
int gslstatus = gsl_spline2d_eval_e(powerspl,lkmin,a,NULL,NULL,&lpk_kmin);

Expand Down Expand Up @@ -1669,7 +1675,7 @@ double ccl_sigmaR(ccl_cosmology *cosmo,double R, int *status)
F.function=&sigmaR_integrand;
F.params=&par;
double sigma_R;
int gslstatus = gsl_integration_cquad(&F, log10(ccl_splines->K_MIN_DEFAULT), log10(ccl_splines->K_MAX),
int gslstatus = gsl_integration_cquad(&F, log10(ccl_splines->K_MIN), log10(ccl_splines->K_MAX),
0.0, ccl_gsl->INTEGRATION_SIGMAR_EPSREL,
workspace,&sigma_R,NULL,NULL);
if(gslstatus != GSL_SUCCESS) {
Expand Down

0 comments on commit 3b5055c

Please sign in to comment.