[nikosarcevic]

this PR investigates a numerical instability in the Limber qag_quad intergation that can produce an isolated discontinuous outlier while GSL QAG reports successful convergence

For the reproducing case, the QAG result was wrong by approximately 4% despite returning GSL_SUCCESS, so the existing CQUAD fallback was not triggered.

This PR:

• adds cquad as an explicit Libmer integration option for comparison
• adds a diagnostic repro and independent numerical valdiation
• tests a twoway split QAG consistency check, falling back to CQUAD when the full and split QAG evaluations disagree
• includes a small performanc benchmark for the additional consistency check.

for the pathological point, the QAG local residual decreases from approx 4.0e-2 to 7.2e-6.
the corrected result agrees with CQUAD and an independently converged Simpson integration.

The existing pyccl/tests/test_cells.py suite passes (16 passed).

This is intended primarily to document and reproduce the issue and test a possible mitigation; the split-QAG implementation is not necessarily intended as the final solution.

[coveralls]

Coverage Report for CI Build 37088652479

Coverage remained the same at 97.471%

Details

• Coverage remained the same as the base build.
• Patch coverage: No coverable lines changed in this PR.
• No coverage regressions found.

Uncovered Changes

No uncovered changes found.

Coverage Regressions

No coverage regressions found.

---

Coverage Stats

Relevant Lines:	6801
Covered Lines:	6629
Line Coverage:	97.47%
Coverage Strength:	4.87 hits per line

---

💛 - Coveralls

[cunhad]

Following up from the Slack thread: while comparing NumCosmo's angular power spectra with CCL's, we (the NumCosmo side) ran into the same QAG false convergence discussed here, and tested the fixes proposed in this PR. The results are below, in case they're useful.

• In lensing it shows up across ℓ rather than across cosmologies. At CCL's defaults, 148 isolated-ℓ outliers (up to 4e-3) in weak and CMB lensing, ℓ 2–3000, five standard cosmologies, all returned with GSL_SUCCESS.
• The split check catches 80 and misses 68, confirmed on CCL built at this PR's head (8978d57). The 68 are the ℓ where QAG's own pass ended with two panels and a wrong sum. The check splits at the same midpoint, so its halves reproduce those two panels and it agrees with itself. For the same reason, a mandatory two-way subdivision would return the same wrong values there.
• On this PR's own clustering case we reproduce the −4.02%. It is a one-panel acceptance, which the check catches.
• Four equal panels in ln k removed all 148 and this PR's outlier, at about 1.15–1.2× the integrand evaluations. CQUAD also removed them all, with the fewest evaluations. A rough wall-time run on a loaded machine gave CQUAD ≈ 1× and this PR's default ≈ 2× the current qag_quad: indicative only.
• Side note, a different method: spline was fine in lensing for us, but on the ISW auto-spectrum it is off by a factor of several at the default DLOGK_INTEGRATION, in the one cosmology checked so far.

Happy to share more details or the scripts if that would help.