1. A boundary matrix

The least phlite needs is a label type for cells and a matrix oracle.

The label

/// A cell of the volume's cubical complex, by its number (see `Volume`).
#[derive(Clone, Copy, Debug, PartialEq, Eq, PartialOrd, Ord, Hash)]
pub struct Cell(pub u64);

// Mark `Cell` as a label for rows and columns. Its order (derived: by number) breaks ties between equal values.
impl BasisElement for Cell {}

BasisElement marks a type as a row or column label. It must be Ord (the order breaks ties between cells with the same filtration value) and Clone, and it is cloned a lot, so it should be small: here one u64.

The oracle

/// Step 1: the boundary matrix, as an oracle. Nothing is stored: a column is computed when it is asked for.
#[derive(Clone, Copy)]
pub struct Boundary<'a>(pub &'a Volume);

impl MatrixOracle for Boundary<'_> {
    type CoefficientField = Z2;
    type ColT = Cell;
    type RowT = Cell;

    fn column(&self, col: Cell) -> impl Iterator<Item = (Z2, Cell)> {
        self.0.faces(col.0).map(|face| (Z2::one(), Cell(face)))
    }
}

Three associated types say what the matrix is over (Z2, whose one non-zero element is zero-sized) and what labels its columns and rows; column returns an iterator over a column’s entries. The boundary is computed from the cell’s number each time: the oracle holds just a reference to the volume.

Running it

To reduce a matrix, phlite also needs the order on rows and which columns to reduce, in that order. With no filtration, with_trivial_filtration() puts every cell at time 0, so rows are ordered by label alone, and the basis must list the columns in the same order:

/// The columns to reduce: every cell up to dimension `max_dim`, by dimension, each in the order of the rows. With no
/// filtration that is just the order on labels: by number.
fn basis_by_number(volume: &Volume, max_dim: usize) -> SplitBasis<Cell> {
    let mut dims: Vec<Vec<Cell>> = vec![vec![]; max_dim + 1];
    for cell in 0..volume.n_cells() {
        let dim = volume.dimension(cell);
        if dim <= max_dim {
            dims[dim].push(Cell(cell));
        }
    }
    SplitBasis::from_dimensions(dims)
}

    // Step 1: no filtration (every cell at time 0): every finite bar has length zero, so only the essential ones are
    // left, which count the holes of a solid cube: one component, and nothing else
    let d = Boundary(&volume).with_trivial_filtration();
    let d = d.with_basis(basis_by_number(&volume, 3));
    let betti = barcode::homology(&d, Algorithm::Standard(Standard::new()), false);
    println!("No filtration: {}", summary(&volume, &betti));

barcode::homology reduces the boundary with the standard algorithm. With every cell at time 0, every finite bar has length zero and is dropped, which leaves the essential bars: the Betti numbers of a solid cube.

No filtration: 1 bars in H0, 0 in H1 and 0 in H2

The basis must be in the rows’ order. The first version of this example sorted the basis by voxel value while the rows, with no filtration, were ordered by number. The reduction then pairs the wrong cells, and this example reported thousands of loops and cavities in a solid cube. phlite can’t check the order for you without visiting every column, so this is worth a test of its own.

You can now compute the homology of anything you can describe cell by cell.