LSSTDESC / CCL

DESC Core Cosmology Library: cosmology routines with validated numerical accuracy
BSD 3-Clause "New" or "Revised" License
141 stars 64 forks source link

Problem calling nonlinear P(k) at k< (kmin of CLASS) #141

Closed c-d-leonard closed 7 years ago

c-d-leonard commented 7 years ago

(I found this problem while working in the neutrinos branch - I think that it is an issue in the master branch as well, but correct me if I'm wrong.)

Short summary: When you try to get the nonlinear power spectrum from CCL below the minimum k that CLASS allows, you get a garbage answer rather than an error message (or, better, the right answer). The linear power spectrum is fine.

Longer notes:

To get the NL power spectrum spline in CCL, the CLASS function "spectra_pk_nl_at_k_and_z" is called from within "ccl_cosmology_compute_power_class" (in ccl_power.c), at a series of k values, defined by double * x = ccl_linear_spacing(kmin,kmax,nk). The power spectrum values at these k's are then used to create a spline. If the kmin in CCL here is at least the CLASS internal "minimum" k value, the right answer is returned.

But, if kmin is less than this, a garbage number is returned (actually, the value of the power spectrum at the largest k for the previous redshift that was calculated). The reason is that "spectra_pk_nl_at_k_and_z" only computes the k>=(CLASS kmin) case. In principle, there is a check in that function which sends an error messages somewhere within CLASS which I guess solves this problem when called in CLASS, but this doesn't seem to trip anything in CC. None of this applies to the calculation of the linear power spectrum, because in that case the appropriate analytical approximation is applied.

The result is that when setting the spline, there is a discontinuous jump (of several orders of magnitude) when you go from the highest k value before (CLASS kmin) to the lowest k value after (CLASS kmin). The result is that when trying to interpolate between these two points, you get a horribly wrong answer, even if the point at which you are trying to call the spline is actually greater than or equal to (CLASS kmin).

So the upshot is that if you try to call the nonlinear power spectrum in CCL at lower than about k = 1.028e-5 h / Mpc, it will give you nonsense, which is probably something we should at least code in an error-and-exit for.

elisachisari commented 7 years ago

This can be integrated with the power spectrum extrapolation to low k as well. Extrapolating below a certain kmin above the CLASS minimum would solve two problems in one go.

damonge commented 7 years ago

Will solve this in the E&Hu branch as discussed and close this when the PR is merged

elisachisari commented 7 years ago

Can this be closed? @damonge

damonge commented 7 years ago

Yep, closed.