Loading...
Searching...
No Matches
Coxeter triangulation

Classes

class  Gudhi::coxeter_triangulation::Vertex_iterator< Permutahedral_representation >
 Iterator over the vertices of a simplex represented by its permutahedral representation. More...
 
class  Gudhi::coxeter_triangulation::Face_iterator< Permutahedral_representation >
 Iterator over the k-faces of a simplex given by its permutahedral representation. More...
 
class  Gudhi::coxeter_triangulation::Coface_iterator< Permutahedral_representation >
 Iterator over the k-cofaces of a simplex given by its permutahedral representation. More...
 
class  Gudhi::coxeter_triangulation::Coxeter_triangulation< Permutahedral_representation_ >
 A class that stores Coxeter triangulation of type \(\tilde{A}_d\). This triangulation has the greatest simplex quality out of all linear transformations of the Freudenthal-Kuhn triangulation. More...
 
class  Gudhi::coxeter_triangulation::Query_result< Simplex_handle >
 The result of a query by an oracle such as Implicit_manifold_intersection_oracle. More...
 
class  Gudhi::coxeter_triangulation::Freudenthal_triangulation< Permutahedral_representation_ >
 A class that stores any affine transformation of the Freudenthal-Kuhn triangulation. More...
 
class  Gudhi::coxeter_triangulation::Implicit_manifold_intersection_oracle< Function_, Domain_function_ >
 An oracle that supports the intersection query on an implicit manifold. More...
 
class  Gudhi::coxeter_triangulation::Mesh_medit
 Structure to store a mesh that can be output in Medit .mesh file format using the output_meshes_to_medit method. More...
 
class  Gudhi::coxeter_triangulation::Manifold_tracing< Triangulation_ >
 A class that assembles methods for manifold tracing algorithm. More...
 
class  Gudhi::coxeter_triangulation::Permutahedral_representation< Vertex_, Ordered_set_partition_ >
 A class that stores the permutahedral representation of a simplex in a Coxeter triangulation or a Freudenthal-Kuhn triangulation. More...
 
class  Gudhi::coxeter_triangulation::Simplex_comparator< Permutahedral_representation_ >
 A comparator class for Permutahedral_representation. The comparison is in lexicographic order first on vertices and then on ordered partitions with sorted parts. The lexicographic order forces that any face is larger than a coface. More...
 

Functions

template<typename... Functions>
Cartesian_product< Functions... > Gudhi::coxeter_triangulation::make_product_function (const Functions &... functions)
 Static constructor of a Cartesian product function. More...
 
template<class Function_ >
Embed_in_Rd< Function_ > Gudhi::coxeter_triangulation::make_embedding (const Function_ &function, std::size_t d)
 Static constructor of an embedding function. More...
 
template<class Function_ >
Linear_transformation< Function_ > Gudhi::coxeter_triangulation::make_linear_transformation (const Function_ &function, const Eigen::MatrixXd &matrix)
 Static constructor of a linearly transformed function. More...
 
template<class Function_ >
Negation< Function_ > Gudhi::coxeter_triangulation::negation (const Function_ &function)
 Static constructor of the negative function. More...
 
template<class Function_ , class Triangulation_ >
PL_approximation< Function_, Triangulation_ > Gudhi::coxeter_triangulation::make_pl_approximation (const Function_ &function, const Triangulation_ &triangulation)
 Static constructor of the piecewise-linear approximation of a function induced by an ambient triangulation. More...
 
template<class Function_ >
Translate< Function_ > Gudhi::coxeter_triangulation::translate (const Function_ &function, Eigen::VectorXd off)
 Static constructor of a translated function. More...
 
template<class Function_ , class Domain_function_ >
Implicit_manifold_intersection_oracle< Function_, Domain_function_ > Gudhi::coxeter_triangulation::make_oracle (const Function_ &function, const Domain_function_ &domain_function)
 Static constructor of an intersection oracle from a function with a domain. More...
 
template<class Function_ >
Implicit_manifold_intersection_oracle< Function_ > Gudhi::coxeter_triangulation::make_oracle (const Function_ &function)
 Static constructor of an intersection oracle from a function without a domain. More...
 
