Skip to content
Snippets Groups Projects
Commit 6c0841d4 authored by Brad Nelson's avatar Brad Nelson
Browse files

initialze homology sketch

parent fa521b7f
No related branches found
No related tags found
No related merge requests found
......@@ -182,6 +182,11 @@ void SimplicialComplex::sortedOrder() {
[this](int i1, int i2) {\
return (full_function[i1].first==full_function[i2].first) ? full_function[i1].second < full_function[i2].second : full_function[i1].first < full_function[i2].first;});
// fill inverse filtration perm
inv_filtration_perm.resize(full_function.size());
for (size_t i = 0; i < filtration_perm.size(); i++) {
inv_filtration_perm[filtration_perm[i]] = i;
}
}
// // reworked extension function
......
......@@ -25,8 +25,13 @@ class SimplicialComplex{
std::vector<int> ncells;
// holds the permutation on cells that gives filtration order
// j = filtration_perm[i] means cell j is in the i'th location in filtration.
std::vector<size_t> filtration_perm;
// holds inverse permutation on cells for filtration order
// i = inv_filtraiton_perm[j] means cell j is in i'th locaiton in filtration.
std::vector<size_t> inv_filtration_perm;
// holds map to critical cell
// function_map[j] = [...] face that causes cell j to appear
std::vector<std::vector<int>> function_map;
......
#include <torch/extension.h>
#include <iostream>
#include <vector>
#include <limits>
/*
Function to build boundary of simplicial complex once sorted order is determined.
Just do single boundary matrix. No blocking by dimension.
column major order
INPUTS:
X - simplicial complex
IMPORTANT: assumes that X has been initialized, and filtration has been extended
MAXDIM - maximum homology dimension
OUTPUTS: boundary matrix in List of Lists format.
*/
std::vector<Cocyle> sorted_boundary(SimplicialComplex &X, int MAXDIM) {
// empty boundary matrix
std::vector<Cocycle> B;
// build boundary in sorted order using filtration_perm
// should also use filtration_perm to permute nzs in rows of columns
std::vector<int> row_inds; // row indices for column
for (size_t j : X.filtration_perm ) {
row_inds.clear(); // clear out column
// go through non-zeros in boundary
for (auto i : X.bdr[j].cochain) {
row_inds.push_back(X.inv_filtration_perm[i]); // push location of bounday cell in filtration
}
// append sorted row_inds to B
B.emplace_back(Cocycle(j, row_inds));
}
return B;
}
/*
Function to reduce boundary matrix.
*/
/*
Standard reduction algorithm on simplicial complex.
INPUTS:
X - simplicial complex
IMPORTANT: assumes that X has been initialized, and filtration has been extended
MAXDIM - maximum homology dimension
OUTPUTS: vector of tensors - t
t[k] is float32 tensor with barcode for dimension k
*/
std::vector<torch::Tensor> persistence_forward_hom(SimplicialComplex &X, int MAXDIM) {
// produce sort permutation on X
X.sortedOrder();
// initialize reutrn diagram
std::vector<torch::Tensor> diagram(MAXDIM+1); // return array
for (int k = 0; k < MAXDIM+1; k++) {
int Nk = X.numPairs(k); // number of bars in dimension k
// allocate return tensor
diagram[k] = torch::empty({Nk,2},
torch::TensorOptions().dtype(torch::kFloat32).layout(torch::kStrided).requires_grad(true));
// allocate critical indices
// TODO: do this in intialization since number of pairs is deterministic
X.backprop_lookup[k].resize(Nk);
}
// keep track of how many pairs we've put in diagram
std::vector<int> nbars(MAXDIM+1);
for (int k = 0; k < MAXDIM+1; k++) {
nbars[k] = 0;
}
// produce boundary matrix
// run standard reduction algorithm
// fill in diagram
return diagram;
}
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment