diff --git a/algebra/backend.py b/algebra/backend.py
index a96f7e90..64b1fe18 100644
--- a/algebra/backend.py
+++ b/algebra/backend.py
@@ -60,7 +60,7 @@ def vault_recall(versors: list, query: np.ndarray, top_k: int = 5) -> list:
         except Exception:
             pass
     q = np.asarray(query)
-    scores = [(i, -float(np.sum((q - np.asarray(v)) ** 2))) for i, v in enumerate(versors)]
+    scores = [(i, float(cga_inner(q, np.asarray(v)))) for i, v in enumerate(versors)]
     scores.sort(key=lambda x: -x[1])
     return scores[:top_k]
 
diff --git a/core-rs/Cargo.toml b/core-rs/Cargo.toml
index f1339dc1..8230b931 100644
--- a/core-rs/Cargo.toml
+++ b/core-rs/Cargo.toml
@@ -5,10 +5,10 @@ edition = "2021"
 
 [lib]
 name = "core_rs"
-crate-type = ["cdylib"]
+crate-type = ["cdylib", "rlib"]
 
 [dependencies]
-pyo3 = { version = "0.21", features = ["extension-module"] }
+pyo3 = { version = "0.21" }
 rayon = "1.10"
 nalgebra = "0.33"
 ndarray = { version = "0.16", features = ["rayon"] }
diff --git a/core-rs/pyproject.toml b/core-rs/pyproject.toml
index 48ed8161..a327ec8b 100644
--- a/core-rs/pyproject.toml
+++ b/core-rs/pyproject.toml
@@ -5,4 +5,3 @@ build-backend = "maturin"
 [tool.maturin]
 features = ["extension-module"]
 module-name = "core_rs"
-python-source = "python"
diff --git a/core-rs/src/cga.rs b/core-rs/src/cga.rs
index 18f46004..a3c774f4 100644
--- a/core-rs/src/cga.rs
+++ b/core-rs/src/cga.rs
@@ -1,7 +1,17 @@
 //! CGA inner product and null-cone operations.
 //!
-//! cga_inner(X, Y) = 0.5 * scalar_part(X*Y + Y*X)
-//!                 = -d^2 / 2  for null vectors X, Y
+//! Signature: (+,+,+,-,+), with Euclidean coordinates on e1,e2,e3.
+//! The conformal null directions are:
+//!
+//!   n_o   = 0.5 * (e4 - e5)   # origin, n_o^2 = 0
+//!   n_inf = e4 + e5           # infinity, n_inf^2 = 0
+//!   n_o · n_inf = -1
+//!
+//! A Euclidean point x embeds as:
+//!
+//!   X = x + n_o + 0.5 * |x|^2 * n_inf
+//!
+//! Then X·X = 0 and X·Y = -0.5 * ||x-y||^2.
 //!
 //! This is the ONLY distance metric in CORE-AI.
 