template<class Cell_complex >
Mesh_medit Gudhi::coxeter_triangulation::build_mesh_from_cell_complex (const Cell_complex &cell_complex, Configuration i_configuration=Configuration(), Configuration b_configuration=Configuration())
 Builds a Gudhi::coxeter_triangulation::Mesh_medit from a Gudhi::coxeter_triangulation::Cell_complex. More...
 
template<typename... Meshes>
void Gudhi::coxeter_triangulation::output_meshes_to_medit (std::size_t amb_d, std::string file_name, const Meshes &... meshes)
 Outputs a text file with specified meshes that can be visualized in Medit. More...
 
template<class Point_range , class Triangulation , class Intersection_oracle , class Out_simplex_map >
void Gudhi::coxeter_triangulation::manifold_tracing_algorithm (const Point_range &seed_points, const Triangulation &triangulation, const Intersection_oracle &oracle, Out_simplex_map &out_simplex_map)
 Static method for Manifold_tracing<Triangulation_>::manifold_tracing_algorithm that computes the set of k-simplices that intersect a boundaryless implicit manifold given by an intersection oracle, where k is the codimension of the manifold. The computation is based on the seed propagation — it starts at the given seed points and then propagates along the manifold. More...
 
template<class Point_range , class Triangulation , class Intersection_oracle , class Out_simplex_map >
void Gudhi::coxeter_triangulation::manifold_tracing_algorithm (const Point_range &seed_points, const Triangulation &triangulation, const Intersection_oracle &oracle, Out_simplex_map &interior_simplex_map, Out_simplex_map &boundary_simplex_map)
 Static method for Manifold_tracing<Triangulation_>::manifold_tracing_algorithm the dimensional manifold given by an intersection oracle, where k is the codimension of the manifold. The computation is based on the seed propagation — it starts at the given seed points and then propagates along the manifold. More...
 
template<class Vertex , class OrderedSetPartition >
std::ostream & Gudhi::coxeter_triangulation::operator<< (std::ostream &os, const Permutahedral_representation< Vertex, OrderedSetPartition > &simplex)
 Print a permutahedral representation to a stream. More...
 

Detailed Description

Author
Siargey Kachanovich

Module overview

Coxeter triangulation module is designed to provide tools for constructing a piecewise-linear approximation of an \(m\)-dimensional smooth manifold embedded in \( \mathbb{R}^d \) using an ambient triangulation. For a more detailed description of the module see [35].

Manifold tracing algorithm

The central piece of the module is the manifold tracing algorithm represented by the class Manifold_tracing. The manifold tracing algorithm takes as input a manifold of some dimension \(m\) embedded in \(\mathbb{R}^d\) represented by an intersection oracle (see Section Intersection oracle), a point on the manifold and an ambient triangulation (see Section Ambient triangulations). The output consists of one map (or two maps in the case of manifolds with boundary) from the \((d-m)\)-dimensional (and \((d-m+1)\)-dimensional in the case of manifolds with boundary) simplices in the ambient triangulation that intersect the manifold to their intersection points. From this output, it is possible to construct the cell complex of the piecewise-linear approximation of the input manifold.

There are two methods that execute the manifold tracing algorithm: the method Manifold_tracing::manifold_tracing_algorithm(seed_points, triangulation, oracle, out_simplex_map) for manifolds without boundary and Manifold_tracing::manifold_tracing_algorithm(seed_points, triangulation, oracle, interior_simplex_map,boundary_simplex_map) for manifolds with boundary. The algorithm functions as follows. It starts at the specified seed points and inserts a \((d-m)\)-dimensional simplices nearby each seed point that intersect the manifold into the output. Starting from this simplex, the algorithm propagates the search for other \((d-m)\)-dimensional simplices that intersect the manifold by marching from a simplex to neighbouring simplices via their common cofaces.

This class Manifold_tracing has one template parameter Triangulation_ which specifies the ambient triangulation which is used by the algorithm. The template type Triangulation_ has to be a model of the concept TriangulationForManifoldTracing.

The module also provides two static methods: manifold_tracing_algorithm(seed_points, triangulation, oracle, out_simplex_map) for manifolds without boundary and manifold_tracing_algorithm(seed_points, triangulation, oracle, interior_simplex_map, boundary_simplex_map) for manifolds with boundary. For these static methods it is not necessary to specify any template arguments.

