Merge pull request #62 from astrofrog/fix-segfault

Fix segmentation fault

.github/workflows/main.yml

@@ -13,31 +13,35 @@ jobs:
       runs-on: |
         linux: ubuntu-22.04
       envs: |
-        - linux: py37-test-numpy116
-        - linux: py37-test-numpy117
         - linux: py38-test-numpy118
-        - linux: py39-test
-        - linux: py310-test
+        - linux: py39-test-numpy119
+        - linux: py39-test-numpy120
+        - linux: py39-test-numpy121
+        - linux: py39-test-numpy122
+        - linux: py39-test-numpy123
+        - linux: py310-test-numpy124
+        - linux: py310-test-numpy125
+        - linux: py311-test-numpy126
 
-        - macos: py36-test-numpy113
-        - macos: py36-test-numpy114
-        - macos: py36-test-numpy115
-        # The following two jobs are disabled due to an issue that seems
-        # specific to Python 3.7 and not worth investigating.
-        # - macos: py37-test-numpy116
-        # - macos: py37-test-numpy117
         - macos: py38-test-numpy118
-        - macos: py39-test
-        - macos: py310-test
+        - macos: py39-test-numpy119
+        - macos: py39-test-numpy120
+        - macos: py39-test-numpy121
+        - macos: py39-test-numpy122
+        - macos: py39-test-numpy123
+        - macos: py310-test-numpy124
+        - macos: py310-test-numpy125
+        - macos: py311-test-numpy126
 
-        - windows: py36-test-numpy113
-        - windows: py36-test-numpy114
-        - windows: py36-test-numpy115
-        - windows: py37-test-numpy116
-        - windows: py37-test-numpy117
         - windows: py38-test-numpy118
-        - windows: py39-test
-        - windows: py310-test
+        - windows: py39-test-numpy119
+        - windows: py39-test-numpy120
+        - windows: py39-test-numpy121
+        - windows: py39-test-numpy122
+        - windows: py39-test-numpy123
+        - windows: py310-test-numpy124
+        - windows: py310-test-numpy125
+        - windows: py311-test-numpy126
 
   publish:
     uses: OpenAstronomy/github-actions-workflows/.github/workflows/publish.yml@main
@@ -46,13 +50,10 @@ jobs:
       test_command: pytest --pyargs fast_histogram -m "not hypothesis"
       sdist-runs-on: ubuntu-22.04
       targets: |
-        - cp*-manylinux_i686
         - cp*-manylinux_x86_64
         - cp*-manylinux_aarch64
-        # - cp*-musllinux_i686
         - cp*-musllinux_x86_64
         # - cp*-musllinux_aarch64
-        - pp*-manylinux_i686
         - pp*-manylinux_x86_64
         # - pp*-manylinux_aarch64
         - cp*-macosx_x86_64

fast_histogram/_histogram_core.c

@@ -171,6 +171,7 @@ static PyObject *_histogram1d(PyObject *self, PyObject *args) {
 
       if (tx >= xmin && tx < xmax) {
         ix = (tx - xmin) * normx;
+        if(ix == nx) ix -= 1;
         count[ix] += 1.;
       }
 
@@ -322,6 +323,8 @@ static PyObject *_histogram2d(PyObject *self, PyObject *args) {
       if (tx >= xmin && tx < xmax && ty >= ymin && ty < ymax) {
         ix = (tx - xmin) * normx;
         iy = (ty - ymin) * normy;
+        if(ix == nx) ix -= 1;
+        if(iy == ny) iy -= 1;
         count[iy + ny * ix] += 1.;
       }
 
@@ -630,11 +633,12 @@ static PyObject *_histogramdd(PyObject *self, PyObject *args) {
           in_range = 0;
         } else {
           local_bin_idx = (tx - xmin) * norms[i];
+          if(local_bin_idx == dims[i]) local_bin_idx -= 1;
           bin_idx += stride[i] * local_bin_idx;
         }
       }
       if (in_range){
-        count[bin_idx] += 1;
+count[bin_idx] += 1;
       }
     }
 
