diff --git a/rascaline/src/systems/neighbors.rs b/rascaline/src/systems/neighbors.rs index ec412c41e..e72141787 100644 --- a/rascaline/src/systems/neighbors.rs +++ b/rascaline/src/systems/neighbors.rs @@ -132,9 +132,9 @@ impl CellList { // number of cells to search in each direction to make sure all possible // pairs below the cutoff are accounted for. let mut n_search = [ - f64::trunc(cutoff * n_cells[0] / distances_between_faces[0]) as i32, - f64::trunc(cutoff * n_cells[1] / distances_between_faces[1]) as i32, - f64::trunc(cutoff * n_cells[2] / distances_between_faces[2]) as i32, + f64::ceil(cutoff * n_cells[0] / distances_between_faces[0]) as i32, + f64::ceil(cutoff * n_cells[1] / distances_between_faces[1]) as i32, + f64::ceil(cutoff * n_cells[2] / distances_between_faces[2]) as i32, ]; let n_cells = [ @@ -407,8 +407,7 @@ impl NeighborsList { mod tests { use approx::assert_ulps_eq; - use crate::Matrix3; - + use crate::types::Matrix3; use super::*; #[test] @@ -525,4 +524,39 @@ mod tests { assert_eq!(&pair.cell_shift_indices, shifts); } } + + #[test] + fn non_cubic_cell() { + let cell = UnitCell::from(Matrix3::new([ + [4.26, -2.45951215, 0.0], + [2.13, 1.22975607, 0.0], + [0.0, 0.0, 50.0], + ])); + let positions = [ + Vector3D::new(1.42, 0.0, 0.0), + Vector3D::new(2.84, 0.0, 0.0), + Vector3D::new(3.55, -1.22975607, 0.0), + Vector3D::new(4.97, -1.22975607, 0.0), + ]; + let neighbors = NeighborsList::new(&positions, cell, 6.4); + + assert_eq!(neighbors.pairs.len(), 90); + + let previously_missing = [ + (0, 3, [-2, 0, 0]), + (0, 3, [-2, 1, 0]), + (0, 3, [-2, 2, 0]), + ]; + + for missing in previously_missing { + let mut found = false; + for pair in &neighbors.pairs { + if pair.first == missing.0 && pair.second == missing.1 + && pair.cell_shift_indices == missing.2 { + found = true; + } + } + assert!(found, "could not find pair {:?}", missing); + } + } }