CGAL 6.1 - Homological Discrete Vector Fields
Loading...
Searching...
No Matches
CGAL::HDVF::Hdvf_duality< CoefficientType, ComplexType > Class Template Reference

#include <CGAL/HDVF/Hdvf_duality.h>

Inherits from

CGAL::HDVF::Hdvf_core< CoefficientType, ComplexType, OSM::Sparse_chain, OSM::Sub_sparse_matrix >.

Definition

template<typename CoefficientType, typename ComplexType>
class CGAL::HDVF::Hdvf_duality< CoefficientType, ComplexType >

The class Hdvf_duality is the implementation of homological discrete vector fields (HDVF for short) for Alexander duality computation.

Warning
The ring of coefficients provided should be a field.

In dimension \(n\), given a complex \(L\) homeomorphic to \(\mathcal S^n\) and a sub-complex \(K\subseteq L\), Alexander duality states that for all \(q\leqslant n\):

\[\tilde H_q(K) \simeq \tilde H^{n-q-1}(L-K)\]

where \(\tilde H_q\) and \(\tilde H^q\) denote reduced homology and cohomology groups.

In [Gonzalez and al. 2025], the authors prove that, even if \(L-K\) is not a sub-complex, it produces a valid chain complex (we call "co-complex" such a complementary of sub-complex). Hence, its homology/cohomology can be computed and for all \(q\leqslant n\):

\[\tilde H_q(K) \simeq \tilde H^{q+1}(L-K)\]

HDVFs provide a fast and convenient mean to compute this isomorphism. In order to work with convenient finite complexes, the complex \(L\) must be homeomorphic to a ball of dimension \(n\) (thus \(\mathcal S^n\) is actually homeomorphic to \(L\) plus an infinite \(n\)-cell closing its boundary).

Perfect HDVFs are first computed over \(K\) and \(L-K\) (providing corresponding relative homology) respectively and Alexander isomorphism gives rise to a pairing between critical cells in \(K\) and \(L-K\), that is a pairing between homology/cohomology generators in \(K\) and \(L-K\).

The class provides HDVF constuction operations: compute_perfect_hdvf() and compute_rand_perfect_hdvf(), which build perfect HDVFs over \(K\) and \(L-K\) respectively. Then, compute_alexander_pairing() computes Alexander isomorphism (and provides a pairing between homology/cohomology generators in \(K\) and \(L-K\)).

Example of Alexander duality isomorphism. The twirl mesh is a subcomplex K of a larger complex L depicted in yellow, homeomorphic to the ball of dimension 3 (right - sectional view).

Figure 92.1 Example of "homological quartet for the twirl model". 1: Homology generators of the twirl \(H_1(K)\), 2: Cohomology generators of the twirl \(H^1(K)\), 3: Homology generators of the complementary of the twirl \(H_1(L-K)\), 4: Cohomology generators of the complementary of the twirl \(H^1(L-K)\). Alexander isomorphism is represented through colours (paired generators have similar colours).


Hence, each hole in \(K\) gives rise to four generators (called its "homological quarted": its homology and cohomology generators in \(K_q\) and the homology and cohomology generators paired with them in \((L-K)_{q+1}\)).

In order to compute relative homology, a sub chain complex mask is used to partially screen the complex L and thus restrict HDVF computation. This mask is called "current mask" (and can be set over K or L-K).

Is model of
HDVF
Template Parameters
CoefficientTypea model of the Ring concept providing the ring used to compute homology.
ComplexTypea model of the AbstractChainComplex concept, providing the type of abstract chain complex used.

[Gonzalez and al. 2025] Gonzalez-Lorenzo, A., Bac, A. & Gazull, YS. A constructive approach of Alexander duality. J Appl. and Comput. Topology 9, 2 (2025).

Examples
HDVF/main_dual_hdvf.cpp.

Public Member Functions

 Hdvf_duality (const ComplexType &L, Sub_chain_complex_mask< CoefficientType, ComplexType > &K, int hdvf_opt=OPT_FULL)
 Hdvf_duality constructor ( from a complex L and a sub-complex K)
 
virtual PairCell find_pair_A (int q, bool &found) const
 Finds a valid PairCell of dimension q / q+1 for A in the current sub chain complex.
 
