geoit 0.0.2

//! Inductive governance families: parameterized governance generators
//! for VGA(n), CGA(n), PGA(n) and the morphisms between them.
//!
//! A `GovernanceFamily` captures the *pattern* of a governance parameterized by
//! spatial dimension n. `GovernanceRegistry` collects families and inter-family
//! morphisms, and `instantiate(n)` produces a complete `GovernanceCategory`.

use crate::algebra::mv::Mv;
use crate::algebra::signature::Signature;
use crate::error::SignatureError;
use crate::governance::category::{GovernanceCategory, GovernedMorphism};
use crate::governance::composition::Embedding;
use crate::governance::construction::Construction;
use crate::governance::expr::Expr;
use crate::governance::field::FieldOp;
use crate::governance::geom_class::GeomClass;
use crate::governance::geom_class::{inner_product_poly, norm_poly};
use crate::governance::governance::Governance;
use crate::governance::poly::Poly;
use crate::governance::reading::VariableMap;
use crate::scalar::Rat;

/// How signature parameters depend on spatial dimension n.
#[derive(Clone, Debug)]
pub enum SignatureFormula {
    /// Cl(0, 0, n) — vanilla Euclidean
    VGA,
    /// Cl(0, 1, n) — projective geometric algebra
    PGA,
    /// Cl(1, 0, n+1) — conformal geometric algebra
    CGA,
}

impl SignatureFormula {
    pub fn instantiate(&self, n: u8) -> Result<Signature, SignatureError> {
        match self {
            SignatureFormula::VGA => Signature::new(0, 0, n),
            SignatureFormula::PGA => Signature::new(0, 1, n),
            SignatureFormula::CGA => Signature::new(1, 0, n + 1),
        }
    }
}

/// Formula for a derived generator in terms of basis generators.
#[derive(Clone, Debug)]
pub enum DerivedGenFormula {
    /// eₒ = (e_{i0} + e_{h0}) / 2 — CGA origin
    ConformalOrigin,
    /// e∞ = e_{h0} - e_{i0} — CGA infinity
    ConformalInfinity,
}

impl DerivedGenFormula {
    pub fn instantiate(&self, sig: &Signature) -> Mv {
        match self {
            DerivedGenFormula::ConformalOrigin => {
                let i0 = sig.resolve_typed('i', 0).unwrap();
                let h0 = sig.resolve_typed('h', 0).unwrap();
                Mv::from_rat_terms(&[(1u64 << i0, Rat::new(1, 2)), (1u64 << h0, Rat::new(1, 2))])
            }
            DerivedGenFormula::ConformalInfinity => {
                let i0 = sig.resolve_typed('i', 0).unwrap();
                let h0 = sig.resolve_typed('h', 0).unwrap();
                Mv::from_rat_terms(&[(1u64 << i0, Rat::from(-1)), (1u64 << h0, Rat::from(1))])
            }
        }
    }
}

/// How to build the class constraints for a given n.
#[derive(Clone, Debug)]
pub enum ClassFormula {
    /// Grade-1 vectors, no constraints (VGA vector).
    Vector,
    /// CGA Point: grade 1, null norm, e∞ normalization.
    ConformalPoint,
    /// PGA Point: grade n (pseudoscalar of Euclidean subspace), weight ≠ 0.
    ProjectivePoint,
}

/// How to build the construction body for a given n.
#[derive(Clone, Debug)]
pub enum ConstructionFormula {
    /// Σᵢ pᵢ·eᵢ — linear combination of generators (VGA vector).
    LinearGenerators,
    /// CGA conformal embedding: Σ pᵢ·hᵢ + (-½Σpᵢ²)·e∞ + eₒ
    ConformalPoint,
    /// PGA point: weight-normalized (n-1)-blade construction.
    ProjectivePoint,
}

/// How dimension n embeds into dimension n+1 within a family.
#[derive(Clone, Debug)]
pub enum StepRule {
    /// gen(k) → gen(k) for all k < source.n(). VGA and PGA pattern.
    Identity,
    /// CGA pattern: conformal pair stays fixed, Euclidean gens map 1:1.
    ConformalExtend,
}

