--- a/README Tue Feb 03 11:31:10 2009 -0800
+++ b/README Sun Mar 22 11:18:07 2009 -0700
@@ -19,10 +19,6 @@
[rlog]: http://www.arg0.net/rlog
[SYNAPS]: http://www-sop.inria.fr/galaad/synaps/
-## Configuration
- The path to CGAL's Makefile is expected to be set in `$CGAL_MAKEFILE`, the rest
- is just the usual CMake configuration
-
## Building
To build examples, create a directory build (to keep everything in one place),
go to that directory and run cmake and make:
--- a/examples/rips/rips.cpp Tue Feb 03 11:31:10 2009 -0800
+++ b/examples/rips/rips.cpp Sun Mar 22 11:18:07 2009 -0700
@@ -1,6 +1,7 @@
#include <topology/rips.h>
#include <topology/filtration.h>
#include <topology/static-persistence.h>
+#include <topology/dynamic-persistence.h>
#include <topology/persistence-diagram.h>
#include <utilities/containers.h> // for BackInsertFunctor
@@ -26,8 +27,9 @@
typedef Generator::Simplex Smplx;
typedef std::vector<Smplx> SimplexVector;
typedef Filtration<SimplexVector, unsigned> Fltr;
-typedef StaticPersistence<> Persistence;
-typedef PersistenceDiagram<> PDgm;
+//typedef StaticPersistence<> Persistence;
+typedef DynamicPersistenceChains<> Persistence;
+typedef PersistenceDiagram<> PDgm;
int main(int argc, char* argv[])
@@ -51,7 +53,7 @@
SimplexVector complex;
// Generate 2-skeleton of the Rips complex for epsilon = 50
- rips.generate(2, 50, make_push_back_functor(complex));
+ rips.generate(2, 10, make_push_back_functor(complex));
std::sort(complex.begin(), complex.end(), Smplx::VertexComparison()); // unnecessary
rInfo("Generated complex of size: %d", complex.size());
@@ -73,4 +75,39 @@
std::ofstream ofs("rips-diagrams");
boost::archive::binary_oarchive oa(ofs);
oa << dgms;
+
+ // Output cycles
+ for (Persistence::OrderIndex cur = p.begin(); cur != p.end(); ++cur)
+ {
+ Persistence::OrderDescriptor::Cycle& cycle = cur->cycle;
+
+ if (!cur->sign()) // only negative simplices have non-empty cycles
+ {
+ Persistence::OrderIndex birth = cur->pair; // the cycle that cur killed was born when we added birth (another simplex)
+
+ const Smplx& b = f.simplex(f.begin() + (birth - p.begin())); // eventually this will be encapsulated behind an interface
+ const Smplx& d = f.simplex(f.begin() + (cur - p.begin()));
+
+ if (b.dimension() != 1) continue;
+ std::cout << "Pair: (" << size(b) << ", " << size(d) << ")" << std::endl;
+ } else if (cur->pair == cur) // positive could be unpaired
+ {
+ const Smplx& b = f.simplex(f.begin() + (cur - p.begin()));
+ if (b.dimension() != 1) continue;
+
+ std::cout << "Unpaired birth: " << size(b) << std::endl;
+ cycle = cur->chain;
+ }
+
+ // Iterate over the cycle
+ for (Persistence::OrderDescriptor::Cycle::const_iterator si = cycle.begin();
+ si != cycle.end(); ++si)
+ {
+ const Smplx& s = f.simplex(f.begin() + (*si - p.begin()));
+ //std::cout << s.dimension() << std::endl;
+ const Smplx::VertexContainer& vertices = s.vertices(); // std::vector<Vertex> where Vertex = Distances::IndexType
+ AssertMsg(vertices.size() == s.dimension() + 1, "dimension of a simplex is one less than the number of its vertices");
+ std::cout << vertices[0] << " " << vertices[1] << std::endl;
+ }
+ }
}
--- a/include/topology/dynamic-persistence.h Tue Feb 03 11:31:10 2009 -0800
+++ b/include/topology/dynamic-persistence.h Sun Mar 22 11:18:07 2009 -0700
@@ -4,8 +4,11 @@
#include "static-persistence.h"
#include <utilities/types.h>
+#include <boost/progress.hpp>
+
#ifdef COUNTERS
static Counter* cTrailLength = GetCounter("persistence/pair/traillength"); // the size of matrix U in RU decomposition
+static Counter* cChainLength = GetCounter("persistence/pair/chainlength"); // the size of matrix V in R=DV decomposition
#endif // COUNTERS
template<class Data, class OrderDescriptor_, class ConsistencyIndex>
@@ -97,6 +100,7 @@
using Parent::begin;
using Parent::end;
+ using Parent::size;
// Struct: TranspositionVisitor
//
@@ -139,6 +143,140 @@
ConsistencyComparison ccmp_;
};
+template<class Data, class OrderDescriptor_, class ConsistencyIndex>
+struct ChainData_: public Data
+{
+ typedef ChainData_<Data, OrderDescriptor_, ConsistencyIndex> Self;
+
+ typedef typename OrderTraits<typename OrderDescriptor_::
+ template RebindData<Self>::
+ other>::Index OrderIndex;
+ typedef typename OrderDescriptor_::
+ Chains::template rebind<OrderIndex>::other::Chain Chain;
+
+ template<class Comparison>
+ struct ConsistencyComparison: public std::binary_function<const OrderIndex&, const OrderIndex&, bool>
+ {
+ ConsistencyComparison(const Comparison& cmp = Comparison()): cmp_(cmp) {}
+
+ bool operator()(const OrderIndex& a, const OrderIndex& b) const { return cmp_(a->consistency, b->consistency); }
+
+ Comparison cmp_;
+ };
+
+ Chain chain;
+ ConsistencyIndex consistency;
+};
+
+/**
+ * Class: DynamicPersistenceChains
+ *
+ * TODO: below comment is incorrect; nothing dynamic about this yet.
+ * Derives from StaticPersistence and allows one to update persistence
+ * after a transposition of two contiguous simplices in a filtration.
+ * In addition to reduced cycles, it stores each OrderElement's chains,
+ * i.e. in addition to matrix R, it stores matrix V in vineyard notation.
+ *
+ * Template parameters:
+ * Data_ - auxilliary contents to store with each OrderElement
+ * OrderDescriptor_ - class describing how the order is stored; it defaults to <VectorOrderDescriptor>
+ * which serves as a prototypical class
+ */
+template<class Data_ = Empty<>,
+ class OrderDescriptor_ = VectorOrderDescriptor<>,
+ class ConsistencyIndex_ = size_t,
+ class ConsistencyComparison_ = std::less<ConsistencyIndex_> >
+class DynamicPersistenceChains:
+ public StaticPersistence<ChainData_<Data_, OrderDescriptor_, ConsistencyIndex_>,
+ OrderDescriptor_>
+{
+ public:
+ typedef Data_ Data;
+ typedef ChainData_<Data_, OrderDescriptor_, ConsistencyIndex_> ChainData;
+ typedef StaticPersistence<ChainData, OrderDescriptor_> Parent;
+
+ typedef typename Parent::Traits Traits;
+ typedef typename Parent::OrderDescriptor OrderDescriptor;
+ typedef typename Parent::OrderComparison OrderComparison;
+ typedef typename Parent::OrderIndex OrderIndex;
+ typedef ConsistencyIndex_ ConsistencyIndex;
+ typedef ThreeOutcomeCompare<
+ typename ChainData::
+ template ConsistencyComparison<ConsistencyComparison_> > ConsistencyComparison;
+
+ /**
+ * Constructor: DynamicPersistenceChains()
+ * TODO: write a description
+ *
+ * Template parameters:
+ * Filtration - filtration of the complex whose persistence we are computing
+ */
+ template<class Filtration> DynamicPersistenceChains(const Filtration& f,
+ const OrderComparison& ocmp = OrderComparison(),
+ const ConsistencyComparison_& ccmp = ConsistencyComparison_());
+
+ void pair_simplices();
+
+ // Function: transpose(i)
+ // Tranpose i and the next element.
+ // Returns: true iff the pairing switched.
+ // TODO
+ //bool transpose(OrderIndex i) { return transpose(i, TranspositionVisitor()); }
+
+ // TODO: the main missing piece to be dynamic
+ //template<class Visitor>
+ //bool transpose(OrderIndex i, const Visitor& visitor = Visitor());
+
+ using Parent::begin;
+ using Parent::end;
+ using Parent::size;
+
+ // Struct: TranspositionVisitor
+ //
+ // For example, a VineardVisitor could implement this archetype.
+ struct TranspositionVisitor
+ {
+ // Function: transpose(i)
+ // This function is called before transposition is processed
+ // (at the very beginning of <transpose(i, visitor)>). It is meant to update any structures
+ // that may need to be updated, but perhaps it has other uses as well.
+ void transpose(OrderIndex i) const {}
+
+ // Function: switched(i, type)
+ // This function is called after the transposition if the switch in pairing has occured.
+ // `i` is the index of the preceding simplex after the transposition.
+ // `type` indicates the <SwitchType>.
+ void switched(OrderIndex i, SwitchType type) const {}
+ };
+
+ protected:
+ using Parent::order;
+
+ private:
+ void swap(OrderIndex i, OrderIndex j);
+ void pairing_switch(OrderIndex i, OrderIndex j);
+
+ struct PairingChainsVisitor: public Parent::PairVisitor
+ {
+ // TODO: this is specialized for std::vector
+ PairingChainsVisitor(OrderIndex bg, ConsistencyComparison ccmp, unsigned size):
+ bg_(bg), ccmp_(ccmp), show_progress(size) {}
+
+ void init(OrderIndex i) const { i->consistency = i - bg_; i->chain.append(i, ccmp_); }
+ void update(OrderIndex j, OrderIndex i) const { j->chain.add(i->pair->chain, ccmp_); }
+ void finished(OrderIndex i) const { CountBy(cChainLength, i->chain.size()); ++show_progress; }
+
+ OrderIndex bg_;
+ ConsistencyComparison ccmp_;
+
+ mutable boost::progress_display
+ show_progress;
+ };
+
+ ConsistencyComparison ccmp_;
+};
+
+
#include "dynamic-persistence.hpp"
#endif // __DYNAMIC_PERSISTENCE_H__
--- a/include/topology/dynamic-persistence.hpp Tue Feb 03 11:31:10 2009 -0800
+++ b/include/topology/dynamic-persistence.hpp Sun Mar 22 11:18:07 2009 -0700
@@ -20,6 +20,8 @@
#endif // COUNTERS
+/* Trails */
+
template<class D, class OD, class CI, class CC>
template<class Filtration>
DynamicPersistenceTrails<D,OD,CI,CC>::
@@ -256,3 +258,22 @@
j_pair->pair = i;
}
}
+
+
+/* Chains */
+
+template<class D, class OD, class CI, class CC>
+template<class Filtration>
+DynamicPersistenceChains<D,OD,CI,CC>::
+DynamicPersistenceChains(const Filtration& f, const OrderComparison& ocmp, const CC& ccmp):
+ Parent(f, ocmp), ccmp_(ccmp)
+{}
+
+template<class D, class OD, class CI, class CC>
+void
+DynamicPersistenceChains<D,OD,CI,CC>::
+pair_simplices()
+{
+ Parent::pair_simplices(begin(), end(), PairingChainsVisitor(begin(), ccmp_, size()));
+}
+