test(fitstats): require Rust ownership for multigroup M2

CHANGELOG.md

@@ -4,6 +4,7 @@
 
 ### Fixed
 
+- Public multigroup M2 fails closed without the compiled Rust core and delegates target/null projected quadratic forms to native `projected_m2`.
 - Public conditional-Rasch M2 fails closed without the compiled Rust core and delegates every result field to the native `m2_cmle_rasch_stat` entrypoint.
 - Reject ambiguous LLM-judge JSON with duplicate keys or unexpected top-level fields; require the exact mode-specific schema including advisory `accepted`.
 - Require compiled Rust ownership for public `s_x2` and `person_fit`, including prior-mean S-X² dispatch, with fail-closed errors when the core is missing.

crates/fast-mlsirm-py/src/lib.rs

@@ -16,7 +16,10 @@ use mlsirm_core::fitstats::{
     benjamini_hochberg as core_benjamini_hochberg, chi2_sf as core_chi2_sf,
     infit_outfit as core_infit_outfit, leniency_residuals as core_leniency_residuals,
     m2_cmle_rasch as core_m2_cmle_rasch, m2_rmsea2 as core_m2, person_fit as core_person_fit,
-    poly_local_dependence as core_poly_ld, poly_m2 as core_poly_m2, s_x2 as core_s_x2, SX2Config,
+    poly_local_dependence as core_poly_ld, poly_m2 as core_poly_m2,
+    projected_m2 as core_projected_m2,
+    projected_m2_workspace_elements as core_projected_m2_workspace_elements,
+    s_x2 as core_s_x2, SX2Config, PROJECTED_M2_MAX_WORKSPACE_ELEMENTS,
 };
 use mlsirm_core::linking::{
     irt_link as core_irt_link, link_fixed_item_parameters as core_link_fixed_item_parameters,
@@ -7034,6 +7037,69 @@ fn m2_stat(
 }
 
 
+
+/// Projected M2 quadratic form ownership entrypoint (dense residual / Delta / Xi).
+#[pyfunction]
+fn projected_m2(
+    residual: PyReadonlyArray1<'_, f64>,
+    delta: PyReadonlyArray2<'_, f64>,
+    xi: PyReadonlyArray2<'_, f64>,
+    n: f64,
+) -> PyResult<f64> {
+    let residual_arr = residual.as_array();
+    let delta_arr = delta.as_array();
+    let xi_arr = xi.as_array();
+    let s = residual_arr.len();
+    if delta_arr.shape()[0] != s {
+        return Err(PyValueError::new_err(
+            "delta rows must match residual length",
+        ));
+    }
+    if xi_arr.shape() != [s, s] {
+        return Err(PyValueError::new_err("xi must be residual_len x residual_len"));
+    }
+    let p = delta_arr.shape()[1];
+    let workspace_elements = core_projected_m2_workspace_elements(s, p)
+        .map_err(PyValueError::new_err)?;
+    if workspace_elements > PROJECTED_M2_MAX_WORKSPACE_ELEMENTS {
+        return Err(PyValueError::new_err(
+            "projected M2 workspace exceeds supported element budget",
+        ));
+    }
+    if !n.is_finite() {
+        return Err(PyValueError::new_err("n must be finite"));
+    }
+    if residual_arr.iter().any(|value| !value.is_finite()) {
+        return Err(PyValueError::new_err("residual values must be finite"));
+    }
+    if delta_arr.iter().any(|value| !value.is_finite()) {
+        return Err(PyValueError::new_err("delta values must be finite"));
+    }
+    if xi_arr.iter().any(|value| !value.is_finite()) {
+        return Err(PyValueError::new_err("xi values must be finite"));
+    }
+    let residual = residual.as_slice()?;
+    let delta_elements = s
+        .checked_mul(p)
+        .ok_or_else(|| PyValueError::new_err("projected M2 dimensions overflow"))?;
+    let xi_elements = s
+        .checked_mul(s)
+        .ok_or_else(|| PyValueError::new_err("projected M2 dimensions overflow"))?;
+    let mut delta_flat = vec![0.0_f64; delta_elements];
+    for row in 0..s {
+        for col in 0..p {
+            delta_flat[row * p + col] = delta_arr[[row, col]];
+        }
+    }
+    let mut xi_flat = vec![0.0_f64; xi_elements];
+    for row in 0..s {
+        for col in 0..s {
+            xi_flat[row * s + col] = xi_arr[[row, col]];
+        }
+    }
+    core_projected_m2(residual, &delta_flat, xi_flat, s, p, n).map_err(PyValueError::new_err)
+}
+
 /// Conditional-Rasch M2 (CMLE ownership path).
 #[pyfunction]
 #[pyo3(signature = (y, observed, n_persons, item_easiness))]
@@ -9169,6 +9235,7 @@ fn fast_mlsirm_core(m: &Bound<'_, PyModule>) -> PyResult<()> {
     m.add_function(wrap_pyfunction!(leniency_residuals_stat, m)?)?;
     m.add_function(wrap_pyfunction!(m2_stat, m)?)?;
     m.add_function(wrap_pyfunction!(m2_cmle_rasch_stat, m)?)?;
+    m.add_function(wrap_pyfunction!(projected_m2, m)?)?;
     m.add_function(wrap_pyfunction!(poly_m2, m)?)?;
     m.add_function(wrap_pyfunction!(poly_local_dependence, m)?)?;
     m.add_function(wrap_pyfunction!(poly_dif, m)?)?;

crates/mlsirm-core/src/fitstats.rs

@@ -1916,8 +1916,48 @@ fn chol_solve(l: &[f64], n: usize, b: &[f64]) -> Vec<f64> {
     x
 }
 
+/// Maximum supported dense workspace for projected M2, in `f64` elements.
+pub const PROJECTED_M2_MAX_WORKSPACE_ELEMENTS: usize = 4_000_000;
+
+/// Return the peak dense element budget required by projected M2.
+pub fn projected_m2_workspace_elements(
+    n_moments: usize,
+    n_parameters: usize,
+) -> Result<usize, String> {
+    let xi_elements = n_moments
+        .checked_mul(n_moments)
+        .ok_or_else(|| "projected M2 dimensions overflow".to_string())?;
+    let delta_elements = n_moments
+        .checked_mul(n_parameters)
+        .ok_or_else(|| "projected M2 dimensions overflow".to_string())?;
+    let parameter_elements = n_parameters
+        .checked_mul(n_parameters)
+        .ok_or_else(|| "projected M2 dimensions overflow".to_string())?;
+    let two_moments = n_moments
+        .checked_mul(2)
+        .ok_or_else(|| "projected M2 dimensions overflow".to_string())?;
+    let two_parameters = n_parameters
+        .checked_mul(2)
+        .ok_or_else(|| "projected M2 dimensions overflow".to_string())?;
+
+    [
+        xi_elements,
+        delta_elements,
+        delta_elements,
+        parameter_elements,
+        two_moments,
+        two_parameters,
+    ]
+    .into_iter()
+    .try_fold(0_usize, |total, elements| {
+        total
+            .checked_add(elements)
+            .ok_or_else(|| "projected M2 dimensions overflow".to_string())
+    })
+}
+
 /// Evaluate the M2 projected quadratic form without forming either inverse.
-fn projected_m2(
+pub fn projected_m2(
     e: &[f64],
     delta: &[f64],
     mut xi: Vec<f64>,
@@ -1927,6 +1967,39 @@ fn projected_m2(
 ) -> Result<f64, String> {
     let s = n_moments;
     let p = n_parameters;
+    let delta_elements = s
+        .checked_mul(p)
+        .ok_or_else(|| "projected M2 dimensions overflow".to_string())?;
+    let xi_elements = s
+        .checked_mul(s)
+        .ok_or_else(|| "projected M2 dimensions overflow".to_string())?;
+    let workspace_elements = projected_m2_workspace_elements(s, p)?;
+    if workspace_elements > PROJECTED_M2_MAX_WORKSPACE_ELEMENTS {
+        return Err(
+            "projected M2 workspace exceeds supported element budget".to_string(),
+        );
+    }
+    if e.len() != s {
+        return Err("residual length must match n_moments".to_string());
+    }
+    if delta.len() != delta_elements {
+        return Err("delta length must equal n_moments * n_parameters".to_string());
+    }
+    if xi.len() != xi_elements {
+        return Err("xi length must equal n_moments * n_moments".to_string());
+    }
+    if !n.is_finite() {
+        return Err("n must be finite".to_string());
+    }
+    if e.iter().any(|value| !value.is_finite()) {
+        return Err("residual values must be finite".to_string());
+    }
+    if delta.iter().any(|value| !value.is_finite()) {
+        return Err("delta values must be finite".to_string());
+    }
+    if xi.iter().any(|value| !value.is_finite()) {
+        return Err("xi values must be finite".to_string());
+    }
     cholesky_lower(&mut xi, s)?;
     let u = chol_solve(&xi, s, e); // Xi^-1 e
     let mut w = vec![0.0_f64; s * p]; // Xi^-1 Delta

crates/mlsirm-core/tests/projected_m2_validation.rs

@@ -0,0 +1,65 @@
+//! Validation contracts for the public projected-M2 Rust boundary.
+
+use mlsirm_core::fitstats::projected_m2;
+
+#[test]
+fn projected_m2_rejects_short_residual_without_panicking() {
+    let result = projected_m2(
+        &[0.25],
+        &[1.0, 0.5],
+        vec![1.0, 0.0, 0.0, 1.0],
+        2,
+        1,
+        100.0,
+    );
+
+    assert!(result.is_err(), "short residual must return Err");
+}
+
+#[test]
+fn projected_m2_rejects_short_delta_without_panicking() {
+    let result = projected_m2(
+        &[0.25, -0.25],
+        &[1.0],
+        vec![1.0, 0.0, 0.0, 1.0],
+        2,
+        1,
+        100.0,
+    );
+
+    assert!(result.is_err(), "short delta must return Err");
+}
+
+#[test]
+fn projected_m2_rejects_short_covariance_without_panicking() {
+    let result = projected_m2(
+        &[0.25, -0.25],
+        &[1.0, 0.5],
+        vec![1.0, 0.0, 0.0],
+        2,
+        1,
+        100.0,
+    );
+
+    assert!(result.is_err(), "short covariance must return Err");
+}
+
+#[test]
+fn projected_m2_rejects_dimension_overflow_without_allocating() {
+    let result = projected_m2(&[], &[], Vec::new(), usize::MAX, 2, 1.0);
+
+    assert!(
+        matches!(result, Err(message) if message.contains("overflow")),
+        "overflowing dimensions must fail before indexing or allocation"
+    );
+}
+
+#[test]
+fn projected_m2_rejects_nonfinite_sample_size_before_arithmetic() {
+    let result = projected_m2(&[0.25], &[1.0], vec![1.0], 1, 1, f64::INFINITY);
+
+    assert!(
+        matches!(result, Err(message) if message.contains("n must be finite")),
+        "non-finite n must fail before numerical execution"
+    );
+}

docs/changelog.d/818-multigroup-m2-rust-projection.md

@@ -0,0 +1,5 @@
+# Multigroup M2 Rust projection ownership
+
+## Fixed
+
+- Public `m2_multigroup` fails closed without the compiled Rust core and delegates target/null projected M2 quadratic forms to `projected_m2`.

python/fast_mlsirm/fitstats.py

@@ -2737,9 +2737,24 @@ def m2_multigroup(
         null_delta[rows, cols] = root_n * null_d
         null_xi_blocks.append(null_xi)
 
-    m2_value = _projected_m2_numpy(residual, delta, _block_diag(xi_blocks), 1.0)
-    null_m2 = _projected_m2_numpy(
-        null_residual, null_delta, _block_diag(null_xi_blocks), 1.0
+    core = _core_module()
+    if core is None or not hasattr(core, "projected_m2"):
+        raise RuntimeError("fit statistics require the compiled Rust core")
+    m2_value = float(
+        core.projected_m2(
+            np.ascontiguousarray(residual, dtype=np.float64),
+            np.ascontiguousarray(delta, dtype=np.float64),
+            np.ascontiguousarray(_block_diag(xi_blocks), dtype=np.float64),
+            1.0,
+        )
+    )
+    null_m2 = float(
+        core.projected_m2(
+            np.ascontiguousarray(null_residual, dtype=np.float64),
+            np.ascontiguousarray(null_delta, dtype=np.float64),
+            np.ascontiguousarray(_block_diag(null_xi_blocks), dtype=np.float64),
+            1.0,
+        )
     )
     df = float(n_groups * s - p)
     null_df = float(n_groups * s - n_groups * y0.shape[1])