/// A parameterized family of governances.
#[derive(Clone, Debug)]
pub struct GovernanceFamily {
    pub name: String,
    pub signature: SignatureFormula,
    pub derived_gen_formulas: Vec<DerivedGenFormula>,
    pub class_formulas: Vec<(String, ClassFormula)>,
    pub construction_formulas: Vec<(String, ConstructionFormula)>,
    pub step_rule: Option<StepRule>,
}

impl GovernanceFamily {
    /// Instantiate this family at spatial dimension n.
    pub fn instantiate(&self, n: u8) -> Result<Governance, crate::error::Error> {
        let sig = self.signature.instantiate(n)?;
        let derived_gens: Vec<Mv> = self
            .derived_gen_formulas
            .iter()
            .map(|f| f.instantiate(&sig))
            .collect();

        let mut geom_classes = Vec::new();
        let mut constructions = Vec::new();

        for (idx, (name, formula)) in self.class_formulas.iter().enumerate() {
            let class = build_class(formula, &sig, &derived_gens, n);
            geom_classes.push(class);

            // Find matching construction formula
            if let Some((_, cf)) = self.construction_formulas.iter().find(|(cn, _)| cn == name) {
                let body = build_construction(cf, &sig, &derived_gens, n);
                constructions.push(Construction {
                    class_index: idx,
                    arity: construction_arity(cf, n),
                    body,
                });
            }
        }

        Ok(Governance::from_parts(
            sig,
            derived_gens,
            geom_classes,
            constructions,
            None,
        ))
    }

    /// Build the step embedding from dimension n to n+1.
    pub fn step_embedding(&self, n: u8) -> Result<Option<Embedding>, crate::error::Error> {
        let rule = match &self.step_rule {
            Some(r) => r,
            None => return Ok(None),
        };

        let source_sig = self.signature.instantiate(n)?;
        let target_sig = self.signature.instantiate(n + 1)?;

        let map: Vec<u8> = match rule {
            StepRule::Identity => (0..source_sig.n()).collect(),
            StepRule::ConformalExtend => (0..source_sig.n()).collect(),
        };

        let embedding =
            Embedding::new(source_sig, target_sig, map).map_err(crate::error::Error::Embedding)?;
        Ok(Some(embedding))
    }
}

/// How one family embeds into another at the same dimension.
#[derive(Clone, Debug)]
pub struct InterFamilyMorphism {
    pub source_family: String,
    pub target_family: String,
    pub name: String,
    /// How source generators map to target generators.
    pub mapping: InterFamilyMapping,
    /// Class preservation: (source_class_name, target_class_name).
    pub class_links: Vec<(String, String)>,
}

#[derive(Clone, Debug)]
pub enum InterFamilyMapping {
    /// VGA→CGA: VGA gen k (all h-type) maps to CGA gen k+2 (past conformal pair + first h).
    VGAtoCGA,
    /// VGA→PGA: VGA gen k maps to PGA gen k+1 (past the degenerate gen).
    VGAtoPGA,
}

impl InterFamilyMapping {
    pub fn instantiate(&self, source_sig: &Signature, _target_sig: &Signature) -> Vec<u8> {
        match self {
            InterFamilyMapping::VGAtoCGA => {
                // VGA Cl(0,0,n): gens are h0..h_{n-1} at flat indices 0..n-1
                // CGA Cl(1,0,n+1): i0 at 0, h0..h_n at 1..n+1
                // VGA h_k maps to CGA flat index k+2 (i0=0, h0=1, then h1=2, h2=3...)
                // Wait: CGA has i=1, d=0, h=n+1. Flat: gen0=i0, gen1=h0, gen2=h1, ...
                // VGA h_k (flat k) maps to CGA h_{k+1} (flat k+2)
                (0..source_sig.n()).map(|k| k + 2).collect()
            }
            InterFamilyMapping::VGAtoPGA => {
                // VGA Cl(0,0,n): gens at flat 0..n-1, all h-type
                // PGA Cl(0,1,n): d0 at flat 0, h0..h_{n-1} at flat 1..n
                // VGA h_k (flat k) maps to PGA h_k (flat k+1)
                (0..source_sig.n()).map(|k| k + 1).collect()
            }
        }
    }
}

