grid1d 0.5.2

//! Property-based tests for grid1d using proptest.
//!
//! These tests verify mathematical invariants that must hold for all valid inputs,
//! providing high confidence in the correctness of the library.

use grid1d::{
    FindIntervalIdOfPoint, Grid1D, Grid1DTrait, Grid1DUniform, Grid1DUnion, HasCoords1D,
    HasDomain1D, HasIntervalIdRange,
    intervals::{
        GetLowerBoundValue, GetUpperBoundValue, IntervalClosed, IntervalFinitePositiveLengthTrait,
        bounded::IntervalFromBounds,
    },
    scalars::{IntervalId, NumIntervals, PositiveNumPoints1D},
};
use proptest::prelude::*;
use sorted_vec::partial::SortedSet;
use std::ops::Deref;
use try_create::TryNew;

// =============================================================================
// Strategy Generators
// =============================================================================

/// Generate valid number of intervals (at least 1)
fn num_intervals_strategy() -> impl Strategy<Value = usize> {
    1..1000usize
}

/// Generate valid domain bounds (lo < hi)
fn domain_bounds_strategy() -> impl Strategy<Value = (f64, f64)> {
    (-1000.0..1000.0f64, -1000.0..1000.0f64).prop_filter_map(
        "lo must be strictly less than hi",
        |(a, b)| {
            if a < b {
                Some((a, b))
            } else if b < a {
                Some((b, a))
            } else {
                None
            }
        },
    )
}

// =============================================================================
// Grid1DUniform Property Tests
// =============================================================================

proptest! {
    #![proptest_config(ProptestConfig::with_cases(500))]

    /// Property: Uniform grid has exactly n+1 points for n intervals
    #[test]
    fn uniform_grid_point_count(
        (lo, hi) in domain_bounds_strategy(),
        n in num_intervals_strategy()
    ) {
        let domain = IntervalClosed::new(lo, hi);
        let num_intervals = NumIntervals::try_new(n).unwrap();
        let grid = Grid1DUniform::new(domain, num_intervals);

        prop_assert_eq!(*grid.num_points().as_ref(), n + 1);
        prop_assert_eq!(*grid.num_intervals().as_ref(), n);
    }

    /// Property: First coordinate equals domain lower bound
    #[test]
    fn uniform_grid_first_coord_equals_lower_bound(
        (lo, hi) in domain_bounds_strategy(),
        n in num_intervals_strategy()
    ) {
        let domain = IntervalClosed::new(lo, hi);
        let num_intervals = NumIntervals::try_new(n).unwrap();
        let grid = Grid1DUniform::new(domain.clone(), num_intervals);

        prop_assert_eq!(grid.coords().first(), domain.lower_bound_value());
    }

    /// Property: Last coordinate equals domain upper bound
    #[test]
    fn uniform_grid_last_coord_equals_upper_bound(
        (lo, hi) in domain_bounds_strategy(),
        n in num_intervals_strategy()
    ) {
        let domain = IntervalClosed::new(lo, hi);
        let num_intervals = NumIntervals::try_new(n).unwrap();
        let grid = Grid1DUniform::new(domain.clone(), num_intervals);

        prop_assert_eq!(grid.coords().last(), domain.upper_bound_value());
    }

    /// Property: All intervals have equal length in uniform grid
    #[test]
    fn uniform_grid_equal_interval_lengths(
        (lo, hi) in domain_bounds_strategy(),
        n in 1..100usize // Smaller range to avoid floating-point precision issues
    ) {
        let domain = IntervalClosed::new(lo, hi);
        let num_intervals = NumIntervals::try_new(n).unwrap();
        let grid = Grid1DUniform::new(domain.clone(), num_intervals);

        let expected_delta = (hi - lo) / (n as f64);
        let tolerance = expected_delta * 1e-9;

        for i in 0..n {
            let interval_id = IntervalId::new(i);
            let length = *grid.interval_length(&interval_id).as_ref();
            prop_assert!(
                (length - expected_delta).abs() < tolerance,
                "Interval {} has length {}, expected {}",
                i, length, expected_delta
            );
        }
    }

    /// Property: Coordinates are strictly increasing
    #[test]
    fn uniform_grid_coords_strictly_increasing(
        (lo, hi) in domain_bounds_strategy(),
        n in num_intervals_strategy()
    ) {
        let domain = IntervalClosed::new(lo, hi);
        let num_intervals = NumIntervals::try_new(n).unwrap();
        let grid = Grid1DUniform::new(domain, num_intervals);

        let coords = grid.coords().as_ref();
        for i in 1..coords.len() {
            prop_assert!(
                coords[i] > coords[i - 1],
                "Coordinates not strictly increasing at index {}: {} <= {}",
                i, coords[i], coords[i - 1]
            );
        }
    }

    /// Property: Point location returns valid interval ID for any point in domain
    #[test]
    fn uniform_grid_point_location_valid(
        (lo, hi) in domain_bounds_strategy(),
        n in 1..100usize,
        t in 0.0..1.0f64
    ) {
        let domain = IntervalClosed::new(lo, hi);
        let num_intervals = NumIntervals::try_new(n).unwrap();
        let grid = Grid1DUniform::new(domain, num_intervals);

        // Generate point within domain using linear interpolation
        let point = lo + t * (hi - lo);
        let interval_id = grid.find_interval_id_of_point(&point).unwrap();

        prop_assert!(
            *interval_id.as_ref() < n,
            "Interval ID {} >= num_intervals {}",
            interval_id.as_ref(), n
        );
    }

    /// Property: Domain length equals sum of all interval lengths
    #[test]
    fn uniform_grid_length_conservation(
        (lo, hi) in domain_bounds_strategy(),
        n in 1..100usize
    ) {
        let domain = IntervalClosed::new(lo, hi);
        let num_intervals = NumIntervals::try_new(n).unwrap();
        let grid = Grid1DUniform::new(domain.clone(), num_intervals);

        let domain_length = *domain.length().as_ref();
        let sum_intervals: f64 = (0..n)
            .map(|i| *grid.interval_length(&IntervalId::new(i)).as_ref())
            .sum();

        let tolerance = domain_length * 1e-10;
        prop_assert!(
            (sum_intervals - domain_length).abs() < tolerance,
            "Sum of intervals {} != domain length {}",
            sum_intervals, domain_length
        );
    }
}

