--- a/CMakeLists.txt Sun Nov 24 13:35:45 2019 -0800
+++ b/CMakeLists.txt Mon Mar 06 12:31:30 2023 -0800
@@ -1,5 +1,11 @@
project (Dionysus)
-cmake_minimum_required (VERSION 2.4)
+cmake_minimum_required (VERSION 3.1)
+
+# Default to Release
+if (NOT CMAKE_BUILD_TYPE)
+ set (CMAKE_BUILD_TYPE "Release" CACHE STRING "Choose the type of build, options are: Debug Release RelWithDebInfo MinSizeRel." FORCE)
+ set_property (CACHE CMAKE_BUILD_TYPE PROPERTY STRINGS "Debug" "Release" "MinSizeRel" "RelWithDebInfo")
+endif (NOT CMAKE_BUILD_TYPE)
option (logging "Build Dionysus with logging on" OFF)
option (counters "Build Dionysus with counters on" OFF)
--- a/bindings/python/CMakeLists.txt Sun Nov 24 13:35:45 2019 -0800
+++ b/bindings/python/CMakeLists.txt Mon Mar 06 12:31:30 2023 -0800
@@ -47,7 +47,7 @@
dionysus/distances.py
)
-get_target_property (_dionysus_location _dionysus LOCATION)
+set (_dionysus_location $<TARGET_FILE:_dionysus>)
add_custom_target (dionysus-link ALL
${CMAKE_COMMAND} -E create_symlink ${_dionysus_location} ${CMAKE_CURRENT_BINARY_DIR}/dionysus/_dionysus.so
DEPENDS _dionysus dionysus)
--- a/examples/homology-zigzags/CMakeLists.txt Sun Nov 24 13:35:45 2019 -0800
+++ b/examples/homology-zigzags/CMakeLists.txt Mon Mar 06 12:31:30 2023 -0800
@@ -1,5 +1,5 @@
set (targets
- rips-pairwise
+ rips-pairwise-zz
M-ZZ
dM-ZZ
iR-ZZ
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/examples/homology-zigzags/rips-pairwise-zz.cpp Mon Mar 06 12:31:30 2023 -0800
@@ -0,0 +1,222 @@
+/***************************************************************************
+
+ rips-pairwise: Rips filtration + persistent homology
+ Copyright (C) 2009-2012 Dmitriy Morozov
+
+ Adapted from its original version ("rips-pairwise.cpp" in
+ the Dionysus library) by Steve Oudot (2012).
+
+ Changelog (2012-11-22):
+
+ - the barcode is now output on a log scale and in a format that is
+ compatible with the Matlab layer of the PLEX 2.5 library,
+
+ - the barcode representation is now by closed intervals instead of
+ half-open intervals.
+
+ This program is free software: you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program. If not, see <http://www.gnu.org/licenses/>.
+
+***************************************************************************/
+
+#include <topology/rips.h>
+#include <topology/filtration.h>
+#include <topology/static-persistence.h>
+#include <topology/dynamic-persistence.h>
+#include <topology/persistence-diagram.h>
+
+#include <geometry/l2distance.h>
+#include <geometry/distances.h>
+
+#include <utilities/containers.h> // for BackInsertFunctor
+#include <utilities/timer.h>
+
+#include <vector>
+
+#include <boost/program_options.hpp>
+
+
+typedef PairwiseDistances<PointContainer, L2Distance> PairDistances;
+typedef PairDistances::DistanceType DistanceType;
+typedef PairDistances::IndexType Vertex;
+
+typedef Rips<PairDistances> Generator;
+typedef Generator::Simplex Smplx;
+typedef Filtration<Smplx> Fltr;
+// typedef StaticPersistence<> Persistence;
+typedef DynamicPersistenceChains<> Persistence;
+typedef PersistenceDiagram<> PDgm;
+
+typedef std::vector<std::list<std::pair<double, double> > > IntervalsVector;
+
+// Forward declarations of auxilliary functions
+void write_intervals(std::ostream& out, const IntervalsVector& intervals, int skeleton_dimension);
+void program_options(int argc, char* argv[], std::string& infilename, Dimension& skeleton, DistanceType& max_distance, std::string& diagram_name);
+
+int main(int argc, char* argv[])
+{
+ Dimension skeleton;
+ DistanceType max_distance;
+ std::string infilename, diagram_name;
+
+ program_options(argc, argv, infilename, skeleton, max_distance, diagram_name);
+ std::ofstream diagram_out(diagram_name.c_str());
+ std::cout << "Diagram: " << diagram_name << std::endl;
+
+ PointContainer points;
+ read_points(infilename, points);
+
+ PairDistances distances(points);
+ Generator rips(distances);
+ Generator::Evaluator size(distances);
+ Fltr f;
+
+ // Generate 2-skeleton of the Rips complex for epsilon = 50
+ rips.generate(skeleton, max_distance, make_push_back_functor(f));
+ std::cout << "# Generated complex of size: " << f.size() << std::endl;
+
+ // Generate filtration with respect to distance and compute its persistence
+ f.sort(Generator::Comparison(distances));
+
+ Timer persistence_timer; persistence_timer.start();
+ Persistence p(f);
+ p.pair_simplices();
+ persistence_timer.stop();
+
+#if 1
+ // Create intervals DS
+ IntervalsVector intervals(skeleton);
+
+ // Output cycles
+ Persistence::SimplexMap<Fltr> m = p.make_simplex_map(f);
+ for (Persistence::iterator cur = p.begin(); cur != p.end(); ++cur)
+ {
+ const Persistence::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 = m[birth];
+ const Smplx& d = m[cur];
+
+ // if (b.dimension() != 1) continue;
+ // std::cout << "Pair: (" << size(b) << ", " << size(d) << ")" << std::endl;
+ if (b.dimension() >= skeleton) continue;
+ // diagram_out << b.dimension() << " " << size(b) << " " << size(d) << std::endl;
+ intervals[b.dimension()].push_back(std::pair<double,double>(size(b), size(d)));
+ } else if (cur->unpaired()) // positive could be unpaired
+ {
+ const Smplx& b = m[cur];
+ // if (b.dimension() != 1) continue;
+
+ // std::cout << "Unpaired birth: " << size(b) << std::endl;
+ // cycle = cur->chain; // TODO
+ if (b.dimension() >= skeleton) continue;
+ // diagram_out << b.dimension() << " " << size(b) << " inf" << std::endl;
+ intervals[b.dimension()].push_back(std::pair<double,double>(size(b), -1));
+ }
+
+ // Iterate over the cycle
+ // for (Persistence::Cycle::const_iterator si = cycle.begin();
+ // si != cycle.end(); ++si)
+ // {
+ // const Smplx& s = m[*si];
+ // //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;
+ // }
+ }
+#endif
+
+ persistence_timer.check("# Persistence timer");
+
+ std::cout << "Writing intervals...";
+ // Note (hack): use epsilons[vertices.size()-2]/2 as minimal value for the log-scale intervals to avoid intervals starting at -infinity and thus a scaling effect
+ write_intervals(diagram_out, intervals, skeleton);
+ std::cout << " done." << std::endl;
+
+}
+
+
+const double LOG2 = std::log(2.0);
+double log2(double x) {
+ return std::log(x) / LOG2;
+}
+
+void write_intervals(std::ostream& out, const IntervalsVector& intervals, int skeleton_dimension) {
+ out << "I = { ";
+ for (int d = 0; d<skeleton_dimension; d++) {
+ out << "[ ";
+ for (std::list<std::pair<double,double> >::const_iterator pit = intervals[d].begin(); pit != intervals[d].end(); pit++) {
+ if (pit->first == 0)
+ out << "[-Inf;";
+ else
+ out << "[" << log2(pit->first) << ";";
+ if (pit->second >= 0)
+ out << log2(pit->second) << "] ";
+ else
+ out << "Inf] ";
+ }
+ // add dummy interval if intervals list is empty
+ if (intervals[d].empty())
+ out << "[0;0] ";
+ out << "] ,";
+ }
+ out << "} " << std::endl;
+}
+
+
+void program_options(int argc, char* argv[], std::string& infilename, Dimension& skeleton, DistanceType& max_distance, std::string& diagram_name)
+{
+ namespace po = boost::program_options;
+
+ po::options_description hidden("Hidden options");
+ hidden.add_options()
+ ("input-file", po::value<std::string>(&infilename), "Point set whose Rips zigzag we want to compute");
+
+ po::options_description visible("Allowed options", 100);
+ visible.add_options()
+ ("help,h", "produce help message")
+ ("skeleton-dimsnion,s", po::value<Dimension>(&skeleton)->default_value(2), "Dimension of the Rips complex we want to compute")
+ ("max-distance,m", po::value<DistanceType>(&max_distance)->default_value(Infinity), "Maximum value for the Rips complex construction")
+ ("diagram,d", po::value<std::string>(&diagram_name), "Filename where to output the persistence diagram");
+#if LOGGING
+ std::vector<std::string> log_channels;
+ visible.add_options()
+ ("log,l", po::value< std::vector<std::string> >(&log_channels), "log channels to turn on (info, debug, etc)");
+#endif
+
+ po::positional_options_description pos;
+ pos.add("input-file", 1);
+
+ po::options_description all; all.add(visible).add(hidden);
+
+ po::variables_map vm;
+ po::store(po::command_line_parser(argc, argv).
+ options(all).positional(pos).run(), vm);
+ po::notify(vm);
+
+#if LOGGING
+ for (std::vector<std::string>::const_iterator cur = log_channels.begin(); cur != log_channels.end(); ++cur)
+ stderrLog.subscribeTo( RLOG_CHANNEL(cur->c_str()) );
+#endif
+
+ if (vm.count("help") || !vm.count("input-file"))
+ {
+ std::cout << "Usage: " << argv[0] << " [options] input-file" << std::endl;
+ std::cout << visible << std::endl;
+ std::abort();
+ }
+}
--- a/examples/homology-zigzags/rips-pairwise.cpp Sun Nov 24 13:35:45 2019 -0800
+++ /dev/null Thu Jan 01 00:00:00 1970 +0000
@@ -1,222 +0,0 @@
-/***************************************************************************
-
- rips-pairwise: Rips filtration + persistent homology
- Copyright (C) 2009-2012 Dmitriy Morozov
-
- Adapted from its original version ("rips-pairwise.cpp" in
- the Dionysus library) by Steve Oudot (2012).
-
- Changelog (2012-11-22):
-
- - the barcode is now output on a log scale and in a format that is
- compatible with the Matlab layer of the PLEX 2.5 library,
-
- - the barcode representation is now by closed intervals instead of
- half-open intervals.
-
- This program is free software: you can redistribute it and/or modify
- it under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- This program is distributed in the hope that it will be useful,
- but WITHOUT ANY WARRANTY; without even the implied warranty of
- MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
- GNU General Public License for more details.
-
- You should have received a copy of the GNU General Public License
- along with this program. If not, see <http://www.gnu.org/licenses/>.
-
-***************************************************************************/
-
-#include <topology/rips.h>
-#include <topology/filtration.h>
-#include <topology/static-persistence.h>
-#include <topology/dynamic-persistence.h>
-#include <topology/persistence-diagram.h>
-
-#include <geometry/l2distance.h>
-#include <geometry/distances.h>
-
-#include <utilities/containers.h> // for BackInsertFunctor
-#include <utilities/timer.h>
-
-#include <vector>
-
-#include <boost/program_options.hpp>
-
-
-typedef PairwiseDistances<PointContainer, L2Distance> PairDistances;
-typedef PairDistances::DistanceType DistanceType;
-typedef PairDistances::IndexType Vertex;
-
-typedef Rips<PairDistances> Generator;
-typedef Generator::Simplex Smplx;
-typedef Filtration<Smplx> Fltr;
-// typedef StaticPersistence<> Persistence;
-typedef DynamicPersistenceChains<> Persistence;
-typedef PersistenceDiagram<> PDgm;
-
-typedef std::vector<std::list<std::pair<double, double> > > IntervalsVector;
-
-// Forward declarations of auxilliary functions
-void write_intervals(std::ostream& out, const IntervalsVector& intervals, int skeleton_dimension);
-void program_options(int argc, char* argv[], std::string& infilename, Dimension& skeleton, DistanceType& max_distance, std::string& diagram_name);
-
-int main(int argc, char* argv[])
-{
- Dimension skeleton;
- DistanceType max_distance;
- std::string infilename, diagram_name;
-
- program_options(argc, argv, infilename, skeleton, max_distance, diagram_name);
- std::ofstream diagram_out(diagram_name.c_str());
- std::cout << "Diagram: " << diagram_name << std::endl;
-
- PointContainer points;
- read_points(infilename, points);
-
- PairDistances distances(points);
- Generator rips(distances);
- Generator::Evaluator size(distances);
- Fltr f;
-
- // Generate 2-skeleton of the Rips complex for epsilon = 50
- rips.generate(skeleton, max_distance, make_push_back_functor(f));
- std::cout << "# Generated complex of size: " << f.size() << std::endl;
-
- // Generate filtration with respect to distance and compute its persistence
- f.sort(Generator::Comparison(distances));
-
- Timer persistence_timer; persistence_timer.start();
- Persistence p(f);
- p.pair_simplices();
- persistence_timer.stop();
-
-#if 1
- // Create intervals DS
- IntervalsVector intervals(skeleton);
-
- // Output cycles
- Persistence::SimplexMap<Fltr> m = p.make_simplex_map(f);
- for (Persistence::iterator cur = p.begin(); cur != p.end(); ++cur)
- {
- const Persistence::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 = m[birth];
- const Smplx& d = m[cur];
-
- // if (b.dimension() != 1) continue;
- // std::cout << "Pair: (" << size(b) << ", " << size(d) << ")" << std::endl;
- if (b.dimension() >= skeleton) continue;
- // diagram_out << b.dimension() << " " << size(b) << " " << size(d) << std::endl;
- intervals[b.dimension()].push_back(std::pair<double,double>(size(b), size(d)));
- } else if (cur->unpaired()) // positive could be unpaired
- {
- const Smplx& b = m[cur];
- // if (b.dimension() != 1) continue;
-
- // std::cout << "Unpaired birth: " << size(b) << std::endl;
- // cycle = cur->chain; // TODO
- if (b.dimension() >= skeleton) continue;
- // diagram_out << b.dimension() << " " << size(b) << " inf" << std::endl;
- intervals[b.dimension()].push_back(std::pair<double,double>(size(b), -1));
- }
-
- // Iterate over the cycle
- // for (Persistence::Cycle::const_iterator si = cycle.begin();
- // si != cycle.end(); ++si)
- // {
- // const Smplx& s = m[*si];
- // //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;
- // }
- }
-#endif
-
- persistence_timer.check("# Persistence timer");
-
- std::cout << "Writing intervals...";
- // Note (hack): use epsilons[vertices.size()-2]/2 as minimal value for the log-scale intervals to avoid intervals starting at -infinity and thus a scaling effect
- write_intervals(diagram_out, intervals, skeleton);
- std::cout << " done." << std::endl;
-
-}
-
-
-const double LOG2 = std::log(2.0);
-double log2(double x) {
- return std::log(x) / LOG2;
-}
-
-void write_intervals(std::ostream& out, const IntervalsVector& intervals, int skeleton_dimension) {
- out << "I = { ";
- for (int d = 0; d<skeleton_dimension; d++) {
- out << "[ ";
- for (std::list<std::pair<double,double> >::const_iterator pit = intervals[d].begin(); pit != intervals[d].end(); pit++) {
- if (pit->first == 0)
- out << "[-Inf;";
- else
- out << "[" << log2(pit->first) << ";";
- if (pit->second >= 0)
- out << log2(pit->second) << "] ";
- else
- out << "Inf] ";
- }
- // add dummy interval if intervals list is empty
- if (intervals[d].empty())
- out << "[0;0] ";
- out << "] ,";
- }
- out << "} " << std::endl;
-}
-
-
-void program_options(int argc, char* argv[], std::string& infilename, Dimension& skeleton, DistanceType& max_distance, std::string& diagram_name)
-{
- namespace po = boost::program_options;
-
- po::options_description hidden("Hidden options");
- hidden.add_options()
- ("input-file", po::value<std::string>(&infilename), "Point set whose Rips zigzag we want to compute");
-
- po::options_description visible("Allowed options", 100);
- visible.add_options()
- ("help,h", "produce help message")
- ("skeleton-dimsnion,s", po::value<Dimension>(&skeleton)->default_value(2), "Dimension of the Rips complex we want to compute")
- ("max-distance,m", po::value<DistanceType>(&max_distance)->default_value(Infinity), "Maximum value for the Rips complex construction")
- ("diagram,d", po::value<std::string>(&diagram_name), "Filename where to output the persistence diagram");
-#if LOGGING
- std::vector<std::string> log_channels;
- visible.add_options()
- ("log,l", po::value< std::vector<std::string> >(&log_channels), "log channels to turn on (info, debug, etc)");
-#endif
-
- po::positional_options_description pos;
- pos.add("input-file", 1);
-
- po::options_description all; all.add(visible).add(hidden);
-
- po::variables_map vm;
- po::store(po::command_line_parser(argc, argv).
- options(all).positional(pos).run(), vm);
- po::notify(vm);
-
-#if LOGGING
- for (std::vector<std::string>::const_iterator cur = log_channels.begin(); cur != log_channels.end(); ++cur)
- stderrLog.subscribeTo( RLOG_CHANNEL(cur->c_str()) );
-#endif
-
- if (vm.count("help") || !vm.count("input-file"))
- {
- std::cout << "Usage: " << argv[0] << " [options] input-file" << std::endl;
- std::cout << visible << std::endl;
- std::abort();
- }
-}