/// The complete system: families + inter-family morphisms.
#[derive(Clone, Debug)]
pub struct GovernanceRegistry {
    pub families: Vec<GovernanceFamily>,
    pub inter_family_morphisms: Vec<InterFamilyMorphism>,
}

impl GovernanceRegistry {
    /// Instantiate all families at dimension n, producing a GovernanceCategory.
    pub fn instantiate(&self, n: u8) -> Result<GovernanceCategory, crate::error::Error> {
        let mut category = GovernanceCategory::new();

        // 1. Instantiate each family
        for family in &self.families {
            let gov = family.instantiate(n)?;
            category.add_governance(&family.name, gov);
        }

        // 2. Compute intra-family step morphisms (n-1 → n)
        if n > 1 {
            for (fam_idx, family) in self.families.iter().enumerate() {
                if let Some(embedding) = family.step_embedding(n - 1)? {
                    let prev_gov = family.instantiate(n - 1)?;
                    let prev_idx =
                        category.add_governance(format!("{}({})", family.name, n - 1), prev_gov);
                    let num_classes = category.governances[fam_idx].geom_classes.len();
                    category.add_morphism(GovernedMorphism {
                        name: format!("{}({}→{})", family.name, n - 1, n),
                        embedding,
                        source_gov_index: prev_idx,
                        target_gov_index: fam_idx,
                        class_map: (0..num_classes).map(Some).collect(),
                    });
                }
            }
        }

        // 3. Compute inter-family morphisms
        for inter in &self.inter_family_morphisms {
            let source_idx = category
                .names
                .iter()
                .position(|n| n == &inter.source_family);
            let target_idx = category
                .names
                .iter()
                .position(|n| n == &inter.target_family);
            if let (Some(si), Some(ti)) = (source_idx, target_idx) {
                let source_sig = category.governances[si].sig;
                let target_sig = category.governances[ti].sig;
                let map = inter.mapping.instantiate(&source_sig, &target_sig);
                if let Ok(embedding) = Embedding::new(source_sig, target_sig, map) {
                    // Resolve class links
                    let source_classes = &category.governances[si].geom_classes;
                    let class_map: Vec<Option<usize>> = (0..source_classes.len())
                        .map(|_| Some(0)) // simplified: map all to class 0
                        .collect();
                    category.add_morphism(GovernedMorphism {
                        name: inter.name.clone(),
                        embedding,
                        source_gov_index: si,
                        target_gov_index: ti,
                        class_map,
                    });
                }
            }
        }

        Ok(category)
    }

    /// Instantiate a single family at dimension n.
    pub fn instantiate_family(
        &self,
        family_name: &str,
        n: u8,
    ) -> Result<Governance, crate::error::Error> {
        let family = self
            .families
            .iter()
            .find(|f| f.name == family_name)
            .ok_or_else(|| {
                crate::error::Error::NotFound(crate::error::NotFoundError::Class(
                    family_name.into(),
                ))
            })?;
        family.instantiate(n)
    }
}

// ═══════════════════════════════════════════════════════════
// STANDARD REGISTRY
// ═══════════════════════════════════════════════════════════

/// The standard registry: VGA, CGA, PGA families with inter-family morphisms.
pub fn standard_registry() -> GovernanceRegistry {
    GovernanceRegistry {
        families: vec![vga_family(), cga_family(), pga_family()],
        inter_family_morphisms: vec![
            InterFamilyMorphism {
                source_family: "VGA".into(),
                target_family: "CGA".into(),
                name: "VGA→CGA".into(),
                mapping: InterFamilyMapping::VGAtoCGA,
                class_links: vec![("Vector".into(), "Point".into())],
            },
            InterFamilyMorphism {
                source_family: "VGA".into(),
                target_family: "PGA".into(),
                name: "VGA→PGA".into(),
                mapping: InterFamilyMapping::VGAtoPGA,
                class_links: vec![],
            },
        ],
    }
}

fn vga_family() -> GovernanceFamily {
    GovernanceFamily {
        name: "VGA".into(),
        signature: SignatureFormula::VGA,
        derived_gen_formulas: vec![],
        class_formulas: vec![("Vector".into(), ClassFormula::Vector)],
        construction_formulas: vec![("Vector".into(), ConstructionFormula::LinearGenerators)],
        step_rule: Some(StepRule::Identity),
    }
}

