(ml5717) Refactored to make dispersal more explicit

Author: juntyr

TODO

@@ -1,6 +1,11 @@
-x check refactored Gillespie implementation
-x check + add contract for alias method dispersal
-x add Gillespie without self-dispersal
+- implement event transform from GPU to CPU
+- implement event deduplication (global vs local? - needs to be global for speciation)
+- implement simple stop condition for time sliced kernel launching
+- implement declarative event filtering which can be put onto the GPU as well
+- implement declarative kernel generic substitution
+- refactor cfg feature cuda and os cudas
+- implement verification
+
 - describe Gillespie in report
 - describe data parallelism idea
 - describe exponential sampling idea

necsim-core/src/cogs/active_lineage_sampler.rs

@@ -3,7 +3,7 @@ use super::{
     RngCore,
 };
 
-use crate::landscape::Location;
+use crate::landscape::IndexedLocation;
 use crate::simulation::partial::active_lineager_sampler::PartialSimulation;
 
 #[allow(clippy::inline_always, clippy::inline_fn_without_body)]
@@ -41,22 +41,22 @@ pub trait ActiveLineageSampler<
             simulation.lineage_store[reference.clone()].time_of_last_event() == event_time
         },
     }, "updates the time of the last event of the returned lineage to the time of the event")]
-    fn pop_active_lineage_location_event_time(
+    fn pop_active_lineage_indexed_location_event_time(
         &mut self,
         time: f64,
         simulation: &mut PartialSimulation<H, G, D, R, S, C, E>,
         rng: &mut G,
-    ) -> Option<(R, Location, f64)>;
+    ) -> Option<(R, IndexedLocation, f64)>;
 
     #[debug_requires(time >= 0.0_f64, "time is non-negative")]
     /*#[debug_requires(( TODO: How can we assert this only for coherent lineage stores?
         simulation.lineage_store.get_active_lineages_at_location(&location).len() <
         (simulation.habitat.get_habitat_at_location(&location) as usize)
     ), "location has habitat capacity for the lineage")]*/
-    fn push_active_lineage_to_location(
+    fn push_active_lineage_to_indexed_location(
         &mut self,
         lineage_reference: R,
-        location: Location,
+        indexed_location: IndexedLocation,
         time: f64,
         simulation: &mut PartialSimulation<H, G, D, R, S, C, E>,
         rng: &mut G,

necsim-core/src/cogs/coalescence_sampler.rs

@@ -1,5 +1,5 @@
 use crate::cogs::RngCore;
-use crate::landscape::Location;
+use crate::landscape::{IndexedLocation, Location};
 
 use super::{Habitat, LineageReference, LineageStore};
 
@@ -9,12 +9,12 @@ pub trait CoalescenceSampler<H: Habitat, G: RngCore, R: LineageReference<H>, S:
     core::fmt::Debug
 {
     #[must_use]
-    #[debug_requires(habitat.get_habitat_at_location(location) > 0, "location is habitable")]
+    #[debug_requires(habitat.get_habitat_at_location(&location) > 0, "location is habitable")]
     fn sample_optional_coalescence_at_location(
         &self,
-        location: &Location,
+        location: Location,
         habitat: &H,
         lineage_store: &S,
         rng: &mut G,
-    ) -> Option<R>;
+    ) -> (IndexedLocation, Option<R>);
 }

necsim-core/src/cogs/event_sampler.rs

@@ -4,7 +4,7 @@ use super::{
     CoalescenceSampler, DispersalSampler, Habitat, LineageReference, LineageStore, RngCore,
 };
 use crate::event::{Event, EventType};
-use crate::landscape::Location;
+use crate::landscape::IndexedLocation;
 use crate::simulation::partial::event_sampler::PartialSimulation;
 
 #[allow(clippy::inline_always, clippy::inline_fn_without_body)]
@@ -31,10 +31,10 @@ pub trait EventSampler<
         "event occurs for lineage_reference"
     )]
     #[debug_ensures(ret.time() == event_time, "event occurs at event_time")]
-    fn sample_event_for_lineage_at_location_time(
+    fn sample_event_for_lineage_at_indexed_location_time(
         &self,
         lineage_reference: R,
-        location: Location,
+        indexed_location: IndexedLocation,
         event_time: f64,
         simulation: &PartialSimulation<H, G, D, R, S, C>,
         rng: &mut G,

necsim-core/src/cogs/habitat.rs

@@ -1,5 +1,5 @@
 use crate::landscape::LandscapeExtent;
-use crate::landscape::Location;
+use crate::landscape::{IndexedLocation, Location};
 
 #[allow(clippy::inline_always, clippy::inline_fn_without_body)]
 #[contract_trait]
@@ -33,11 +33,13 @@ pub trait Habitat: core::fmt::Debug {
 #[contract_trait]
 pub trait HabitatToU64Injection: Habitat {
     #[must_use]
-    #[debug_requires(self.get_extent().contains(location), "location is inside habitat extent")]
-    #[debug_requires(index_at_location < self.get_habitat_at_location(location) as usize)]
-    fn map_indexed_location_to_u64_injective(
-        &self,
-        location: &Location,
-        index_at_location: usize,
-    ) -> u64;
+    #[debug_requires(
+        self.get_extent().contains(indexed_location.location()),
+        "location is inside habitat extent"
+    )]
+    #[debug_requires(
+        indexed_location.index() < self.get_habitat_at_location(indexed_location.location()),
+        "index is within the location's habitat capacity"
+    )]
+    fn map_indexed_location_to_u64_injective(&self, indexed_location: &IndexedLocation) -> u64;
 }

