Fix #2132

package/CHANGELOG

@@ -17,6 +17,8 @@ mm/dd/yy
   * 0.19.3
 
 Fixes
+  * fixed lib.nsgrid segfault for coordinates very close to the box boundary
+    (Issue #2132, PR #2136)
 
 Deprecations
 

package/MDAnalysis/lib/nsgrid.pyx

@@ -539,7 +539,9 @@ cdef class _NSGrid(object):
     cdef ns_int ncoords  # number of coordinates
     cdef ns_int[DIM] ncells  # individual cells in every dimension
     cdef ns_int[DIM] cell_offsets  # Cell Multipliers
-    cdef real[DIM] cellsize  # cell size in every dimension
+    # cellsize MUST be double precision, otherwise coord2cellid() may fail for
+    # coordinates very close to the upper box boundaries! See Issue #2132
+    cdef dreal[DIM] cellsize  # cell size in every dimension
     cdef ns_int nbeads_per_cell  # maximum beads
     cdef ns_int *nbeads  # size (Number of beads in every cell)
     cdef ns_int *beadids  # size * nbeads_per_cell (Beadids in every cell)
@@ -572,14 +574,14 @@ cdef class _NSGrid(object):
 
         # Calculate best cutoff
         self.cutoff = cutoff
-        # First, we add a small margin to the cell size so that we can safely
-        # use the condition d <= cutoff (instead of d < cutoff) for neighbor
-        # search.
-        relative_cutoff_margin = 1.0e-8
-        while self.cutoff == cutoff:
-            self.cutoff = cutoff * (1.0 + relative_cutoff_margin)
-            relative_cutoff_margin *= 10.0
         if not force:
+            # First, we add a small margin to the cell size so that we can safely
+            # use the condition d <= cutoff (instead of d < cutoff) for neighbor
+            # search.
+            relative_cutoff_margin = 1.0e-8
+            while self.cutoff == cutoff:
+                self.cutoff = cutoff * (1.0 + relative_cutoff_margin)
+                relative_cutoff_margin *= 10.0
             bbox_vol = box.c_pbcbox.box[XX][XX] * box.c_pbcbox.box[YY][YY] * box.c_pbcbox.box[YY][YY]
             while bbox_vol/self.cutoff**3 > max_size:
                 self.cutoff *= 1.2

testsuite/MDAnalysisTests/lib/test_nsgrid.py

@@ -47,7 +47,6 @@ def run_grid_search(u, ref_id, cutoff=3):
 
     return searcher.search(searchcoords)
 
-
 def test_pbc_box():
     """Check that PBC box accepts only well-formated boxes"""
     pbc = True
@@ -214,3 +213,20 @@ def test_nsgrid_selfsearch(box, result):
         searchresults = searcher.self_search()
     pairs = searchresults.get_pairs()
     assert_equal(len(pairs)//2, result)
+
+def test_nsgrid_probe_close_to_box_boundary():
+    # FastNS.search used to segfault with this box, cutoff and reference
+    # coordinate prior to PR #2136, so we ensure that this remains fixed.
+    # See Issue #2132 for further information.
+    ref = np.array([[55.783722, 44.190044, -54.16671]], dtype=np.float32)
+    box = np.array([53.785854, 43.951054, 57.17597, 90., 90., 90.], dtype=np.float32)
+    cutoff = 3.0
+    # search within a configuration where we know the expected outcome:
+    conf = np.ones((1, 3), dtype=np.float32)
+    searcher = nsgrid.FastNS(cutoff, conf, box)
+    results = searcher.search(ref)
+    # check if results are as expected:
+    expected_pairs = np.zeros((1, 2), dtype=np.int64)
+    expected_dists = np.array([2.3689647], dtype=np.float64)
+    assert_equal(results.get_pairs(), expected_pairs)
+    assert_allclose(results.get_pair_distances(), expected_dists, rtol=1.e-6)