fn cga_family() -> GovernanceFamily {
    GovernanceFamily {
        name: "CGA".into(),
        signature: SignatureFormula::CGA,
        derived_gen_formulas: vec![
            DerivedGenFormula::ConformalOrigin,
            DerivedGenFormula::ConformalInfinity,
        ],
        class_formulas: vec![("Point".into(), ClassFormula::ConformalPoint)],
        construction_formulas: vec![("Point".into(), ConstructionFormula::ConformalPoint)],
        step_rule: Some(StepRule::ConformalExtend),
    }
}

fn pga_family() -> GovernanceFamily {
    GovernanceFamily {
        name: "PGA".into(),
        signature: SignatureFormula::PGA,
        derived_gen_formulas: vec![],
        class_formulas: vec![("Point".into(), ClassFormula::ProjectivePoint)],
        construction_formulas: vec![("Point".into(), ConstructionFormula::ProjectivePoint)],
        step_rule: Some(StepRule::Identity),
    }
}

// ═══════════════════════════════════════════════════════════
// CLASS AND CONSTRUCTION BUILDERS
// ═══════════════════════════════════════════════════════════

fn build_class(formula: &ClassFormula, sig: &Signature, derived_gens: &[Mv], n: u8) -> GeomClass {
    match formula {
        ClassFormula::Vector => GeomClass::grades_only(&[1]),
        ClassFormula::ConformalPoint => {
            let grade_mask = 0b10u64; // grade 1
            let vm = VariableMap::for_grade_mask(sig, grade_mask);
            let null_eq = norm_poly(sig, grade_mask, vm.num_vars, &vm.mask_to_var);
            let einf = &derived_gens[1];
            let ip_eq = inner_product_poly(
                einf,
                sig,
                grade_mask,
                Rat::ONE,
                vm.num_vars,
                &vm.mask_to_var,
            );
            GeomClass {
                grade_mask,
                equations: vec![null_eq, ip_eq],
                inequalities: vec![],
                field_op: FieldOp::ScalarProduct,
                expected_profile: None,
            }
        }
        ClassFormula::ProjectivePoint => {
            let grade = n; // PGA points live at grade n
            let grade_mask = 1u64 << grade;
            let vm = VariableMap::for_grade_mask(sig, grade_mask);
            // Weight inequality: the purely-Euclidean component must be nonzero
            // For PGA(n), the Euclidean pseudoscalar blade is e_{1}e_{2}...e_{n}
            let euclidean_mask: u64 = ((1u64 << n) - 1) << 1; // bits 1..n (skip d0 at bit 0)
            if let Some(&weight_var) = vm.mask_to_var.get(&euclidean_mask) {
                let weight_poly = Poly::variable(weight_var, vm.num_vars);
                GeomClass {
                    grade_mask,
                    equations: vec![],
                    inequalities: vec![weight_poly],
                    field_op: FieldOp::LeftContraction,
                    expected_profile: None,
                }
            } else {
                GeomClass {
                    grade_mask,
                    equations: vec![],
                    inequalities: vec![],
                    field_op: FieldOp::LeftContraction,
                    expected_profile: None,
                }
            }
        }
    }
}

fn construction_arity(formula: &ConstructionFormula, n: u8) -> usize {
    match formula {
        ConstructionFormula::LinearGenerators => n as usize,
        ConstructionFormula::ConformalPoint => n as usize,
        ConstructionFormula::ProjectivePoint => n as usize,
    }
}