// =============================================================================
// Grid1DNonUniform Property Tests
// =============================================================================

proptest! {
    #![proptest_config(ProptestConfig::with_cases(300))]

    /// Property: Non-uniform grid preserves coordinate count
    #[test]
    fn non_uniform_grid_coord_count(
        (lo, hi) in domain_bounds_strategy(),
        extra_points in 0..50usize
    ) {
        let mut coords_vec = vec![lo, hi];
        // Add some random points in between
        for i in 1..=extra_points {
            let t = i as f64 / (extra_points + 1) as f64;
            coords_vec.push(lo + t * (hi - lo));
        }

        let sorted = SortedSet::from_unsorted(coords_vec.clone());
        let expected_len = sorted.len();

        let grid = Grid1D::<IntervalClosed<f64>>::try_from_sorted(sorted).unwrap();

        prop_assert_eq!(grid.coords().len(), expected_len);
    }

    /// Property: Non-uniform grid coordinates are sorted
    #[test]
    fn non_uniform_grid_coords_sorted(
        (lo, hi) in domain_bounds_strategy(),
        extra_points in 0..50usize
    ) {
        let mut coords_vec = vec![lo, hi];
        for i in 1..=extra_points {
            let t = i as f64 / (extra_points + 1) as f64;
            coords_vec.push(lo + t * (hi - lo));
        }

        let grid = Grid1D::<IntervalClosed<f64>>::try_from_sorted(
            SortedSet::from_unsorted(coords_vec)
        ).unwrap();

        let coords = grid.coords().as_ref();
        for i in 1..coords.len() {
            prop_assert!(coords[i] > coords[i - 1]);
        }
    }
}

// =============================================================================
// Grid Refinement Property Tests
// =============================================================================

