geographdb_core/algorithms/

persistence.rs

//! Temporal persistence barcode for 4D graphs.
//!
//! Sweeps a temporal window across the time axis and tracks how the strongly
//! connected component landscape evolves. The result is a persistence barcode:
//! each component has a (birth_t, death_t) lifetime, forming the discrete
//! analogue of 0-dimensional persistent homology with time as the filtration
//! parameter.

use crate::algorithms::four_d::{
    strongly_connected_components_4d, GraphNode4D, SpatialRegion, TemporalWindow,
    TraversalContext4D,
};
use glam::Vec3;

/// One measurement point in the temporal persistence sweep.
#[derive(Debug, Clone)]
pub struct TemporalPersistencePoint {
    /// Centre of the temporal window at this measurement.
    pub t: u64,
    /// Nodes active under the (sphere, time) constraints.
    pub active: usize,
    /// Number of strongly connected components.
    pub n_components: usize,
    /// Size of the largest SCC.
    pub largest_size: usize,
    /// `largest_size / active`, or 0.0 when `active == 0`.
    pub fraction_largest: f32,
}

/// One bar in the persistence barcode.
///
/// Represents a connected component that was first observed at `birth` and
/// last observed at `death` (or still alive at sweep end when `death` is
/// `None`).
#[derive(Debug, Clone, PartialEq)]
pub struct TemporalBarcode {
    /// Time slice where this component first appeared.
    pub birth: u64,
    /// Time slice where this component was last seen, or `None` if it
    /// survived to the end of the sweep.
    pub death: Option<u64>,
    /// Largest size this component reached across its lifetime.
    pub peak_size: usize,
}

/// Sweep the temporal window across `time_slices`, computing the SCC landscape
/// at each step.
///
/// The spatial query is a sphere of `sphere_radius` centred at `sphere_center`.
/// `time_half_width` forms the window `[t - half, t + half]` around each slice.
pub fn temporal_persistence_sweep(
    nodes: &[GraphNode4D],
    sphere_center: Vec3,
    sphere_radius: f32,
    time_slices: &[u64],
    time_half_width: u64,
) -> Vec<TemporalPersistencePoint> {
    time_slices
        .iter()
        .map(|&t| {
            let ctx = TraversalContext4D {
                spatial_region: Some(SpatialRegion::Sphere {
                    center: sphere_center,
                    radius: sphere_radius,
                }),
                time_window: Some(TemporalWindow {
                    start: t.saturating_sub(time_half_width),
                    end: t.saturating_add(time_half_width),
                }),
                ..TraversalContext4D::default()
            };

            let sccs = strongly_connected_components_4d(nodes, &ctx);
            let active: usize = sccs.iter().map(|c| c.len()).sum();
            let n_components = sccs.len();
            let largest_size = sccs.iter().map(|c| c.len()).max().unwrap_or(0);
            let fraction_largest = if active == 0 {
                0.0
            } else {
                largest_size as f32 / active as f32
            };

            TemporalPersistencePoint {
                t,
                active,
                n_components,
                largest_size,
                fraction_largest,
            }
        })
        .collect()
}

/// Derive a persistence barcode from a sweep.
///
/// Tracks the number of components across time slices. Each time the count
/// rises, new components are "born"; each time it falls, existing ones "die"
/// (merge or disappear). Returns one `TemporalBarcode` per birth event.
///
/// The `points` must be ordered by increasing `t` (as produced by
/// `temporal_persistence_sweep`).
pub fn compute_temporal_barcode(points: &[TemporalPersistencePoint]) -> Vec<TemporalBarcode> {
    // Track open bars: one entry per "live" component count slot.
    // Each birth event opens a bar; each death event closes the oldest open bar.
    let mut open: Vec<(u64, usize)> = Vec::new(); // (birth_t, peak_size)
    let mut bars: Vec<TemporalBarcode> = Vec::new();

    let mut prev_count = 0usize;

    for pt in points {
        let curr = pt.n_components;

        // New components born this step.
        if curr > prev_count {
            for _ in 0..(curr - prev_count) {
                open.push((pt.t, pt.largest_size));
            }
        }

        // Update peak sizes for all open bars at each step.
        for ob in &mut open {
            ob.1 = ob.1.max(pt.largest_size);
        }

        // Components that died this step — close the most recently opened bars.
        if curr < prev_count {
            let died = prev_count - curr;
            let close_from = open.len().saturating_sub(died);
            for (birth, peak) in open.drain(close_from..) {
                bars.push(TemporalBarcode {
                    birth,
                    death: Some(pt.t),
                    peak_size: peak,
                });
            }
        }

        prev_count = curr;
    }

    // Remaining open bars survive to end of sweep.
    for (birth, peak) in open {
        bars.push(TemporalBarcode {
            birth,
            death: None,
            peak_size: peak,
        });
    }

    bars.sort_by_key(|b| b.birth);
    bars
}