fn build_construction(
    formula: &ConstructionFormula,
    sig: &Signature,
    _derived_gens: &[Mv],
    n: u8,
) -> Expr {
    match formula {
        ConstructionFormula::LinearGenerators => {
            // Σᵢ pᵢ·gen(i) for i in 0..n (all generators)
            build_linear_sum(0, n as usize, 0)
        }
        ConstructionFormula::ConformalPoint => {
            // P(x₁,...,xₙ) = Σ xᵢ·h_{i+1} + (-½Σxᵢ²)·e∞ + eₒ
            // h generators start at flat index 2 in CGA (i0=0, h0=1, h1=2, ...)
            // Actually: gen layout is i0, h0, h1, ..., h_n.
            // The Euclidean gens are h1..h_n (flat 2..n+1).
            let euclidean_start = sig.i() + sig.d() + 1; // skip i0 and h0 (conformal pair)

            // Euclidean part: Σ pᵢ · gen(euclidean_start + i)
            let euclidean = build_linear_sum(0, n as usize, euclidean_start as usize);

            // r² = Σ pᵢ²
            let r_squared = build_sum_of_squares(0, n as usize);

            // (-1/2) * r² * e∞
            let neg_half = Expr::Literal(crate::scalar::Scalar::Rat(Rat::new(-1, 2)));
            let conformal = Expr::Mul(
                Box::new(Expr::Mul(Box::new(neg_half), Box::new(r_squared))),
                Expr::dgen(1), // e∞ is derived gen 1
            );

            // euclidean + conformal + eₒ
            Expr::Add(
                Box::new(Expr::Add(Box::new(euclidean), Box::new(conformal))),
                Expr::dgen(0), // eₒ is derived gen 0
            )
        }
        ConstructionFormula::ProjectivePoint => {
            // PGA point construction: weight-normalized n-blade
            // For PGA(n) with sig Cl(0,1,n):
            // Generators: d0 (flat 0), h0..h_{n-1} (flat 1..n)
            // Point = Σ pᵢ · (appropriate (n-1)-blade involving d0) + e_{1..n}
            // This is a simplified construction that builds the blade directly
            build_pga_point_construction(n)
        }
    }
}

/// Build Σ param(i) * gen(gen_start + i) for i in 0..count.
fn build_linear_sum(param_start: usize, count: usize, gen_start: usize) -> Expr {
    if count == 0 {
        return Expr::Literal(crate::scalar::Scalar::from(0i64));
    }
    let mut result = *Expr::mul(Expr::param(param_start), Expr::gen(gen_start as u8));
    for i in 1..count {
        result = Expr::Add(
            Box::new(result),
            Expr::mul(
                Expr::param(param_start + i),
                Expr::gen((gen_start + i) as u8),
            ),
        );
    }
    result
}

/// Build Σ param(i)² for i in start..start+count.
fn build_sum_of_squares(start: usize, count: usize) -> Expr {
    let mut result = *Expr::mul(Expr::param(start), Expr::param(start));
    for i in 1..count {
        result = Expr::Add(
            Box::new(result),
            Expr::mul(Expr::param(start + i), Expr::param(start + i)),
        );
    }
    result
}

/// Build PGA point construction for PGA(n).
///
/// PGA(n) sig = Cl(0, 1, n). Generators: d0 (flat 0), h0..h_{n-1} (flat 1..n).
/// A point is a grade-n element with the pattern:
///   Σ_{k=0..n} (-1)^k · param(k) · (d0 ∧ all h-gens except h_k) + (h0 ∧ ... ∧ h_{n-1})
///
/// Example PGA(3): p0·e023 − p1·e013 + p2·e012 + e123
fn build_pga_point_construction(n: u8) -> Expr {
    // gen layout: d0=flat(0), h0=flat(1), h1=flat(2), ..., h_{n-1}=flat(n)
    let n = n as usize;

    // Build the Euclidean pseudoscalar: h0 ∧ h1 ∧ ... ∧ h_{n-1} = gen(1) * gen(2) * ... * gen(n)
    let euclidean_ps = build_gen_product((1..=n).map(|g| g as u8).collect());

    // Build the parametric blades: for each k in 0..n,
    // blade_k = d0 ∧ (all h-gens except h_k) with sign (-1)^k
    let mut terms: Vec<Box<Expr>> = Vec::new();
    for k in 0..n {
        // Generators for this blade: d0 (flat 0) plus all h-gens except h_k
        // h_k is at flat index k+1
        let mut gens: Vec<u8> = vec![0]; // d0
        for j in 0..n {
            if j != k {
                gens.push((j + 1) as u8); // h_j at flat j+1
            }
        }
        let blade = build_gen_product(gens);
        let param_times_blade = Expr::mul(Expr::param(k), Box::new(blade));
        if k % 2 == 1 {
            terms.push(Expr::neg(param_times_blade));
        } else {
            terms.push(param_times_blade);
        }
    }

    // Sum all parametric terms + Euclidean pseudoscalar
    let mut result = *terms.remove(0);
    for t in terms {
        result = Expr::Add(Box::new(result), t);
    }
    Expr::Add(Box::new(result), Box::new(euclidean_ps))
}

