From 0dcebad6a0c9d9f7968b582522fbbf548758d9c0 Mon Sep 17 00:00:00 2001 From: mcarrier Date: Mon, 26 Feb 2018 12:53:46 +0000 Subject: git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/Nerve_GIC@3257 636b058d-ea47-450e-bf9e-a15bfbe3eedb Former-commit-id: 25ee894f17d60c05f961faec7bfa41f510d91827 --- src/Nerve_GIC/include/gudhi/GIC.h | 56 ++++++++++++++++++++++++++++++++------- 1 file changed, 47 insertions(+), 9 deletions(-) (limited to 'src/Nerve_GIC/include/gudhi/GIC.h') diff --git a/src/Nerve_GIC/include/gudhi/GIC.h b/src/Nerve_GIC/include/gudhi/GIC.h index 642b88b0..c5fd4d22 100644 --- a/src/Nerve_GIC/include/gudhi/GIC.h +++ b/src/Nerve_GIC/include/gudhi/GIC.h @@ -115,7 +115,7 @@ class Cover_complex { std::vector distribution; std::map > - cover; // function associating to each data point its vectors of cover elements to which it belongs. + cover; // function associating to each data point the vector of cover elements to which it belongs. std::map > cover_back; // inverse of cover, in order to get the data points associated to a specific cover element. std::map cover_std; // standard function (induced by func) used to compute the extended persistence @@ -1044,13 +1044,31 @@ class Cover_complex { minf = std::min(minf, it->second); } + /*int magic[] = {-2}; + st.insert_simplex(magic, -3); + for (auto const& simplex : simplices) { - // Add a simplex and a cone on it - std::vector splx = simplex; - splx.push_back(-2); - st.insert_simplex_and_subfaces(splx); + std::vector splx = simplex; splx.push_back(-2); + st.insert_simplex_and_subfaces(splx, -3); } + for (auto const& simplex : simplices) { + if(simplex.size() == 1){ + st.insert_simplex(it->first, -2 + (it->second - minf)/(maxf - minf)); + int[] cone_edge = {-2,it->first}; + st.insert_simplex(cone_edge, 2 - (it->second - minf)/(maxf - minf)); + } + else{ + st.insert_simplex(simplex, -2.5); + std::vector splx = simplex; splx.push_back(-2); + st.insert_simplex(splx, 0); + } + } + + st.make_filtration_non_decreasing();*/ + + + // Build filtration for (auto simplex : st.complex_simplex_range()) { double filta = std::numeric_limits::lowest(); @@ -1062,7 +1080,7 @@ class Cover_complex { continue; } filta = std::max(-2 + (cover_std[vertex] - minf) / (maxf - minf), filta); - filts = std::max(2 - (cover_std[vertex] - minf) / (maxf - minf), filts); + filts = std::max( 2 - (cover_std[vertex] - minf) / (maxf - minf), filts); } if (ascending) st.assign_filtration(simplex, filta); @@ -1071,11 +1089,10 @@ class Cover_complex { } int magic[] = {-2}; st.assign_filtration(st.find(magic), -3); + st.initialize_filtration(); // Compute PD - st.initialize_filtration(); - Gudhi::persistent_cohomology::Persistent_cohomology pcoh(st); - pcoh.init_coefficients(2); + Gudhi::persistent_cohomology::Persistent_cohomology pcoh(st); pcoh.init_coefficients(2); pcoh.compute_persistent_cohomology(); // Output PD @@ -1201,6 +1218,27 @@ class Cover_complex { } } + private: + std::vector > subfaces(std::vector simplex){ + if (simplex.size() == 1){ + std::vector > dummy; dummy.clear(); + std::vector empty; empty.clear(); + dummy.push_bakc(empty); dummy.push_back(simplex); return dummy; + } + else{ + int popped_vertex = simplex[simplex.size()-1]; + std::vector popped_simplex = simplex; popped_simplex.pop_back(); + std::vector > subf1 = subfaces(popped_simplex); std::vector > subf2; + for (int i = 0; i < subf1.size(); i++){ + std::vector face = subf1[i]; + face.push_back(popped_vertex); + subf2.push_back(face); + } + subf1.insert(subf1.end(), subf2.begin(), subf2.end() ); + return subf1; + } + } + public: /** \brief Computes the simplices of the simplicial complex. */ -- cgit v1.2.3 From 6141868f2336a2650e0f9ac755735aee30d56356 Mon Sep 17 00:00:00 2001 From: mcarrier Date: Tue, 27 Feb 2018 15:32:07 +0000 Subject: added parallelization with TBB + some map changed to vector + construction of PD git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/Nerve_GIC@3258 636b058d-ea47-450e-bf9e-a15bfbe3eedb Former-commit-id: b1990a9f139015217e5a59cc87e1f6b5296be6b3 --- src/Nerve_GIC/include/gudhi/GIC.h | 287 ++++++++++++++++++---------------- src/cmake/modules/GUDHI_modules.cmake | 2 +- 2 files changed, 154 insertions(+), 135 deletions(-) (limited to 'src/Nerve_GIC/include/gudhi/GIC.h') diff --git a/src/Nerve_GIC/include/gudhi/GIC.h b/src/Nerve_GIC/include/gudhi/GIC.h index c5fd4d22..1805959b 100644 --- a/src/Nerve_GIC/include/gudhi/GIC.h +++ b/src/Nerve_GIC/include/gudhi/GIC.h @@ -23,6 +23,12 @@ #ifndef GIC_H_ #define GIC_H_ +#ifdef GUDHI_USE_TBB +#include +#include +#include +#endif + #include #include #include @@ -99,9 +105,8 @@ class Cover_complex { int data_dimension; // dimension of input data. int n; // number of points. - std::map func; // function used to compute the output simplicial complex. - std::map - func_color; // function used to compute the colors of the nodes of the output simplicial complex. + std::vector func; // function used to compute the output simplicial complex. + std::vector func_color; // function used to compute the colors of the nodes of the output simplicial complex. bool functional_cover = false; // whether we use a cover with preimages of a function or not. Graph one_skeleton_OFF; // one-skeleton given by the input OFF file (if it exists). @@ -114,7 +119,7 @@ class Cover_complex { Persistence_diagram PD; std::vector distribution; - std::map > + std::vector > cover; // function associating to each data point the vector of cover elements to which it belongs. std::map > cover_back; // inverse of cover, in order to get the data points associated to a specific cover element. @@ -140,8 +145,8 @@ class Cover_complex { // Point comparator struct Less { - Less(std::map func) { Fct = func; } - std::map Fct; + Less(std::vector func) { Fct = func; } + std::vector Fct; bool operator()(int a, int b) { if (Fct[a] == Fct[b]) return a < b; @@ -276,6 +281,7 @@ class Cover_complex { point_cloud.emplace_back(point.begin(), point.begin() + data_dimension); boost::add_vertex(one_skeleton_OFF); vertices.push_back(boost::add_vertex(one_skeleton)); + std::vector dummy; dummy.clear(); cover.push_back(dummy); i++; } } @@ -431,17 +437,29 @@ class Cover_complex { if (distances.size() == 0) compute_pairwise_distances(distance); - // #pragma omp parallel for - for (int i = 0; i < N; i++) { - SampleWithoutReplacement(n, m, samples); - double hausdorff_dist = 0; - for (int j = 0; j < n; j++) { - double mj = distances[j][samples[0]]; - for (int k = 1; k < m; k++) mj = std::min(mj, distances[j][samples[k]]); - hausdorff_dist = std::max(hausdorff_dist, mj); + #ifdef GUDHI_USE_TBB + tbb::parallel_for(0, N, [&](int i){ + SampleWithoutReplacement(n, m, samples); + double hausdorff_dist = 0; + for (int j = 0; j < n; j++) { + double mj = distances[j][samples[0]]; + for (int k = 1; k < m; k++) mj = std::min(mj, distances[j][samples[k]]); + hausdorff_dist = std::max(hausdorff_dist, mj); + } + delta += hausdorff_dist / N; + }); + #else + for (int i = 0; i < N; i++) { + SampleWithoutReplacement(n, m, samples); + double hausdorff_dist = 0; + for (int j = 0; j < n; j++) { + double mj = distances[j][samples[0]]; + for (int k = 1; k < m; k++) mj = std::min(mj, distances[j][samples[k]]); + hausdorff_dist = std::max(hausdorff_dist, mj); + } + delta += hausdorff_dist / N; } - delta += hausdorff_dist / N; - } + #endif if (verbose) std::cout << "delta = " << delta << std::endl; set_graph_from_rips(delta, distance); @@ -466,7 +484,7 @@ class Cover_complex { while (std::getline(input, line)) { std::stringstream stream(line); stream >> f; - func.emplace(i, f); + func.push_back(f); i++; } functional_cover = true; @@ -480,7 +498,7 @@ class Cover_complex { * */ void set_function_from_coordinate(int k) { - for (int i = 0; i < n; i++) func.emplace(i, point_cloud[i][k]); + for (int i = 0; i < n; i++) func.push_back(point_cloud[i][k]); char coordinate[100]; sprintf(coordinate, "coordinate %d", k); functional_cover = true; @@ -495,7 +513,7 @@ class Cover_complex { */ template void set_function_from_range(InputRange const& function) { - for (int i = 0; i < n; i++) func.emplace(i, function[i]); + for (int i = 0; i < n; i++) func.push_back(function[i]); functional_cover = true; } @@ -713,37 +731,70 @@ class Cover_complex { funcstd[i] = 0.5 * (u + v); } - if (verbose) std::cout << "Computing connected components..." << std::endl; - // #pragma omp parallel for - for (int i = 0; i < res; i++) { - // Compute connected components - Graph G = one_skeleton.create_subgraph(); - int num = preimages[i].size(); - std::vector component(num); - for (int j = 0; j < num; j++) boost::add_vertex(index[vertices[preimages[i][j]]], G); - boost::connected_components(G, &component[0]); - int max = 0; - - // For each point in preimage - for (int j = 0; j < num; j++) { - // Update number of components in preimage - if (component[j] > max) max = component[j]; - - // Identify component with Cantor polynomial N^2 -> N - int identifier = (std::pow(i + component[j], 2) + 3 * i + component[j]) / 2; - - // Update covers - cover[preimages[i][j]].push_back(identifier); - cover_back[identifier].push_back(preimages[i][j]); - cover_fct[identifier] = i; - cover_std[identifier] = funcstd[i]; - cover_color[identifier].second += func_color[preimages[i][j]]; - cover_color[identifier].first += 1; - } + #ifdef GUDHI_USE_TBB + tbb::task_scheduler_init init(4); + if (verbose) std::cout << "Computing connected components (parallelized)..." << std::endl; + tbb::parallel_for(0, res, [&](int i){ + // Compute connected components + Graph G = one_skeleton.create_subgraph(); + int num = preimages[i].size(); + std::vector component(num); + for (int j = 0; j < num; j++) boost::add_vertex(index[vertices[preimages[i][j]]], G); + boost::connected_components(G, &component[0]); + int max = 0; + + // For each point in preimage + for (int j = 0; j < num; j++) { + // Update number of components in preimage + if (component[j] > max) max = component[j]; + + // Identify component with Cantor polynomial N^2 -> N + int identifier = ((i + component[j])*(i + component[j]) + 3 * i + component[j]) / 2; + + // Update covers + cover[preimages[i][j]].push_back(identifier); + cover_back[identifier].push_back(preimages[i][j]); + cover_fct[identifier] = i; + cover_std[identifier] = funcstd[i]; + cover_color[identifier].second += func_color[preimages[i][j]]; + cover_color[identifier].first += 1; + } - // Maximal dimension is total number of connected components - id += max + 1; - } + // Maximal dimension is total number of connected components + id += max + 1; + }); + #else + if (verbose) std::cout << "Computing connected components..." << std::endl; + for (int i = 0; i < res; i++) { + // Compute connected components + Graph G = one_skeleton.create_subgraph(); + int num = preimages[i].size(); + std::vector component(num); + for (int j = 0; j < num; j++) boost::add_vertex(index[vertices[preimages[i][j]]], G); + boost::connected_components(G, &component[0]); + int max = 0; + + // For each point in preimage + for (int j = 0; j < num; j++) { + // Update number of components in preimage + if (component[j] > max) max = component[j]; + + // Identify component with Cantor polynomial N^2 -> N + int identifier = (std::pow(i + component[j], 2) + 3 * i + component[j]) / 2; + + // Update covers + cover[preimages[i][j]].push_back(identifier); + cover_back[identifier].push_back(preimages[i][j]); + cover_fct[identifier] = i; + cover_std[identifier] = funcstd[i]; + cover_color[identifier].second += func_color[preimages[i][j]]; + cover_color[identifier].first += 1; + } + + // Maximal dimension is total number of connected components + id += max + 1; + } + #endif maximal_dim = id - 1; for (std::map >::iterator iit = cover_color.begin(); iit != cover_color.end(); iit++) @@ -805,25 +856,48 @@ class Cover_complex { std::vector mindist(n); for (int j = 0; j < n; j++) mindist[j] = std::numeric_limits::max(); + // Compute the geodesic distances to subsamples with Dijkstra - // #pragma omp parallel for - for (int i = 0; i < m; i++) { - if (verbose) std::cout << "Computing geodesic distances to seed " << i << "..." << std::endl; - int seed = voronoi_subsamples[i]; - std::vector dmap(n); - boost::dijkstra_shortest_paths( - one_skeleton, vertices[seed], - boost::weight_map(weight).distance_map(boost::make_iterator_property_map(dmap.begin(), index))); - - for (int j = 0; j < n; j++) - if (mindist[j] > dmap[j]) { - mindist[j] = dmap[j]; - if (cover[j].size() == 0) - cover[j].push_back(i); - else - cover[j][0] = i; - } - } + #ifdef GUDHI_USE_TBB + if (verbose) std::cout << "Computing geodesic distances (parallelized)..." << std::endl; + tbb::mutex coverMutex; tbb::mutex mindistMutex; + tbb::parallel_for(0, m, [&](int i){ + int seed = voronoi_subsamples[i]; + std::vector dmap(n); + boost::dijkstra_shortest_paths( + one_skeleton, vertices[seed], + boost::weight_map(weight).distance_map(boost::make_iterator_property_map(dmap.begin(), index))); + + coverMutex.lock(); mindistMutex.lock(); + for (int j = 0; j < n; j++) + if (mindist[j] > dmap[j]) { + mindist[j] = dmap[j]; + if (cover[j].size() == 0) + cover[j].push_back(i); + else + cover[j][0] = i; + } + coverMutex.unlock(); mindistMutex.unlock(); + }); + #else + for (int i = 0; i < m; i++) { + if (verbose) std::cout << "Computing geodesic distances to seed " << i << "..." << std::endl; + int seed = voronoi_subsamples[i]; + std::vector dmap(n); + boost::dijkstra_shortest_paths( + one_skeleton, vertices[seed], + boost::weight_map(weight).distance_map(boost::make_iterator_property_map(dmap.begin(), index))); + + for (int j = 0; j < n; j++) + if (mindist[j] > dmap[j]) { + mindist[j] = dmap[j]; + if (cover[j].size() == 0) + cover[j].push_back(i); + else + cover[j][0] = i; + } + } + #endif for (int i = 0; i < n; i++) { cover_back[cover[i][0]].push_back(i); @@ -863,7 +937,7 @@ class Cover_complex { while (std::getline(input, line)) { std::stringstream stream(line); stream >> f; - func_color.emplace(i, f); + func_color.push_back(f); i++; } color_name = color_file_name; @@ -876,7 +950,7 @@ class Cover_complex { * */ void set_color_from_coordinate(int k = 0) { - for (int i = 0; i < n; i++) func_color[i] = point_cloud[i][k]; + for (int i = 0; i < n; i++) func_color.push_back(point_cloud[i][k]); color_name = "coordinate "; color_name.append(std::to_string(k)); } @@ -888,7 +962,7 @@ class Cover_complex { * */ void set_color_from_vector(std::vector color) { - for (unsigned int i = 0; i < color.size(); i++) func_color[i] = color[i]; + for (unsigned int i = 0; i < color.size(); i++) func_color.push_back(color[i]); } public: // Create a .dot file that can be compiled with neato to produce a .pdf file. @@ -1044,52 +1118,19 @@ class Cover_complex { minf = std::min(minf, it->second); } - /*int magic[] = {-2}; - st.insert_simplex(magic, -3); - + // Build filtration for (auto const& simplex : simplices) { std::vector splx = simplex; splx.push_back(-2); st.insert_simplex_and_subfaces(splx, -3); } - for (auto const& simplex : simplices) { - if(simplex.size() == 1){ - st.insert_simplex(it->first, -2 + (it->second - minf)/(maxf - minf)); - int[] cone_edge = {-2,it->first}; - st.insert_simplex(cone_edge, 2 - (it->second - minf)/(maxf - minf)); - } - else{ - st.insert_simplex(simplex, -2.5); - std::vector splx = simplex; splx.push_back(-2); - st.insert_simplex(splx, 0); - } - } - - st.make_filtration_non_decreasing();*/ - - - - // Build filtration - for (auto simplex : st.complex_simplex_range()) { - double filta = std::numeric_limits::lowest(); - double filts = filta; - bool ascending = true; - for (auto vertex : st.simplex_vertex_range(simplex)) { - if (vertex == -2) { - ascending = false; - continue; - } - filta = std::max(-2 + (cover_std[vertex] - minf) / (maxf - minf), filta); - filts = std::max( 2 - (cover_std[vertex] - minf) / (maxf - minf), filts); - } - if (ascending) - st.assign_filtration(simplex, filta); - else - st.assign_filtration(simplex, filts); + for (std::map::iterator it = cover_std.begin(); it != cover_std.end(); it++) { + int vertex = it->first; float val = it->second; + int vert[] = {vertex}; int edge[] = {vertex, -2}; + st.assign_filtration(st.find(vert), -2 + (val - minf)/(maxf - minf)); + st.assign_filtration(st.find(edge), 2 - (val - minf)/(maxf - minf)); } - int magic[] = {-2}; - st.assign_filtration(st.find(magic), -3); - st.initialize_filtration(); + st.make_filtration_non_decreasing(); // Compute PD Gudhi::persistent_cohomology::Persistent_cohomology pcoh(st); pcoh.init_coefficients(2); @@ -1099,7 +1140,7 @@ class Cover_complex { int max_dim = st.dimension(); for (int i = 0; i < max_dim; i++) { std::vector > bars = pcoh.intervals_in_dimension(i); - int num_bars = bars.size(); + int num_bars = bars.size(); if(i == 0) num_bars -= 1; if(verbose) std::cout << num_bars << " interval(s) in dimension " << i << ":" << std::endl; for (int j = 0; j < num_bars; j++) { double birth = bars[j].first; @@ -1218,27 +1259,6 @@ class Cover_complex { } } - private: - std::vector > subfaces(std::vector simplex){ - if (simplex.size() == 1){ - std::vector > dummy; dummy.clear(); - std::vector empty; empty.clear(); - dummy.push_bakc(empty); dummy.push_back(simplex); return dummy; - } - else{ - int popped_vertex = simplex[simplex.size()-1]; - std::vector popped_simplex = simplex; popped_simplex.pop_back(); - std::vector > subf1 = subfaces(popped_simplex); std::vector > subf2; - for (int i = 0; i < subf1.size(); i++){ - std::vector face = subf1[i]; - face.push_back(popped_vertex); - subf2.push_back(face); - } - subf1.insert(subf1.end(), subf2.begin(), subf2.end() ); - return subf1; - } - } - public: /** \brief Computes the simplices of the simplicial complex. */ @@ -1249,8 +1269,7 @@ class Cover_complex { } if (type == "Nerve") { - for(auto& simplex : cover) - simplices.push_back(simplex.second); + for(int i = 0; i < n; i++) simplices.push_back(cover[i]); std::sort(simplices.begin(), simplices.end()); std::vector >::iterator it = std::unique(simplices.begin(), simplices.end()); simplices.resize(std::distance(simplices.begin(), it)); diff --git a/src/cmake/modules/GUDHI_modules.cmake b/src/cmake/modules/GUDHI_modules.cmake index f95d0c34..276fb2cc 100644 --- a/src/cmake/modules/GUDHI_modules.cmake +++ b/src/cmake/modules/GUDHI_modules.cmake @@ -17,7 +17,7 @@ function(add_gudhi_module file_path) endfunction(add_gudhi_module) option(WITH_GUDHI_BENCHMARK "Activate/desactivate benchmark compilation" OFF) -option(WITH_GUDHI_EXAMPLE "Activate/desactivate examples compilation and installation" OFF) +option(WITH_GUDHI_EXAMPLE "Activate/desactivate examples compilation and installation" ON) option(WITH_GUDHI_PYTHON "Activate/desactivate python module compilation and installation" ON) option(WITH_GUDHI_TEST "Activate/desactivate examples compilation and installation" ON) option(WITH_GUDHI_UTILITIES "Activate/desactivate utilities compilation and installation" ON) -- cgit v1.2.3 From 2e60e1e4737729a1f772a9cbdb83583fbb166f71 Mon Sep 17 00:00:00 2001 From: mcarrier Date: Wed, 28 Mar 2018 09:41:38 +0000 Subject: git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/Nerve_GIC@3309 636b058d-ea47-450e-bf9e-a15bfbe3eedb Former-commit-id: 1af850057913b691b237bde016ef11de9d65b91b --- src/Nerve_GIC/include/gudhi/GIC.h | 1 - src/Nerve_GIC/test/test_GIC.cpp | 3 +++ 2 files changed, 3 insertions(+), 1 deletion(-) (limited to 'src/Nerve_GIC/include/gudhi/GIC.h') diff --git a/src/Nerve_GIC/include/gudhi/GIC.h b/src/Nerve_GIC/include/gudhi/GIC.h index 1805959b..80d4c024 100644 --- a/src/Nerve_GIC/include/gudhi/GIC.h +++ b/src/Nerve_GIC/include/gudhi/GIC.h @@ -732,7 +732,6 @@ class Cover_complex { } #ifdef GUDHI_USE_TBB - tbb::task_scheduler_init init(4); if (verbose) std::cout << "Computing connected components (parallelized)..." << std::endl; tbb::parallel_for(0, res, [&](int i){ // Compute connected components diff --git a/src/Nerve_GIC/test/test_GIC.cpp b/src/Nerve_GIC/test/test_GIC.cpp index d633753c..e3067d35 100644 --- a/src/Nerve_GIC/test/test_GIC.cpp +++ b/src/Nerve_GIC/test/test_GIC.cpp @@ -39,6 +39,7 @@ BOOST_AUTO_TEST_CASE(check_nerve) { N.set_type("Nerve"); std::string cloud_file_name("data/cloud"); N.read_point_cloud(cloud_file_name); + N.set_color_from_coordinate(); std::string graph_file_name("data/graph"); N.set_graph_from_file(graph_file_name); std::string cover_file_name("data/cover"); @@ -58,6 +59,7 @@ BOOST_AUTO_TEST_CASE(check_GIC) { GIC.set_type("GIC"); std::string cloud_file_name("data/cloud"); GIC.read_point_cloud(cloud_file_name); + GIC.set_color_from_coordinate(); std::string graph_file_name("data/graph"); GIC.set_graph_from_file(graph_file_name); std::string cover_file_name("data/cover"); @@ -77,6 +79,7 @@ BOOST_AUTO_TEST_CASE(check_voronoiGIC) { GIC.set_type("GIC"); std::string cloud_file_name("data/cloud"); GIC.read_point_cloud(cloud_file_name); + GIC.set_color_from_coordinate(); std::string graph_file_name("data/graph"); GIC.set_graph_from_file(graph_file_name); GIC.set_cover_from_Voronoi(Gudhi::Euclidean_distance(), 2); -- cgit v1.2.3 From eb1a95034622365303f30a6987c15f9fe80bd8ab Mon Sep 17 00:00:00 2001 From: mcarrier Date: Thu, 29 Mar 2018 15:38:20 +0000 Subject: git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/Nerve_GIC@3320 636b058d-ea47-450e-bf9e-a15bfbe3eedb Former-commit-id: 21f1574da074759d38553d2a7c87f6c979f77934 --- src/Nerve_GIC/include/gudhi/GIC.h | 20 ++++---------------- 1 file changed, 4 insertions(+), 16 deletions(-) (limited to 'src/Nerve_GIC/include/gudhi/GIC.h') diff --git a/src/Nerve_GIC/include/gudhi/GIC.h b/src/Nerve_GIC/include/gudhi/GIC.h index 80d4c024..bd15225f 100644 --- a/src/Nerve_GIC/include/gudhi/GIC.h +++ b/src/Nerve_GIC/include/gudhi/GIC.h @@ -143,18 +143,6 @@ class Cover_complex { std::string point_cloud_name; std::string color_name; - // Point comparator - struct Less { - Less(std::vector func) { Fct = func; } - std::vector Fct; - bool operator()(int a, int b) { - if (Fct[a] == Fct[b]) - return a < b; - else - return Fct[a] < Fct[b]; - } - }; - // Remove all edges of a graph. void remove_edges(Graph& G) { boost::graph_traits::edge_iterator ei, ei_end; @@ -280,8 +268,7 @@ class Cover_complex { point.assign(std::istream_iterator(iss), std::istream_iterator()); point_cloud.emplace_back(point.begin(), point.begin() + data_dimension); boost::add_vertex(one_skeleton_OFF); - vertices.push_back(boost::add_vertex(one_skeleton)); - std::vector dummy; dummy.clear(); cover.push_back(dummy); + vertices.push_back(boost::add_vertex(one_skeleton)); cover.emplace_back(); i++; } } @@ -675,7 +662,7 @@ class Cover_complex { // Sort points according to function values std::vector points(n); for (int i = 0; i < n; i++) points[i] = i; - std::sort(points.begin(), points.end(), Less(this->func)); + std::sort(points.begin(), points.end(), [=](const int & p1, const int & p2){return (this->func[p1] < this->func[p2]);}); int id = 0; int pos = 0; @@ -1179,8 +1166,9 @@ class Cover_complex { int id = std::floor(u * (this->n)); boot[j] = id; Cboot.point_cloud[j] = this->point_cloud[id]; - Cboot.func.emplace(j, this->func[id]); + Cboot.func.push_back(this->func[id]); } + Cboot.set_color_from_vector(Cboot.func); for (int j = 0; j < n; j++) { std::vector dist(n); for (int k = 0; k < n; k++) dist[k] = distances[boot[j]][boot[k]]; -- cgit v1.2.3 From 78335c71e46bd3b77d1595edef63cedbe6cf006c Mon Sep 17 00:00:00 2001 From: mcarrier Date: Fri, 30 Mar 2018 08:37:43 +0000 Subject: changed functions for bootstrap git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/Nerve_GIC@3322 636b058d-ea47-450e-bf9e-a15bfbe3eedb Former-commit-id: b26ac991e90767fb0de3d91d01ffe1ef25c572fe --- src/Nerve_GIC/example/CoordGIC.cpp | 3 +++ src/Nerve_GIC/include/gudhi/GIC.h | 45 +++++++++++++++++++++----------------- 2 files changed, 28 insertions(+), 20 deletions(-) (limited to 'src/Nerve_GIC/include/gudhi/GIC.h') diff --git a/src/Nerve_GIC/example/CoordGIC.cpp b/src/Nerve_GIC/example/CoordGIC.cpp index c03fcbb3..7e595382 100644 --- a/src/Nerve_GIC/example/CoordGIC.cpp +++ b/src/Nerve_GIC/example/CoordGIC.cpp @@ -66,6 +66,9 @@ int main(int argc, char **argv) { GIC.find_simplices(); + GIC.compute_distribution(10); + GIC.compute_p_value(); + GIC.plot_DOT(); Gudhi::Simplex_tree<> stree; diff --git a/src/Nerve_GIC/include/gudhi/GIC.h b/src/Nerve_GIC/include/gudhi/GIC.h index bd15225f..e6c508fc 100644 --- a/src/Nerve_GIC/include/gudhi/GIC.h +++ b/src/Nerve_GIC/include/gudhi/GIC.h @@ -415,9 +415,7 @@ class Cover_complex { template double set_graph_from_automatic_rips(Distance distance, int N = 100) { int m = floor(n / std::exp((1 + rate_power) * std::log(std::log(n) / std::log(rate_constant)))); - m = std::min(m, n - 1); - std::vector samples(m); - double delta = 0; + m = std::min(m, n - 1); std::vector samples(m); double delta = 0; if (verbose) std::cout << n << " points in R^" << data_dimension << std::endl; if (verbose) std::cout << "Subsampling " << m << " points" << std::endl; @@ -1152,23 +1150,25 @@ class Cover_complex { * @param[in] N number of bootstrap iterations. * */ - template - void compute_distribution(int N = 100) { - if (distribution.size() >= N) { + void compute_distribution(unsigned int N = 100) { + unsigned int sz = distribution.size(); + if (sz >= N) { std::cout << "Already done!" << std::endl; } else { - for (int i = 0; i < N - distribution.size(); i++) { - Cover_complex Cboot; - Cboot.n = this->n; + for (unsigned int i = 0; i < N - sz; i++) { + if (verbose) std::cout << "Computing " << i << "th bootstrap, bottleneck distance = "; + + Cover_complex Cboot; Cboot.n = this->n; Cboot.data_dimension = this->data_dimension; Cboot.type = this->type; Cboot.functional_cover = true; + std::vector boot(this->n); for (int j = 0; j < this->n; j++) { double u = GetUniform(); - int id = std::floor(u * (this->n)); - boot[j] = id; - Cboot.point_cloud[j] = this->point_cloud[id]; - Cboot.func.push_back(this->func[id]); + int id = std::floor(u * (this->n)); boot[j] = id; + Cboot.point_cloud.push_back(this->point_cloud[id]); Cboot.cover.emplace_back(); Cboot.func.push_back(this->func[id]); + boost::add_vertex(Cboot.one_skeleton_OFF); Cboot.vertices.push_back(boost::add_vertex(Cboot.one_skeleton)); } Cboot.set_color_from_vector(Cboot.func); + for (int j = 0; j < n; j++) { std::vector dist(n); for (int k = 0; k < n; k++) dist[k] = distances[boot[j]][boot[k]]; @@ -1181,8 +1181,9 @@ class Cover_complex { Cboot.set_cover_from_function(); Cboot.find_simplices(); Cboot.compute_PD(); - - distribution.push_back(Gudhi::persistence_diagram::bottleneck_distance(this->PD, Cboot.PD)); + double db = Gudhi::persistence_diagram::bottleneck_distance(this->PD, Cboot.PD); + if (verbose) std::cout << db << std::endl; + distribution.push_back(db); } std::sort(distribution.begin(), distribution.end()); @@ -1196,7 +1197,7 @@ class Cover_complex { * */ double compute_distance_from_confidence_level(double alpha) { - int N = distribution.size(); + unsigned int N = distribution.size(); return distribution[std::floor(alpha * N)]; } @@ -1207,9 +1208,11 @@ class Cover_complex { * */ double compute_confidence_level_from_distance(double d) { - int N = distribution.size(); - for (int i = 0; i < N; i++) - if (distribution[i] > d) return i * 1.0 / N; + unsigned int N = distribution.size(); + double level = 1; + for (unsigned int i = 0; i < N; i++) + if (distribution[i] > d){ level = i * 1.0 / N; break; } + return level; } public: @@ -1221,7 +1224,9 @@ class Cover_complex { double distancemin = -std::numeric_limits::lowest(); int N = PD.size(); for (int i = 0; i < N; i++) distancemin = std::min(distancemin, 0.5 * (PD[i].second - PD[i].first)); - return 1 - compute_confidence_level_from_distance(distancemin); + double p_value = 1 - compute_confidence_level_from_distance(distancemin); + if (verbose) std::cout << "p value = " << p_value << std::endl; + return p_value; } // ******************************************************************************************************************* -- cgit v1.2.3 From 2528ddff5d820020374ece89228006409c224e78 Mon Sep 17 00:00:00 2001 From: mcarrier Date: Mon, 28 May 2018 02:39:20 +0000 Subject: added mutex git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/Nerve_GIC@3468 636b058d-ea47-450e-bf9e-a15bfbe3eedb Former-commit-id: 0527349c20f3e44c4f3db196df13ccdd4c265bbb --- src/Nerve_GIC/include/gudhi/GIC.h | 16 ++++++++++++---- 1 file changed, 12 insertions(+), 4 deletions(-) (limited to 'src/Nerve_GIC/include/gudhi/GIC.h') diff --git a/src/Nerve_GIC/include/gudhi/GIC.h b/src/Nerve_GIC/include/gudhi/GIC.h index 8834858c..4bd2c849 100644 --- a/src/Nerve_GIC/include/gudhi/GIC.h +++ b/src/Nerve_GIC/include/gudhi/GIC.h @@ -416,7 +416,7 @@ class Cover_complex { template double set_graph_from_automatic_rips(Distance distance, int N = 100) { int m = floor(n / std::exp((1 + rate_power) * std::log(std::log(n) / std::log(rate_constant)))); - m = std::min(m, n - 1); std::vector samples(m); double delta = 0; + m = std::min(m, n - 1); double delta = 0; if (verbose) std::cout << n << " points in R^" << data_dimension << std::endl; if (verbose) std::cout << "Subsampling " << m << " points" << std::endl; @@ -424,7 +424,9 @@ class Cover_complex { if (distances.size() == 0) compute_pairwise_distances(distance); #ifdef GUDHI_USE_TBB + tbb::mutex deltamutex; tbb::parallel_for(0, N, [&](int i){ + std::vector samples(m); SampleWithoutReplacement(n, m, samples); double hausdorff_dist = 0; for (int j = 0; j < n; j++) { @@ -432,10 +434,13 @@ class Cover_complex { for (int k = 1; k < m; k++) mj = std::min(mj, distances[j][samples[k]]); hausdorff_dist = std::max(hausdorff_dist, mj); } + deltamutex.lock(); delta += hausdorff_dist / N; + deltamutex.unlock(); }); #else for (int i = 0; i < N; i++) { + std::vector samples(m); SampleWithoutReplacement(n, m, samples); double hausdorff_dist = 0; for (int j = 0; j < n; j++) { @@ -718,12 +723,11 @@ class Cover_complex { } #ifdef GUDHI_USE_TBB - if (verbose) std::cout << "Computing connected components (parallelized)..." << std::endl; + if (verbose) std::cout << "Computing connected components (parallelized)..." << std::endl; tbb::mutex covermutex, idmutex; tbb::parallel_for(0, res, [&](int i){ // Compute connected components Graph G = one_skeleton.create_subgraph(); - int num = preimages[i].size(); - std::vector component(num); + int num = preimages[i].size(); std::vector component(num); for (int j = 0; j < num; j++) boost::add_vertex(index[vertices[preimages[i][j]]], G); boost::connected_components(G, &component[0]); int max = 0; @@ -737,16 +741,20 @@ class Cover_complex { int identifier = ((i + component[j])*(i + component[j]) + 3 * i + component[j]) / 2; // Update covers + covermutex.lock(); cover[preimages[i][j]].push_back(identifier); cover_back[identifier].push_back(preimages[i][j]); cover_fct[identifier] = i; cover_std[identifier] = funcstd[i]; cover_color[identifier].second += func_color[preimages[i][j]]; cover_color[identifier].first += 1; + covermutex.unlock(); } // Maximal dimension is total number of connected components + idmutex.lock(); id += max + 1; + idmutex.unlock(); }); #else if (verbose) std::cout << "Computing connected components..." << std::endl; -- cgit v1.2.3 From 68b93015eaf44d45a3a85747b4f3c53bb755e8af Mon Sep 17 00:00:00 2001 From: mcarrier Date: Mon, 11 Jun 2018 02:50:07 +0000 Subject: small change in bootstrap code git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/Nerve_GIC@3577 636b058d-ea47-450e-bf9e-a15bfbe3eedb Former-commit-id: da8a421316d545b400f82eca3f426912b4767a8d --- src/Nerve_GIC/include/gudhi/GIC.h | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) (limited to 'src/Nerve_GIC/include/gudhi/GIC.h') diff --git a/src/Nerve_GIC/include/gudhi/GIC.h b/src/Nerve_GIC/include/gudhi/GIC.h index 4bd2c849..2a50acd7 100644 --- a/src/Nerve_GIC/include/gudhi/GIC.h +++ b/src/Nerve_GIC/include/gudhi/GIC.h @@ -1184,8 +1184,8 @@ class Cover_complex { } Cboot.set_graph_from_automatic_rips(Gudhi::Euclidean_distance()); - Cboot.set_automatic_resolution(); Cboot.set_gain(); + Cboot.set_automatic_resolution(); Cboot.set_cover_from_function(); Cboot.find_simplices(); Cboot.compute_PD(); @@ -1206,7 +1206,9 @@ class Cover_complex { */ double compute_distance_from_confidence_level(double alpha) { unsigned int N = distribution.size(); - return distribution[std::floor(alpha * N)]; + double d = distribution[std::floor(alpha * N)]; + if (verbose) std::cout << "Distance corresponding to confidence " << alpha << " is " << d << std::endl; + return d; } public: @@ -1220,6 +1222,7 @@ class Cover_complex { double level = 1; for (unsigned int i = 0; i < N; i++) if (distribution[i] > d){ level = i * 1.0 / N; break; } + if (verbose) std::cout << "Confidence level of distance " << d << " is " << level << std::endl; return level; } @@ -1231,7 +1234,7 @@ class Cover_complex { double compute_p_value() { double distancemin = -std::numeric_limits::lowest(); int N = PD.size(); - for (int i = 0; i < N; i++) distancemin = std::min(distancemin, 0.5 * (PD[i].second - PD[i].first)); + for (int i = 0; i < N; i++) distancemin = std::min(distancemin, 0.5 * std::abs(PD[i].second - PD[i].first)); double p_value = 1 - compute_confidence_level_from_distance(distancemin); if (verbose) std::cout << "p value = " << p_value << std::endl; return p_value; -- cgit v1.2.3 From 5f8c3b95f41f623be2f999ce38640c39f4c448da Mon Sep 17 00:00:00 2001 From: mcarrier Date: Mon, 11 Jun 2018 08:53:27 +0000 Subject: minor change git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/Nerve_GIC@3578 636b058d-ea47-450e-bf9e-a15bfbe3eedb Former-commit-id: 262dc4f45eba9f5ada23b078d52be54b4eb60a41 --- src/Nerve_GIC/include/gudhi/GIC.h | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) (limited to 'src/Nerve_GIC/include/gudhi/GIC.h') diff --git a/src/Nerve_GIC/include/gudhi/GIC.h b/src/Nerve_GIC/include/gudhi/GIC.h index 2a50acd7..6b8df6dc 100644 --- a/src/Nerve_GIC/include/gudhi/GIC.h +++ b/src/Nerve_GIC/include/gudhi/GIC.h @@ -1232,8 +1232,7 @@ class Cover_complex { * */ double compute_p_value() { - double distancemin = -std::numeric_limits::lowest(); - int N = PD.size(); + double distancemin = 0; int N = PD.size(); for (int i = 0; i < N; i++) distancemin = std::min(distancemin, 0.5 * std::abs(PD[i].second - PD[i].first)); double p_value = 1 - compute_confidence_level_from_distance(distancemin); if (verbose) std::cout << "p value = " << p_value << std::endl; -- cgit v1.2.3 From d34354883a5021ab9beaa03f7739378caf3e5683 Mon Sep 17 00:00:00 2001 From: mcarrier Date: Mon, 11 Jun 2018 14:43:39 +0000 Subject: git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/Nerve_GIC@3579 636b058d-ea47-450e-bf9e-a15bfbe3eedb Former-commit-id: 4fe9f0185421f9c2be697f0db323c042f7c8e7ca --- src/Nerve_GIC/include/gudhi/GIC.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) (limited to 'src/Nerve_GIC/include/gudhi/GIC.h') diff --git a/src/Nerve_GIC/include/gudhi/GIC.h b/src/Nerve_GIC/include/gudhi/GIC.h index 6b8df6dc..8cd7bdbf 100644 --- a/src/Nerve_GIC/include/gudhi/GIC.h +++ b/src/Nerve_GIC/include/gudhi/GIC.h @@ -1232,7 +1232,7 @@ class Cover_complex { * */ double compute_p_value() { - double distancemin = 0; int N = PD.size(); + double distancemin = std::numeric_limits::max(); int N = PD.size(); for (int i = 0; i < N; i++) distancemin = std::min(distancemin, 0.5 * std::abs(PD[i].second - PD[i].first)); double p_value = 1 - compute_confidence_level_from_distance(distancemin); if (verbose) std::cout << "p value = " << p_value << std::endl; -- cgit v1.2.3 From 24683e8d0134e7ede7033f05b83c44dd72ff7e3a Mon Sep 17 00:00:00 2001 From: mcarrier Date: Tue, 18 Sep 2018 22:01:20 +0000 Subject: git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/Nerve_GIC@3897 636b058d-ea47-450e-bf9e-a15bfbe3eedb Former-commit-id: f5fe42a6f42443ca7217c69f8a18780cdb85669e --- src/Nerve_GIC/include/gudhi/GIC.h | 18 ++++++++++++++++++ src/cython/cython/nerve_gic.pyx | 9 +++++++++ 2 files changed, 27 insertions(+) (limited to 'src/Nerve_GIC/include/gudhi/GIC.h') diff --git a/src/Nerve_GIC/include/gudhi/GIC.h b/src/Nerve_GIC/include/gudhi/GIC.h index 2c37dfae..242ecf7d 100644 --- a/src/Nerve_GIC/include/gudhi/GIC.h +++ b/src/Nerve_GIC/include/gudhi/GIC.h @@ -389,6 +389,24 @@ class Cover_complex { distances[index[boost::source(*ei, one_skeleton)]][index[boost::target(*ei, one_skeleton)]]); } + public: + /** \brief Reads and stores the distance matrices from vector stored in memory. + * + * @param[in] distance_matrix input vector representing the distance matrix. + * + */ + template + void set_distances_from_range(InputRange const & distance_matrix) { + n = distance_matrix.size(); data_dimension = 0; point_cloud_name = "matrix"; + for(int i = 0; i < n; i++){ + point_cloud.emplace_back(); + boost::add_vertex(one_skeleton_OFF); + vertices.push_back(boost::add_vertex(one_skeleton)); + cover.emplace_back(); + } + distances = distance_matrix; + } + public: // Pairwise distances. /** \private \brief Computes all pairwise distances. */ diff --git a/src/cython/cython/nerve_gic.pyx b/src/cython/cython/nerve_gic.pyx index 01dd0a4b..5f01b379 100644 --- a/src/cython/cython/nerve_gic.pyx +++ b/src/cython/cython/nerve_gic.pyx @@ -68,6 +68,7 @@ cdef extern from "Nerve_gic_interface.h" namespace "Gudhi": void plot_DOT() void plot_OFF() void set_point_cloud_from_range(vector[vector[double]] cloud) + void set_distances_from_range(vector[vector[double]] distance_matrix) # CoverComplex python interface cdef class CoverComplex: @@ -111,6 +112,14 @@ cdef class CoverComplex: """ return self.thisptr.set_point_cloud_from_range(cloud) + def set_distances_from_range(self, distance_matrix): + """ Reads and stores the input distance matrix from a vector stored in memory. + + :param distance_matrix: Input vector containing the distance matrix. + :type distance_matrix: vector[vector[double]] + """ + return self.thisptr.set_distances_from_range(distance_matrix) + def compute_confidence_level_from_distance(self, distance): """Computes the confidence level of a specific bottleneck distance threshold. -- cgit v1.2.3 From cae907acde52894ea1a20f3a0c0c22538c1b4128 Mon Sep 17 00:00:00 2001 From: mcarrier Date: Thu, 20 Sep 2018 21:44:01 +0000 Subject: fixed small bug on persistence diagrams git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/Nerve_GIC@3901 636b058d-ea47-450e-bf9e-a15bfbe3eedb Former-commit-id: cbeecd2f928b4c4df13acf5a7a2083882fdc98d2 --- src/Nerve_GIC/include/gudhi/GIC.h | 55 +++++++++++++++++++++++++++------------ 1 file changed, 39 insertions(+), 16 deletions(-) (limited to 'src/Nerve_GIC/include/gudhi/GIC.h') diff --git a/src/Nerve_GIC/include/gudhi/GIC.h b/src/Nerve_GIC/include/gudhi/GIC.h index 242ecf7d..30f89d65 100644 --- a/src/Nerve_GIC/include/gudhi/GIC.h +++ b/src/Nerve_GIC/include/gudhi/GIC.h @@ -397,12 +397,16 @@ class Cover_complex { */ template void set_distances_from_range(InputRange const & distance_matrix) { - n = distance_matrix.size(); data_dimension = 0; point_cloud_name = "matrix"; - for(int i = 0; i < n; i++){ - point_cloud.emplace_back(); - boost::add_vertex(one_skeleton_OFF); - vertices.push_back(boost::add_vertex(one_skeleton)); - cover.emplace_back(); + if(point_cloud.size() == 0){ + n = distance_matrix.size(); + point_cloud_name = "matrix"; + data_dimension = 0; + for(int i = 0; i < n; i++){ + point_cloud.emplace_back(); + boost::add_vertex(one_skeleton_OFF); + vertices.push_back(boost::add_vertex(one_skeleton)); + cover.emplace_back(); + } } distances = distance_matrix; } @@ -536,9 +540,17 @@ class Cover_complex { * */ void set_function_from_coordinate(int k) { - for (int i = 0; i < n; i++) func.push_back(point_cloud[i][k]); - functional_cover = true; - cover_name = "coordinate " + std::to_string(k); + if(point_cloud[0].size() > 0){ + for (int i = 0; i < n; i++) func.push_back(point_cloud[i][k]); + functional_cover = true; + cover_name = "coordinate " + std::to_string(k); + } + else{ + std::cout << "Only pairwise distances provided---cannot access " << k << "th coordinate; returning null vector instead" << std::endl; + for (int i = 0; i < n; i++) func.push_back(0.0); + functional_cover = true; + cover_name = "null"; + } } public: // Set function from vector. @@ -989,9 +1001,17 @@ class Cover_complex { * */ void set_color_from_coordinate(int k = 0) { - for (int i = 0; i < n; i++) func_color.push_back(point_cloud[i][k]); - color_name = "coordinate "; - color_name.append(std::to_string(k)); + if(point_cloud[0].size() > 0){ + for (int i = 0; i < n; i++) func_color.push_back(point_cloud[i][k]); + color_name = "coordinate "; + color_name.append(std::to_string(k)); + } + else{ + std::cout << "Only pairwise distances provided---cannot access " << k << "th coordinate; returning null vector instead" << std::endl; + for (int i = 0; i < n; i++) func.push_back(0.0); + functional_cover = true; + cover_name = "null"; + } } public: // Set color from vector. @@ -1157,15 +1177,18 @@ class Cover_complex { // Build filtration for (auto const& simplex : simplices) { - std::vector splx = simplex; splx.push_back(-2); + std::vector splx = simplex; + splx.push_back(-2); st.insert_simplex_and_subfaces(splx, -3); } for (std::map::iterator it = cover_std.begin(); it != cover_std.end(); it++) { int vertex = it->first; float val = it->second; int vert[] = {vertex}; int edge[] = {vertex, -2}; - st.assign_filtration(st.find(vert), -2 + (val - minf)/(maxf - minf)); - st.assign_filtration(st.find(edge), 2 - (val - minf)/(maxf - minf)); + if(st.find(vert) != st.null_simplex()){ + st.assign_filtration(st.find(vert), -2 + (val - minf)/(maxf - minf)); + st.assign_filtration(st.find(edge), 2 - (val - minf)/(maxf - minf)); + } } st.make_filtration_non_decreasing(); @@ -1267,7 +1290,7 @@ class Cover_complex { unsigned int N = distribution.size(); double level = 1; for (unsigned int i = 0; i < N; i++) - if (distribution[i] > d){ level = i * 1.0 / N; break; } + if (distribution[i] >= d){ level = i * 1.0 / N; break; } if (verbose) std::cout << "Confidence level of distance " << d << " is " << level << std::endl; return level; } -- cgit v1.2.3