summaryrefslogtreecommitdiff
path: root/src
diff options
context:
space:
mode:
authorROUVREAU Vincent <vincent.rouvreau@inria.fr>2019-06-27 10:48:03 +0200
committerROUVREAU Vincent <vincent.rouvreau@inria.fr>2019-06-27 10:48:03 +0200
commit7ba3590fe76357367ee0ea1eaaef95e6343dc222 (patch)
tree9b8e8d44e81a66ac923dd82788dbbe0d9b2541ff /src
parentdf294a362cfb6d57a881650cb0f0914812afc3aa (diff)
parent9f249e001089af45186aaf0d196d6300705eee31 (diff)
Merge master modifications
Diffstat (limited to 'src')
-rw-r--r--src/Persistent_cohomology/doc/Intro_persistent_cohomology.h69
-rw-r--r--src/Persistent_cohomology/example/plain_homology.cpp27
-rw-r--r--src/Persistent_cohomology/include/gudhi/Persistent_cohomology.h18
3 files changed, 59 insertions, 55 deletions
diff --git a/src/Persistent_cohomology/doc/Intro_persistent_cohomology.h b/src/Persistent_cohomology/doc/Intro_persistent_cohomology.h
index a16591ce..6fd19706 100644
--- a/src/Persistent_cohomology/doc/Intro_persistent_cohomology.h
+++ b/src/Persistent_cohomology/doc/Intro_persistent_cohomology.h
@@ -144,6 +144,8 @@ diagram.
3 1 0.104347 inf
3 2 0.138335 inf \endcode
+More details on the <a href="../../ripscomplex/">Rips complex utilities</a> dedicated page.
+
\li <a href="_persistent_cohomology_2rips_multifield_persistence_8cpp-example.html">
Persistent_cohomology/rips_multifield_persistence.cpp</a> computes the Rips complex of a point cloud and outputs its
persistence diagram with a family of field coefficients.
@@ -156,6 +158,8 @@ The file should contain square or lower triangular distance matrix with semicolo
The code do not check if it is dealing with a distance matrix. It is the user responsibility to provide a valid input.
Please refer to data/distance_matrix/lower_triangular_distance_matrix.csv for an example of a file.
+More details on the <a href="../../ripscomplex/">Rips complex utilities</a> dedicated page.
+
\li <a href="_rips_complex_2rips_correlation_matrix_persistence_8cpp-example.html">
Rips_complex/rips_correlation_matrix_persistence.cpp</a>
computes the Rips complex of a correlation matrix and outputs its persistence diagram.
@@ -165,6 +169,8 @@ It is the user responsibility to ensure that this is the case. The input is to b
triangular matrix.
Please refer to data/correlation_matrix/lower_triangular_correlation_matrix.csv for an example of a file.
+More details on the <a href="../../ripscomplex/">Rips complex utilities</a> dedicated page.
+
\li <a href="_alpha_complex_2alpha_complex_3d_persistence_8cpp-example.html">
Alpha_complex/alpha_complex_3d_persistence.cpp</a> computes the persistent homology with
\f$\mathbb{Z}/2\mathbb{Z}\f$ coefficients of the alpha complex on points sampling from an OFF file.
@@ -175,48 +181,33 @@ Alpha_complex/alpha_complex_3d_persistence.cpp</a> computes the persistent homol
2 1 0.0934117 1.00003
2 2 0.56444 1.03938 \endcode
-\li <a href="_alpha_complex_2exact_alpha_complex_3d_persistence_8cpp-example.html">
-Alpha_complex/exact_alpha_complex_3d_persistence.cpp</a> computes the persistent homology with
-\f$\mathbb{Z}/2\mathbb{Z}\f$ coefficients of the alpha complex on points sampling from an OFF file.
-Here, as CGAL computes the exact values, it is slower, but it is necessary when points are on a grid
-for instance.
-\code $> ./exact_alpha_complex_3d_persistence ../../data/points/sphere3D_pts_on_grid.off -p 2 -m 0.1 \endcode
+More details on the <a href="../../alphacomplex/">Alpha complex utilities</a> dedicated page.
+
+CGAL can be forced to compute the exact values, it is slower, but it is necessary when points are on a grid
+for instance (the fast version `--fast` would give incorrect values).
+\code $> ./alpha_complex_3d_persistence ../../data/points/sphere3D_pts_on_grid.off --exact -p 2 -m 0.1 \endcode
\code Simplex_tree dim: 3
2 0 0 inf
2 2 0.0002 0.2028 \endcode
-\li <a href="_alpha_complex_2weighted_alpha_complex_3d_persistence_8cpp-example.html">
-Alpha_complex/weighted_alpha_complex_3d_persistence.cpp</a> computes the persistent homology with
+It can also compute the persistent homology with
\f$\mathbb{Z}/2\mathbb{Z}\f$ coefficients of the weighted alpha complex on points sampling from an OFF file
and a weights file.
-\code $> ./weighted_alpha_complex_3d_persistence ../../data/points/tore3D_300.off
-../../data/points/tore3D_300.weights -p 2 -m 0.45 \endcode
+\code $> ./alpha_complex_3d_persistence ../../data/points/tore3D_300.off
+--weight-file ../../data/points/tore3D_300.weights -p 2 -m 0.45 \endcode
\code Simplex_tree dim: 3
2 0 -1 inf
2 1 -0.931784 0.000103311
2 1 -0.906588 2.60165e-05
2 2 -0.43556 0.0393798 \endcode
-\li <a href="_alpha_complex_2alpha_complex_persistence_8cpp-example.html">
-Alpha_complex/alpha_complex_persistence.cpp</a> computes the persistent homology with
-\f$\mathbb{Z}/p\mathbb{Z}\f$ coefficients of the alpha complex on points sampling from an OFF file.
-\code $> ./alpha_complex_persistence -r 32 -p 2 -m 0.45 ../../data/points/tore3D_300.off \endcode
-\code Alpha complex is of dimension 3 - 9273 simplices - 300 vertices.
-Simplex_tree dim: 3
-2 0 0 inf
-2 1 0.0682162 1.0001
-2 1 0.0934117 1.00003
-2 2 0.56444 1.03938 \endcode
-
-\li <a href="_alpha_complex_2periodic_alpha_complex_3d_persistence_8cpp-example.html">
-Alpha_complex/periodic_alpha_complex_3d_persistence.cpp</a> computes the persistent homology with
+One can also compute the persistent homology with
\f$\mathbb{Z}/2\mathbb{Z}\f$ coefficients of the periodic alpha complex on points sampling from an OFF file.
The second parameter is a \ref FileFormatsIsoCuboid file with coordinates of the periodic cuboid.
Note that the lengths of the sides of the periodic cuboid have to be the same.
-\code $> ./periodic_alpha_complex_3d_persistence ../../data/points/grid_10_10_10_in_0_1.off
-../../data/points/iso_cuboid_3_in_0_1.txt -p 3 -m 1.0 \endcode
-\code Periodic Delaunay computed.
-Simplex_tree dim: 3
+\code $> ./alpha_complex_3d_persistence ../../data/points/grid_10_10_10_in_0_1.off
+--cuboid-file ../../data/points/iso_cuboid_3_in_0_1.txt -p 3 -m 1.0 \endcode
+\code Simplex_tree dim: 3
3 0 0 inf
3 1 0.0025 inf
3 1 0.0025 inf
@@ -226,18 +217,17 @@ Simplex_tree dim: 3
3 2 0.005 inf
3 3 0.0075 inf \endcode
-\li <a href="_persistent_cohomology_2weighted_periodic_alpha_complex_3d_persistence_8cpp-example.html">
-Persistent_cohomology/weighted_periodic_alpha_complex_3d_persistence.cpp</a> computes the persistent homology with
+In order to compute the persistent homology with
\f$\mathbb{Z}/2\mathbb{Z}\f$ coefficients of the periodic alpha complex on weighted points from an OFF file. The
additional parameters of this program are:<br>
(a) The file with the weights of points. The file consist of a sequence of numbers (as many as points).
Note that the weight of each single point have to be bounded by 1/64 times the square of the cuboid edge length.<br>
(b) A \ref FileFormatsIsoCuboid file with coordinates of the periodic cuboid.
Note that the lengths of the sides of the periodic cuboid have to be the same.<br>
-\code $> ./weighted_periodic_alpha_complex_3d_persistence ../../data/points/shifted_sphere.off
-../../data/points/shifted_sphere.weights ../../data/points/iso_cuboid_3_in_0_10.txt 3 1.0 \endcode
-\code Weighted Periodic Delaunay computed.
-Simplex_tree dim: 3
+\code $> ./alpha_complex_3d_persistence ../../data/points/shifted_sphere.off
+--weight-file ../../data/points/shifted_sphere.weights
+--cuboid-file ../../data/points/iso_cuboid_3_in_0_10.txt -p 3 -m 1.0 \endcode
+\code Simplex_tree dim: 3
3 0 -0.0001 inf
3 1 16.0264 inf
3 1 16.0273 inf
@@ -247,6 +237,19 @@ Simplex_tree dim: 3
3 2 36.8838 inf
3 3 58.6783 inf \endcode
+\li <a href="_alpha_complex_2alpha_complex_persistence_8cpp-example.html">
+Alpha_complex/alpha_complex_persistence.cpp</a> computes the persistent homology with
+\f$\mathbb{Z}/p\mathbb{Z}\f$ coefficients of the alpha complex on points sampling from an OFF file.
+\code $> ./alpha_complex_persistence -r 32 -p 2 -m 0.45 ../../data/points/tore3D_300.off \endcode
+\code Alpha complex is of dimension 3 - 9273 simplices - 300 vertices.
+Simplex_tree dim: 3
+2 0 0 inf
+2 1 0.0682162 1.0001
+2 1 0.0934117 1.00003
+2 2 0.56444 1.03938 \endcode
+
+More details on the <a href="../../alphacomplex/">Alpha complex utilities</a> dedicated page.
+
\li <a href="_persistent_cohomology_2plain_homology_8cpp-example.html">
Persistent_cohomology/plain_homology.cpp</a> computes the plain homology of a simple simplicial complex without
filtration values.
diff --git a/src/Persistent_cohomology/example/plain_homology.cpp b/src/Persistent_cohomology/example/plain_homology.cpp
index b8c8b1f9..fbb25cea 100644
--- a/src/Persistent_cohomology/example/plain_homology.cpp
+++ b/src/Persistent_cohomology/example/plain_homology.cpp
@@ -40,26 +40,34 @@ int main() {
ST st;
/* Complex to build. */
- /* 1 3 */
- /* o---o */
- /* /X\ / */
+ /* 1 3 5 */
+ /* o---o---o */
+ /* / \ / */
/* o---o o */
/* 2 0 4 */
- const short triangle012[] = {0, 1, 2};
+ const short edge01[] = {0, 1};
+ const short edge02[] = {0, 2};
+ const short edge12[] = {1, 2};
const short edge03[] = {0, 3};
const short edge13[] = {1, 3};
+ const short edge35[] = {3, 5};
const short vertex4[] = {4};
- st.insert_simplex_and_subfaces(triangle012);
+ st.insert_simplex_and_subfaces(edge01);
+ st.insert_simplex_and_subfaces(edge02);
+ st.insert_simplex_and_subfaces(edge12);
st.insert_simplex_and_subfaces(edge03);
st.insert_simplex(edge13);
+ st.insert_simplex_and_subfaces(edge35);
st.insert_simplex(vertex4);
// Sort the simplices in the order of the filtration
st.initialize_filtration();
// Class for homology computation
- Persistent_cohomology pcoh(st);
+ // By default, since the complex has dimension 1, only 0-dimensional homology would be computed.
+ // Here we also want persistent homology to be computed for the maximal dimension in the complex (persistence_dim_max = true)
+ Persistent_cohomology pcoh(st, true);
// Initialize the coefficient field Z/2Z for homology
pcoh.init_coefficients(2);
@@ -72,13 +80,14 @@ int main() {
// 2 0 0 inf
// 2 0 0 inf
// 2 1 0 inf
- // means that in Z/2Z-homology, the Betti numbers are b0=2 and b1=1.
+ // 2 1 0 inf
+ // means that in Z/2Z-homology, the Betti numbers are b0=2 and b1=2.
pcoh.output_diagram();
- // Print the Betti numbers are b0=2 and b1=1.
+ // Print the Betti numbers are b0=2 and b1=2.
std::cout << std::endl;
std::cout << "The Betti numbers are : ";
- for (int i = 0; i < st.dimension(); i++)
+ for (int i = 0; i < 3; i++)
std::cout << "b" << i << " = " << pcoh.betti_number(i) << " ; ";
std::cout << std::endl;
}
diff --git a/src/Persistent_cohomology/include/gudhi/Persistent_cohomology.h b/src/Persistent_cohomology/include/gudhi/Persistent_cohomology.h
index e62ccbc8..bee54ded 100644
--- a/src/Persistent_cohomology/include/gudhi/Persistent_cohomology.h
+++ b/src/Persistent_cohomology/include/gudhi/Persistent_cohomology.h
@@ -88,9 +88,13 @@ class Persistent_cohomology {
*
* @param[in] cpx Complex for which the persistent homology is computed.
* cpx is a model of FilteredComplex
+ *
+ * @param[in] persistence_dim_max if true, the persistent homology for the maximal dimension in the
+ * complex is computed. If false, it is ignored. Default is false.
+ *
* @exception std::out_of_range In case the number of simplices is more than Simplex_key type numeric limit.
*/
- explicit Persistent_cohomology(FilteredComplex& cpx)
+ explicit Persistent_cohomology(FilteredComplex& cpx, bool persistence_dim_max = false)
: cpx_(&cpx),
dim_max_(cpx.dimension()), // upper bound on the dimension of the simplices
coeff_field_(), // initialize the field coefficient structure.
@@ -116,18 +120,6 @@ class Persistent_cohomology {
++idx_fil;
dsets_.make_set(cpx_->key(sh));
}
- }
-
- /** \brief Initializes the Persistent_cohomology class.
- *
- * @param[in] cpx Complex for which the persistent homology is compiuted.
- * cpx is a model of FilteredComplex
- *
- * @param[in] persistence_dim_max if true, the persistent homology for the maximal dimension in the
- * complex is computed. If false, it is ignored. Default is false.
- */
- Persistent_cohomology(FilteredComplex& cpx, bool persistence_dim_max)
- : Persistent_cohomology(cpx) {
if (persistence_dim_max) {
++dim_max_;
}