necsim-core/src/cogs/lineage_store/contract.rs

@@ -1,3 +1,5 @@
+use core::convert::TryFrom;
+
 use crate::cogs::{CoherentLineageStore, Habitat, LineageReference};
 use crate::landscape::Location;
 
@@ -18,19 +20,14 @@ pub fn explicit_lineage_store_lineage_at_location_contract<
         None => return false,
     };
 
-    let location = match lineage.location() {
-        Some(location) => location,
+    let indexed_location = match lineage.indexed_location() {
+        Some(indexed_location) => indexed_location,
         None => return false,
     };
 
-    let lineages_at_location = &store.get_active_lineages_at_location(location);
-
-    let index_at_location = match lineage.index_at_location() {
-        Some(index_at_location) => index_at_location,
-        None => return false,
-    };
+    let lineages_at_location = &store.get_active_lineages_at_location(indexed_location.location());
 
-    match lineages_at_location.get(index_at_location) {
+    match lineages_at_location.get(indexed_location.index() as usize) {
         Some(reference_at_location) => reference_at_location == &input_reference,
         None => false,
     }
@@ -50,8 +47,15 @@ pub(super) fn explicit_lineage_store_invariant_contract<
     lineages_at_location
         .iter()
         .enumerate()
-        .all(|(i, reference)| {
-            store[reference.clone()].location() == Some(location)
-                && store[reference.clone()].index_at_location() == Some(i)
+        .all(|(index, reference)| {
+            match (
+                u32::try_from(index),
+                store[reference.clone()].indexed_location(),
+            ) {
+                (Ok(index), Some(indexed_location)) => {
+                    indexed_location.location() == location && indexed_location.index() == index
+                }
+                _ => false,
+            }
         })
 }

necsim-core/src/cogs/lineage_store/mod.rs

@@ -1,7 +1,7 @@
 use core::ops::Index;
 
 use super::{Habitat, LineageReference};
-use crate::landscape::Location;
+use crate::landscape::{IndexedLocation, Location};
 use crate::lineage::Lineage;
 
 mod contract;
@@ -58,8 +58,9 @@ pub trait CoherentLineageStore<H: Habitat, R: LineageReference<H>>: LineageStore
         contract::explicit_lineage_store_invariant_contract(self, &location),
         "invariant of lineage-location bijection holds"
     )]
+    #[debug_ensures(self[old(reference.clone())].is_active(), "lineage was activated")]
     #[debug_ensures(
-        self[old(reference.clone())].location() == Some(&old(location.clone())),
+        self[old(reference.clone())].indexed_location().map(IndexedLocation::location) == Some(&old(location.clone())),
         "lineage was added to location"
     )]
     #[debug_ensures(
@@ -82,24 +83,25 @@ pub trait CoherentLineageStore<H: Habitat, R: LineageReference<H>>: LineageStore
     )]
     #[debug_requires(
         contract::explicit_lineage_store_invariant_contract(
-            self, self[reference.clone()].location().unwrap()
+            self, self[reference.clone()].indexed_location().unwrap().location()
         ), "invariant of lineage-location bijection holds"
     )]
+    #[debug_ensures(!self[old(reference.clone())].is_active(), "lineage was deactivated")]
     #[debug_ensures(
         !contract::explicit_lineage_store_lineage_at_location_contract(
             self, old(reference.clone())
         ), "lineage was removed from the location and index it references"
     )]
     #[debug_ensures(
         contract::explicit_lineage_store_invariant_contract(
-            self, &old(self[reference.clone()].location().unwrap().clone())
+            self, &old(self[reference.clone()].indexed_location().unwrap().location().clone())
         ), "maintains invariant of lineage-location bijection"
     )]
     #[debug_ensures(
-        ret == old(self[reference.clone()].location().unwrap().clone()),
-        "returns the individual's prior location"
+        ret == old(self[reference.clone()].indexed_location().unwrap().clone()),
+        "returns the individual's prior indexed_location"
     )]
-    fn pop_lineage_from_its_location(&mut self, reference: R) -> Location;
+    fn pop_lineage_from_its_location(&mut self, reference: R) -> IndexedLocation;
 }
 
 #[allow(clippy::inline_always, clippy::inline_fn_without_body)]