Ambient triangulations

The ambient triangulations supported by the manifold tracing algorithm have to be models of the concept TriangulationForManifoldTracing. This module offers two such models: the class Freudenthal_triangulation and the derived class Coxeter_triangulation.

Both these classes encode affine transformations of the so-called Freudenthal-Kuhn triangulation of \(\mathbb{R}^d\). The Freudenthal-Kuhn triangulation of \(\mathbb{R}^d\) is defined as the simplicial subdivision of the unit cubic partition of \(\mathbb{R}^d\). Each simplex is encoded using the permutahedral representation, which consists of an integer-valued vector \(y\) that positions the simplex in a specific cube in the cubical partition and an ordered partition \(\omega\) of the set \(\{1,\ldots,d+1\}\), which positions the simplex in the simplicial subdivision of the cube. The default constructor Freudenthal_triangulation(d) the Freudenthal-Kuhn triangulation of \(\mathbb{R}^d\). The class Freudenthal_triangulation can also encode any affine transformation of the Freudenthal-Kuhn triangulation of \(\mathbb{R}^d\) using an invertible matrix \(\Lambda\) and an offset vector \(b\) that can be specified in the constructor and which can be changed using the methods change_matrix and change_offset. The class Coxeter_triangulation is derived from Freudenthal_triangulation and its default constructor Coxeter_triangulation(d) builds a Coxeter triangulation of type \(\tilde{A}_d\), which has the best simplex quality of all linear transformations of the Freudenthal-Kuhn triangulation of \(\mathbb{R}^d\).

Coxeter (on the left) and Freudenthal-Kuhn triangulation (on the right)

Intersection oracle

The input \(m\)-dimensional manifold in \(\mathbb{R}^d\) needs to be given via the intersection oracle that answers the following query: given a \((d-m)\)-dimensional simplex, does it intersect the manifold? The concept IntersectionOracle describes all requirements for an intersection oracle class to be compatible with the class Manifold_tracing. This module offers one model of the concept IntersectionOracle, which is the class Implicit_manifold_intersection_oracle. This class represents a manifold given as the zero-set of a specified function \(F: \mathbb{R}^d \rightarrow \mathbb{R}^{d-m}\). The function \(F\) is given by a class which is a model of the concept FunctionForImplicitManifold. There are multiple function classes that are already implemented in this module.

The base function classes above can be composed or modified into new functions using the following classes and methods:

It is also possible to implement your own function as detailed in this Example with a custom function.

Cell complex construction

The output of the manifold tracing algorithm can be transformed into the Hasse diagram of a cell complex that approximates the input manifold using the class Cell_complex. The type of the cells in the Hasse diagram is Hasse_cell<int, double, bool> provided by the module Hasse diagram. The cells in the cell complex given by an object of the class Cell_complex are accessed through several maps that are accessed through the following methods.

The use and interfaces of this Cell_complex is limited to the Coxeter_triangulation implementation.

Examples

Basic example without boundaries

#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;
}
Data structure to store a cell in a Hasse diagram.
Definition: Hasse_diagram_cell.h:49
A class that constructs the cell complex from the output provided by the class Gudhi::coxeter_triangu...
Definition: Cell_complex.h:38
A class that stores Coxeter triangulation of type . This triangulation has the greatest simplex quali...
Definition: Coxeter_triangulation.h:45
A class that assembles methods for manifold tracing algorithm.
Definition: Manifold_tracing.h:39
void manifold_tracing_algorithm(const Point_range &seed_points, const Triangulation &triangulation, const Intersection_oracle &oracle, Out_simplex_map &out_simplex_map)
Static method for Manifold_tracing<Triangulation_>::manifold_tracing_algorithm that computes the set ...
Definition: Manifold_tracing.h:223
Implicit_manifold_intersection_oracle< Function_, Domain_function_ > make_oracle(const Function_ &function, const Domain_function_ &domain_function)
Static constructor of an intersection oracle from a function with a domain.
Definition: Implicit_manifold_intersection_oracle.h:240
A class for the function that defines an m-dimensional implicit sphere embedded in the d-dimensional ...
Definition: Function_Sm_in_Rd.h:27

The program output is:

Vertices:
0 : (-0.680375, 0.523483)
1 : (0.147642, 0.887879)
2 : (-0.847996, 0.30801)
3 : (-0.881369, 0.0951903)
4 : (0.638494, -0.550215)
5 : (0.415344, 0.843848)
6 : (0.812453, -0.0815816)
7 : (0.319625, -0.7709)
8 : (0.319625, 0.889605)
9 : (0.579487, 0.638553)
10 : (-0.680375, -0.461325)
11 : (-0.364269, -0.760962)
Edges:
3 2
3 10
10 11
11 7
7 4
2 0
0 1
6 9
6 4
1 8
8 5
5 9

Example with boundaries

Here is an example of constructing a piecewise-linear approximation of a flat torus embedded in \(\mathbb{R}^4\), rotated by a random rotation in \(\mathbb{R}^4\) and cut by a hyperplane.

// workaround for the annoying boost message in boost 1.69
#define BOOST_PENDING_INTEGER_LOG2_HPP
#include <boost/integer/integer_log2.hpp>
// end workaround
#include <iostream>
#include <gudhi/Coxeter_triangulation.h>
#include <gudhi/Functions/Function_affine_plane_in_Rd.h>
#include <gudhi/Functions/Function_Sm_in_Rd.h>
#include <gudhi/Functions/Cartesian_product.h>
#include <gudhi/Functions/Linear_transformation.h>
#include <gudhi/Implicit_manifold_intersection_oracle.h>
#include <gudhi/Manifold_tracing.h>
#include <gudhi/Coxeter_triangulation/Cell_complex/Cell_complex.h>
#include <gudhi/Functions/random_orthogonal_matrix.h> // requires CGAL
#include <gudhi/IO/build_mesh_from_cell_complex.h>
#include <gudhi/IO/output_meshes_to_medit.h>
using namespace Gudhi::coxeter_triangulation;
int main(int argc, char** argv) {
// Creating a circle S1 in R2 of specified radius
double radius = 1.0;
Function_Sm_in_Rd fun_circle(radius, 1);
// Creating a flat torus S1xS1 in R4 from two circle functions
auto fun_flat_torus = make_product_function(fun_circle, fun_circle);
// Apply a random rotation in R4
auto matrix = random_orthogonal_matrix(4);
auto fun_flat_torus_rotated = make_linear_transformation(fun_flat_torus, matrix);
// Computing the seed of the function fun_flat_torus
Eigen::VectorXd seed = fun_flat_torus_rotated.seed();
// Defining a domain function that defines the boundary, which is a hyperplane passing by the origin and orthogonal to
// x.
Eigen::MatrixXd normal_matrix = Eigen::MatrixXd::Zero(4, 1);
for (std::size_t i = 0; i < 4; i++) normal_matrix(i, 0) = -seed(i);
Function_affine_plane_in_Rd fun_bound(normal_matrix, -seed / 2);
// Defining the intersection oracle
auto oracle = make_oracle(fun_flat_torus_rotated, fun_bound);
// Define a Coxeter triangulation scaled by a factor lambda.
// The triangulation is translated by a random vector to avoid violating the genericity hypothesis.
double lambda = 0.2;
Coxeter_triangulation<> cox_tr(oracle.amb_d());
cox_tr.change_offset(Eigen::VectorXd::Random(oracle.amb_d()));
cox_tr.change_matrix(lambda * cox_tr.matrix());
// Manifold tracing algorithm
using Out_simplex_map = typename MT::Out_simplex_map;
std::vector<Eigen::VectorXd> seed_points(1, seed);
Out_simplex_map interior_simplex_map, boundary_simplex_map;
manifold_tracing_algorithm(seed_points, cox_tr, oracle, interior_simplex_map, boundary_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, boundary_simplex_map);
// Output the cell complex to a file readable by medit
output_meshes_to_medit(3, "flat_torus_with_boundary",
build_mesh_from_cell_complex(cell_complex, Configuration(true, true, true, 1, 5, 3),
Configuration(true, true, true, 2, 13, 14)));
return 0;
}
Linear_transformation< Function_ > make_linear_transformation(const Function_ &function, const Eigen::MatrixXd &matrix)
Static constructor of a linearly transformed function.
Definition: Linear_transformation.h:80
void output_meshes_to_medit(std::size_t amb_d, std::string file_name, const Meshes &... meshes)
Outputs a text file with specified meshes that can be visualized in Medit.
Definition: output_meshes_to_medit.h:81
Mesh_medit build_mesh_from_cell_complex(const Cell_complex &cell_complex, Configuration i_configuration=Configuration(), Configuration b_configuration=Configuration())
Builds a Gudhi::coxeter_triangulation::Mesh_medit from a Gudhi::coxeter_triangulation::Cell_complex.
Definition: build_mesh_from_cell_complex.h:122
Cartesian_product< Functions... > make_product_function(const Functions &... functions)
Static constructor of a Cartesian product function.
Definition: Cartesian_product.h:149
A class for the function that defines an m-dimensional implicit affine plane embedded in d-dimensiona...
Definition: Function_affine_plane_in_Rd.h:27

The output in medit is:

Output from the example of a flat torus with boundary

Example with a custom function

In the following more complex example, we define a custom function for the implicit manifold.

#include <iostream>
#include <gudhi/Coxeter_triangulation.h>
#include <gudhi/Functions/Function_Sm_in_Rd.h>
#include <gudhi/Implicit_manifold_intersection_oracle.h>
#include <gudhi/Manifold_tracing.h>
#include <gudhi/Coxeter_triangulation/Cell_complex/Cell_complex.h>
#include <gudhi/Functions/Linear_transformation.h>
#include <gudhi/IO/build_mesh_from_cell_complex.h>
#include <gudhi/IO/output_meshes_to_medit.h>
using namespace Gudhi::coxeter_triangulation;
/* A definition of a function that defines a 2d surface embedded in R^4, but that normally
* lives on a complex projective plane.
* In terms of harmonic coordinates [x:y:z] of points on the complex projective plane,
* the equation of the manifold is x^3*y + y^3*z + z^3*x = 0.
* The embedding consists of restricting the manifold to the affine subspace z = 1.
*/
struct Function_surface_on_CP2_in_R4 {
Eigen::VectorXd operator()(const Eigen::VectorXd& p) const {
// The real and imaginary parts of the variables x and y
double xr = p(0), xi = p(1), yr = p(2), yi = p(3);
Eigen::VectorXd result(cod_d());
// Squares and cubes of real and imaginary parts used in the computations
double xr2 = xr * xr, xi2 = xi * xi, yr2 = yr * yr, yi2 = yi * yi, xr3 = xr2 * xr, xi3 = xi2 * xi, yr3 = yr2 * yr,
yi3 = yi2 * yi;
// The first coordinate of the output is Re(x^3*y + y^3 + x)
result(0) = xr3 * yr - 3 * xr * xi2 * yr - 3 * xr2 * xi * yi + xi3 * yi + yr3 - 3 * yr * yi2 + xr;
// The second coordinate of the output is Im(x^3*y + y^3 + x)
result(1) = 3 * xr2 * xi * yr + xr3 * yi - 3 * xr * xi2 * yi - xi3 * yr + 3 * yr2 * yi - yi3 + xi;
return result;
}
std::size_t amb_d() const { return 4; };
std::size_t cod_d() const { return 2; };
Eigen::VectorXd seed() const {
Eigen::VectorXd result = Eigen::VectorXd::Zero(4);
return result;
}
Function_surface_on_CP2_in_R4() {}
};
int main(int argc, char** argv) {
// The function for the (non-compact) manifold
Function_surface_on_CP2_in_R4 fun;
// Seed of the function
Eigen::VectorXd seed = fun.seed();
// Creating the function that defines the boundary of a compact region on the manifold
double radius = 3.0;
Function_Sm_in_Rd fun_sph(radius, 3, seed);
// Defining the intersection oracle
auto oracle = make_oracle(fun, fun_sph);
// Define a Coxeter triangulation scaled by a factor lambda.
// The triangulation is translated by a random vector to avoid violating the genericity hypothesis.
double lambda = 0.2;
Coxeter_triangulation<> cox_tr(oracle.amb_d());
cox_tr.change_offset(Eigen::VectorXd::Random(oracle.amb_d()));
cox_tr.change_matrix(lambda * cox_tr.matrix());
// Manifold tracing algorithm
using Out_simplex_map = typename MT::Out_simplex_map;
std::vector<Eigen::VectorXd> seed_points(1, seed);
Out_simplex_map interior_simplex_map, boundary_simplex_map;
manifold_tracing_algorithm(seed_points, cox_tr, oracle, interior_simplex_map, boundary_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, boundary_simplex_map);
// Output the cell complex to a file readable by medit
output_meshes_to_medit(3, "manifold_on_CP2_with_boundary",
build_mesh_from_cell_complex(cell_complex, Configuration(true, true, true, 1, 5, 3),
Configuration(true, true, true, 2, 13, 14)));
return 0;
}