/// Build the geometric product of a list of generators: gen(gs[0]) * gen(gs[1]) * ...
fn build_gen_product(gs: Vec<u8>) -> Expr {
    assert!(!gs.is_empty());
    let mut result = Expr::Generator(gs[0]);
    for &g in &gs[1..] {
        result = Expr::Mul(Box::new(result), Expr::gen(g));
    }
    result
}

// ═══════════════════════════════════════════════════════════
// TESTS
// ═══════════════════════════════════════════════════════════

#[cfg(test)]
mod tests {
    use super::*;
    use crate::algebra::ops;
    use crate::governance;
    use crate::scalar::Scalar;

    #[test]
    fn vga_family_instantiate() {
        let reg = standard_registry();
        let gov = reg.instantiate_family("VGA", 3).unwrap();
        assert_eq!(gov.sig.n(), 3);
        assert_eq!(gov.sig.i(), 0);
        assert_eq!(gov.sig.d(), 0);
        assert_eq!(gov.sig.h(), 3);
        assert_eq!(gov.geom_classes.len(), 1);
        assert_eq!(gov.constructions.len(), 1);
        assert_eq!(gov.constructions[0].arity, 3);
    }

    #[test]
    fn vga_construct_and_govern() {
        let reg = standard_registry();
        let gov = reg.instantiate_family("VGA", 3).unwrap();
        let params = vec![Scalar::from(3i64), Scalar::from(4i64), Scalar::from(5i64)];
        let mv = gov.construct(0, &params).unwrap();
        let geoit = governance::govern(&mv, &gov, 0).unwrap();
        assert!(geoit.is_satisfied());
        assert_eq!(geoit.read_all().unwrap(), params);
    }

    #[test]
    fn vga_scales_to_5d() {
        let reg = standard_registry();
        let gov = reg.instantiate_family("VGA", 5).unwrap();
        assert_eq!(gov.sig.n(), 5);
        let params: Vec<Scalar> = (1..=5).map(|n| Scalar::from(n as i64)).collect();
        let mv = gov.construct(0, &params).unwrap();
        let geoit = governance::govern(&mv, &gov, 0).unwrap();
        assert!(geoit.is_satisfied());
    }

    #[test]
    fn cga_family_instantiate() {
        let reg = standard_registry();
        let gov = reg.instantiate_family("CGA", 2).unwrap();
        assert_eq!(gov.sig, Signature::new(1, 0, 3).unwrap());
        assert_eq!(gov.derived_gens.len(), 2);
        // Verify null vectors
        let eo = &gov.derived_gens[0];
        let einf = &gov.derived_gens[1];
        assert!(ops::norm_squared(eo, &gov.sig).is_zero());
        assert!(ops::norm_squared(einf, &gov.sig).is_zero());
    }

    #[test]
    fn cga_construct_point() {
        let reg = standard_registry();
        let gov = reg.instantiate_family("CGA", 2).unwrap();
        let params = vec![Scalar::from(1i64), Scalar::from(2i64)];
        let mv = gov.construct(0, &params).unwrap();
        let geoit = governance::govern(&mv, &gov, 0).unwrap();
        assert!(geoit.is_satisfied());
    }

    #[test]
    fn cga_scales_to_3d() {
        let reg = standard_registry();
        let gov = reg.instantiate_family("CGA", 3).unwrap();
        assert_eq!(gov.sig, Signature::new(1, 0, 4).unwrap());
        assert_eq!(gov.derived_gens.len(), 2);
        let params = vec![Scalar::from(1i64), Scalar::from(2i64), Scalar::from(3i64)];
        let mv = gov.construct(0, &params).unwrap();
        let geoit = governance::govern(&mv, &gov, 0).unwrap();
        assert!(geoit.is_satisfied());
    }