@@ -108,59 +110,24 @@ pub trait CoherentLineageStore<H: Habitat, R: LineageReference<H>>: LineageStore
 pub trait IncoherentLineageStore<H: Habitat, R: LineageReference<H>>: LineageStore<H, R> {
     #[debug_requires(self.get(reference.clone()).is_some(), "lineage reference is valid")]
     #[debug_requires(!self[reference.clone()].is_active(), "lineage is inactive")]
-    /*#[debug_requires(
-        !contract::explicit_lineage_store_lineage_at_location_contract(self, reference.clone()),
-        "lineage is not at the location and index it references"
-    )]
-    #[debug_requires(
-        contract::explicit_lineage_store_invariant_contract(self, &location),
-        "invariant of lineage-location bijection holds"
-    )]
+    #[debug_ensures(self[old(reference.clone())].is_active(), "lineage was activated")]
     #[debug_ensures(
-        self[old(reference.clone())].location() == Some(&old(location.clone())),
-        "lineage was added to location"
-    )]
-    #[debug_ensures(
-        contract::explicit_lineage_store_lineage_at_location_contract(
-            self, old(reference.clone())
-        ), "lineage is at the location and index it references"
+        self[old(reference.clone())].indexed_location() == Some(&old(indexed_location.clone())),
+        "lineage was added to indexed_location"
     )]
-    #[debug_ensures(
-        contract::explicit_lineage_store_invariant_contract(self, &old(location.clone())),
-        "maintains invariant of lineage-location bijection"
-    )]*/
-    fn insert_lineage_to_location_at_index(
+    fn insert_lineage_to_indexed_location(
         &mut self,
         reference: R,
-        location: Location,
-        index_at_location: usize,
+        indexed_location: IndexedLocation,
     );
 
     #[must_use]
     #[debug_requires(self.get(reference.clone()).is_some(), "lineage reference is valid")]
     #[debug_requires(self[reference.clone()].is_active(), "lineage is active")]
-    /*#[debug_requires(
-        contract::explicit_lineage_store_lineage_at_location_contract(self, reference.clone()),
-        "lineage is at the location and index it references"
-    )]
-    #[debug_requires(
-        contract::explicit_lineage_store_invariant_contract(
-            self, self[reference.clone()].location().unwrap()
-        ), "invariant of lineage-location bijection holds"
-    )]
-    #[debug_ensures(
-        !contract::explicit_lineage_store_lineage_at_location_contract(
-            self, old(reference.clone())
-        ), "lineage was removed from the location and index it references"
-    )]
-    #[debug_ensures(
-        contract::explicit_lineage_store_invariant_contract(
-            self, &old(self[reference.clone()].location().unwrap().clone())
-        ), "maintains invariant of lineage-location bijection"
-    )]*/
+    #[debug_ensures(!self[old(reference.clone())].is_active(), "lineage was deactivated")]
     #[debug_ensures(
-        ret == old(self[reference.clone()].location().unwrap().clone()),
+        ret == old(self[reference.clone()].indexed_location().unwrap().clone()),
         "returns the individual's prior location"
     )]
-    fn extract_lineage_from_its_location(&mut self, reference: R) -> Location;
+    fn extract_lineage_from_its_location(&mut self, reference: R) -> IndexedLocation;
 }

necsim-core/src/event.rs

@@ -1,7 +1,7 @@
 use core::marker::PhantomData;
 
 use crate::cogs::{Habitat, LineageReference};
-use crate::landscape::Location;
+use crate::landscape::IndexedLocation;
 
 pub struct Event<H: Habitat, R: LineageReference<H>> {
     time: f64,
@@ -46,8 +46,8 @@ impl<H: Habitat, R: LineageReference<H>> Event<H, R> {
 pub enum EventType<H: Habitat, R: LineageReference<H>> {
     Speciation,
     Dispersal {
-        origin: Location,
-        target: Location,
+        origin: IndexedLocation,
+        target: IndexedLocation,
         coalescence: Option<R>,
         _marker: PhantomData<H>,
     },

necsim-core/src/landscape/location.rs

@@ -23,3 +23,39 @@ impl Location {
         self.y
     }
 }
+
+impl From<IndexedLocation> for Location {
+    fn from(indexed_location: IndexedLocation) -> Location {
+        indexed_location.location
+    }
+}
+
+#[derive(Eq, PartialEq, Clone, Hash, Debug)]
+#[cfg_attr(feature = "cuda", derive(DeviceCopy))]
+#[allow(clippy::module_name_repetitions)]
+pub struct IndexedLocation {
+    location: Location,
+    pub(crate) index: u32,
+}
+
+impl IndexedLocation {
+    #[must_use]
+    #[debug_ensures(
+        ret.location.x() == old(location.x()) && ret.location.y() == old(location.y()),
+        "stores location"
+    )]
+    #[debug_ensures(ret.index() == index, "stores index")]
+    pub fn new(location: Location, index: u32) -> Self {
+        Self { location, index }
+    }
+
+    #[must_use]
+    pub fn location(&self) -> &Location {
+        &self.location
+    }
+
+    #[must_use]
+    pub fn index(&self) -> u32 {
+        self.index
+    }
+}

necsim-core/src/landscape/mod.rs

@@ -2,4 +2,4 @@ mod extent;
 mod location;
 
 pub use extent::LandscapeExtent;
-pub use location::Location;
+pub use location::{IndexedLocation, Location};