6. Leaving cells out

This is where a lazy oracle pulls ahead of a stored matrix. A cell in a zero-length apparent pair needs no reduction and gives no bar, and the reduction can recognise it again as an owner when it needs it (through first_column_in_row). So it needn’t be in the basis at all:

/// Step 6: whether a cell is in a zero-length apparent pair, as either its birth or its death. Such cells needn't be
/// in the basis at all: the reduction finds them again as owners (through `first_column_in_row`) when it needs them.
fn in_zero_apparent_pair(volume: &Volume, cell: u64) -> bool {
    let value = volume.value(cell);
    // A birth: its oldest coface has the same value, and the cell is that coface's youngest face
    let birth = oldest_coface(volume, cell)
        .is_some_and(|coface| volume.value(coface) == value && youngest_face(volume, coface) == Some(cell));
    // A death: its youngest face has the same value, and the cell is that face's oldest coface
    let death = youngest_face(volume, cell)
        .is_some_and(|face| volume.value(face) == value && oldest_coface(volume, face) == Some(cell));
    birth || death
}

/// The columns to reduce, up to dimension `max_dim`, without the cells of zero-length apparent pairs. Cells are
/// streamed (never all held at once), so memory goes on the cells that are kept.
fn reduction_basis(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 && !in_zero_apparent_pair(volume, cell) {
            dims[dim].push(Cell(cell));
        }
    }
    for cells in &mut dims {
        cells.sort_by_key(|&cell| key(volume, cell.0));
    }
    SplitBasis::from_dimensions(dims)
}

A cell is in a zero-length apparent pair as a birth if its oldest coface has the same value and has it as youngest face, or as a death if its youngest face has the same value and has it as oldest coface. Both checks are local: a few voxel lookups. Cells are streamed, so only those kept are ever stored.

    // Step 6: leave the cells of zero-length apparent pairs out of the basis altogether
    let reduced = timed("build the reduction basis", || reduction_basis(&volume, 2));
    let reduced_size = basis_size(&reduced);
    let delta = Coboundary(&volume).with_basis(reduced);
    let lazy = timed("cohomology, clearing, reduction basis", || {
        barcode::cohomology(&delta, Algorithm::Clearing(algorithm.clone()))
    });
    println!(
        "    basis: {reduced_size} cells instead of {full_size} ({:.1}%)",
        100.0 * reduced_size as f64 / full_size as f64
    );
    if let (Some(full), Some(lazy)) = (full_peak, peak_memory_mb()) {
        println!("    peak memory: {lazy:.0} MB instead of {full:.0} MB");
    }

On a 96³ volume (7 million cells, 6.1 million of them up to dimension 2):

build every cell up to dimension 2                    2965.6 ms
cohomology, clearing, hooks                           3017.3 ms
build the reduction basis                             2211.4 ms
cohomology, clearing, reduction basis                 2020.2 ms
    basis: 338183 cells instead of 6110496 (5.5%)
    peak memory: 78 MB instead of 326 MB

The basis is 18× smaller, the reduction 1.5× faster, and peak memory a quarter. Building the basis is now the largest cost, and it is embarrassingly parallel (phlite_rips builds its basis on all cores).

For comparison, storing this complex explicitly, as PHAT or the phlite Python package do, takes about 36 bytes per cell and 8 per boundary entry, roughly 400 MB here, before any reduction.

Trade-offs

• Only the clearing algorithms can leave cells out: the standard algorithm needs every column.
• The test must be exactly right: leave out a cell that isn’t in a zero-length apparent pair and its bar is lost.
• It pays when most pairs are zero-length apparent pairs (lower-star filtrations of smooth data, Vietoris–Rips). On noise, where few pairs have length zero, the test costs more than it saves.