#[cfg(test)]
mod tests {
    use super::*;
    use crate::algorithms::four_d::{GraphProperties, TemporalEdge};

    fn node(id: u64, x: f32, successors: Vec<TemporalEdge>) -> GraphNode4D {
        GraphNode4D {
            id,
            x,
            y: 0.0,
            z: 0.0,
            begin_ts: 0,
            end_ts: 100,
            properties: GraphProperties::new(),
            successors,
        }
    }

    fn edge(dst: u64) -> TemporalEdge {
        TemporalEdge {
            dst,
            weight: 1.0,
            begin_ts: 0,
            end_ts: 100,
        }
    }

    fn big_sphere() -> Vec3 {
        Vec3::ZERO
    }

    // ── temporal_persistence_sweep ───────────────────────────────────────────

    #[test]
    fn test_empty_nodes_zero_components() {
        let pts = temporal_persistence_sweep(&[], big_sphere(), 1000.0, &[50], 10);
        assert_eq!(pts.len(), 1);
        assert_eq!(pts[0].active, 0);
        assert_eq!(pts[0].n_components, 0);
    }

    #[test]
    fn test_disconnected_nodes_n_components_equals_n_active() {
        // 4 isolated nodes → 4 SCCs (each node is its own component).
        let nodes: Vec<GraphNode4D> = (0..4).map(|i| node(i, i as f32, vec![])).collect();
        let pts = temporal_persistence_sweep(&nodes, big_sphere(), 1000.0, &[50], 10);
        assert_eq!(pts[0].n_components, 4);
        assert_eq!(pts[0].active, 4);
        assert_eq!(pts[0].largest_size, 1);
    }

    #[test]
    fn test_fully_connected_cycle_is_one_component() {
        // 0 → 1 → 2 → 0 forms a single SCC.
        let nodes = vec![
            node(0, 0.0, vec![edge(1)]),
            node(1, 1.0, vec![edge(2)]),
            node(2, 2.0, vec![edge(0)]),
        ];
        let pts = temporal_persistence_sweep(&nodes, big_sphere(), 1000.0, &[50], 10);
        assert_eq!(pts[0].n_components, 1);
        assert_eq!(pts[0].largest_size, 3);
        assert!((pts[0].fraction_largest - 1.0).abs() < 1e-6);
    }

    #[test]
    fn test_time_filter_outside_node_validity_gives_zero() {
        // Nodes valid in [0, 100]; window centred at t=500 → nothing active.
        let nodes: Vec<GraphNode4D> = (0..3).map(|i| node(i, i as f32, vec![])).collect();
        let pts = temporal_persistence_sweep(&nodes, big_sphere(), 1000.0, &[500], 10);
        assert_eq!(pts[0].active, 0);
        assert_eq!(pts[0].n_components, 0);
        assert_eq!(pts[0].fraction_largest, 0.0);
    }

    // ── compute_temporal_barcode ─────────────────────────────────────────────

    #[test]
    fn test_barcode_birth_and_death() {
        // Simulate: at t=10 one component appears, at t=50 it merges away.
        let points = vec![
            TemporalPersistencePoint {
                t: 0,
                active: 0,
                n_components: 0,
                largest_size: 0,
                fraction_largest: 0.0,
            },
            TemporalPersistencePoint {
                t: 10,
                active: 3,
                n_components: 1,
                largest_size: 3,
                fraction_largest: 1.0,
            },
            TemporalPersistencePoint {
                t: 30,
                active: 3,
                n_components: 1,
                largest_size: 3,
                fraction_largest: 1.0,
            },
            TemporalPersistencePoint {
                t: 50,
                active: 0,
                n_components: 0,
                largest_size: 0,
                fraction_largest: 0.0,
            },
        ];
        let bars = compute_temporal_barcode(&points);
        assert_eq!(bars.len(), 1);
        assert_eq!(bars[0].birth, 10);
        assert_eq!(bars[0].death, Some(50));
        assert_eq!(bars[0].peak_size, 3);
    }

    #[test]
    fn test_barcode_surviving_component_has_no_death() {
        // Component born at t=10, never dies → death = None.
        let points = vec![
            TemporalPersistencePoint {
                t: 0,
                active: 0,
                n_components: 0,
                largest_size: 0,
                fraction_largest: 0.0,
            },
            TemporalPersistencePoint {
                t: 10,
                active: 2,
                n_components: 1,
                largest_size: 2,
                fraction_largest: 1.0,
            },
            TemporalPersistencePoint {
                t: 20,
                active: 2,
                n_components: 1,
                largest_size: 2,
                fraction_largest: 1.0,
            },
        ];
        let bars = compute_temporal_barcode(&points);
        assert_eq!(bars.len(), 1);
        assert_eq!(bars[0].birth, 10);
        assert_eq!(bars[0].death, None);
    }
}