Commit 11513a94 authored by Goran Jelic-Cizmek's avatar Goran Jelic-Cizmek
Browse files

Fixed GSL error handler bailing

With just density, the l=2,4 multipoles are extremely small, and the GSL
integrator would abort the computation, so now I just turn it off for both
`compute_corrfunc` and `compute_multipole`.
parent 3b2aa9aa
......@@ -5,6 +5,16 @@
# remember to add it inside this cdef.
cdef extern from "<gsl/gsl_errno.h>":
cdef struct gsl_error_handler_t:
pass
cdef gsl_error_handler_t *gsl_set_error_handler_off()
cdef void gsl_set_error_handler(gsl_error_handler_t *)
cdef extern from "<gsl/gsl_spline.h>":
cdef struct gsl_spline:
......
......@@ -1356,7 +1356,10 @@ cdef class Coffe:
_check_parameter('r', r, (int, float), 0, 1500)
_check_parameter('mu', mu, (int, float), -1, 1)
return ccoffe.coffe_integrate(
cdef ccoffe.gsl_error_handler_t *default_handler = \
ccoffe.gsl_set_error_handler_off()
result = ccoffe.coffe_integrate(
&self._parameters,
&self._background,
&self._integral,
......@@ -1378,6 +1381,10 @@ cdef class Coffe:
ccoffe.DOUBLE_INTEGRATED, ccoffe.CORRFUNC
)
ccoffe.gsl_set_error_handler(default_handler)
return result
def compute_multipole(
self, *,
......@@ -1413,7 +1420,10 @@ cdef class Coffe:
_check_parameter('r', r, (int, float), 0, 1500)
_check_parameter('l', l, (int,), 0, 10)
return ccoffe.coffe_integrate(
cdef ccoffe.gsl_error_handler_t *default_handler = \
ccoffe.gsl_set_error_handler_off()
result = ccoffe.coffe_integrate(
&self._parameters,
&self._background,
&self._integral,
......@@ -1435,6 +1445,10 @@ cdef class Coffe:
ccoffe.DOUBLE_INTEGRATED, ccoffe.MULTIPOLES
)
ccoffe.gsl_set_error_handler(default_handler)
return result
def compute_multipoles_bulk(self, recompute : bool = False):
"""
......
......@@ -243,6 +243,20 @@ class TestCoffe:
assert covariance_matrix(result, rstep=111).size == 0
def test_error_handler(self):
"""
Checks that the GSL error handler doesn't abort the computation
"""
# initialize COFFE with some cosmology
cosmo = coffe.Coffe()
cosmo.set_parameters(
sep=np.linspace(10, 300, 100),
has_rsd=False,
omega_m=0.31,
)
cosmo.compute_multipole(l=4, r=10, z=1.5)
class TestRepresentation:
"""
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment