blob: dfaaffa8ea63e6f8cb5f6f0659b141d43613f6f8 (
plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
|
#include <iostream>
#include <gudhi/Coxeter_triangulation.h>
#include <gudhi/Implicit_manifold_intersection_oracle.h> // for Gudhi::coxeter_triangulation::make_oracle
#include <gudhi/Manifold_tracing.h>
#include <gudhi/Coxeter_triangulation/Cell_complex/Cell_complex.h>
#include <gudhi/Functions/Function_Sm_in_Rd.h>
using namespace Gudhi::coxeter_triangulation;
int main(int argc, char** argv) {
// Oracle is a circle of radius 1
double radius = 1.;
auto oracle = make_oracle(Function_Sm_in_Rd(radius, 1));
// Define a Coxeter triangulation.
Coxeter_triangulation<> cox_tr(oracle.amb_d());
// Theory forbids that a vertex of the triangulation lies exactly on the circle.
// Add some offset to avoid algorithm degeneracies.
cox_tr.change_offset(-Eigen::VectorXd::Random(oracle.amb_d()));
// For a better manifold approximation, one can change the circle radius value or change the linear transformation
// matrix.
// The number of points and edges will increase with a better resolution.
//cox_tr.change_matrix(0.5 * cox_tr.matrix());
// Manifold tracing algorithm
using Out_simplex_map = typename Manifold_tracing<Coxeter_triangulation<> >::Out_simplex_map;
std::vector<Eigen::VectorXd> seed_points(1, oracle.seed());
Out_simplex_map interior_simplex_map;
manifold_tracing_algorithm(seed_points, cox_tr, oracle, interior_simplex_map);
// Constructing the cell complex
std::size_t intr_d = oracle.amb_d() - oracle.cod_d();
Cell_complex<Out_simplex_map> cell_complex(intr_d);
cell_complex.construct_complex(interior_simplex_map);
// List of Hasse_cell pointers to retrieve vertices values from edges
std::map<Cell_complex<Out_simplex_map>::Hasse_cell*, std::size_t> vi_map;
std::size_t index = 0;
std::clog << "Vertices:" << std::endl;
for (const auto& cp_pair : cell_complex.cell_point_map()) {
std::clog << index << " : (" << cp_pair.second(0) << ", " << cp_pair.second(1) << ")" << std::endl;
vi_map.emplace(cp_pair.first, index++);
}
std::clog << "Edges:" << std::endl;
for (const auto& sc_pair : cell_complex.interior_simplex_cell_map(1)) {
Cell_complex<Out_simplex_map>::Hasse_cell* edge_cell = sc_pair.second;
for (const auto& vi_pair : edge_cell->get_boundary()) std::clog << vi_map[vi_pair.first] << " ";
std::clog << std::endl;
}
return 0;
}
|