The output in medit looks as follows:

Output from the example with a custom function

Iterator types for Permutahedral_representation

Function Documentation

◆ build_mesh_from_cell_complex()

template<class Cell_complex >
Mesh_medit Gudhi::coxeter_triangulation::build_mesh_from_cell_complex ( const Cell_complex cell_complex,
Configuration  i_configuration = Configuration(),
Configuration  b_configuration = Configuration() 
)

◆ make_embedding()

template<class Function_ >
Embed_in_Rd< Function_ > Gudhi::coxeter_triangulation::make_embedding ( const Function_ &  function,
std::size_t  d 
)

Static constructor of an embedding function.

Parameters
[in]functionThe function to be embedded in higher dimension.
[in]dEmbedding dimension.
Template Parameters
Function_The function template parameter. Should be a model of the concept FunctionForImplicitManifold.

◆ make_linear_transformation()

template<class Function_ >
Linear_transformation< Function_ > Gudhi::coxeter_triangulation::make_linear_transformation ( const Function_ &  function,
const Eigen::MatrixXd &  matrix 
)

Static constructor of a linearly transformed function.

Parameters
[in]functionThe function to be linearly transformed.
[in]matrixThe transformation matrix. Its dimension should be d*d, where d is the domain (ambient) dimension of 'function'.
Template Parameters
Function_The function template parameter. Should be a model of the concept FunctionForImplicitManifold.

◆ make_oracle() [1/2]

template<class Function_ >
Implicit_manifold_intersection_oracle< Function_ > Gudhi::coxeter_triangulation::make_oracle ( const Function_ &  function)

Static constructor of an intersection oracle from a function without a domain.

Parameters
functionThe input function that represents the implicit manifold without boundary.

◆ make_oracle() [2/2]

template<class Function_ , class Domain_function_ >
Implicit_manifold_intersection_oracle< Function_, Domain_function_ > Gudhi::coxeter_triangulation::make_oracle ( const Function_ &  function,
const Domain_function_ &  domain_function 
)

Static constructor of an intersection oracle from a function with a domain.

Parameters
functionThe input function that represents the implicit manifold before the restriction with the domain.
domain_functionThe input domain function that can be used to define an implicit manifold with boundary.

◆ make_pl_approximation()

template<class Function_ , class Triangulation_ >
PL_approximation< Function_, Triangulation_ > Gudhi::coxeter_triangulation::make_pl_approximation ( const Function_ &  function,
const Triangulation_ &  triangulation 
)

Static constructor of the piecewise-linear approximation of a function induced by an ambient triangulation.

Parameters
[in]functionThe function.
[in]triangulationThe ambient triangulation.
Template Parameters
Function_The function template parameter. Should be a model of the concept FunctionForImplicitManifold.

◆ make_product_function()

template<typename... Functions>
Cartesian_product< Functions... > Gudhi::coxeter_triangulation::make_product_function ( const Functions &...  functions)

Static constructor of a Cartesian product function.

Parameters
[in]functionsThe functions the zero-sets of which are factors in the Cartesian product of the resulting function.
Template Parameters
FunctionsA pack template parameter for functions. All functions should be models of the concept FunctionForImplicitManifold.

◆ manifold_tracing_algorithm() [1/2]

template<class Point_range , class Triangulation , class Intersection_oracle , class Out_simplex_map >
void Gudhi::coxeter_triangulation::manifold_tracing_algorithm ( const Point_range &  seed_points,
const Triangulation &  triangulation,
const Intersection_oracle &  oracle,
Out_simplex_map &  interior_simplex_map,
Out_simplex_map &  boundary_simplex_map 
)