    #[test]
    fn registry_instantiate_category() {
        let reg = standard_registry();
        let cat = reg.instantiate(3).unwrap();
        // Should have VGA(3), CGA(3), PGA(3) plus step-morphism sources
        assert!(cat.governances.len() >= 3);
        // Should have inter-family morphisms
        assert!(cat.morphism_by_name("VGA→CGA").is_some());
    }

    #[test]
    fn pga_family_instantiate() {
        let reg = standard_registry();
        let gov = reg.instantiate_family("PGA", 3).unwrap();
        assert_eq!(gov.sig, Signature::new(0, 1, 3).unwrap());
        assert_eq!(gov.geom_classes.len(), 1);
    }

    #[test]
    fn step_embedding_vga() {
        let fam = vga_family();
        let emb = fam.step_embedding(2).unwrap().unwrap();
        assert_eq!(emb.source_sig, Signature::new(0, 0, 2).unwrap());
        assert_eq!(emb.target_sig, Signature::new(0, 0, 3).unwrap());
        assert_eq!(emb.generator_map, vec![0, 1]);
    }

    #[test]
    fn step_embedding_cga() {
        let fam = cga_family();
        let emb = fam.step_embedding(2).unwrap().unwrap();
        // CGA(2) = Cl(1,0,3), CGA(3) = Cl(1,0,4)
        assert_eq!(emb.source_sig, Signature::new(1, 0, 3).unwrap());
        assert_eq!(emb.target_sig, Signature::new(1, 0, 4).unwrap());
        // All 4 source gens map to same flat index in target
        assert_eq!(emb.generator_map, vec![0, 1, 2, 3]);
    }

    #[test]
    fn pga_construct_and_govern() {
        let reg = standard_registry();
        let gov = reg.instantiate_family("PGA", 3).unwrap();
        let params = vec![Scalar::from(2i64), Scalar::from(3i64), Scalar::from(4i64)];
        let mv = gov.construct(0, &params).unwrap();
        let geoit = governance::govern(&mv, &gov, 0).unwrap();
        assert!(geoit.is_satisfied());
        // Verify the Mv is grade 3 (point in PGA(3))
        for (blade, coeff) in geoit.mv().blades_bk() {
            if !coeff.is_zero() {
                assert_eq!(blade.grade(), 3, "PGA point should be pure grade 3");
            }
        }
        // Verify parameters round-trip
        let extracted = geoit.read_all().unwrap();
        assert_eq!(extracted, params);
    }

    #[test]
    fn pga_scales_to_4d() {
        let reg = standard_registry();
        let gov = reg.instantiate_family("PGA", 4).unwrap();
        assert_eq!(gov.sig, Signature::new(0, 1, 4).unwrap());
        assert_eq!(gov.constructions[0].arity, 4);
        let params: Vec<Scalar> = (1..=4).map(|n| Scalar::from(n as i64)).collect();
        let mv = gov.construct(0, &params).unwrap();
        let geoit = governance::govern(&mv, &gov, 0).unwrap();
        assert!(geoit.is_satisfied());
        // Grade 4 in PGA(4)
        for (blade, coeff) in geoit.mv().blades_bk() {
            if !coeff.is_zero() {
                assert_eq!(blade.grade(), 4);
            }
        }
    }

    #[test]
    fn pga_matches_hand_written() {
        // Verify the family-generated PGA(3) produces the same Mv as circuit_pga's manual construction
        let reg = standard_registry();
        let gov = reg.instantiate_family("PGA", 3).unwrap();
        let params = vec![Scalar::from(2i64), Scalar::from(3i64), Scalar::from(4i64)];
        let mv = gov.construct(0, &params).unwrap();
        // PGA(3) point should have the e123 component (weight = 1)
        // e123 = gen(1)*gen(2)*gen(3) = mask 0b1110
        assert_eq!(mv.coefficient(0b1110), Scalar::from(1i64));
        // e023 component = p0 = 2 (mask 0b1101)
        assert_eq!(mv.coefficient(0b1101), Scalar::from(2i64));
        // e013 component = -p1 = -3 (mask 0b1011)
        assert_eq!(mv.coefficient(0b1011), Scalar::from(-3i64));
        // e012 component = p2 = 4 (mask 0b0111)
        assert_eq!(mv.coefficient(0b0111), Scalar::from(4i64));
    }
}