diff --git a/.github/workflows/rust.yml b/.github/workflows/rust.yml index 14470bd..29c8e50 100644 --- a/.github/workflows/rust.yml +++ b/.github/workflows/rust.yml @@ -12,9 +12,17 @@ jobs: override: true components: rustfmt, clippy profile: minimal + + - uses: actions/setup-python@v2 + with: + python-version: 3.9 + architecture: x64 + - name: Install deps run: | sudo apt-get -y update sudo apt-get -y install build-essential libssl-dev + pip install numpy + - run: cargo build --features geos/static --verbose - run: cargo test --features geos/static --verbose diff --git a/Cargo.toml b/Cargo.toml index eabebdd..d8c5885 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -20,16 +20,21 @@ path = "src/devel/make_bitmap.rs" [dependencies] geos = { version = "8" } lazy_static = "1.4.0" -numpy = "0.14" +numpy = { version = "0.14", features = [ "rayon" ] } pyo3 = { version = "0.14" } roaring = "0.7.0" rust-embed = "5.9.0" xz2 = "0.1.6" -openssl = { version = "0.10", features = ["vendored"] } +openssl = { version = "0.10" } +ndarray = { version = "0.15.3", features = [ "rayon" ] } + +[dev-dependencies] +rayon = "1" [patch.crates-io] geos-sys = { version = "2", git = "https://github.com/gauteh/geos-sys", branch = "static-build" } geos = { version = "8", git = "https://github.com/gauteh/geos", branch = "static" } +#geos-sys = { path = "../../dev/geos-sys" } [build-dependencies] reqwest = { version = "0.11", features = [ "blocking" ] } diff --git a/pyproject.toml b/pyproject.toml index f6575d4..f634445 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -3,4 +3,4 @@ requires = ["maturin>=0.11,<0.12"] build-backend = "maturin" [tool.maturin] -cargo-extra-args = "--features extension-module,geos/static" +cargo-extra-args = "--features extension-module,geos/static,openssl/vendored" diff --git a/src/lib.rs b/src/lib.rs index de1f2d4..84979d1 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -127,6 +127,20 @@ impl RoaringLandmask { ) .to_owned() } + + pub fn contains_many_par( + &self, + py: Python, + x: PyReadonlyArrayDyn, + y: PyReadonlyArrayDyn, + ) -> Py> { + let x = x.as_array(); + let y = y.as_array(); + + use ndarray::Zip; + let contains = Zip::from(&x).and(&y).par_map_collect(|x, y| self.contains(*x, *y)); + PyArray::from_owned_array(py, contains).to_owned() + } } #[cfg(test)] @@ -157,4 +171,72 @@ mod tests { b.iter(|| mask.contains(5., 65.6)) } + + #[bench] + fn test_contains_many(b: &mut Bencher) { + let mask = RoaringLandmask::new().unwrap(); + + let (x, y): (Vec, Vec) = (0..360 * 2) + .map(|v| v as f64 * 0.5 - 180.) + .map(|x| { + (0..180 * 2) + .map(|y| y as f64 * 0.5 - 90.) + .map(move |y| (x, y)) + }) + .flatten() + .unzip(); + + pyo3::prepare_freethreaded_python(); + pyo3::Python::with_gil(|py| { + + let x = PyArray::from_vec(py, x); + let y = PyArray::from_vec(py, y); + + println!("testing {} points..", x.len()); + + b.iter(|| { + let len = x.len(); + + let x = x.to_dyn().readonly(); + let y = y.to_dyn().readonly(); + + let onland = mask.contains_many(py, x, y); + assert!(onland.as_ref(py).len() == len); + }) + }) + } + + #[bench] + fn test_contains_many_par(b: &mut Bencher) { + let mask = RoaringLandmask::new().unwrap(); + + let (x, y): (Vec, Vec) = (0..360 * 2) + .map(|v| v as f64 * 0.5 - 180.) + .map(|x| { + (0..180 * 2) + .map(|y| y as f64 * 0.5 - 90.) + .map(move |y| (x, y)) + }) + .flatten() + .unzip(); + + pyo3::prepare_freethreaded_python(); + pyo3::Python::with_gil(|py| { + + let x = PyArray::from_vec(py, x); + let y = PyArray::from_vec(py, y); + + println!("testing {} points..", x.len()); + + b.iter(|| { + let len = x.len(); + + let x = x.to_dyn().readonly(); + let y = y.to_dyn().readonly(); + + let onland = mask.contains_many_par(py, x, y); + assert!(onland.as_ref(py).len() == len); + }) + }) + } } diff --git a/src/shapes.rs b/src/shapes.rs index a1d1776..f13d397 100644 --- a/src/shapes.rs +++ b/src/shapes.rs @@ -23,14 +23,24 @@ impl Drop for Gshhg { } } -// PreparedGeometry is Send, Geometry is Send. *mut Geometry is never modified. +// PreparedGeometry is Send+Sync, Geometry is Send+Sync. *mut Geometry is never modified. unsafe impl Send for Gshhg {} +unsafe impl Sync for Gshhg {} + +// `PreparededGeometry::contains` needs a call to `contains` before it is thread-safe: +// https://github.com/georust/geos/issues/95 +fn warmup_prepped(prepped: &PreparedGeometry<'_>) { + let point = CoordSeq::new_from_vec(&[&[0.0, 0.0]]).unwrap(); + let point = Geometry::create_point(point).unwrap(); + prepped.contains(&point).unwrap(); +} impl Clone for Gshhg { fn clone(&self) -> Self { let gptr = self.geom.clone(); debug_assert!(gptr != self.geom); let prepped = unsafe { (&*gptr).to_prepared_geom().unwrap() }; + warmup_prepped(&prepped); Gshhg { geom: gptr, @@ -45,6 +55,7 @@ impl Gshhg { let gptr = Box::into_raw(bxd); let prepped = unsafe { (&*gptr).to_prepared_geom() } .map_err(|_| io::Error::new(io::ErrorKind::InvalidData, "cannot prepare geomtry"))?; + warmup_prepped(&prepped); Ok(Gshhg { geom: gptr, @@ -53,15 +64,19 @@ impl Gshhg { } pub fn from_compressed>(path: P) -> io::Result { + let g = Gshhg::get_geometry_from_compressed(path)?; + + Gshhg::from_geom(g) + } + + pub fn get_geometry_from_compressed>(path: P) -> io::Result> { let fd = File::open(path)?; let fd = io::BufReader::new(fd); let mut fd = xz2::bufread::XzDecoder::new(fd); let mut buf = Vec::new(); fd.read_to_end(&mut buf)?; - let g = geos::Geometry::new_from_wkb(&buf).unwrap(); - - Gshhg::from_geom(g) + Ok(geos::Geometry::new_from_wkb(&buf).unwrap()) } } diff --git a/tests/test_geos_par_prepped.rs b/tests/test_geos_par_prepped.rs new file mode 100644 index 0000000..4460fa8 --- /dev/null +++ b/tests/test_geos_par_prepped.rs @@ -0,0 +1,54 @@ +use roaring_landmask::Gshhg; +use geos::{CoordSeq, Geom, Geometry}; + +#[ignore] +#[test] +fn test_par_prepped_no_warmup() { + use rayon::prelude::*; + + // let s = Gshhg::new().unwrap(); + // let prepped = &s.prepped; + + let g = Gshhg::get_geometry_from_compressed( + "gshhs/gshhs_f_-180.000000E-90.000000N180.000000E90.000000N.wkb.xz", + ).unwrap(); + let prepped = g.to_prepared_geom().unwrap(); + + (0..10000).into_par_iter().for_each(|k| { + + let x = k % 180; + let y = (k / 180) % 90; + + + let point = CoordSeq::new_from_vec(&[&[x as f64, y as f64]]).unwrap(); + let point = Geometry::create_point(point).unwrap(); + prepped.contains(&point).unwrap(); + }); +} + +#[test] +fn test_par_prepped_with_warmup() { + use rayon::prelude::*; + + // let s = Gshhg::new().unwrap(); + // let prepped = &s.prepped; + + let g = Gshhg::get_geometry_from_compressed( + "gshhs/gshhs_f_-180.000000E-90.000000N180.000000E90.000000N.wkb.xz", + ).unwrap(); + let prepped = g.to_prepared_geom().unwrap(); + + let point = CoordSeq::new_from_vec(&[&[10., 50.]]).unwrap(); + let point = Geometry::create_point(point).unwrap(); + prepped.contains(&point).unwrap(); + + (0..10000).into_par_iter().for_each(|k| { + let x = k % 180; + let y = (k / 180) % 90; + + + let point = CoordSeq::new_from_vec(&[&[x as f64, y as f64]]).unwrap(); + let point = Geometry::create_point(point).unwrap(); + prepped.contains(&point).unwrap(); + }); +} diff --git a/tests/test_roaring.py b/tests/test_roaring.py index e2ec849..2e49ee3 100644 --- a/tests/test_roaring.py +++ b/tests/test_roaring.py @@ -21,3 +21,37 @@ def test_landmask_many(benchmark): print ("points:", len(xx.ravel())) benchmark(l.contains_many, xx.ravel(), yy.ravel()) + +def test_landmask_many_few(benchmark): + l = RoaringLandmask.new() + + x = np.linspace(-180, 180, 10) + y = np.linspace(-90, 90, 5) + + xx, yy = np.meshgrid(x,y) + + print ("points:", len(xx.ravel())) + benchmark(l.contains_many, xx.ravel(), yy.ravel()) + +def test_landmask_many_par(benchmark): + l = RoaringLandmask.new() + + x = np.arange(-180, 180, .5) + y = np.arange(-90, 90, .5) + + xx, yy = np.meshgrid(x,y) + + print ("points:", len(xx.ravel())) + benchmark(l.contains_many_par, xx.ravel(), yy.ravel()) + +def test_landmask_many_par_few(benchmark): + l = RoaringLandmask.new() + + x = np.linspace(-180, 180, 10) + y = np.linspace(-90, 90, 5) + + xx, yy = np.meshgrid(x,y) + + print ("points:", len(xx.ravel())) + benchmark(l.contains_many_par, xx.ravel(), yy.ravel()) +