tawasuyu/hifas

1//! `cohomology` — la única cuenta que importa (§6.4, §7).
2//!
3//! `rank δ`, `dim H⁰`, `dim H¹` y una base del cokernel para localizar el nudo.
4
5use crate::error::Result;
6use crate::linalg::{F64Backend, LinAlg, Scalar};
7use crate::sheaf::Sheaf;
8
9/// Resultado del cálculo cohomológico sobre un haz.
10#[derive(Debug, Clone, PartialEq)]
11pub struct Cohomology {
12    /// `dim C⁰ = Σ dim F(v)`.
13    pub dim_c0: usize,
14    /// `dim C¹ = Σ dim F(e)`.
15    pub dim_c1: usize,
16    /// `rank δ`.
17    pub rank_delta: usize,
18    /// `dim H⁰ = dim_c0 − rank_delta` (secciones globales coherentes).
19    pub dim_h0: usize,
20    /// `dim H¹ = dim_c1 − rank_delta` (obstrucciones).
21    pub dim_h1: usize,
22    /// Si `dim_h1 > 0`, una base del cokernel: cada vector señala aristas
23    /// que participan en la obstrucción (para localizar el nudo).
24    pub obstruction_basis: Vec<Vec<Scalar>>,
25}
26
27/// Calcula la cohomología del haz con el backend `f64` por defecto (§7.1).
28pub fn compute(sheaf: &Sheaf) -> Result<Cohomology> {
29    compute_with(sheaf, &F64Backend::new())
30}
31
32/// Calcula la cohomología del haz con un backend de álgebra lineal explícito.
33///
34/// Pasos (§7): ensamblar `δ`, `rank(δ)`, `dim H⁰ = dimC⁰ − rank`,
35/// `dim H¹ = dimC¹ − rank`, y si `H¹ ≠ 0`, una base del cokernel.
36pub fn compute_with<L: LinAlg>(sheaf: &Sheaf, backend: &L) -> Result<Cohomology> {
37    let delta = sheaf.coboundary()?;
38    let dim_c0 = sheaf.dim_c0();
39    let dim_c1 = sheaf.dim_c1();
40    let rank_delta = backend.rank(&delta);
41
42    // rank ≤ min(dimC⁰, dimC¹), así que las restas nunca subdesbordan.
43    let dim_h0 = dim_c0 - rank_delta;
44    let dim_h1 = dim_c1 - rank_delta;
45
46    let obstruction_basis = if dim_h1 > 0 {
47        backend.cokernel_basis(&delta)
48    } else {
49        Vec::new()
50    };
51
52    Ok(Cohomology {
53        dim_c0,
54        dim_c1,
55        rank_delta,
56        dim_h0,
57        dim_h1,
58        obstruction_basis,
59    })
60}
61
62#[cfg(test)]
63mod tests {
64    use super::*;
65    use crate::gf2::Gf2Backend;
66
67    /// Prueba de cordura del haz constante (§7.2): para stalk de dimensión `d`,
68    /// `dim H⁰ = (#componentes)·d` y `dim H¹ = (primer nº de Betti)·d`.
69    ///
70    /// Aquí verificamos con `d = 1`, donde H⁰ = #componentes y H¹ = b₁. Además
71    /// cross-validamos que el backend exacto `GF(2)` coincide con `f64` en las
72    /// dimensiones (sin epsilon): si difirieran, uno de los dos estaría mal.
73    fn coho_constante(nv: usize, edges: &[(usize, usize)]) -> Cohomology {
74        let sheaf = Sheaf::constant(nv, edges, 1);
75        let f64_coho = compute(&sheaf).expect("haz constante bien formado");
76        let gf2_coho = compute_with(&sheaf, &Gf2Backend).expect("haz constante");
77        assert_eq!(
78            (f64_coho.dim_h0, f64_coho.dim_h1),
79            (gf2_coho.dim_h0, gf2_coho.dim_h1),
80            "f64 y GF(2) deben coincidir en las dimensiones cohomológicas"
81        );
82        f64_coho
83    }
84
85    #[test]
86    fn arbol_no_tiene_obstruccion() {
87        // Camino de 4 vértices (árbol conexo): H⁰ = 1, H¹ = 0.
88        let coho = coho_constante(4, &[(0, 1), (1, 2), (2, 3)]);
89        assert_eq!(coho.dim_h0, 1, "un árbol es conexo");
90        assert_eq!(coho.dim_h1, 0, "un árbol no tiene ciclos");
91        assert!(coho.obstruction_basis.is_empty());
92    }
93
94    #[test]
95    fn triangulo_tiene_un_ciclo() {
96        // Triángulo: conexo con un ciclo. H⁰ = 1, H¹ = 1.
97        let coho = coho_constante(3, &[(0, 1), (1, 2), (2, 0)]);
98        assert_eq!(coho.dim_h0, 1);
99        assert_eq!(coho.dim_h1, 1, "b₁ del triángulo es 1");
100        assert_eq!(coho.obstruction_basis.len(), 1);
101        // La obstrucción toca las 3 aristas del ciclo (soporte no nulo en todas).
102        let v = &coho.obstruction_basis[0];
103        assert_eq!(v.len(), 3);
104        assert!(v.iter().all(|x| x.abs() > 1e-9), "el nudo es el ciclo entero");
105    }
106
107    #[test]
108    fn dos_triangulos_pegados_tienen_dos_ciclos() {
109        // Dos triángulos que comparten la arista (1,2): V=4, E=5, C=1.
110        // b₁ = E − V + C = 5 − 4 + 1 = 2.
111        let edges = [(0, 1), (0, 2), (1, 2), (1, 3), (2, 3)];
112        let coho = coho_constante(4, &edges);
113        assert_eq!(coho.dim_h0, 1);
114        assert_eq!(coho.dim_h1, 2, "dos ciclos independientes");
115        assert_eq!(coho.obstruction_basis.len(), 2);
116    }
117
118    #[test]
119    fn bosque_desconexo_cuenta_componentes() {
120        // Dos aristas disjuntas: V=4, E=2, C=2, sin ciclos.
121        let coho = coho_constante(4, &[(0, 1), (2, 3)]);
122        assert_eq!(coho.dim_h0, 2, "dos componentes conexas");
123        assert_eq!(coho.dim_h1, 0);
124    }
125
126    #[test]
127    fn stalk_de_dimension_mayor_escala_betti() {
128        // Mismo triángulo pero con stalk d=2: H⁰ = 1·2, H¹ = 1·2.
129        let sheaf = Sheaf::constant(3, &[(0, 1), (1, 2), (2, 0)], 2);
130        let coho = compute(&sheaf).unwrap();
131        assert_eq!(coho.dim_h0, 2);
132        assert_eq!(coho.dim_h1, 2);
133    }
134}