proptest! {
    #![proptest_config(ProptestConfig::with_cases(200))]

    /// Property: Uniform refinement doubles intervals when adding 1 point per interval
    #[test]
    fn uniform_refinement_doubles_intervals(
        (lo, hi) in domain_bounds_strategy(),
        n in 1..50usize
    ) {
        let domain = IntervalClosed::new(lo, hi);
        let grid = Grid1D::uniform(domain, NumIntervals::try_new(n).unwrap());

        let refinement = grid.refine_uniform(&PositiveNumPoints1D::try_new(1).unwrap());
        let refined = refinement.refined_grid();

        prop_assert_eq!(
            *refined.num_intervals().as_ref(),
            n * 2,
            "Refined intervals should be 2x original"
        );
    }

    /// Property: Refinement preserves domain
    #[test]
    fn refinement_preserves_domain(
        (lo, hi) in domain_bounds_strategy(),
        n in 1..50usize,
        extra_points in 1..5usize
    ) {
        let domain = IntervalClosed::new(lo, hi);
        let grid = Grid1D::uniform(domain.clone(), NumIntervals::try_new(n).unwrap());

        let refinement = grid.refine_uniform(&PositiveNumPoints1D::try_new(extra_points).unwrap());
        let refined = refinement.refined_grid();

        prop_assert_eq!(refined.domain().lower_bound_value(), &lo);
        prop_assert_eq!(refined.domain().upper_bound_value(), &hi);
    }

    /// Property: Refinement mapping is valid - each refined interval maps to a valid original
    #[test]
    fn refinement_mapping_valid(
        (lo, hi) in domain_bounds_strategy(),
        n in 1..50usize
    ) {
        let domain = IntervalClosed::new(lo, hi);
        let grid = Grid1D::uniform(domain, NumIntervals::try_new(n).unwrap());

        let refinement = grid.refine_uniform(&PositiveNumPoints1D::try_new(1).unwrap());

        for (refined_id, original_id) in refinement.iter_refined_with_mapping() {
            prop_assert!(
                *original_id.as_ref() < n,
                "Original interval ID {} out of range [0, {})",
                original_id.as_ref(), n
            );
            prop_assert!(
                *refined_id.as_ref() < *refinement.refined_grid().num_intervals().as_ref(),
                "Refined interval ID out of range"
            );
        }
    }
}

// =============================================================================
// Grid Union Property Tests
// =============================================================================

proptest! {
    #![proptest_config(ProptestConfig::with_cases(200))]

    /// Property: Union of identical grids has same coordinate count
    #[test]
    fn union_identical_grids_same_coords(
        (lo, hi) in domain_bounds_strategy(),
        n in 1..50usize
    ) {
        let domain = IntervalClosed::new(lo, hi);
        let grid = Grid1D::uniform(domain, NumIntervals::try_new(n).unwrap());

        let union = Grid1DUnion::try_new(&grid, &grid).unwrap();

        prop_assert_eq!(
            union.coords().len(),
            grid.coords().len(),
            "Union of identical grids should have same coord count"
        );
    }

    /// Property: Union preserves domain
    #[test]
    fn union_preserves_domain(
        (lo, hi) in domain_bounds_strategy(),
        n1 in 1..50usize,
        n2 in 1..50usize
    ) {
        let domain = IntervalClosed::new(lo, hi);
        let grid1 = Grid1D::uniform(domain.clone(), NumIntervals::try_new(n1).unwrap());
        let grid2 = Grid1D::uniform(domain, NumIntervals::try_new(n2).unwrap());

        let union = Grid1DUnion::try_new(&grid1, &grid2).unwrap();

        prop_assert_eq!(union.domain().lower_bound_value(), &lo);
        prop_assert_eq!(union.domain().upper_bound_value(), &hi);
    }

    /// Property: Union has at least as many points as either input grid
    #[test]
    fn union_has_superset_coords(
        (lo, hi) in domain_bounds_strategy(),
        n1 in 1..50usize,
        n2 in 1..50usize
    ) {
        let domain = IntervalClosed::new(lo, hi);
        let grid1 = Grid1D::uniform(domain.clone(), NumIntervals::try_new(n1).unwrap());
        let grid2 = Grid1D::uniform(domain, NumIntervals::try_new(n2).unwrap());

        let union = Grid1DUnion::try_new(&grid1, &grid2).unwrap();

        prop_assert!(
            union.coords().len() >= grid1.coords().len(),
            "Union should have >= coords than grid1"
        );
        prop_assert!(
            union.coords().len() >= grid2.coords().len(),
            "Union should have >= coords than grid2"
        );
    }

    /// Property: Union mapping is valid
    #[test]
    fn union_mapping_valid(
        (lo, hi) in domain_bounds_strategy(),
        n1 in 1..30usize,
        n2 in 1..30usize
    ) {
        let domain = IntervalClosed::new(lo, hi);
        let grid1 = Grid1D::uniform(domain.clone(), NumIntervals::try_new(n1).unwrap());
        let grid2 = Grid1D::uniform(domain, NumIntervals::try_new(n2).unwrap());

        let union = Grid1DUnion::try_new(&grid1, &grid2).unwrap();

        for (unified_id, a_id, b_id) in union.iter_interval_mappings() {
            prop_assert!(
                *a_id.as_ref() < n1,
                "Grid A interval ID {} out of range", a_id.as_ref()
            );
            prop_assert!(
                *b_id.as_ref() < n2,
                "Grid B interval ID {} out of range", b_id.as_ref()
            );
            prop_assert!(
                *unified_id.as_ref() < *union.num_refined_intervals().as_ref(),
                "Unified interval ID out of range"
            );
        }
    }
}