virtual PairCell find_pair_A (int q, bool &found, size_t gamma) const
 Finds a valid PairCell for A containing gamma in the current sub chain complex (a cell of dimension q)
 
virtual std::vector< PairCellfind_pairs_A (int q, bool &found) const
 Finds all valid PairCell of dimension q / q+1 in the current sub chain complex for A.
 
virtual std::vector< PairCellfind_pairs_A (int q, bool &found, size_t gamma) const
 Finds all valid PairCell for A containing gamma in the current sub chain complex (a cell of dimension q)
 
void set_mask_K ()
 Sets the current sub chain complex masks over K.
 
void set_mask_L_K ()
 Sets the current sub chain complex masks over L-K.
 
Sub_chain_complex_mask< CoefficientType, ComplexType > get_current_mask ()
 Returns the value of the current sub chain complex mask.
 
std::vector< PairCellcompute_perfect_hdvf (bool verbose=false)
 Computes a perfect HDVF over the current sub chain complex.
 
std::vector< PairCellcompute_rand_perfect_hdvf (bool verbose=false)
 Computes a random perfect HDVF over the current sub chain complex.
 
std::vector< PairCellcompute_pairing_hdvf ()
 Computes a "pairing" HDVF between K and L-K.
 
std::vector< PairCellcompute_rand_pairing_hdvf ()
 Computes a random "pairing" HDVF between K and L-K.
 
vector< vector< size_t > > get_flag (FlagType flag) const
 Gets cells with a given flag in any dimension in the current sub chain complex.
 
vector< size_t > get_flag_dim (FlagType flag, int q) const
 Gets cells with a given flag in dimension q in the current sub chain complex.
 
virtual ostream & print_reduction (ostream &out=cout)
 Prints the homology and cohomology reduction information for K and L-K.
 
virtual ostream & print_reduction_sub (ostream &out=cout)
 Prints the homology and cohomology reduction information for the current such chain complex.
 
ostream & print_bnd_pairing (ostream &out=cout)
 Prints the reduced boundary over critical cells of K and L-K.
 
virtual vector< vector< int > > get_psc_labels () const
 Exports primary/secondary/critical labels of the current sub chain complex for vtk export.
 
virtual CChain get_homology_chain (size_t cell, int dim) const
 Exports homology generators of the current sub chain complex associated to cell (critical cell) of dimension q (used by vtk export).
 
virtual CChain get_cohomology_chain (size_t cell, int dim) const
 Exports cohomology generators of the current sub chain complex associated to cell (critical cell) of dimension q (used by vtk export).
 
- Public Member Functions inherited from CGAL::HDVF::Hdvf_core< CoefficientType, ComplexType, OSM::Sparse_chain, OSM::Sub_sparse_matrix >
 Hdvf_core (const ComplexType &K, int hdvf_opt=OPT_FULL)
 Default constructor.
 
 Hdvf_core (const Hdvf_core &hdvf)
 
 ~Hdvf_core ()
 
virtual PairCell find_pair_A (int q, bool &found) const
 Finds a valid PairCell of dimension q / q+1 for A.
 
virtual PairCell find_pair_A (int q, bool &found, size_t gamma) const
 Finds a valid PairCell for A containing gamma (a cell of dimension q)
 
virtual std::vector< PairCellfind_pairs_A (int q, bool &found) const
 Finds all valid PairCell of dimension q / q+1 for A.
 
virtual std::vector< PairCellfind_pairs_A (int q, bool &found, size_t gamma) const
 Finds all valid PairCell for A containing gamma (a cell of dimension q)
 
void A (size_t gamma1, size_t gamma2, int q)
 A operation: pairs critical cells.
 
std::vector< PairCellcompute_perfect_hdvf (bool verbose=false)
 Computes a perfect HDVF.
 
std::vector< PairCellcompute_rand_perfect_hdvf (bool verbose=false)
 Computes a random perfect HDVF.
 
bool is_perfect_hdvf ()
 Tests if a HDVF is perfect.
 
virtual std::vector< std::vector< size_t > > get_flag (FlagType flag) const
 Gets cells with a given flag in any dimension.
 
virtual std::vector< size_t > get_flag_dim (FlagType flag, int q) const
 Gets cells with a given flag in dimension q.
 
FlagType get_cell_flag (int q, size_t tau) const
 Gets the flag of the cell tau in dimension q.
 
