From e55ea02b21fbb14eeb92a559c96b77bab59ad76c Mon Sep 17 00:00:00 2001 From: Steven Murray Date: Thu, 27 Jun 2024 12:11:08 +0200 Subject: [PATCH] maint: add warnings tests ci --- .github/workflows/warnings-tests.yaml | 53 +++++++++++++++++++++ development/growth_factor_integration.ipynb | 8 ++-- src/hmf/density_field/halofit.py | 4 +- src/hmf/density_field/transfer.py | 2 +- src/hmf/halos/mass_definitions.py | 19 ++++---- src/hmf/helpers/sample.py | 17 ++++--- src/hmf/mass_function/integrate_hmf.py | 28 +++++------ tests/test_fcoll.py | 8 +--- tests/test_filters.py | 9 +--- tests/test_hmf.py | 5 +- tests/test_integrate_hmf.py | 7 +-- tests/test_sample.py | 15 ++---- 12 files changed, 109 insertions(+), 66 deletions(-) create mode 100644 .github/workflows/warnings-tests.yaml diff --git a/.github/workflows/warnings-tests.yaml b/.github/workflows/warnings-tests.yaml new file mode 100644 index 0000000..cb4acbf --- /dev/null +++ b/.github/workflows/warnings-tests.yaml @@ -0,0 +1,53 @@ +name: Warnings Tests + +# Test on all pushes, except when the push is literally just a tag (because we +# tag automatically via CI, and therefore there's no extra code in that push). +# Also, only test on pull requests into master/dev. +on: + push: + tags-ignore: + - 'v*' + pull_request: + branches: + - 'master' + - 'dev' + +jobs: + tests: + if: "!contains(github.event.pull_request.labels.*.name, 'auto-pr')" + env: + ENV_NAME: tests + PYTHON: ${{ matrix.python-version }} + OS: ${{ matrix.os }} + name: Testing + # needs: [linter] + runs-on: ${{ matrix.os }} + strategy: + fail-fast: false + matrix: + os: [ubuntu-latest] + python-version: ["3.11"] + steps: + - uses: actions/checkout@master + with: + fetch-depth: 1 + - uses: actions/setup-python@v4 + with: + python-version: ${{ matrix.python-version }} + - uses: fortran-lang/setup-fortran@v1 + id: setup-fortran + with: + compiler: gcc + version: 13 + - name: Install + run: pip install .[dev] + - name: Run Tests + run: | + python -m pytest -Werror --cov=hmf --cov-config=.coveragerc --cov-report xml:./coverage.xml --durations=25 + - uses: codecov/codecov-action@v3 + if: matrix.os == 'ubuntu-latest' && success() && !contains(github.event.pull_request.labels.*.name, 'auto-pr') + with: + file: ./coverage.xml # optional + fail_ci_if_error: true + verbose: true + token: ${{ secrets.CODECOV_TOKEN }} diff --git a/development/growth_factor_integration.ipynb b/development/growth_factor_integration.ipynb index d935179..b686b12 100644 --- a/development/growth_factor_integration.ipynb +++ b/development/growth_factor_integration.ipynb @@ -93,10 +93,10 @@ " integrand = 1.0 / (np.exp(lna) * cosmo.efunc(_zvec)) ** 3\n", "\n", " if not getvec:\n", - " integral = intg.simpson(np.exp(lna) * integrand, lna,even=\"avg\")\n", + " integral = intg.simpson(np.exp(lna) * integrand, x=lna)\n", " dplus = 5.0 * cosmo.Om0 * cosmo.efunc(z) * integral / 2.0\n", " else:\n", - " integral = intg.cumtrapz(np.exp(lna) * integrand, dx=dlna, initial=0.0)\n", + " integral = intg.cumulative_trapezoid(np.exp(lna) * integrand, dx=dlna, initial=0.0)\n", " dplus = 5.0 * cosmo.Om0 * cosmo.efunc(_zvec) * integral / 2.0\n", "\n", " return dplus" @@ -245,10 +245,10 @@ " integrand = 1.0 / (np.exp(lna) * cosmo.efunc(_zvec)) ** 3\n", "\n", " if not getvec:\n", - " integral = intg.simpson(integrand, np.exp(lna),even=\"avg\")\n", + " integral = intg.simpson(integrand, x=np.exp(lna))\n", " dplus = 5.0 * cosmo.Om0 * cosmo.efunc(z) * integral / 2.0\n", " else:\n", - " integral = intg.cumtrapz(np.exp(lna) * integrand, dx=dlna, initial=0.0)\n", + " integral = intg.cumulative_trapezoid(np.exp(lna) * integrand, dx=dlna, initial=0.0)\n", " dplus = 5.0 * cosmo.Om0 * cosmo.efunc(_zvec) * integral / 2.0\n", "\n", " return dplus" diff --git a/src/hmf/density_field/halofit.py b/src/hmf/density_field/halofit.py index ff05d52..cc6cb7c 100644 --- a/src/hmf/density_field/halofit.py +++ b/src/hmf/density_field/halofit.py @@ -48,7 +48,7 @@ def _get_spec( def get_log_sigma2(lnr): R = np.exp(lnr) integrand = delta_k * np.exp(-((k * R) ** 2)) - return np.log(_simps(integrand, np.log(k))) + return np.log(_simps(integrand, x=np.log(k))) def get_sigma_abs(lnr): return np.abs(get_log_sigma2(lnr)) @@ -77,7 +77,7 @@ def get_sigma_abs(lnr): return knl, n_eff, n_curv -def halofit(k, delta_k, sigma_8=None, z=0, cosmo=None, takahashi=True): +def halofit(k, delta_k, *, sigma_8=None, z=0, cosmo=None, takahashi=True): """ Implementation of HALOFIT (Smith+2003). diff --git a/src/hmf/density_field/transfer.py b/src/hmf/density_field/transfer.py index 6010b68..6ff83ba 100644 --- a/src/hmf/density_field/transfer.py +++ b/src/hmf/density_field/transfer.py @@ -323,5 +323,5 @@ def nonlinear_delta_k(self): .. math:: \Delta_k = \frac{k^3 P_{\rm nl}(k)}{2\pi^2} """ return _hfit( - self.k, self.delta_k, self.sigma_8, self.z, self.cosmo, self.takahashi + self.k, self.delta_k, z=self.z, cosmo=self.cosmo, takahashi=self.takahashi ) diff --git a/src/hmf/halos/mass_definitions.py b/src/hmf/halos/mass_definitions.py index 44765cf..2a11b00 100644 --- a/src/hmf/halos/mass_definitions.py +++ b/src/hmf/halos/mass_definitions.py @@ -180,6 +180,8 @@ def change_definition( if not hasattr(rhos, "__len__"): rhos = [rhos] + + if not hasattr(c, "__len__"): c = [c] c_new = np.array( @@ -356,7 +358,7 @@ def _find_new_concentration(rho_s, halo_density, h=None, x_guess=5.0): # provide lower and upper limits for the root finder. To balance stability and # performance, we do so iteratively: if there is no result within relatively # aggressive limits, we try again with more conservative limits. - args = rho_s, halo_density + target = halo_density / rho_s x = None i = 0 XDELTA_GUESS_FACTORS = [5.0, 10.0, 20.0, 100.0, 10000.0] @@ -366,22 +368,23 @@ def _find_new_concentration(rho_s, halo_density, h=None, x_guess=5.0): def h(x): return np.log(1.0 + x) - x / (1.0 + x) - fnc = ( - lambda x, rhos, density_threshold: rhos * h(x) * 3.0 / x**3 - density_threshold - ) + def fnc(x): + return h(x) * 3.0 / x**3 - target while x is None and i < len(XDELTA_GUESS_FACTORS): try: xmin = x_guess / XDELTA_GUESS_FACTORS[i] xmax = x_guess * XDELTA_GUESS_FACTORS[i] - x = sp.optimize.brentq(fnc, xmin, xmax, args) - except Exception: + x = sp.optimize.brentq(fnc, xmin, xmax) + except Exception as e: + warnings.warn(f"raised following error: {e}") i += 1 if x is None: raise OptimizationException( - "Could not determine x where the density threshold %.2f is satisfied." - % halo_density + "Could not determine x where the density threshold " + f"{halo_density / rho_s:.2f} is satisfied. Largest x-range tried was " + f"x={xmin} -- {xmax} which brackets {fnc(xmin)} -- {fnc(xmax)}" ) return x diff --git a/src/hmf/helpers/sample.py b/src/hmf/helpers/sample.py index 98919e2..ebd7700 100644 --- a/src/hmf/helpers/sample.py +++ b/src/hmf/helpers/sample.py @@ -13,20 +13,24 @@ def _prepare_mf(log_mmin, **mf_kwargs): h = hmf.MassFunction(Mmin=log_mmin, **mf_kwargs) mask = h.ngtm > 0 - icdf = _spline((h.ngtm[mask] / h.ngtm[0])[::-1], np.log10(h.m[mask][::-1]), k=3) + _icdf = np.log10(h.ngtm[mask] / h.ngtm[0]) + icdf = _spline(_icdf[::-1], np.log10(h.m[mask][::-1]), k=3) return icdf, h -def _choose_halo_masses_num(N, icdf, xmin=0): +def _choose_halo_masses_num(N, icdf, xmin=0, rng=None): # Generate random variates from 0 to maxcum - x = np.random.uniform(low=xmin, high=1.0, size=int(N)) + if rng is None: + rng = np.random.default_rng() + x = rng.uniform(low=xmin, high=1.0, size=int(N)) + logm = icdf(np.log10(x)) # Generate halo masses from mf distribution - return 10 ** icdf(x) + return 10**logm -def sample_mf(N, log_mmin, sort=False, **mf_kwargs): +def sample_mf(N, log_mmin, sort=False, rng=None, **mf_kwargs): """ Create a sample of halo masses from a theoretical mass function. @@ -60,8 +64,7 @@ def sample_mf(N, log_mmin, sort=False, **mf_kwargs): >>> m,hmf = sample_mf(1e6,10.0,hmf_model="PS",Mmax=17) """ icdf, h = _prepare_mf(log_mmin, **mf_kwargs) - - m = _choose_halo_masses_num(N, icdf, xmin=h.ngtm.min() / h.ngtm[0]) + m = _choose_halo_masses_num(N, icdf, xmin=h.ngtm.min() / h.ngtm[0], rng=rng) if sort: m.sort() diff --git a/src/hmf/mass_function/integrate_hmf.py b/src/hmf/mass_function/integrate_hmf.py index bb771ce..5fe6d13 100644 --- a/src/hmf/mass_function/integrate_hmf.py +++ b/src/hmf/mass_function/integrate_hmf.py @@ -13,7 +13,7 @@ class NaNException(Exception): pass -def hmf_integral_gtm(M, dndm, mass_density=False): +def hmf_integral_gtm(m, dndm, mass_density=False): """ Cumulatively integrate dn/dm. @@ -52,15 +52,17 @@ def hmf_integral_gtm(M, dndm, mass_density=False): True """ + n = len(m) + # Eliminate NaN's - m = M[np.logical_not(np.isnan(dndm))] - dndm = dndm[np.logical_not(np.isnan(dndm))] + mask = np.isfinite(dndm) + m = m[mask] + dndm = dndm[mask] dndlnm = m * dndm if len(m) < 4: raise NaNException( - "There are too few real numbers in dndm: len(dndm) = %s, #NaN's = %s" - % (len(M), len(M) - len(dndm)) + f"There are too few real numbers in dndm: len(dndm) = {n}, #NaN's = {n - len(m)}" ) # Calculate the mass function (and its integral) from the highest M up to 10**18 @@ -70,13 +72,9 @@ def hmf_integral_gtm(M, dndm, mass_density=False): mf = mf_func(m_upper) if not mass_density: - int_upper = intg.simpson( - np.exp(mf), dx=m_upper[2] - m_upper[1], even="first" - ) + int_upper = intg.simpson(np.exp(mf), dx=m_upper[2] - m_upper[1]) else: - int_upper = intg.simpson( - np.exp(m_upper + mf), dx=m_upper[2] - m_upper[1], even="first" - ) + int_upper = intg.simpson(np.exp(m_upper + mf), dx=m_upper[2] - m_upper[1]) else: int_upper = 0 @@ -84,16 +82,16 @@ def hmf_integral_gtm(M, dndm, mass_density=False): if not mass_density: ngtm = np.concatenate( ( - intg.cumtrapz(dndlnm[::-1], dx=np.log(m[1]) - np.log(m[0]))[::-1], + intg.cumulative_trapezoid(dndlnm[::-1], dx=np.log(m[1] / m[0]))[::-1], np.zeros(1), ) ) else: ngtm = np.concatenate( ( - intg.cumtrapz(m[::-1] * dndlnm[::-1], dx=np.log(m[1]) - np.log(m[0]))[ - ::-1 - ], + intg.cumulative_trapezoid( + m[::-1] * dndlnm[::-1], dx=np.log(m[1]) - np.log(m[0]) + )[::-1], np.zeros(1), ) ) diff --git a/tests/test_fcoll.py b/tests/test_fcoll.py index 968f694..8233d0e 100644 --- a/tests/test_fcoll.py +++ b/tests/test_fcoll.py @@ -80,12 +80,8 @@ def test_ranges_cut(self, peacock, Mmin, Mmax): anl = fcoll_Peacock(np.sqrt(peacock.nu)) num = peacock.rho_gtm / peacock.mean_density0 - err = np.abs((num - anl) / anl)[ - np.logical_and(peacock.m > 10**10, peacock.m < 10**15) - ] - err = err[np.logical_not(np.isnan(err))] - print(np.max(err)) - assert np.max(err) < 0.4 + mask = np.logical_and(peacock.m > 10**10, peacock.m < 10**15) + np.testing.assert_allclose(num[mask], anl[mask], rtol=0.4) @pytest.fixture(scope="class") def ps(self): diff --git a/tests/test_filters.py b/tests/test_filters.py index c3c768c..fa20e2a 100644 --- a/tests/test_filters.py +++ b/tests/test_filters.py @@ -137,12 +137,10 @@ def test_sigma_Rhalf(self, cls): t = 2 + 2 + 1 true = 1.0 / (2 * pi**2 * t * thisr**t) - with warnings.catch_warnings(record=True) as w: + with pytest.warns(UserWarning): # should also raise a warning R = 0.5 s2 = cls.sigma(R)[0] ** 2 - assert w - print(s2, true) assert np.isclose(s2, true) def test_sigma1_Rhalf(self, cls): @@ -151,13 +149,10 @@ def test_sigma1_Rhalf(self, cls): t = 4 + 2 + 1 true = 1.0 / (2 * pi**2 * t * thisr**t) - with warnings.catch_warnings(record=True) as w: + with pytest.warns(UserWarning): # should also raise a warning R = 0.5 s2 = cls.sigma(R, 1)[0] ** 2 - assert w - - print(s2, true) assert np.isclose(s2, true) def test_dlnssdlnr_Rhalf(self, cls): diff --git a/tests/test_hmf.py b/tests/test_hmf.py index 3a8f1e8..689f6b9 100644 --- a/tests/test_hmf.py +++ b/tests/test_hmf.py @@ -1,9 +1,9 @@ """Tests of HMF.""" +import pytest from pytest import raises import numpy as np -import warnings from hmf import MassFunction @@ -47,6 +47,5 @@ def test_str_filter(): def test_mass_nonlinear_outside_range(): h = MassFunction(Mmin=8, Mmax=9) - with warnings.catch_warnings(record=True) as w: + with pytest.warns(UserWarning): assert h.mass_nonlinear > 0 - assert len(w) diff --git a/tests/test_integrate_hmf.py b/tests/test_integrate_hmf.py index 1944d68..5274dbd 100644 --- a/tests/test_integrate_hmf.py +++ b/tests/test_integrate_hmf.py @@ -48,9 +48,10 @@ def test_high_z(self): m = np.logspace(10, 18, 500) dndm = self.tggd(m, 9.0, -1.93, 0.4) ngtm = self.anl_int(m, 9.0, -1.93, 0.4) - - print(ngtm / hmf_integral_gtm(m, dndm)) - assert np.allclose(ngtm, hmf_integral_gtm(m, dndm), rtol=0.03) + mask = m < 10**15 + np.testing.assert_allclose( + ngtm[mask], hmf_integral_gtm(m, dndm)[mask], rtol=0.03, atol=1e-8 + ) # def test_low_mmax_z0(self): # m = np.logspace(10,15,500) diff --git a/tests/test_sample.py b/tests/test_sample.py index e9c1d6a..6a218f1 100644 --- a/tests/test_sample.py +++ b/tests/test_sample.py @@ -1,5 +1,3 @@ -from pytest import raises - import numpy as np from scipy.interpolate import InterpolatedUnivariateSpline as spline @@ -7,8 +5,8 @@ def test_circular(): - np.random.seed(1234) - m, h = sample_mf(1e5, 11, transfer_model="EH") + rng = np.random.default_rng(1234) + m, h = sample_mf(1e5, 11, transfer_model="EH", rng=rng) centres, hist = dndm_from_sample(m, 1e5 / h.ngtm[0]) s = spline(np.log10(h.m), np.log10(h.dndm)) @@ -18,12 +16,9 @@ def test_circular(): def test_mmax_big(): - # raises ValueError because ngtm=0 exactly at m=18 - # due to hard limit of integration in integrate_hmf. - np.random.seed(12345) + rng = np.random.default_rng(12345) - m, h = sample_mf(1e5, 11, transfer_model="EH", Mmax=18) + m, h = sample_mf(1e5, 11, transfer_model="EH", Mmax=18, rng=rng) # centres,hist = dndm_from_sample(m,1e5/h.ngtm[0]) # print centres,hist - with raises(ValueError): - dndm_from_sample(m, 1e5 / h.ngtm[0]) + dndm_from_sample(m, 1e5 / h.ngtm[0])