// =============================================================================
// Interval Property Tests
// =============================================================================

proptest! {
    #![proptest_config(ProptestConfig::with_cases(500))]

    /// Property: Interval length is always positive for valid intervals
    #[test]
    fn interval_length_positive(
        (lo, hi) in domain_bounds_strategy()
    ) {
        let interval = IntervalClosed::new(lo, hi);
        let length = *interval.length().as_ref();

        prop_assert!(length > 0.0, "Interval length should be positive");
        prop_assert_eq!(length, hi - lo, "Length should equal hi - lo");
    }

    /// Property: Midpoint is within interval
    #[test]
    fn interval_midpoint_in_bounds(
        (lo, hi) in domain_bounds_strategy()
    ) {
        let interval = IntervalClosed::new(lo, hi);
        let mid = interval.midpoint();

        prop_assert!(mid >= lo, "Midpoint {} < lower bound {}", mid, lo);
        prop_assert!(mid <= hi, "Midpoint {} > upper bound {}", mid, hi);
    }

    /// Property: Midpoint is equidistant from bounds
    #[test]
    fn interval_midpoint_equidistant(
        (lo, hi) in domain_bounds_strategy()
    ) {
        let interval = IntervalClosed::new(lo, hi);
        let mid = interval.midpoint();

        let dist_lo = mid - lo;
        let dist_hi = hi - mid;
        let tolerance = (hi - lo) * 1e-10;

        prop_assert!(
            (dist_lo - dist_hi).abs() < tolerance,
            "Midpoint not equidistant: {} from lo, {} from hi",
            dist_lo, dist_hi
        );
    }
}

// =============================================================================
// Serialization Roundtrip Property Tests
// =============================================================================

/// Check if two f64 values are approximately equal (relative tolerance)
fn approx_eq(a: f64, b: f64) -> bool {
    if a == b {
        return true;
    }
    let abs_diff = (a - b).abs();
    let max_val = a.abs().max(b.abs());
    // Use relative tolerance for larger values, absolute for smaller
    abs_diff <= max_val * 1e-14 || abs_diff <= 1e-15
}

/// Check if two f64 slices are approximately equal
fn slices_approx_eq(a: &[f64], b: &[f64]) -> bool {
    a.len() == b.len() && a.iter().zip(b.iter()).all(|(x, y)| approx_eq(*x, *y))
}

proptest! {
    #![proptest_config(ProptestConfig::with_cases(100))]

    /// Property: Grid survives JSON serialization roundtrip (within floating-point tolerance)
    #[test]
    fn grid_json_roundtrip(
        (lo, hi) in domain_bounds_strategy(),
        n in 1..50usize
    ) {
        let domain = IntervalClosed::new(lo, hi);
        let grid = Grid1D::uniform(domain, NumIntervals::try_new(n).unwrap());

        let json = serde_json::to_string(&grid).unwrap();
        let deserialized: Grid1D<IntervalClosed<f64>> = serde_json::from_str(&json).unwrap();

        // Use approximate equality due to JSON floating-point serialization precision
        prop_assert!(
            slices_approx_eq(grid.coords().deref(), deserialized.coords().deref()),
            "Coordinates don't match after roundtrip"
        );
        prop_assert_eq!(grid.num_intervals(), deserialized.num_intervals());
    }

    /// Property: Interval survives JSON serialization roundtrip (within floating-point tolerance)
    #[test]
    fn interval_json_roundtrip(
        (lo, hi) in domain_bounds_strategy()
    ) {
        let interval = IntervalClosed::new(lo, hi);

        let json = serde_json::to_string(&interval).unwrap();
        let deserialized: IntervalClosed<f64> = serde_json::from_str(&json).unwrap();

        // Use approximate equality due to JSON floating-point serialization precision
        prop_assert!(
            approx_eq(*interval.lower_bound_value(), *deserialized.lower_bound_value()),
            "Lower bounds don't match after roundtrip"
        );
        prop_assert!(
            approx_eq(*interval.upper_bound_value(), *deserialized.upper_bound_value()),
            "Upper bounds don't match after roundtrip"
        );
    }
}