@@ -788,6 +792,7 @@ static PyObject *_histogram1d_weighted(PyObject *self, PyObject *args) {
 
       if (tx >= xmin && tx < xmax) {
         ix = (tx - xmin) * normx;
+        if(ix == nx) ix -= 1;
         count[ix] += tw;
       }
 
@@ -951,6 +956,8 @@ static PyObject *_histogram2d_weighted(PyObject *self, PyObject *args) {
       if (tx >= xmin && tx < xmax && ty >= ymin && ty < ymax) {
         ix = (tx - xmin) * normx;
         iy = (ty - ymin) * normy;
+        if(ix == nx) ix -= 1;
+        if(iy == ny) iy -= 1;
         count[iy + ny * ix] += tw;
       }
 
@@ -1264,6 +1271,7 @@ static PyObject *_histogramdd_weighted(PyObject *self, PyObject *args) {
           in_range = 0;
         } else {
           local_bin_idx = (tx - xmin) * norms[i];
+          if(local_bin_idx == dims[i]) local_bin_idx -= 1;
           bin_idx += stride[i] * local_bin_idx;
         }
       }

fast_histogram/tests/test_histogram.py

@@ -9,44 +9,46 @@
 
 from ..histogram import histogram1d, histogram2d, histogramdd
 
-# NOTE: for now we don't test the full range of floating-point values in the
-# tests below, because Numpy's behavior isn't always deterministic in some
-# of the extreme regimes. We should add manual (non-hypothesis and not
-# comparing to Numpy) test cases.
 
-@given(values=arrays(dtype='<f8', shape=st.integers(0, 200),
-                     elements=st.floats(-1000, 1000), unique=True),
+# NOTE: For now we randomly generate values ourselves when comparing to Numpy -
+# ideally we would make use of hypothesis to do this but when we do this we run
+# into issues with the numpy implementation too. We have a separate test to
+# make sure our implementation works with arbitrary hypothesis-generated
+# arrays.
+
+
+
+@given(size=st.integers(0, 50),
        nx=st.integers(1, 10),
        xmin=st.floats(-1e10, 1e10),
        xmax=st.floats(-1e10, 1e10),
        weights=st.booleans(),
        dtype=st.sampled_from(['>f4', '<f4', '>f8', '<f8']))
-@settings(max_examples=500)
-def test_1d_compare_with_numpy(values, nx, xmin, xmax, weights, dtype):
+@settings(max_examples=1000)
+def test_1d_compare_with_numpy(size, nx, xmin, xmax, weights, dtype):
 
-    if xmax <= xmin:
-        return
+    # For now we randomly generate values ourselves - ideally we would make use
+    # of hypothesis to do this but when we do this we run into issues with the
+    # numpy implementation too. We have a separate test to make sure our
+    # implementation works with arbitrary hypothesis-generated arrays.
+    values = np.random.uniform(-1000, 1000, size * 2).astype(dtype)
 
-    values = values.astype(dtype)
+    # Numpy will automatically cast the bounds to the dtype so we should do this
+    # here to get consistent results
+    xmin = float(np.array(xmin, dtype=dtype))
+    xmax = float(np.array(xmax, dtype=dtype))
 
-    size = len(values) // 2
+    if xmax <= xmin:
+        assume(False)
 
     if weights:
         w = values[:size]
     else:
         w = None
 
-    x = values[size:size * 2]
+    x = values[size:]
 
-    try:
-        reference = np.histogram(x, bins=nx, weights=w, range=(xmin, xmax))[0]
-    except ValueError:
-        if 'f4' in str(x.dtype):
-            # Numpy has a bug in certain corner cases
-            # https://github.com/numpy/numpy/issues/11586
-            return
-        else:
-            raise
+    reference = np.histogram(x, bins=nx, weights=w, range=(xmin, xmax))[0]
 
     # First, check the Numpy result because it sometimes doesn't make sense. See
     # bug report https://github.com/numpy/numpy/issues/9435
@@ -74,38 +76,42 @@ def test_1d_compare_with_numpy(values, nx, xmin, xmax, weights, dtype):
     fastdd = histogramdd((x,), bins=nx, weights=w, range=[(xmin, xmax)])
     np.testing.assert_array_equal(fast, fastdd)
 
-@given(values=arrays(dtype='<f8', shape=st.integers(0, 300),
-                     elements=st.floats(-1000, 1000), unique=True),
+@given(size=st.integers(0, 50),
        nx=st.integers(1, 10),
        xmin=st.floats(-1e10, 1e10), xmax=st.floats(-1e10, 1e10),
        ny=st.integers(1, 10),
        ymin=st.floats(-1e10, 1e10), ymax=st.floats(-1e10, 1e10),
        weights=st.booleans(),
        dtype=st.sampled_from(['>f4', '<f4', '>f8', '<f8']))
-@settings(max_examples=500)
-def test_2d_compare_with_numpy(values, nx, xmin, xmax, ny, ymin, ymax, weights, dtype):
+@settings(max_examples=1000)
+def test_2d_compare_with_numpy(size, nx, xmin, xmax, ny, ymin, ymax, weights, dtype):
+
+    # For now we randomly generate values ourselves - ideally we would make use
+    # of hypothesis to do this but when we do this we run into issues with the
+    # numpy implementation too. We have a separate test to make sure our
+    # implementation works with arbitrary hypothesis-generated arrays.
+    values = np.random.uniform(-1000, 1000, size * 3).astype(dtype)
+
+    # Numpy will automatically cast the bounds to the dtype so we should do this
+    # here to get consistent results
+    xmin = float(np.array(xmin, dtype=dtype))
+    xmax = float(np.array(xmax, dtype=dtype))
+    ymin = float(np.array(ymin, dtype=dtype))
+    ymax = float(np.array(ymax, dtype=dtype))
 
     if xmax <= xmin or ymax <= ymin:
         return
 
-    values = values.astype(dtype)
-
-    size = len(values) // 3
-
     if weights:
         w = values[:size]
     else:
         w = None
 
     x = values[size:size * 2]
-    y = values[size * 2:size * 3]
+    y = values[size * 2:]
 
-    try:
-        reference = np.histogram2d(x, y, bins=(nx, ny), weights=w,
-                                   range=((xmin, xmax), (ymin, ymax)))[0]
-    except Exception:
-        # If Numpy fails, we skip the comparison since this isn't our fault
-        return
+    reference = np.histogram2d(x, y, bins=(nx, ny), weights=w,
+                                range=((xmin, xmax), (ymin, ymax)))[0]
 
     # First, check the Numpy result because it sometimes doesn't make sense. See
     # bug report https://github.com/numpy/numpy/issues/9435.
@@ -131,49 +137,43 @@ def test_2d_compare_with_numpy(values, nx, xmin, xmax, ny, ymin, ymax, weights,
                          range=((xmin, xmax), (ymin, ymax)))
     np.testing.assert_array_equal(fast, fastdd)
 
-@given(values=arrays(dtype='<f8', shape=st.integers(0, 1000),
-                     elements=st.floats(-1000, 1000), unique=True),
+@given(size=st.integers(0, 50),
        hist_size=st.integers(1, 1e5),
-       bins=arrays(elements=st.integers(1, 10), shape=(10,), dtype=np.int32),
+       bins=arrays(elements=st.integers(1, 10), shape=st.integers(1, 5), dtype=np.int32),
        ranges=arrays(elements=st.floats(1e-10, 1e5), dtype='<f8',
                      shape=(10,), unique=True),
        weights=st.booleans(),
        dtype=st.sampled_from(['>f4', '<f4', '>f8', '<f8']))
-@settings(max_examples=200, suppress_health_check=[HealthCheck.too_slow], deadline=None)
-def test_dd_compare_with_numpy(values, hist_size, bins, ranges, weights, dtype):
-
-    # To avoid generating huge histograms that take a long time, we only take
-    # as many dimensions as we can such that the total hist_size is still within the
-    # limit. If `hist_size = 1`, we will take all the leading ones in `bins`.
-    _bins = []
-    accum_size = 1
-    for i in range(10):
-        if bins[i] * accum_size > hist_size:
-            break
-        _bins.append(bins[i])
-        accum_size *= bins[i]
-    ndim = len(_bins)
-    values = values.astype(dtype)
+@settings(max_examples=1000)
+def test_dd_compare_with_numpy(size, hist_size, bins, ranges, weights, dtype):
+
+    ndim = len(bins)
+
+    # For now we randomly generate values ourselves - ideally we would make use
+    # of hypothesis to do this but when we do this we run into issues with the
+    # numpy implementation too. We have a separate test to make sure our
+    # implementation works with arbitrary hypothesis-generated arrays.
+    values = np.random.uniform(-1000, 1000, size * (ndim + 1)).astype(dtype)
+
     ranges = ranges.astype(dtype)
     ranges = ranges[:ndim]
+
     # Ranges are symmetric because otherwise the probability of samples falling inside
     # is just too small and we would just be testing a bunch of empty histograms.
     ranges = np.vstack((-ranges, ranges)).T
 
-    size = len(values) // (ndim + 1)
+    # Numpy will automatically cast the bounds to the dtype so we should do this
+    # here to get consistent results
+    ranges = ranges.astype(dtype)
 
     if weights:
         w = values[:size]
     else:
         w = None
 
     sample = tuple(values[size*(i+1):size*(i+2)] for i in range(ndim))
-    # for simplicity using the same range in all dimensions
-    try:
-        reference = np.histogramdd(sample, bins=_bins, weights=w, range=ranges)[0]
-    except Exception:
-        # If Numpy fails, we skip the comparison since this isn't our fault
-        return
+
+    reference = np.histogramdd(sample, bins=bins, weights=w, range=ranges)[0]
 
     # First, check the Numpy result because it sometimes doesn't make sense. See
     # bug report https://github.com/numpy/numpy/issues/9435.
@@ -188,7 +188,7 @@ def test_dd_compare_with_numpy(values, hist_size, bins, ranges, weights, dtype):
         n_inside = np.sum(inside)
         assume(n_inside == np.sum(reference))
 
-    fast = histogramdd(sample, bins=_bins, weights=w, range=ranges)
+    fast = histogramdd(sample, bins=bins, weights=w, range=ranges)
 
     if sample[0].dtype.kind == 'f' and sample[0].dtype.itemsize == 4:
         rtol = 1e-7

pyproject.toml

@@ -2,10 +2,12 @@
 requires = ["setuptools",
             "setuptools_scm",
             "wheel",
-            "oldest-supported-numpy"]
+            "oldest-supported-numpy;python_version<'3.9'",
+            "numpy>=1.25;python_version>='3.9'"]
 build-backend = 'setuptools.build_meta'
 
 [tool.cibuildwheel.linux]
+skip = "pp310* pp311* pp312*"
 archs = ["auto", "aarch64"]
 
 [tool.cibuildwheel.macos]

tox.ini

@@ -26,11 +26,19 @@ deps =
     numpy116: numpy==1.16.*
     numpy117: numpy==1.17.*
     numpy118: numpy==1.18.*
+    numpy119: numpy==1.19.*
+    numpy120: numpy==1.20.*
+    numpy121: numpy==1.21.*
+    numpy122: numpy==1.22.*
+    numpy123: numpy==1.23.*
+    numpy124: numpy==1.24.*
+    numpy125: numpy==1.25.*
+    numpy126: numpy==1.26.*
 extras =
     test
 commands =
     pip freeze
-    pytest --pyargs fast_histogram {posargs}
+    pytest --pyargs fast_histogram --hypothesis-show-statistics {posargs}
 
 [testenv:style]
 skip_install = true