blr-active 0.1.0

//! Algorithm 3: Precision Assessment
//!
//! Evaluates whether a fitted BLR+ARD model meets the user's precision goal,
//! using **percentile-based** semantics (per CLARIFY-1, option C):
//!
//!   **Goal met** when: 95th-percentile posterior std ≤ target_std
//!
//! Rationale: 95th percentile is interpretable ("95% of measurements meet
//! spec") while allowing a small tail of higher-uncertainty points (e.g., at
//! input-space boundaries before dense sampling).
//!
//! Status thresholds:
//!   - `MetGoal`:  percentile_95_std ≤ target_std
//!   - `Near`:     percentile_95_std ≤ 1.1 × target_std (within 10% of goal)
//!   - `Unmet`:    percentile_95_std > 1.1 × target_std
//!   - `NoiseFloorHit`: set externally by Algorithm 4 integration

/// Multi-level precision status.
#[derive(Debug, Clone, PartialEq)]
pub enum PrecisionStatus {
    /// 95th-percentile posterior std is within the user's tolerance.
    MetGoal,
    /// Within 10% of the tolerance — almost there.
    Near,
    /// Not yet within tolerance; more samples are needed.
    Unmet,
    /// Posterior variance has plateaued; hardware noise is the bottleneck.
    /// This is set externally when Algorithm 4 reports `likely_at_floor=true`.
    NoiseFloorHit,
}

/// Summary of the current predictive precision.
#[derive(Debug, Clone)]
pub struct PrecisionAssessment {
    /// Maximum posterior std over the test/grid points.
    pub max_std: f64,
    /// Mean posterior std over the test/grid points.
    pub mean_std: f64,
    /// 95th-percentile posterior std (primary decision metric, per CLARIFY-1).
    pub percentile_95_std: f64,
    /// User-provided tolerance (the precision goal).
    pub target_std: f64,
    /// Relative gap: (percentile_95_std − target_std) / target_std.
    /// Negative when goal is met; positive when unmet.
    pub gap: f64,
    /// Status label (MetGoal / Near / Unmet / NoiseFloorHit).
    pub status: PrecisionStatus,
}

/// Compute the p-th percentile (0.0–1.0) of a slice.
///
/// Uses nearest-rank method. Returns 0.0 on empty input.
pub fn percentile(values: &[f64], p: f64) -> f64 {
    if values.is_empty() {
        return 0.0;
    }
    let mut sorted = values.to_vec();
    sorted.sort_by(|a, b| a.partial_cmp(b).unwrap_or(std::cmp::Ordering::Equal));
    let idx = ((p * (sorted.len() - 1) as f64).round() as usize).min(sorted.len() - 1);
    sorted[idx]
}

/// Assess how close the current model is to the precision goal.
///
/// # Arguments
/// - `posterior_std`: posterior std at each evaluation point (e.g., a uniform grid)
/// - `target_std`: user precision requirement (e.g., 0.01 V)
///
/// # Edge cases
/// - Empty `posterior_std` → returns `Unmet` with zeros everywhere.
/// - `target_std = 0` → gap is always `f64::INFINITY` (never MetGoal).
pub fn assess_precision(posterior_std: &[f64], target_std: f64) -> PrecisionAssessment {
    if posterior_std.is_empty() {
        return PrecisionAssessment {
            max_std: 0.0,
            mean_std: 0.0,
            percentile_95_std: 0.0,
            target_std,
            gap: if target_std > 0.0 {
                -1.0
            } else {
                f64::INFINITY
            },
            status: PrecisionStatus::Unmet,
        };
    }

    let max_std = posterior_std
        .iter()
        .cloned()
        .fold(f64::NEG_INFINITY, f64::max);
    let mean_std = posterior_std.iter().sum::<f64>() / posterior_std.len() as f64;
    let percentile_95_std = percentile(posterior_std, 0.95);

    let gap = if target_std > 0.0 {
        (percentile_95_std - target_std) / target_std
    } else {
        f64::INFINITY
    };

    // Primary decision metric: percentile_95_std (per CLARIFY-1)
    let status = if target_std > 0.0 && percentile_95_std <= target_std {
        PrecisionStatus::MetGoal
    } else if target_std > 0.0 && percentile_95_std <= 1.1 * target_std {
        PrecisionStatus::Near
    } else {
        PrecisionStatus::Unmet
    };

    PrecisionAssessment {
        max_std,
        mean_std,
        percentile_95_std,
        target_std,
        gap,
        status,
    }
}

