blr-active 0.1.0

/// Integration tests for Algorithm 1: Posterior Variance Computation.
///
/// Tests interact with the blr_core public API (posterior_variance, posterior_std_grid)
/// and verify correctness against manually computed expected values.
use blr_active::active_learning::variance::{
    posterior_std, posterior_std_grid, posterior_variance,
};

// ─── Manual calculation baseline ─────────────────────────────────────────────
//
// For Σ = I (identity), φ = e_1 = (1, 0, ..., 0), β = 1:
//   σ²_* = 1/β + φ^T Σ φ = 1 + 1 = 2  →  σ_* = √2 ≈ 1.41421356
//
// For Σ = 2·I, φ = (1, 1), β = 2:
//   σ²_* = 1/2 + (1,1) · 2·I · (1,1)^T = 0.5 + 2·(1+1) = 0.5 + 4 = 4.5
//   σ_* = √4.5 ≈ 2.12132034

#[test]
fn test_variance_identity_cov_e1() {
    // Σ = I₄, φ = e₁ = (1,0,0,0), β = 1 → var = 1 + 1 = 2
    let d = 4usize;
    let sigma_cov: Vec<f64> = (0..d * d)
        .map(|k| if k % (d + 1) == 0 { 1.0 } else { 0.0 })
        .collect();
    let phi = vec![1.0, 0.0, 0.0, 0.0];
    let var = posterior_variance(1.0, &sigma_cov, d, &phi, 1);
    assert!((var[0] - 2.0).abs() < 1e-10, "expected 2.0, got {}", var[0]);
}

#[test]
fn test_variance_scaled_identity_two_points() {
    // Σ = 2·I₂, β = 2
    // φ₁ = (1,1) → var₁ = 0.5 + (1,1)·2·(1,1)^T = 0.5 + 4 = 4.5
    // φ₂ = (0,1) → var₂ = 0.5 + (0,1)·2·(0,1)^T = 0.5 + 2 = 2.5
    let d = 2usize;
    let sigma_cov = vec![2.0, 0.0, 0.0, 2.0];
    let phi_test = vec![1.0, 1.0, 0.0, 1.0]; // two test points, row-major
    let var = posterior_variance(2.0, &sigma_cov, d, &phi_test, 2);
    assert!((var[0] - 4.5).abs() < 1e-10, "expected 4.5, got {}", var[0]);
    assert!((var[1] - 2.5).abs() < 1e-10, "expected 2.5, got {}", var[1]);
}

#[test]
fn test_std_is_sqrt_of_variance() {
    let d = 3usize;
    let sigma_cov = vec![0.5, 0.1, 0.0, 0.1, 0.3, 0.0, 0.0, 0.0, 0.2];
    let phi = vec![1.0, -0.5, 0.3];
    let var = posterior_variance(1.5, &sigma_cov, d, &phi, 1);
    let std = posterior_std(1.5, &sigma_cov, d, &phi, 1);
    assert!(
        (std[0] - var[0].sqrt()).abs() < 1e-12,
        "std should be sqrt(var)"
    );
}

#[test]
fn test_grid_length_matches_resolution() {
    let d = 2usize;
    let sigma_cov = vec![1.0, 0.0, 0.0, 1.0];
    let feature_fn = |x: f64| vec![1.0, x];
    let (grid, stds) = posterior_std_grid(1.0, &sigma_cov, d, 0.0, 10.0, 100, &feature_fn);
    assert_eq!(grid.len(), 100);
    assert_eq!(stds.len(), 100);
}

#[test]
fn test_grid_all_finite_nonneg() {
    let d = 3usize;
    let sigma_cov = vec![0.3, 0.0, 0.0, 0.0, 0.5, 0.1, 0.0, 0.1, 0.2];
    let feature_fn = |x: f64| vec![1.0, x, x * x];
    let (_, stds) = posterior_std_grid(2.0, &sigma_cov, d, -10.0, 10.0, 200, &feature_fn);
    for (i, s) in stds.iter().enumerate() {
        assert!(s.is_finite(), "grid[{i}] std is not finite: {s}");
        assert!(*s >= 0.0, "grid[{i}] std is negative: {s}");
    }
}

#[test]
fn test_performance_100_points_fast() {
    // Verify 100-point × 8-feature grid completes well under 1 ms
    // (checked as a wall-clock guard, not as a hard performance assertion)
    let d = 8usize;
    let sigma_cov: Vec<f64> = (0..d * d)
        .map(|k| if k % (d + 1) == 0 { 0.1 } else { 0.0 })
        .collect();
    let feature_fn = |x: f64| {
        vec![
            1.0,
            x,
            x * x,
            x.sin(),
            x.cos(),
            x.exp().min(1e6),
            x.ln().abs(),
            x * x * x,
        ]
    };
    let start = std::time::Instant::now();
    let (_, _) = posterior_std_grid(1.0, &sigma_cov, d, 0.0, 1.0, 100, &feature_fn);
    let elapsed = start.elapsed();
    // Allow generous 100ms — should be <1ms in practice; this just guards regressions
    assert!(
        elapsed.as_millis() < 100,
        "100-point grid took {:?}",
        elapsed
    );
}

#[test]
fn test_degenerate_single_point_no_panic() {
    let d = 2usize;
    let sigma_cov = vec![1.0, 0.0, 0.0, 1.0];
    let feature_fn = |_: f64| vec![0.0, 0.0];
    let (grid, stds) = posterior_std_grid(1.0, &sigma_cov, d, 5.0, 5.0, 2, &feature_fn);
    // resolution clamped to 2, so grid has 2 points (both equal to 5.0)
    assert!(grid.len() >= 2);
    for s in &stds {
        assert!(s.is_finite());
    }
}