int get_hdvf_opts () const
 Gets HDVF computation option.
 
const RMatrixget_f (int q) const
 Gets the row-major matrix of \(f\) (from the reduction associated to the HDVF).
 
const CMatrixget_g (int q) const
 Gets the column-major matrix of \(g\) (from the reduction associated to the HDVF).
 
const CMatrixget_h (int q) const
 Gets the column-major matrix of \(h\) (from the reduction associated to the HDVF).
 
const CMatrixget_dd (int q) const
 Gets the column-major matrix of \(\partial'\), reduced boundary operator (from the reduction associated to the HDVF).
 
std::ostream & print_matrices (std::ostream &out=std::cout) const
 Prints the matrices of the reduction.
 
std::ostream & print_reduction (std::ostream &out=std::cout) const
 Prints the homology and cohomology reduction information.
 
virtual std::vector< std::vector< int > > get_psc_labels () const
 Exports primary/secondary/critical labels (in particular for vtk export)
 
virtual CChain get_homology_chain (size_t cell, int q) const
 Gets homology generators associated to cell (critical cell) of dimension q (used by vtk export).
 
virtual CChain get_cohomology_chain (size_t cell, int dim) const
 Gets cohomology generators associated to cell (critical cell) of dimension q (used by vtk export).
 
std::ostream & save_hdvf_reduction (std::ostream &out)
 Saves a HDVF together with the associated reduction (f, g, h, d matrices)
 
std::istream & load_hdvf_reduction (std::istream &out)
 Loads a HDVF together with the associated reduction (f, g, h, d matrices)
 

Additional Inherited Members

- Public Types inherited from CGAL::HDVF::Hdvf_core< CoefficientType, ComplexType, OSM::Sparse_chain, OSM::Sub_sparse_matrix >
typedef OSM::Sparse_chain< CoefficientType, OSM::COLUMN > CChain
 Type of column-major chains.
 
typedef OSM::Sparse_chain< CoefficientType, OSM::ROW > RChain
 Type of row-major chains.
 
typedef OSM::Sub_sparse_matrix< CoefficientType, OSM::COLUMN > CMatrix
 Type of column-major sparse matrices.
 
typedef OSM::Sub_sparse_matrix< CoefficientType, OSM::ROW > RMatrix
 Type of row-major sparse matrices.
 
- Protected Member Functions inherited from CGAL::HDVF::Hdvf_core< CoefficientType, ComplexType, OSM::Sparse_chain, OSM::Sub_sparse_matrix >
OSM::Sparse_chain< CoefficientType, ChainTypeFlag > projection (const OSM::Sparse_chain< CoefficientType, ChainTypeFlag > &chain, FlagType flag, int q) const
 
void progress_bar (size_t i, size_t n)
 
- Protected Attributes inherited from CGAL::HDVF::Hdvf_core< CoefficientType, ComplexType, OSM::Sparse_chain, OSM::Sub_sparse_matrix >
std::vector< std::vector< FlagType > > _flag
 
std::vector< size_t > _nb_P
 
std::vector< size_t > _nb_S
 
std::vector< size_t > _nb_C
 
std::vector< RMatrix_F_row
 
std::vector< CMatrix_G_col
 
std::vector< CMatrix_H_col
 
std::vector< CMatrix_DD_col
 
const ComplexType & _K
 
int _hdvf_opt
 

Constructor & Destructor Documentation

◆ Hdvf_duality()

template<typename CoefficientType , typename ComplexType >
CGAL::HDVF::Hdvf_duality< CoefficientType, ComplexType >::Hdvf_duality ( const ComplexType &  L,
Sub_chain_complex_mask< CoefficientType, ComplexType > &  K,
int  hdvf_opt = OPT_FULL 
)

Hdvf_duality constructor ( from a complex L and a sub-complex K)

L is a complex of a given dimension \(n\) homeomorphic to \(\mathcal B^n\) and K is a sub-complex of L described by a bitboard (cells of K have a bit set to 1, cells of K have a bit set to 0).

Initially, the sub chain complex mask is set to K.

Parameters
[in]LA complex of a given dimension \(n\) homeomorphic to \(\mathcal B^n\).
[in]KA sub complex of L encoded through a bitboard.
[in]hdvf_optOption for HDVF computation (OPT_BND, OPT_F, OPT_G or OPT_FULL).

Member Function Documentation

◆ compute_pairing_hdvf()

template<typename CoefficientType , typename ComplexType >
std::vector< PairCell > CGAL::HDVF::Hdvf_duality< CoefficientType, ComplexType >::compute_pairing_hdvf

Computes a "pairing" HDVF between K and L-K.

Warning
Run compute_perfect_hdvf() first (to build perfect HDVFs over K and L-K respectively).

The function computes a perfect HDVF over remaining critical cells. Each pair of cells inserted with the A() operation maps corresponding homology/cohomology generators in the Alexander isomorphism.

Returns
The vector of paired critical cells (encoding Alexander isomorphism).
Examples
HDVF/main_dual_hdvf.cpp.

◆ compute_perfect_hdvf()

template<typename CoefficientType , typename ComplexType >
std::vector< PairCell > CGAL::HDVF::Hdvf_duality< CoefficientType, ComplexType >::compute_perfect_hdvf ( bool  verbose = false)

Computes a perfect HDVF over the current sub chain complex.

As long as valid pairs for A exist in the current sub chain complex, the function selects the first available pair (returned by find_pair_A()) and applies the corresponding A() operation. If the Ring of coefficients is a field, this operation always produces a perfect HDVF (ie. the reduced boundary is null and the reduction provides homology and cohomology information). Otherwise the operation produces a maximal HDVF with a residual boundary matrix over critical cells.

If the HDVF is initially not trivial (some cells have already been paired), the function completes it into a perfect HDVF.

Parameters
[in]verboseIf this parameter is true, all intermediate reductions are printed out.
Returns
The vector of all PairCell paired with A.

◆ compute_rand_pairing_hdvf()

template<typename CoefficientType , typename ComplexType >
std::vector< PairCell > CGAL::HDVF::Hdvf_duality< CoefficientType, ComplexType >::compute_rand_pairing_hdvf

Computes a random "pairing" HDVF between K and L-K.

Warning
Run compute_perfect_hdvf() first (to build perfect HDVFs over K and L-K respectively).

The function computes a random perfect HDVF over remaining critical cells. Each pair of cells inserted with the `A() operation maps corresponding homology/cohomology generators in the Alexander isomorphism.

Returns
The vector of paired critical cells (encoding Alexander isomorphism).

◆ compute_rand_perfect_hdvf()

template<typename CoefficientType , typename ComplexType >
std::vector< PairCell > CGAL::HDVF::Hdvf_duality< CoefficientType, ComplexType >::compute_rand_perfect_hdvf ( bool  verbose = false)

Computes a random perfect HDVF over the current sub chain complex.

As long as valid pairs for A exist in the current sub chain complex, the function selects a random pair (among pairs returned by find_pairs_A()) and applies the corresponding A() operation. If the Ring of coefficients is a field, this operation always produces a perfect HDVF (that is the reduced boundary is null and the reduction provides homology and cohomology information).

If the HDVF is initially not trivial (some cells have already been paired), the function randomly completes it into a perfect HDVF.

Warning
This method is slower that compute_perfect_hdvf() (finding out all possible valid pairs requires additional time).
Parameters
[in]verboseIf this parameter is true, all intermediate reductions are printed out.
Returns
The vector of all pairs of cells used for apply A.

◆ find_pair_A() [1/2]

template<typename CoefficientType , typename ComplexType >
PairCell CGAL::HDVF::Hdvf_duality< CoefficientType, ComplexType >::find_pair_A ( int  q,
bool &  found 
) const
virtual

Finds a valid PairCell of dimension q / q+1 for A in the current sub chain complex.

The function searches a pair of critical cells, in the current sub chain complex, \((\gamma_1, \gamma2)\) of dimension q / q+1, valid for A (ie. such that \(\langle \partial_{q+1}(\gamma_2), \gamma_1 \rangle\) invertible). It returns the first valid pair found by iterators.

Parameters
[in]qLower dimension of the pair.
[in]foundReference to a Boolean variable. The method sets found to true if a valid pair is found, false otherwise.

Reimplemented from CGAL::HDVF::Hdvf_core< CoefficientType, ComplexType, OSM::Sparse_chain, OSM::Sub_sparse_matrix >.

◆ find_pair_A() [2/2]

template<typename CoefficientType , typename ComplexType >
PairCell CGAL::HDVF::Hdvf_duality< CoefficientType, ComplexType >::find_pair_A ( int  q,
bool &  found,
size_t  gamma 
) const
virtual

Finds a valid PairCell for A containing gamma in the current sub chain complex (a cell of dimension q)

The function searches a cell \(\gamma'\) in the current sub chain complex such that one of the following conditions holds:

  • \(\gamma'\) has dimension q+1 and \((\gamma, \gamma')\) is valid for A (ie. such that \(\langle \partial_{q+1}(\gamma'), \gamma \rangle\) invertible),
  • \(\gamma'\) has dimension q-1 and \((\gamma', \gamma)\) is valid for A (ie. such that \(\langle \partial_{q}(\gamma), \gamma' \rangle\) invertible).
Parameters
[in]qDimension of the cell gamma.
[in]foundReference to a Boolean variable. The method sets found to true if a valid pair is found, false otherwise.
[in]gammaIndex of a cell to pair.

Reimplemented from CGAL::HDVF::Hdvf_core< CoefficientType, ComplexType, OSM::Sparse_chain, OSM::Sub_sparse_matrix >.

◆ find_pairs_A() [1/2]

template<typename CoefficientType , typename ComplexType >
std::vector< PairCell > CGAL::HDVF::Hdvf_duality< CoefficientType, ComplexType >::find_pairs_A ( int  q,
bool &  found 
) const
virtual

Finds all valid PairCell of dimension q / q+1 in the current sub chain complex for A.

The function searches all pairs of critical cells \((\gamma_1, \gamma2)\) in the current sub chain complex of dimension q / q+1, valid for A (ie. such that \(\langle \partial_{q+1}(\gamma_2), \gamma_1 \rangle\) invertible). It returns a vector of such pairs.

Parameters
[in]qLower dimension of the pair.
[in]foundReference to a Boolean variable. The method sets found to true if a valid pair is found, false otherwise.

Reimplemented from CGAL::HDVF::Hdvf_core< CoefficientType, ComplexType, OSM::Sparse_chain, OSM::Sub_sparse_matrix >.

◆ find_pairs_A() [2/2]

template<typename CoefficientType , typename ComplexType >
std::vector< PairCell > CGAL::HDVF::Hdvf_duality< CoefficientType, ComplexType >::find_pairs_A ( int  q,
bool &  found,
size_t  gamma 
) const
virtual

Finds all valid PairCell for A containing gamma in the current sub chain complex (a cell of dimension q)

The function searches all CRITICAL cells \(\gamma'\) in the current sub chain complex such that one of the following conditions holds:

  • \(\gamma'\) has dimension q+1 and \((\gamma, \gamma')\) is valid for A (ie. such that \(\langle \partial_{q+1}(\gamma'), \gamma \rangle\) invertible),
  • \(\gamma'\) has dimension q-1 and \((\gamma', \gamma)\) is valid for A (ie. such that \(\langle \partial_{q}(\gamma), \gamma' \rangle\) invertible). It returns a vector of such pairs.
Parameters
[in]qDimension of the cell gamma.
[in]foundReference to a Boolean variable. The method sets found to true if a valid pair is found, false otherwise.
[in]gammaIndex of a cell to pair.

Reimplemented from CGAL::HDVF::Hdvf_core< CoefficientType, ComplexType, OSM::Sparse_chain, OSM::Sub_sparse_matrix >.

◆ get_cohomology_chain()

template<typename CoefficientType , typename ComplexType >
virtual CChain CGAL::HDVF::Hdvf_duality< CoefficientType, ComplexType >::get_cohomology_chain ( size_t  cell,
int  dim 
) const
virtual

Exports cohomology generators of the current sub chain complex associated to cell (critical cell) of dimension q (used by vtk export).

The method exports the chain \(f^\star(\sigma)\) for \(\sigma\) the cell of index cell and dimension q.

Returns
A column-major chain.

Reimplemented from CGAL::HDVF::Hdvf_core< CoefficientType, ComplexType, OSM::Sparse_chain, OSM::Sub_sparse_matrix >.

◆ get_flag()

template<typename CoefficientType , typename ComplexType >
vector< vector< size_t > > CGAL::HDVF::Hdvf_duality< CoefficientType, ComplexType >::get_flag ( FlagType  flag) const
virtual

Gets cells with a given flag in any dimension in the current sub chain complex.

The function returns a vector containing, for each dimension, the vector of cells with a given flag.

Parameters
[in]flagFlag to select.

Reimplemented from CGAL::HDVF::Hdvf_core< CoefficientType, ComplexType, OSM::Sparse_chain, OSM::Sub_sparse_matrix >.

◆ get_flag_dim()

template<typename CoefficientType , typename ComplexType >
vector< size_t > CGAL::HDVF::Hdvf_duality< CoefficientType, ComplexType >::get_flag_dim ( FlagType  flag,
int  q 
) const
virtual

Gets cells with a given flag in dimension q in the current sub chain complex.

The function returns the vector of cells of dimension q with a given flag.

Parameters
[in]flagFlag to select.
[in]qDimension visited.

Reimplemented from CGAL::HDVF::Hdvf_core< CoefficientType, ComplexType, OSM::Sparse_chain, OSM::Sub_sparse_matrix >.

◆ get_homology_chain()

template<typename CoefficientType , typename ComplexType >
virtual CChain CGAL::HDVF::Hdvf_duality< CoefficientType, ComplexType >::get_homology_chain ( size_t  cell,
int  dim 
) const
virtual

Exports homology generators of the current sub chain complex associated to cell (critical cell) of dimension q (used by vtk export).

The method exports the chain \(g(\sigma)\) for \(\sigma\) the cell of index cell and dimension q.

Returns
A column-major chain.

Reimplemented from CGAL::HDVF::Hdvf_core< CoefficientType, ComplexType, OSM::Sparse_chain, OSM::Sub_sparse_matrix >.

◆ get_psc_labels()

template<typename CoefficientType , typename ComplexType >
virtual vector< vector< int > > CGAL::HDVF::Hdvf_duality< CoefficientType, ComplexType >::get_psc_labels ( ) const
virtual

Exports primary/secondary/critical labels of the current sub chain complex for vtk export.

The method exports the labels of every cells in each dimension.

Returns
A vector containing, for each dimension, the vector of labels by cell index.

Reimplemented from CGAL::HDVF::Hdvf_core< CoefficientType, ComplexType, OSM::Sparse_chain, OSM::Sub_sparse_matrix >.

◆ print_bnd_pairing()

template<typename CoefficientType , typename ComplexType >
ostream & CGAL::HDVF::Hdvf_duality< CoefficientType, ComplexType >::print_bnd_pairing ( ostream &  out = cout)

Prints the reduced boundary over critical cells of K and L-K.

The method prints out the reduced boundary matrix in each dimension, restricted to critical cells of K and L-K (ie. the matrix used to compute Alexander pairing).

Warning
Call this method after compute_perfect_hdvf().

By default, outputs the complex to std::cout.

◆ print_reduction()

template<typename CoefficientType , typename ComplexType >
virtual ostream & CGAL::HDVF::Hdvf_duality< CoefficientType, ComplexType >::print_reduction ( ostream &  out = cout)
virtual

Prints the homology and cohomology reduction information for K and L-K.

Prints \(f^*\), \(g\) \(\partial'\) the reduced boundary over each critical cell.

By default, outputs the complex to std::cout.

◆ print_reduction_sub()

template<typename CoefficientType , typename ComplexType >
virtual ostream & CGAL::HDVF::Hdvf_duality< CoefficientType, ComplexType >::print_reduction_sub ( ostream &  out = cout)
virtual

Prints the homology and cohomology reduction information for the current such chain complex.

Prints \(f^*\), \(g\) \(\partial'\) the reduced boundary over each critical cell.

By default, outputs the complex to std::cout.

◆ set_mask_K()

template<typename CoefficientType , typename ComplexType >
void CGAL::HDVF::Hdvf_duality< CoefficientType, ComplexType >::set_mask_K ( )

Sets the current sub chain complex masks over K.

Further HDVF computations will be restricted to K (ie. computation of reduced homology).

◆ set_mask_L_K()

template<typename CoefficientType , typename ComplexType >
void CGAL::HDVF::Hdvf_duality< CoefficientType, ComplexType >::set_mask_L_K ( )

Sets the current sub chain complex masks over L-K.

Further HDVF computations will be restricted to L-K (ie. computation of reduced homology).