Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Make tracer creation more robust #928

Merged
merged 22 commits into from
Apr 28, 2022
Merged

Make tracer creation more robust #928

merged 22 commits into from
Apr 28, 2022

Conversation

tilmantroester
Copy link
Contributor

@tilmantroester tilmantroester commented Apr 14, 2022

This is an attempt to make CCL more robust to bad n(z). See e.g. #927, #919.

  • Added spline integration for the n(z) normalisation integral and the lensing kernel integral
  • Added flags in gsl_params to switch between spline and GSL quadrature integration for the two cases. Since spline integration is more robust and seems to be faster in the case where the n(z) is not super smooth or not very finely sampled, this is chosen as the default.
  • When using GSL quadrature for the lensing integral, the integration error is compared to the peak of the lensing kernel, instead of just the relative accuracy of the integral. This solves issues with numerical convergence of the integrals in the tails of the kernel where the absolute values are small but the integration method struggles.
  • Added option to choose lower number of samples for the lensing kernel than samples in the n(z). This is done using the n_samples argument, like in CMBLensingTracer. Right now there is an inconsitency in the API, in that the tSZ, CIB, and ISW tracers use n_chi for the same thing. This should be unified but as this involves a change to existing API, I've not included it here. @damonge Let me know what you think.
  • Added the option to call Tracer.get_kernels without specifying chi. In this case, the internal spline representation is returned. This is mostly useful for debugging purposes.
  • Enabled warnings by default (can be turned off with pyccl.debug_mode(False). Many of the warnings are there for good reason and warning the user that internal computations might produce garbage results is probably preferable to silent failures.
  • Added a warning when the status message gets overwritten. Right now there's no error stack or traceback. If multiple errors get written to cosmo.status_message, earlier messages get lost, which makes debugging difficult. Partially addresses ENH clean out the C layer warnings #672. But this should really be handled cleaner (let's just switch to C++ and get exceptions and smart pointers!)

@coveralls
Copy link

coveralls commented Apr 14, 2022

Pull Request Test Coverage Report for Build 2236361717

  • 27 of 27 (100.0%) changed or added relevant lines in 1 file are covered.
  • No unchanged relevant lines lost coverage.
  • Overall coverage increased (+0.02%) to 97.077%

Totals Coverage Status
Change from base Build 2219940793: 0.02%
Covered Lines: 4451
Relevant Lines: 4585

💛 - Coveralls

@tilmantroester tilmantroester marked this pull request as ready for review April 20, 2022 15:03
@tilmantroester tilmantroester marked this pull request as draft April 21, 2022 22:18
@tilmantroester tilmantroester marked this pull request as ready for review April 25, 2022 12:03
Copy link
Collaborator

@damonge damonge left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks a lot @tilmantroester . A few comments below

pyccl/tracers.py Outdated Show resolved Hide resolved
pyccl/tracers.py Outdated Show resolved Hide resolved
src/ccl_tracers.c Show resolved Hide resolved
src/ccl_tracers.c Show resolved Hide resolved
pyccl/tests/test_tracers.py Outdated Show resolved Hide resolved
@damonge
Copy link
Collaborator

damonge commented Apr 25, 2022

On cleaning up warnings: I agree... we should try to sort that out before V3
On nchi vs. nsamples: I also agree, but let's not break the API yet. Perhaps this could be added to #889

@nikfilippas
Copy link
Contributor

nikfilippas commented Apr 25, 2022

Unrelated to the title of this PR but definitely related to C-level warnings: currently, GSL raises an error when a value outside of the interpolation range is requested. For some reason, this is caught and converted to a warning in CCL, and the code ends up failing with an unexplained error CCLError: Error 1: . Would we consider promoting GSL out-of-bounds warnings to errors (perhaps in this PR)?

Try, for example, ccl.CosmologyVanillaLCDM().sigmaM(M=1e5, a=1).

As the vast majority (citation needed?) of GSL failures in CCL are caused by the queried values being out of bounds, maybe it's worth either catching this case specifically and promoting it to an error (GSL_EDOM=-1 here) or just promoting every GSL failure to an error.

@tilmantroester
Copy link
Contributor Author

tilmantroester commented Apr 26, 2022

Ok, I think I've addressed the comments. Still need to think about how to handle the case of sparsely sampled n(z).

Re: error handling. The more I look into this, the more jank I find. Probably best to completely rewrite the ccl_set_status and ccl_raise_* with something more robust. If someone knows a C library that does this nicely, that would be good to know. The requirement of being python friendly but C-only makes this awkard though.

@damonge
Copy link
Collaborator

damonge commented Apr 26, 2022

OK, a more complete solution to error handling is for another PR. @tilmantroester see my suggestion above for the poorly-sampled case. Not sure if it's clear. I can provide a diagram if not.

@tilmantroester
Copy link
Contributor Author

@damonge The poorly-sampled case should be addressed now.
@nikfilippas I added a warning if the n(z) is poorly sampled that points used to changing the integration method, which is set by a flag in gsl_params. Once your PR is in that warning needs to be adjusted accordingly.

Copy link
Collaborator

@damonge damonge left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

LGTM!

@tilmantroester tilmantroester merged commit fd54a21 into master Apr 28, 2022
@damonge damonge deleted the tracer_integrals branch April 28, 2022 12:40
@nikfilippas
Copy link
Contributor

nikfilippas commented May 5, 2022

@tilmantroester I noticed while running some tests that this PR introduces a bunch of unnecessary warnings in the python unit tests and in the benchmarks. It would be nice to get rid of them in a separate PR to keep the output as clean as possible.

UPDATE: I addressed that in PR889.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

4 participants