Static method for Manifold_tracing<Triangulation_>::manifold_tracing_algorithm the dimensional manifold given by an intersection oracle, where k is the codimension of the manifold. The computation is based on the seed propagation — it starts at the given seed points and then propagates along the manifold.

Template Parameters
Point_rangeRange of points of type Eigen::VectorXd.
Triangulation_The type of the ambient triangulation. Needs to be a model of the concept TriangulationForManifoldTracing.
Intersection_oracleIntersection oracle that represents the manifold. Needs to be a model of the concept IntersectionOracle.
Out_simplex_mapNeeds to be Manifold_tracing<Triangulation_>::Out_simplex_map.
Parameters
[in]seed_pointsThe range of points on the manifold from which the computation begins.
[in]triangulationThe ambient triangulation.
[in]oracleThe intersection oracle for the manifold. The ambient dimension needs to match the dimension of the triangulation.
[out]interior_simplex_mapThe output map, where the keys are k-simplices in the input triangulation that intersect the relative interior of the input manifold and the mapped values are the intersection points.
[out]boundary_simplex_mapThe output map, where the keys are k-simplices in the input triangulation that intersect the boundary of the input manifold and the mapped values are the intersection points.

◆ manifold_tracing_algorithm() [2/2]

template<class Point_range , class Triangulation , class Intersection_oracle , class Out_simplex_map >
void Gudhi::coxeter_triangulation::manifold_tracing_algorithm ( const Point_range &  seed_points,
const Triangulation &  triangulation,
const Intersection_oracle &  oracle,
Out_simplex_map &  out_simplex_map 
)

Static method for Manifold_tracing<Triangulation_>::manifold_tracing_algorithm that computes the set of k-simplices that intersect a boundaryless implicit manifold given by an intersection oracle, where k is the codimension of the manifold. The computation is based on the seed propagation — it starts at the given seed points and then propagates along the manifold.

Template Parameters
Point_rangeRange of points of type Eigen::VectorXd.
Triangulation_The type of the ambient triangulation. Needs to be a model of the concept TriangulationForManifoldTracing.
Intersection_oracleIntersection oracle that represents the manifold. Needs to be a model of the concept IntersectionOracle.
Out_simplex_mapNeeds to be Manifold_tracing<Triangulation_>::Out_simplex_map.
Parameters
[in]seed_pointsThe range of points on the manifold from which the computation begins.
[in]triangulationThe ambient triangulation.
[in]oracleThe intersection oracle for the manifold. The ambient dimension needs to match the dimension of the triangulation.
[out]out_simplex_mapThe output map, where the keys are k-simplices in the input triangulation that intersect the input manifold and the mapped values are the intersection points.

◆ negation()

template<class Function_ >
Negation< Function_ > Gudhi::coxeter_triangulation::negation ( const Function_ &  function)

Static constructor of the negative function.

Parameters
[in]functionThe function to be translated. domain (ambient) dimension of 'function'.
Template Parameters
Function_The function template parameter. Should be a model of the concept FunctionForImplicitManifold.

◆ operator<<()

template<class Vertex , class OrderedSetPartition >
std::ostream & Gudhi::coxeter_triangulation::operator<< ( std::ostream &  os,
const Permutahedral_representation< Vertex, OrderedSetPartition > &  simplex 
)

Print a permutahedral representation to a stream.

Parameters
[in]osThe output stream.
[in]simplexA simplex represented by its permutahedral representation.

◆ output_meshes_to_medit()

template<typename... Meshes>
void Gudhi::coxeter_triangulation::output_meshes_to_medit ( std::size_t  amb_d,
std::string  file_name,
const Meshes &...  meshes 
)

Outputs a text file with specified meshes that can be visualized in Medit.

Parameters
[in]amb_dAmbient dimension. Can be 2 or 3.
[in]file_nameThe name of the output file.
[in]meshesA pack of meshes to be specified separated by commas.

◆ translate()

template<class Function_ >
Translate< Function_ > Gudhi::coxeter_triangulation::translate ( const Function_ &  function,
Eigen::VectorXd  off 
)

Static constructor of a translated function.

Parameters
[in]functionThe function to be translated.
[in]offThe offset vector. The dimension should correspond to the domain (ambient) dimension of 'function'.
Template Parameters
Function_The function template parameter. Should be a model of the concept FunctionForImplicitManifold.