Module ateams.arithmetic

Uses a variant of the Chen/Kerber (2011) and PHAT (2017) twist_reduce algorithm to compute the persistent homology of the complex specified by the flat boundary matrix and the filtration over the field $\mathbb{Z}/p\mathbb{Z}$.

Args

boundary

Full boundary matrix (i.e. Cubical.matrices.full).

filtration

A list of indices specifying the order in which the cells in the complex will be added. The assumed order is adding cells in order of dimension.

homology

The homology group for which persistence is computed.

p

Characteristic of the field $\mathbb{Z}/p\mathbb{Z}$.

breaks

Index locations for cells of each degree (i.e. Cubical.breaks).

Computes the persistence pairs of the complex corresponding to the provided boundary matrix and filtration.

Args

boundary : np.array

Flattened boundary matrix given by Complex.matrices.full.

filtration : np.array

Permutation on the order of the columns of the boundary matrix. Here, we assume that only cells of dimension homology are being permuted. Shuffling the order of cells of different dimensions will result in incorrect computations.

homology : int

Homology group we're interested in; corresponds to the dimension of permuted cells.

breaks : np.array

Index ranges for cells by dimension, given by Complex.breaks.

Returns

A list of [birth, death] pairs.

Returns a basis for the kernel of A.

Args

M : MatrixReduction

A MatrixReduction object.

A : np.array

A 2D NumPy array to be reduced.

Returns

A 2D NumPy array representing the reduced basis of $\mathrm{ker}(A)$.

Draws a uniform random sample from the kernel of A using the MatrixReduction object M.

Args

M : MatrixReduction

Contains the RREF form of a matrix.

A : np.ndarray

2D NumPy array whose kernel is computed and sampled from.

Returns

A vector sampled at uniform random from the kernel of the matrix stored in M.

Uses the LinBox (Block?) Lanczos black-box solver to get a uniform random element of the kernel of the coboundary.

Args

coboundary

Memory-contiguous one-dimensional coboundary array, like that of the Cubical Complex object.

zeros

Memory-contiguous array of row indices that will be included in the coboundary submatrix.

faces

The number of faces per plaquette.

columns

Integer representing the number of columns in the coboundary submatrix; should be the total number of faces in the complex.

p

Characteristic of the field $\mathbb Z/p\mathbb Z$.

maxTries

How many times do we try to get a nonzero solution from the black-box solver? Default is 8.

Returns

A C++ std::vector with spin assignments.

Uses the LinBox (Block?) Lanczos black-box solver to get a uniform random element of the kernel of the coboundary.

Args

coboundary

Memory-contiguous one-dimensional coboundary array, like that of the Cubical Complex object.

zeros

Memory-contiguous array of row indices that will be included in the coboundary submatrix.

p

Characteristic of the field $\mathbb Z/p\mathbb Z$.

maxTries

How many times do we try to get a nonzero solution from the black-box solver? Default is 8.

Returns

A C++ std::vector with spin assignments.

A class for putting matrices over finite fields into reduced row echelon form (RREF); works faster than comparable method in the galois library, accelerated by sparse matrix methods and cached finite field arithmetic.

Args

characteristic : int

Characteristic of the finite field.

parallel : bool=False

Should computations be done in parallel? Defaults to False.

minBlockSize : int=32

Smallest number of columns reduced by a single thread at one time.

maxBlockSize : int=64

Largest number of columns reduced by a single thread at one time.

cores : int=2

Number of cores/threads available for use at run-time.

Suppose we need to compute a basis for the kernel of an $m \times n$ matrix $A$ with entries in the finite field $\mathbb F_3$. The standard algorithm for kernel computation is as follows:

1. create an augmented matrix $C$ by adjoining the identity matrix $I_n$ on the right of $A^\top$, so $C= \left[ A^\top \mid I_n \right]$;
2. put $C$ in RREF, so $C = \left[ P \mid K \right]$;
3. the rows of $K$ corresponding to zero rows of $P$ form a basis for the kernel of $A$.