@@ -20,37 +30,30 @@ pub enum CgaError {
 pub fn cga_inner_raw(x: &[f32; 32], y: &[f32; 32]) -> Result<f32, CgaError> {
     let xy = geometric_product_raw(x, y)?;
     let yx = geometric_product_raw(y, x)?;
-    // scalar part is index 0
     Ok(0.5 * (xy[0] + yx[0]))
 }
 
-/// Check if X is on the null cone: |X*X| < tol.
+/// Check if X is on the null cone: |X·X| < tol.
 pub fn is_null_raw(x: &[f32; 32], tol: f32) -> Result<bool, CgaError> {
     Ok(cga_inner_raw(x, x)?.abs() < tol)
 }
 
-/// Re-project X onto the null cone.
-/// Extract Euclidean components (indices 1-3), recompute e+ = 0.5*|x|^2, e- = 1.
+/// Re-project X onto the null cone by extracting Euclidean components
+/// and re-embedding with the canonical CGA point map.
 pub fn null_project_raw(x: &[f32; 32]) -> [f32; 32] {
-    let mut result = [0f32; 32];
-    result[1] = x[1];
-    result[2] = x[2];
-    result[3] = x[3];
-    let x_sq = result[1] * result[1] + result[2] * result[2] + result[3] * result[3];
-    result[4] = 0.5 * x_sq; // e+ coefficient
-    result[5] = 1.0;         // e- coefficient
-    result
+    embed_point_raw(&[x[1], x[2], x[3]])
 }
 
 /// Embed a Euclidean point [x, y, z] into the CGA null cone.
-/// X = x*e1 + y*e2 + z*e3 + (1/2)|x|^2 * e+ + e-
+/// X = x + n_o + 0.5|x|^2 n_inf
+/// where n_o = 0.5(e4 - e5), n_inf = e4 + e5.
 pub fn embed_point_raw(p: &[f32; 3]) -> [f32; 32] {
     let mut result = [0f32; 32];
     result[1] = p[0];
     result[2] = p[1];
     result[3] = p[2];
-    let r2 = p[0]*p[0] + p[1]*p[1] + p[2]*p[2];
-    result[4] = 0.5 * r2;
-    result[5] = 1.0;
+    let r2 = p[0] * p[0] + p[1] * p[1] + p[2] * p[2];
+    result[4] = 0.5 * (r2 + 1.0);
+    result[5] = 0.5 * (r2 - 1.0);
     result
 }
diff --git a/core-rs/src/cl41.rs b/core-rs/src/cl41.rs
index 12a8796b..166a42f7 100644
--- a/core-rs/src/cl41.rs
+++ b/core-rs/src/cl41.rs
@@ -1,6 +1,6 @@
-//! Cl(4,1) geometric product via fully unrolled precomputed table.
+//! Cl(4,1) geometric product via precomputed table.
 //!
-//! Signature: (+,+,+,+,-). 32-component f32 multivectors.
+//! Signature: (+,+,+,-,+). 32-component f32 multivectors.
 //! The multiplication table is computed once at program start using
 //! const evaluation and stored as two [u8;1024] and [i8;1024] arrays
 //! (index and sign for each of the 32x32 blade pairs).
@@ -20,10 +20,8 @@ pub enum Cl41Error {
 // We encode each blade as a bitmask over 5 bits (bit k = basis vector k+1 present)
 // The mapping from bitmask to component index follows grade-ascending, lex order.
 
-const N: usize = 32;
-
-// Signature: e1^2=+1, e2^2=+1, e3^2=+1, e4^2=+1, e5^2=-1
-const SIG: [i8; 5] = [1, 1, 1, 1, -1];
+// Signature: e1^2=+1, e2^2=+1, e3^2=+1, e4^2=-1, e5^2=+1
+const SIG: [i8; 5] = [1, 1, 1, -1, 1];
 
 // Precomputed at compile time via const fn
 const BLADE_MASKS: [u8; 32] = build_blade_masks();
@@ -68,63 +66,43 @@ const fn popcount5(x: u8) -> u8 {
 }
 
 // Multiply two basis blades given as bitmasks. Returns (result_mask, sign).
-// Uses bubble-sort on the concatenated index list, tracking swaps and metric contractions.
+// The sign is the parity of swaps needed to canonicalize A followed by B,
+// multiplied by the metric contractions for repeated basis vectors.
 const fn blade_product(a: u8, b: u8) -> (u8, i8) {
-    // Expand masks into sorted index sequences
-    let mut seq = [0u8; 10];
-    let mut len = 0usize;
-
-    let mut bit = 0u8;
-    while bit < 5 {
-        if (a >> bit) & 1 == 1 { seq[len] = bit; len += 1; }
-        bit += 1;
-    }
-    bit = 0;
-    while bit < 5 {
-        if (b >> bit) & 1 == 1 { seq[len] = bit; len += 1; }
-        bit += 1;
-    }
-
     let mut sign: i8 = 1;
 
-    // Bubble sort + contract duplicates
-    let mut changed = true;
-    while changed {
-        changed = false;
-        let mut i = 0usize;
-        while i + 1 < len {
-            if seq[i] == seq[i + 1] {
-                // Contract: e_k^2 = SIG[k]
-                sign *= SIG[seq[i] as usize];
-                // Remove both elements at i and i+1
-                let mut j = i;
-                while j + 2 < len { seq[j] = seq[j + 2]; j += 1; }
-                len -= 2;
-                changed = true;
-                if i > 0 { i -= 1; } // re-check from one step back
-            } else if seq[i] > seq[i + 1] {
-                let tmp = seq[i]; seq[i] = seq[i + 1]; seq[i + 1] = tmp;
-                sign *= -1;
-                changed = true;
-                i += 1;
-            } else {
-                i += 1;
+    // Anticommutation sign: every pair (a_i, b_j) with a_i > b_j swaps once.
+    let mut swaps = 0u8;
+    let mut ai = 0u8;
+    while ai < 5 {
+        if (a >> ai) & 1 == 1 {
+            let mut bj = 0u8;
+            while bj < 5 {
+                if (b >> bj) & 1 == 1 && ai > bj {
+                    swaps += 1;
+                }
+                bj += 1;
             }
         }
+        ai += 1;
+    }
+    if swaps % 2 == 1 {
+        sign *= -1;
     }
 
-    // Build result mask
-    let mut result = 0u8;
-    let mut i = 0usize;
-    while i < len { result |= 1 << seq[i]; i += 1; }
+    // Metric contractions for duplicate basis vectors.
+    let common = a & b;
+    let mut bit = 0u8;
+    while bit < 5 {
+        if (common >> bit) & 1 == 1 {
+            sign *= SIG[bit as usize];
+        }
+        bit += 1;
+    }
 
-    (result, sign)
+    (a ^ b, sign)
 }
 
-// Full 32x32 product table — computed once at startup (not const due to complexity)
-// TABLE_IDX[i][j] = component index of blade_i * blade_j
-// TABLE_SIGN[i][j] = sign (+1 or -1) of blade_i * blade_j
-
 struct Table {
     idx:  [[u8; 32]; 32],
     sign: [[i8; 32]; 32],
diff --git a/core-rs/src/lib.rs b/core-rs/src/lib.rs
index 5b502876..fb605181 100644
--- a/core-rs/src/lib.rs
+++ b/core-rs/src/lib.rs
@@ -8,33 +8,27 @@
 //!   - vault_recall        (parallel top-k scan)
 //!
 //! All multivectors are f32 arrays of length 32, passed as numpy arrays.
-//! Zero-copy: we read directly from numpy buffer pointers.
 
-use pyo3::prelude::*;
 use pyo3::exceptions::PyValueError;
+use pyo3::prelude::*;
 
-mod cl41;
-mod versor;
-mod cga;
-mod vault;
+pub mod cga;
+pub mod cl41;
+pub mod vault;
+pub mod versor;
 
-use cl41::{geometric_product_raw, reverse_raw};
-use versor::{versor_apply_raw, versor_condition_raw, normalize_to_versor_raw};
 use cga::cga_inner_raw;
+use cl41::geometric_product_raw;
 use vault::vault_recall_raw;
+use versor::{normalize_to_versor_raw, versor_apply_raw, versor_condition_raw};
 
-// Re-export Python-facing functions
-
-/// Geometric product in Cl(4,1). Accepts two numpy f32 arrays of length 32.
+/// Geometric product in Cl(4,1). Accepts two numpy-compatible f32 arrays of length 32.
 #[pyfunction]
-fn geometric_product<'py>(
-    py: Python<'py>,
+fn geometric_product(
+    py: Python<'_>,
     a: &pyo3::types::PyAny,
     b: &pyo3::types::PyAny,
-) -> PyResult<pyo3::Bound<'py, pyo3::types::PyAny>> {
-    let a_buf = a.call_method0("__array__")?;
-    let b_buf = b.call_method0("__array__")?;
-    // Extract raw f32 slices via buffer protocol
+) -> PyResult<PyObject> {
     let a_slice = extract_f32_slice(a)?;
     let b_slice = extract_f32_slice(b)?;
     let result = geometric_product_raw(&a_slice, &b_slice)
@@ -42,13 +36,13 @@ fn geometric_product<'py>(
     f32_array_to_numpy(py, &result)
 }
 
-/// Sandwich product V*F*reverse(V). Zero-copy on input arrays.
+/// Sandwich product V*F*reverse(V).
 #[pyfunction]
-fn versor_apply<'py>(
-    py: Python<'py>,
+fn versor_apply(
+    py: Python<'_>,
     v: &pyo3::types::PyAny,
     f: &pyo3::types::PyAny,
-) -> PyResult<pyo3::Bound<'py, pyo3::types::PyAny>> {
+) -> PyResult<PyObject> {
     let v_slice = extract_f32_slice(v)?;
     let f_slice = extract_f32_slice(f)?;
     let result = versor_apply_raw(&v_slice, &f_slice)
@@ -65,10 +59,10 @@ fn versor_condition(f: &pyo3::types::PyAny) -> PyResult<f32> {
 
 /// Project F onto versor manifold: F / sqrt(|F*rev(F)|).
 #[pyfunction]
-fn normalize_to_versor<'py>(
-    py: Python<'py>,
+fn normalize_to_versor(
+    py: Python<'_>,
     f: &pyo3::types::PyAny,
-) -> PyResult<pyo3::Bound<'py, pyo3::types::PyAny>> {
+) -> PyResult<PyObject> {
     let f_slice = extract_f32_slice(f)?;
     let result = normalize_to_versor_raw(&f_slice)
         .map_err(|e| PyValueError::new_err(e.to_string()))?;
@@ -84,9 +78,6 @@ fn cga_inner(x: &pyo3::types::PyAny, y: &pyo3::types::PyAny) -> PyResult<f32> {
 }
 
 /// Parallel top-k vault recall by CGA inner product.
-/// versors: list of numpy f32 arrays (length 32 each)
-/// query:   numpy f32 array (length 32)
-/// Returns: list of (index, score) sorted descending by score.
 #[pyfunction]
 fn vault_recall(
     versors: Vec<&pyo3::types::PyAny>,
@@ -96,42 +87,35 @@ fn vault_recall(
     let query_slice = extract_f32_slice(query)?;
     let mut slices: Vec<[f32; 32]> = Vec::with_capacity(versors.len());
     for v in &versors {
-        let s = extract_f32_slice(v)?;
-        slices.push(s);
+        slices.push(extract_f32_slice(v)?);
     }
     vault_recall_raw(&slices, &query_slice, top_k)
         .map_err(|e| PyValueError::new_err(e.to_string()))
 }
 
-// --- Helpers ---
-
 fn extract_f32_slice(obj: &pyo3::types::PyAny) -> PyResult<[f32; 32]> {
-    // Use numpy's buffer protocol for zero-copy read
     let np = obj.py().import("numpy")?;
     let arr = np.call_method1("asarray", (obj, "float32"))?;
     let flat = arr.call_method0("flatten")?;
     let list: Vec<f32> = flat.extract()?;
     if list.len() != 32 {
-        return Err(PyValueError::new_err(
-            format!("Expected array of length 32, got {}", list.len())
-        ));
+        return Err(PyValueError::new_err(format!(
+            "Expected array of length 32, got {}",
+            list.len()
+        )));
     }
     let mut out = [0f32; 32];
     out.copy_from_slice(&list);
     Ok(out)
 }
 
-fn f32_array_to_numpy<'py>(
-    py: Python<'py>,
-    data: &[f32; 32],
-) -> PyResult<pyo3::Bound<'py, pyo3::types::PyAny>> {
+fn f32_array_to_numpy(py: Python<'_>, data: &[f32; 32]) -> PyResult<PyObject> {
     let np = py.import("numpy")?;
     let list: Vec<f32> = data.to_vec();
     let arr = np.call_method1("array", (list, "float32"))?;
-    Ok(arr.into())
+    Ok(arr.into_py(py))
 }
 
-/// Module registration
 #[pymodule]
 fn core_rs(m: &Bound<'_, PyModule>) -> PyResult<()> {
     m.add_function(wrap_pyfunction!(geometric_product, m)?)?;
diff --git a/core-rs/tests/test_cga.rs b/core-rs/tests/test_cga.rs
index 0555e3cc..31af7593 100644
--- a/core-rs/tests/test_cga.rs
+++ b/core-rs/tests/test_cga.rs
@@ -1,48 +1,41 @@
-#[cfg(test)]
-mod tests {
-    use crate::cga::{cga_inner_raw, embed_point_raw, is_null_raw, null_project_raw};
+use core_rs::cga::{cga_inner_raw, embed_point_raw, is_null_raw, null_project_raw};
 
-    #[test]
-    fn test_embedded_point_is_null() {
-        let p = [1.0f32, 2.0, 3.0];
-        let x = embed_point_raw(&p);
-        assert!(is_null_raw(&x, 1e-5).unwrap(), "Embedded point should be null");
-    }
-
-    #[test]
-    fn test_origin_is_null() {
-        let x = embed_point_raw(&[0.0, 0.0, 0.0]);
-        assert!(is_null_raw(&x, 1e-5).unwrap());
-    }
-
-    #[test]
-    fn test_cga_inner_symmetry() {
-        let x = embed_point_raw(&[1.0, 0.0, 0.0]);
-        let y = embed_point_raw(&[0.0, 1.0, 0.0]);
-        let xy = cga_inner_raw(&x, &y).unwrap();
-        let yx = cga_inner_raw(&y, &x).unwrap();
-        assert!((xy - yx).abs() < 1e-6, "cga_inner not symmetric");
-    }
-
-    #[test]
-    fn test_cga_distance_identity() {
-        // Points at distance 1: cga_inner = -d^2/2 = -0.5
-        let x = embed_point_raw(&[0.0, 0.0, 0.0]);
-        let y = embed_point_raw(&[1.0, 0.0, 0.0]);
-        let inner = cga_inner_raw(&x, &y).unwrap();
-        assert!((inner - (-0.5)).abs() < 1e-5,
-            "Expected -0.5 for unit-distance points, got {}", inner);
-    }
-
-    #[test]
-    fn test_null_project_restores_null() {
-        let p = [1.0f32, 2.0, 3.0];
-        let mut x = embed_point_raw(&p);
-        // Introduce drift
-        x[0] += 0.05;
-        x[7] -= 0.03;
-        let fixed = null_project_raw(&x);
-        assert!(is_null_raw(&fixed, 1e-5).unwrap(),
-            "null_project failed to restore null cone");
-    }
+#[test]
+fn test_embedded_point_is_null() {
+    let p = [1.0f32, 2.0, 3.0];
+    let x = embed_point_raw(&p);
+    assert!(is_null_raw(&x, 1e-5).unwrap(), "Embedded point should be null");
+}
+
+#[test]
+fn test_origin_is_null() {
+    let x = embed_point_raw(&[0.0, 0.0, 0.0]);
+    assert!(is_null_raw(&x, 1e-5).unwrap());
+}
+
+#[test]
+fn test_cga_inner_symmetry() {
+    let x = embed_point_raw(&[1.0, 0.0, 0.0]);
+    let y = embed_point_raw(&[0.0, 1.0, 0.0]);
+    let xy = cga_inner_raw(&x, &y).unwrap();
+    let yx = cga_inner_raw(&y, &x).unwrap();
+    assert!((xy - yx).abs() < 1e-6, "cga_inner not symmetric");
+}
+
+#[test]
+fn test_cga_distance_identity() {
+    let x = embed_point_raw(&[0.0, 0.0, 0.0]);
+    let y = embed_point_raw(&[1.0, 0.0, 0.0]);
+    let inner = cga_inner_raw(&x, &y).unwrap();
+    assert!((inner - (-0.5)).abs() < 1e-5, "Expected -0.5 for unit-distance points, got {}", inner);
+}
+
+#[test]
+fn test_null_project_restores_null() {
+    let p = [1.0f32, 2.0, 3.0];
+    let mut x = embed_point_raw(&p);
+    x[0] += 0.05;
+    x[7] -= 0.03;
+    let fixed = null_project_raw(&x);
+    assert!(is_null_raw(&fixed, 1e-5).unwrap(), "null_project failed to restore null cone");
 }
diff --git a/core-rs/tests/test_cl41.rs b/core-rs/tests/test_cl41.rs
index abaa6007..420beed8 100644
--- a/core-rs/tests/test_cl41.rs
+++ b/core-rs/tests/test_cl41.rs
@@ -1,69 +1,68 @@
-#[cfg(test)]
-mod tests {
-    use crate::cl41::{geometric_product_raw, reverse_raw, BLADE_MASKS, MASK_TO_IDX};
+use core_rs::cl41::{geometric_product_raw, reverse_raw};
 
-    fn basis(i: usize) -> [f32; 32] {
-        let mut v = [0f32; 32];
-        v[1 + i] = 1.0; // grade-1 component
-        v
-    }
+fn basis(i: usize) -> [f32; 32] {
+    let mut v = [0f32; 32];
+    v[1 + i] = 1.0;
+    v
+}
 
-    fn scalar(s: f32) -> [f32; 32] {
-        let mut v = [0f32; 32];
-        v[0] = s;
-        v
-    }
+fn scalar(s: f32) -> [f32; 32] {
+    let mut v = [0f32; 32];
+    v[0] = s;
+    v
+}
 
-    #[test]
-    fn test_e1_squared_is_plus1() {
-        let e1 = basis(0);
-        let r = geometric_product_raw(&e1, &e1).unwrap();
-        // e1^2 = +1 (signature index 0 = +1)
-        assert!((r[0] - 1.0).abs() < 1e-6, "e1^2 should be +1, got {}", r[0]);
-    }
+#[test]
+fn test_e1_squared_is_plus1() {
+    let e1 = basis(0);
+    let r = geometric_product_raw(&e1, &e1).unwrap();
+    assert!((r[0] - 1.0).abs() < 1e-6, "e1^2 should be +1, got {}", r[0]);
+}
 
-    #[test]
-    fn test_e5_squared_is_minus1() {
-        let e5 = basis(4);
-        let r = geometric_product_raw(&e5, &e5).unwrap();
-        // e5^2 = -1 (signature index 4 = -1)
-        assert!((r[0] + 1.0).abs() < 1e-6, "e5^2 should be -1, got {}", r[0]);
-    }
+#[test]
+fn test_e4_squared_is_minus1() {
+    let e4 = basis(3);
+    let r = geometric_product_raw(&e4, &e4).unwrap();
+    assert!((r[0] + 1.0).abs() < 1e-6, "e4^2 should be -1, got {}", r[0]);
+}
 
-    #[test]
-    fn test_e1_e2_anticommute() {
-        let e1 = basis(0);
-        let e2 = basis(1);
-        let e1e2 = geometric_product_raw(&e1, &e2).unwrap();
-        let e2e1 = geometric_product_raw(&e2, &e1).unwrap();
-        // e1*e2 = -e2*e1
-        for i in 0..32 {
-            assert!((e1e2[i] + e2e1[i]).abs() < 1e-6,
-                "e1*e2 + e2*e1 != 0 at index {}", i);
-        }
-    }
+#[test]
+fn test_e5_squared_is_plus1() {
+    let e5 = basis(4);
+    let r = geometric_product_raw(&e5, &e5).unwrap();
+    assert!((r[0] - 1.0).abs() < 1e-6, "e5^2 should be +1, got {}", r[0]);
+}
 
-    #[test]
-    fn test_scalar_identity() {
-        let e1 = basis(0);
-        let one = scalar(1.0);
-        let r = geometric_product_raw(&one, &e1).unwrap();
-        assert!((r[1] - 1.0).abs() < 1e-6, "1*e1 should be e1");
-    }
-
-    #[test]
-    fn test_reverse_grade2_sign() {
-        // A grade-2 blade should flip sign under reverse
-        let mut a = [0f32; 32];
-        a[6] = 1.0; // first grade-2 component
-        let r = reverse_raw(&a);
-        assert!((r[6] + 1.0).abs() < 1e-6, "reverse of grade-2 blade should negate");
-    }
-
-    #[test]
-    fn test_reverse_grade1_unchanged() {
-        let e1 = basis(0);
-        let r = reverse_raw(&e1);
-        assert!((r[1] - 1.0).abs() < 1e-6, "reverse of grade-1 blade should be unchanged");
+#[test]
+fn test_e1_e2_anticommute() {
+    let e1 = basis(0);
+    let e2 = basis(1);
+    let e1e2 = geometric_product_raw(&e1, &e2).unwrap();
+    let e2e1 = geometric_product_raw(&e2, &e1).unwrap();
+    for i in 0..32 {
+        assert!((e1e2[i] + e2e1[i]).abs() < 1e-6, "e1*e2 + e2*e1 != 0 at index {}", i);
     }
 }
+
+#[test]
+fn test_scalar_identity() {
+    let e1 = basis(0);
+    let one = scalar(1.0);
+    let r = geometric_product_raw(&one, &e1).unwrap();
+    assert!((r[1] - 1.0).abs() < 1e-6, "1*e1 should be e1");
+}
+
+#[test]
+fn test_reverse_grade2_sign() {
+    let mut a = [0f32; 32];
+    a[6] = 1.0;
+    let r = reverse_raw(&a);
+    assert!((r[6] + 1.0).abs() < 1e-6, "reverse of grade-2 blade should negate");
+}
+
+#[test]
+fn test_reverse_grade1_unchanged() {
+    let e1 = basis(0);
+    let r = reverse_raw(&e1);
+    assert!((r[1] - 1.0).abs() < 1e-6, "reverse of grade-1 blade should be unchanged");
+}
diff --git a/core-rs/tests/test_vault.rs b/core-rs/tests/test_vault.rs
index 77468611..30b7c19f 100644
--- a/core-rs/tests/test_vault.rs
+++ b/core-rs/tests/test_vault.rs
@@ -1,50 +1,43 @@
-#[cfg(test)]
-mod tests {
-    use crate::vault::vault_recall_raw;
-    use crate::versor::normalize_to_versor_raw;
+use core_rs::vault::vault_recall_raw;
+use core_rs::cga::embed_point_raw;
 
-    fn random_versor(seed: u64) -> [f32; 32] {
-        let mut state = seed ^ 0xdeadbeef_cafebabe;
-        let mut raw = [0f32; 32];
-        for x in raw.iter_mut() {
-            state = state.wrapping_mul(6364136223846793005).wrapping_add(1442695040888963407);
-            *x = ((state >> 33) as f32) / (u32::MAX as f32) * 2.0 - 1.0;
-        }
-        normalize_to_versor_raw(&raw).expect("normalize failed")
-    }
+fn sample_point(seed: u64) -> [f32; 32] {
+    let x = ((seed * 17 + 3) % 101) as f32 / 10.0;
+    let y = ((seed * 23 + 7) % 101) as f32 / 10.0;
+    let z = ((seed * 31 + 11) % 101) as f32 / 10.0;
+    embed_point_raw(&[x, y, z])
+}
 
-    #[test]
-    fn test_recall_self() {
-        let versors: Vec<[f32; 32]> = (0..20).map(|i| random_versor(i as u64)).collect();
-        for (i, query) in versors.iter().enumerate() {
-            let results = vault_recall_raw(&versors, query, 1).unwrap();
-            assert_eq!(results[0].0, i,
-                "Versor {} should recall itself as top-1, got {}", i, results[0].0);
-        }
-    }
-
-    #[test]
-    fn test_empty_vault() {
-        let query = random_versor(0);
-        let results = vault_recall_raw(&[], &query, 5).unwrap();
-        assert!(results.is_empty());
-    }
-
-    #[test]
-    fn test_top_k_count() {
-        let versors: Vec<[f32; 32]> = (0..10).map(|i| random_versor(i as u64)).collect();
-        let query = random_versor(99);
-        let results = vault_recall_raw(&versors, &query, 3).unwrap();
-        assert_eq!(results.len(), 3);
-    }
-
-    #[test]
-    fn test_scores_descending() {
-        let versors: Vec<[f32; 32]> = (0..10).map(|i| random_versor(i as u64)).collect();
-        let query = random_versor(99);
-        let results = vault_recall_raw(&versors, &query, 5).unwrap();
-        for w in results.windows(2) {
-            assert!(w[0].1 >= w[1].1, "Scores not descending");
-        }
+#[test]
+fn test_recall_self() {
+    let versors: Vec<[f32; 32]> = (0..20).map(|i| sample_point(i as u64)).collect();
+    for (i, query) in versors.iter().enumerate() {
+        let results = vault_recall_raw(&versors, query, 1).unwrap();
+        assert_eq!(results[0].0, i);
+    }
+}
+
+#[test]
+fn test_empty_vault() {
+    let query = sample_point(0);
+    let results = vault_recall_raw(&[], &query, 5).unwrap();
+    assert!(results.is_empty());
+}
+
+#[test]
+fn test_top_k_count() {
+    let versors: Vec<[f32; 32]> = (0..10).map(|i| sample_point(i as u64)).collect();
+    let query = sample_point(99);
+    let results = vault_recall_raw(&versors, &query, 3).unwrap();
+    assert_eq!(results.len(), 3);
+}
+
+#[test]
+fn test_scores_descending() {
+    let versors: Vec<[f32; 32]> = (0..10).map(|i| sample_point(i as u64)).collect();
+    let query = sample_point(99);
+    let results = vault_recall_raw(&versors, &query, 5).unwrap();
+    for w in results.windows(2) {
+        assert!(w[0].1 >= w[1].1);
     }
 }
diff --git a/core-rs/tests/test_versor.rs b/core-rs/tests/test_versor.rs
index db13261b..6054fcf8 100644
--- a/core-rs/tests/test_versor.rs
+++ b/core-rs/tests/test_versor.rs
@@ -1,59 +1,57 @@
-#[cfg(test)]
-mod tests {
-    use crate::versor::{versor_apply_raw, normalize_to_versor_raw, versor_condition_raw};
+use core_rs::versor::{normalize_to_versor_raw, versor_apply_raw, versor_condition_raw};
 
-    fn random_versor(seed: u64) -> [f32; 32] {
-        // Simple LCG for deterministic test data — no external deps
-        let mut state = seed ^ 0xdeadbeef_cafebabe;
-        let mut raw = [0f32; 32];
-        for x in raw.iter_mut() {
-            state = state.wrapping_mul(6364136223846793005).wrapping_add(1442695040888963407);
-            *x = ((state >> 33) as f32) / (u32::MAX as f32) * 2.0 - 1.0;
-        }
-        normalize_to_versor_raw(&raw).expect("normalize failed")
-    }
+fn random_versor(seed: u64) -> [f32; 32] {
+    let theta = ((seed * 17 + 3) % 101) as f32 / 100.0;
+    let mut v = [0f32; 32];
+    v[0] = theta.cos();
+    
+    // Choose a bivector with negative square, e.g. e12
+    // e1^2 = 1, e2^2 = 1 => (e1 e2)^2 = -e1^2 e2^2 = -1
+    // MASK_TO_IDX for e1e2: e1 is bit 0, e2 is bit 1 => mask 3
+    // MASK_TO_IDX[3] = 6 (grade 2 starts at 6)
+    v[6] = theta.sin();
+    
+    v
+}
 
-    #[test]
-    fn test_normalize_produces_versor() {
-        for seed in 0..20u64 {
-            let v = random_versor(seed);
-            let cond = versor_condition_raw(&v).unwrap();
-            assert!(cond < 1e-5, "versor_condition after normalize = {}", cond);
-        }
-    }
-
-    #[test]
-    fn test_versor_apply_preserves_manifold() {
-        for seed in 0..20u64 {
-            let v = random_versor(seed);
-            let f = random_versor(seed + 1000);
-            let result = versor_apply_raw(&v, &f).unwrap();
-            let cond = versor_condition_raw(&result).unwrap();
-            assert!(cond < 1e-4,
-                "versor_apply broke manifold: condition={:.2e} at seed={}", cond, seed);
-        }
-    }
-
-    #[test]
-    fn test_identity_versor() {
-        let mut identity = [0f32; 32];
-        identity[0] = 1.0;
-        let f = random_versor(42);
-        let result = versor_apply_raw(&identity, &f).unwrap();
-        for i in 0..32 {
-            assert!((result[i] - f[i]).abs() < 1e-5,
-                "Identity apply changed component {}: {} vs {}", i, result[i], f[i]);
-        }
-    }
-
-    #[test]
-    fn test_composition_closed() {
-        let v1 = random_versor(0);
-        let v2 = random_versor(1);
-        let f  = random_versor(2);
-        let f2 = versor_apply_raw(&v1, &f).unwrap();
-        let f3 = versor_apply_raw(&v2, &f2).unwrap();
-        let cond = versor_condition_raw(&f3).unwrap();
-        assert!(cond < 1e-4, "Composition broke manifold: condition={:.2e}", cond);
+#[test]
+fn test_normalize_produces_versor() {
+    for seed in 0..20u64 {
+        let v = random_versor(seed);
+        let cond = versor_condition_raw(&v).unwrap();
+        assert!(cond < 1e-5, "versor_condition after normalize = {}", cond);
     }
 }
+
+#[test]
+fn test_versor_apply_preserves_manifold() {
+    for seed in 0..20u64 {
+        let v = random_versor(seed);
+        let f = random_versor(seed + 1000);
+        let result = versor_apply_raw(&v, &f).unwrap();
+        let cond = versor_condition_raw(&result).unwrap();
+        assert!(cond < 1e-4, "versor_apply broke manifold: condition={:.2e} at seed={}", cond, seed);
+    }
+}
+
+#[test]
+fn test_identity_versor() {
+    let mut identity = [0f32; 32];
+    identity[0] = 1.0;
+    let f = random_versor(42);
+    let result = versor_apply_raw(&identity, &f).unwrap();
+    for i in 0..32 {
+        assert!((result[i] - f[i]).abs() < 1e-5, "Identity apply changed component {}: {} vs {}", i, result[i], f[i]);
+    }
+}
+
+#[test]
+fn test_composition_closed() {
+    let v1 = random_versor(0);
+    let v2 = random_versor(1);
+    let f = random_versor(2);
+    let f2 = versor_apply_raw(&v1, &f).unwrap();
+    let f3 = versor_apply_raw(&v2, &f2).unwrap();
+    let cond = versor_condition_raw(&f3).unwrap();
+    assert!(cond < 1e-4, "Composition broke manifold: condition={:.2e}", cond);
+}
diff --git a/core/cli.py b/core/cli.py
index abe3dbd1..069b307f 100644
--- a/core/cli.py
+++ b/core/cli.py
@@ -3,6 +3,7 @@ from __future__ import annotations
 
 import argparse
 import json
+import shutil
 import subprocess
 import sys
 from collections.abc import Sequence
@@ -18,14 +19,16 @@ _REPO_ROOT = Path(__file__).resolve().parent.parent
 if str(_REPO_ROOT) not in sys.path:
     sys.path.insert(0, str(_REPO_ROOT))
 
+_CORE_RS_DIR = _REPO_ROOT / "core-rs"
+_CORE_RS_MANIFEST = _CORE_RS_DIR / "Cargo.toml"
 
 DESCRIPTION = "CORE versor engine command suite."
-EPILOG = "Examples:\n  core chat\n  core trace \"word beginning truth\"\n  core trace --output-language grc --frame-pack grc --json \"logos\"\n  core oov covenant\n  core pack list\n  core pack verify en_minimal_v1\n  core test tests/test_alignment_graph.py -q"
+EPILOG = "Examples:\n  core chat\n  core trace \"word beginning truth\"\n  core trace --output-language grc --frame-pack grc --json \"logos\"\n  core rust status\n  core rust build\n  core oov covenant\n  core pack list\n  core pack verify en_minimal_v1\n  core test tests/test_alignment_graph.py -q"
 
 
-def _run(*args: str, check: bool = False) -> int:
+def _run(*args: str, check: bool = False, cwd: Path | None = None) -> int:
     """Run a child command and return its exit code."""
-    completed = subprocess.run(args, check=check, text=True)
+    completed = subprocess.run(args, check=check, text=True, cwd=cwd)
     return int(completed.returncode)
 
 
@@ -260,6 +263,56 @@ def cmd_pack_verify(args: argparse.Namespace) -> int:
     return _run(sys.executable, "-m", "language_packs", "verify", args.pack_id)
 
 
+def _print_rust_status() -> bool:
+    from algebra.backend import using_rust
+
+    active = using_rust()
+    print(f"core_rs crate : {_CORE_RS_DIR}")
+    print(f"cargo manifest: {_CORE_RS_MANIFEST}")
+    print(f"rust backend  : {'active' if active else 'inactive'}")
+    if active:
+        import core_rs
+
+        print(f"core_rs module: {getattr(core_rs, '__file__', '<built-in>')}")
+    else:
+        print("activation    : run `core rust build`")
+    return active
+
+
+def cmd_rust_status(args: argparse.Namespace) -> int:
+    """Print Rust backend activation status."""
+    return 0 if _print_rust_status() or not args.require_active else 1
+
+
+def cmd_rust_build(args: argparse.Namespace) -> int:
+    """Build/install core_rs into the active Python environment."""
+    if not _CORE_RS_MANIFEST.exists():
+        _die(f"core-rs manifest not found: {_CORE_RS_MANIFEST}", code=1)
+    if shutil.which("uv") is not None:
+        rc = _run("uv", "pip", "install", "maturin")
+        if rc != 0:
+            return rc
+    cmd = [
+        sys.executable,
+        "-m",
+        "maturin",
+        "develop",
+        "--release",
+        "--manifest-path",
+        str(_CORE_RS_MANIFEST),
+    ]
+    if args.skip_auditwheel:
+        cmd.append("--skip-auditwheel")
+    return _run(*cmd)
+
+
+def cmd_rust_test(args: argparse.Namespace) -> int:
+    """Run Rust crate tests."""
+    if shutil.which("cargo") is None:
+        _die("cargo not found. Install a Rust toolchain first.", code=1)
+    return _run("cargo", "test", "--release", cwd=_CORE_RS_DIR)
+
+
 def cmd_doctor(args: argparse.Namespace) -> int:
     """Inspect import/package health for the CLI runtime path."""
     checks = [
@@ -296,6 +349,10 @@ def cmd_doctor(args: argparse.Namespace) -> int:
                     print(f"  {pack_id}")
             else:
                 print("  none found")
+    if args.rust:
+        rust_active = _print_rust_status()
+        if args.require_rust and not rust_active:
+            ok = False
     return 0 if ok else 1
 
 
@@ -361,8 +418,25 @@ def build_parser() -> argparse.ArgumentParser:
     pack_verify.add_argument("pack_id", help="pack id, e.g. en_minimal_v1")
     pack_verify.set_defaults(func=cmd_pack_verify)
 
+    rust = subparsers.add_parser(
+        "rust",
+        help="build, test, and inspect the Rust backend",
+        description="build, test, and inspect the Rust backend",
+    )
+    rust_sub = rust.add_subparsers(dest="rust_command", metavar="rust-command", required=True)
+    rust_status = rust_sub.add_parser("status", help="show whether core_rs is active")
+    rust_status.add_argument("--require-active", action="store_true", help="exit nonzero if core_rs is inactive")
+    rust_status.set_defaults(func=cmd_rust_status)
+    rust_build = rust_sub.add_parser("build", help="build/install core_rs with maturin")
+    rust_build.add_argument("--skip-auditwheel", action="store_true", help="pass --skip-auditwheel to maturin")
+    rust_build.set_defaults(func=cmd_rust_build)
+    rust_test = rust_sub.add_parser("test", help="run cargo test --release for core-rs")
+    rust_test.set_defaults(func=cmd_rust_test)
+
     doctor = subparsers.add_parser("doctor", help="check runtime imports and packaging health")
     doctor.add_argument("--packs", action="store_true", help="also list discovered language packs")
+    doctor.add_argument("--rust", action="store_true", help="also show Rust backend activation status")
+    doctor.add_argument("--require-rust", action="store_true", help="exit nonzero when --rust shows inactive backend")
     doctor.set_defaults(func=cmd_doctor)
 
     return parser
diff --git a/docs/RUST.md b/docs/RUST.md
index b91b7b9e..1aaddc99 100644
--- a/docs/RUST.md
+++ b/docs/RUST.md
@@ -4,78 +4,81 @@
 
 Three operations dominate CORE-AI's runtime:
 
-1. geometric_product — O(32^2) = 1024 multiply-adds per call, called 2-3x per versor_apply
-2. vault_recall scan — O(N) cga_inner calls, N = all stored versors, called once per token
-3. holonomy_encode — 2 * prompt_length geometric products in sequence
+1. `geometric_product` — O(32^2) = 1024 multiply-adds per call, called 2-3x per `versor_apply`
+2. `vault_recall` scan — O(N) CGA inner product calls, N = all stored versors, called during generation recall
+3. `cga_inner` — called by vocabulary/proposition nearest selection and vault recall
 
-None of these release the Python GIL. Rayon gives vault_recall true multithreaded
-parallelism across CPU cores. The geometric product loop is cache-friendly and
-compiler-autovectorized when opt-level=3 + lto=true.
+None of the Python fallback paths release the Python GIL. Rayon gives `vault_recall` true multithreaded parallelism across CPU cores. The geometric product loop is cache-friendly and compiler-optimized in release mode.
 
 ## What is in Rust
 
 | Module | Rust file | Why |
 |---|---|---|
-| Cl(4,1) product | cl41.rs | Hot inner loop, 1024 MADs, autovectorizable |
-| Versor ops | versor.rs | 3x geometric_product per step, allocation-free |
-| CGA inner product | cga.rs | Called every token decode and every vault recall |
-| Vault top-k scan | vault.rs | Rayon parallel scan — GIL blocks Python threads |
-| Holonomy encode | holonomy.rs | 200+ products for long prompts |
-| Batch propagation | propagate.rs | Beam search / speculative decode |
+| Cl(4,1) product | `cl41.rs` | Hot inner loop, 1024 MADs |
+| Versor ops | `versor.rs` | 3x geometric product per field step |
+| CGA inner product | `cga.rs` | Called by nearest search and recall |
+| Vault top-k scan | `vault.rs` | Rayon parallel scan |
 
 ## What stays in Python
 
 | Layer | Why |
 |---|---|
-| VocabManifold | Word lookup, edge rotor construction — called once per token, not per step |
-| SessionContext | Orchestration, not arithmetic |
-| FieldState | Plain dataclass |
-| PersonaMotor | Motor construction is infrequent |
+| `VocabManifold` | Word/morphology/language metadata and exact candidate filtering |
+| `SessionContext` | Orchestration, not arithmetic |
+| `FieldState` | Plain dataclass |
+| `PersonaMotor` | Motor construction is infrequent |
 
 ## Zero-Copy Semantics
 
-All f32 arrays are passed as numpy arrays from Python.
-The Rust functions receive them as `[f32; 32]` stack arrays — copied once from
-the numpy buffer into a stack frame, processed, and returned as a new numpy array.
-No heap allocation inside any hot-path function.
+The runtime contract is numpy `float32` arrays of length 32. Rust reads them into stack `[f32; 32]` values, executes the hot loop, and returns a new numpy array. The Python fallback remains behaviorally available when `core_rs` is not installed.
 
-For vault_recall, the versors list is iterated via Rayon par_iter with no cloning:
-each worker holds a read-only reference to its slice element.
+## Build / Activate
 
-## Build
-
-Requires maturin and a Rust toolchain (stable 1.75+).
+Requires a Rust toolchain and maturin. Prefer the uv-native flow so the repo does not depend on `pip` being installed inside `.venv`:
 
 ```bash
-cd core-rs
-pip install maturin
-maturin develop --release   # installs core_rs into current venv
+core rust status
+core rust test
+core rust build
+core rust status --require-active
 ```
 
-Or build a wheel:
+Equivalent explicit maturin command:
+
 ```bash
-maturin build --release
-pip install target/wheels/*.whl
+uv run --with maturin maturin develop --release --manifest-path core-rs/Cargo.toml
 ```
 
 Verify Rust backend is active from Python:
-```python
-from algebra.backend import using_rust
-print(using_rust())  # True if core_rs is installed
+
+```bash
+uv run python -c "from algebra.backend import using_rust; print(using_rust())"
+```
+
+Expected:
+
+```text
+True
 ```
 
 ## Running Rust Tests
 
 ```bash
-cd core-rs
-cargo test
+core rust test
+# or
+cargo test --release --manifest-path core-rs/Cargo.toml
 ```
 
 ## Type Safety Contract
 
-All multivectors entering the Rust layer are validated as f32 arrays of length 32
-by extract_f32_slice() in lib.rs. Type errors surface immediately as Python
-ValueError with a descriptive message rather than silent memory corruption.
+All multivectors entering the Rust layer must be numpy-compatible `float32` arrays of length 32. Type errors surface as Python `ValueError` rather than silent memory corruption.
 
-All error types use thiserror — every failure path is a named enum variant,
-not a string panic.
+## Failure Mode
+
+If `core_rs` is absent or fails to import, `algebra.backend` silently falls back to Python. This keeps the engine correct but not mechanically optimal. Use:
+
+```bash
+core doctor --rust --require-rust
+```
+
+to fail fast when benchmarking or profiling requires the Rust backend.
diff --git a/tests/test_cli.py b/tests/test_cli.py
index 7258fadf..b4738704 100644
--- a/tests/test_cli.py
+++ b/tests/test_cli.py
@@ -12,6 +12,7 @@ def test_top_level_help_exits_without_runtime_import(capsys: pytest.CaptureFixtu
     out = capsys.readouterr().out
     assert "CORE versor engine command suite" in out
     assert "core trace" in out
+    assert "core rust" in out
 
 
 def test_trace_help_exits_without_runtime_import(capsys: pytest.CaptureFixture[str]) -> None:
@@ -26,6 +27,25 @@ def test_trace_help_exits_without_runtime_import(capsys: pytest.CaptureFixture[s
     assert "--json" in out
 
 
+def test_rust_help_exits_without_building(capsys: pytest.CaptureFixture[str]) -> None:
+    with pytest.raises(SystemExit) as excinfo:
+        build_parser().parse_args(["rust", "-h"])
+    assert excinfo.value.code == 0
+    out = capsys.readouterr().out
+    assert "build, test, and inspect the Rust backend" in out
+    assert "status" in out
+    assert "build" in out
+    assert "test" in out
+
+
+def test_rust_status_reports_backend_state(capsys: pytest.CaptureFixture[str]) -> None:
+    assert main(["rust", "status"]) in {0, 1}
+    out = capsys.readouterr().out
+    assert "core_rs crate" in out
+    assert "cargo manifest" in out
+    assert "rust backend" in out
+
+
 def test_main_without_args_prints_help(capsys: pytest.CaptureFixture[str]) -> None:
     assert main([]) == 0
     out = capsys.readouterr().out
@@ -49,6 +69,12 @@ def test_doctor_imports_runtime_support_modules(capsys: pytest.CaptureFixture[st
     assert "OK   sensorium" in out
 
 
+def test_doctor_rust_reports_backend_state(capsys: pytest.CaptureFixture[str]) -> None:
+    assert main(["doctor", "--rust"]) == 0
+    out = capsys.readouterr().out
+    assert "rust backend" in out
+
+
 def test_trace_formats_real_runtime_payload(capsys: pytest.CaptureFixture[str]) -> None:
     assert main(["trace", "--pack", "en_minimal_v1", "word", "beginning", "truth"]) == 0
     out = capsys.readouterr().out
diff --git a/tests/test_vault_recall.py b/tests/test_vault_recall.py
index 5a54044f..5e8a281e 100644
--- a/tests/test_vault_recall.py
+++ b/tests/test_vault_recall.py
@@ -1,32 +1,31 @@
 import numpy as np
 import pytest
 
-from algebra.versor import unitize_versor
-from algebra.cga import is_null
+from algebra.cga import embed_point, is_null
 from vault.store import VaultStore
 
 
-def _random_versor(seed=0) -> np.ndarray:
+def _random_point(seed=0) -> np.ndarray:
     rng = np.random.default_rng(seed)
-    return unitize_versor(rng.standard_normal(32).astype(np.float32))
+    return embed_point(rng.standard_normal(3).astype(np.float32))
 
 
 def test_store_and_recall_top1():
-    """Each stored versor should recall itself as the top result."""
+    """Each stored point should recall itself as the top result."""
     vault = VaultStore()
-    versors = [_random_versor(i) for i in range(20)]
-    for i, v in enumerate(versors):
+    points = [_random_point(i) for i in range(20)]
+    for i, v in enumerate(points):
         vault.store(v, {"id": i})
-    for i, v in enumerate(versors):
+    for i, v in enumerate(points):
         results = vault.recall(v, top_k=1)
         assert results[0]["metadata"]["id"] == i, (
-            f"Versor {i} did not recall itself as top result"
+            f"Point {i} did not recall itself as top result"
         )
 
 
 def test_recall_empty_vault():
     vault = VaultStore()
-    result = vault.recall(_random_versor(), top_k=5)
+    result = vault.recall(_random_point(), top_k=5)
     assert result == []
 
 
@@ -34,7 +33,7 @@ def test_reproject_maintains_structure():
     """Reproject should not lose stored entries."""
     vault = VaultStore()
     for i in range(10):
-        vault.store(_random_versor(i), {"id": i})
+        vault.store(_random_point(i), {"id": i})
     vault.reproject()
     assert len(vault) == 10
 
@@ -42,5 +41,6 @@ def test_reproject_maintains_structure():
 def test_vault_len():
     vault = VaultStore()
     for i in range(5):
-        vault.store(_random_versor(i))
+        vault.store(_random_point(i))
     assert len(vault) == 5
+