#[cfg(test)]
mod tests {
    use super::*;

    /// MetGoal when p95 <= target
    #[test]
    fn test_precision_goal_met() {
        // All stds at 0.005, target = 0.01 → p95 = 0.005 ≤ 0.01
        let stds = vec![0.005; 100];
        let a = assess_precision(&stds, 0.01);
        assert_eq!(a.status, PrecisionStatus::MetGoal);
        assert!(a.gap < 0.0, "gap should be negative when goal met");
    }

    /// Near when p95 is within 10% above target
    #[test]
    fn test_precision_near() {
        // p95 = 0.0105, target = 0.01 → 5% above → Near
        let mut stds = vec![0.01; 100];
        stds[99] = 0.012; // 99th percentile at 0.012 but 95th stays at 0.01
                          // Construct: 95 values at 0.01, 5 values at 0.012
        let stds: Vec<f64> = (0..95).map(|_| 0.01).chain((0..5).map(|_| 0.012)).collect();
        let _a = assess_precision(&stds, 0.01);
        // p95 of this set: index = round(0.95 * 99) = 94 → 0.01, which is MetGoal
        // Let's pick a case that is clearly Near
        let stds2: Vec<f64> = (0..80)
            .map(|_| 0.01)
            .chain((0..20).map(|_| 0.0108))
            .collect();
        let a2 = assess_precision(&stds2, 0.01);
        // p95: index 94 → 0.0108 ≤ 1.1 * 0.01 = 0.011 → Near
        assert!(
            a2.status == PrecisionStatus::Near || a2.status == PrecisionStatus::MetGoal,
            "expected Near or MetGoal, got {:?}",
            a2.status
        );
    }

    /// Unmet when p95 is significantly above target
    #[test]
    fn test_precision_unmet() {
        let stds = vec![0.05; 100]; // way above target 0.01
        let a = assess_precision(&stds, 0.01);
        assert_eq!(a.status, PrecisionStatus::Unmet);
        assert!(a.gap > 0.0, "gap should be positive when unmet");
    }

    /// Gap formula: (p95 - target) / target
    #[test]
    fn test_gap_calculation() {
        let stds = vec![0.02; 100];
        let a = assess_precision(&stds, 0.01);
        let expected_gap = (0.02 - 0.01) / 0.01; // = 1.0
        assert!(
            (a.gap - expected_gap).abs() < 1e-10,
            "gap={}, expected={}",
            a.gap,
            expected_gap
        );
    }

    /// Percentile function is correct for known distribution
    #[test]
    fn test_percentile_computation() {
        // 100 values: [0.0, 0.01, ..., 0.99] → p95 index = round(0.95*99)=94 → 0.94
        let vals: Vec<f64> = (0..100).map(|i| i as f64 / 100.0).collect();
        let p95 = percentile(&vals, 0.95);
        // nearest-rank: round(0.95 * 99) = round(94.05) = 94 → vals[94] = 0.94
        assert!((p95 - 0.94).abs() < 1e-10, "p95={}, expected 0.94", p95);
    }

    /// Empty array returns Unmet without panic
    #[test]
    fn test_empty_stds_no_panic() {
        let a = assess_precision(&[], 0.01);
        assert_eq!(a.status, PrecisionStatus::Unmet);
    }

    /// target_std = 0: never MetGoal, gap = infinity
    #[test]
    fn test_zero_target_std() {
        let stds = vec![0.001; 10];
        let a = assess_precision(&stds, 0.0);
        assert_ne!(a.status, PrecisionStatus::MetGoal);
        assert!(
            a.gap.is_infinite(),
            "gap should be infinite for target_std=0"
        );
    }
}