Code to perform these computations looks like

    import numpy as np
    from ateams.arithmetic import MatrixReduction, FINT, Kernel

    F = 3
    m, n = 50, 60
    A = np.random.randint(0, F, size=(m,n))
    I = np.identity(n)
    C = np.concatenate([A.T, I], axis=1, dtype=FINT)

    M = MatrixReduction(F)
    M.RREF(C, m)

    high = M.HighestZeroRow(m)
    reduced = M.ToArray()[:,:m]
    basis = reduced[high:]

But can be shortened to

    import numpy as np
    from ateams.arithmetic import MatrixReduction, FINT, Kernel

    F = 3
    m, n = 50, 60
    A = np.random.randint(0, F, size=(m,n), dtype=FINT)

    M = MatrixReduction(F)
    basis = Kernel(M, A)

Using the Kernel() convenience method.

Computes the persistent homology of a complex with coefficients in an arbitrary finite field.

Args

characteristic : int

Characteristic of the finite field over which we do computations.

flattened : list

List-of-lists representing the original flattened boundary matrix for the complex.

homology : int=-1

Degree of homology observed for homological percolation events; if not specified, the homology of the codimension-1 skeleton is reported.

The Persistence.ComputePercolationEvents() method computes the birth times of giant cycles given a complete filtration — that is, a filtration where the terminal element is the entire cubical complex. For example, the code below finds the birth times of giant cycles on the default $3 \times 3$ cubical torus: the filtration stored in filtration adds each cube to the cubical complex in order of construction. At completion, events contains the times {17, 21} which correspond to crossings of the meridian and equator of the torus, respectively.

    from ateams.arithmetic import Persistence
    from ateams.structures import Complex

    L = Complex().fromCorners([3,3], field=3)
    P = Persistence(L.field.characteristic, L.flattened, homology=1)

    filtration = np.arange(len(L.flattened))
    events = P.ComputePercolationEvents(filtration)

The Persistence.ComputeBettiNumbers() method provides functionality to compute the Betti numbers for an entire subcomplex. The subcomplex is specified by a list of indices corresponding to the cells included in the subcomplex. Initializing the Persistence object as before, we can get the Betti numbers of the $2$-torus encoded by the Complex above by passing a list of all cells' indices in the flattened boundary matrix: since there are $36$ cells (including vertices, edges, and squares), the subcomplex is just the range of integers $0-36$:

    from ateams.arithmetic import Persistence
    from ateams.structures import Complex

    L = Complex().fromCorners([3,3], field=3)
    P = Persistence(L.field.characteristic, L.flattened)

    subcomplex = np.arange(len(L.flattened))

    bettis = P.ComputeBettiNumbers(subcomplex)

The result in bettis is the list [1,2,1]. By deleting an edge, we reduce the number of squares by two, so we should expect the homology of a one-edge-deleted $2$-torus to have betti numbers [1,2,0]. To access the indices of the flattened boundary matrix corresponding to cells of each dimension, we can use the Complex.tranches property which maps an integer dimension to a (start, stop) pair: for example, L.tranches[1] is [9, 27] in the complex above, so we know that all entries in the flattened boundary matrix between indices $9$ and $27$ (right-exclusive) correspond to edges. To test our theory, we construct a subcomplex by deleting the first edge added:

    from ateams.arithmetic import Persistence
    from ateams.structures import Complex

    L = Complex().fromCorners([3,3], field=3)
    P = Persistence(L.field.characteristic, L.flattened)

    subcomplex = np.arange(len(L.flattened))
    firstEdge = L.tranches[1][0]
    subcomplex = np.concatenate([subcomplex[:firstEdge], subcomplex[firstEdge+1:]])

    bettis = P.ComputeBettiNumbers(subcomplex)

The result in bettis is the list [1,2,0], as expected. Note that we only need to pass a homology parameter to the Persistence constructor if we call Persistence.ComputePercolationEvents().