chain_rep_cycles.h

Go to the documentation of this file.

/*    This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
 *    See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
 *    Author(s):       Hannah Schreiber
 *
 *    Copyright (C) 2022 Inria
 *
 *    Modification(s):
 *      - YYYY/MM Author: Description of the modification
 */

 
#ifndef PM_CHAIN_REP_CYCLES_H
#define PM_CHAIN_REP_CYCLES_H
 
#include <vector>
 
#ifdef GUDHI_USE_TBB
#include <tbb/parallel_for.h>
#endif
 
#include <gudhi/Debug_utils.h>
#include <gudhi/persistence_matrix_options.h>
 
namespace Gudhi {
namespace persistence_matrix {

struct Dummy_chain_representative_cycles {
  friend void swap([[maybe_unused]] Dummy_chain_representative_cycles& d1,
                   [[maybe_unused]] Dummy_chain_representative_cycles& d2) noexcept
  {}
};

template <class Master_matrix>
class Chain_representative_cycles
{
 public:
  using Bar = typename Master_matrix::Bar;                           
  using Cycle = typename Master_matrix::Cycle;                       
  using Column_container = typename Master_matrix::Column_container; 
  using Index = typename Master_matrix::Index;                       
  using ID_index = typename Master_matrix::ID_index;                 
  using Dimension = typename Master_matrix::Dimension;               

  Chain_representative_cycles() = default;

  void update_all_representative_cycles(Dimension dim = Master_matrix::template get_null_value<Dimension>());

  void update_representative_cycle(const Bar& bar);

  const std::vector<Cycle>& get_all_representative_cycles() const;
  const Cycle& get_representative_cycle(const Bar& bar) const;

  friend void swap(Chain_representative_cycles& base1, Chain_representative_cycles& base2) noexcept
  {
    base1.representativeCycles_.swap(base2.representativeCycles_);
    base1.birthToCycle_.swap(base2.birthToCycle_);
  }
 
 protected:
  void _reset();
 
 private:
  using Master_chain_matrix = typename Master_matrix::Master_chain_matrix;
 
  std::vector<Cycle> representativeCycles_; 
  std::vector<Index> birthToCycle_;         
 
  // access to inheriting Chain_matrix class
  constexpr Master_chain_matrix* _matrix() { return static_cast<Master_chain_matrix*>(this); }
 
  constexpr const Master_chain_matrix* _matrix() const { return static_cast<const Master_chain_matrix*>(this); }
};
 
template <class Master_matrix>
inline void Chain_representative_cycles<Master_matrix>::update_all_representative_cycles(Dimension dim)
{
  Index nberColumns = _matrix()->get_number_of_columns();
  Index nullValue = Master_matrix::template get_null_value<Index>();
 
  auto get_position = [&](ID_index pivot) {
    if constexpr (Master_matrix::Option_list::has_vine_update) {
      if constexpr (Master_matrix::Option_list::has_map_column_container) {
        return _matrix()->map_.at(pivot);
      } else {
        return _matrix()->map_[pivot];
      }
    } else {
      return pivot;
    }
  };
  auto get_col_index = [&](ID_index pivot) -> Index {
    if constexpr (Master_matrix::Option_list::has_map_column_container) {
      auto it = _matrix()->pivotToColumnIndex_.find(pivot);
      return it == _matrix()->pivotToColumnIndex_.end() ? nullValue : it->second;
    } else {
      return _matrix()->pivotToColumnIndex_[pivot];
    }
  };
 
  birthToCycle_.clear();
  birthToCycle_.resize(nberColumns, nullValue);
  representativeCycles_.clear();
 
#ifdef GUDHI_USE_TBB
  Index c = 0;
  for (Index i = 0; i < nberColumns; i++) {
    Index colIdx = get_col_index(i);
    if (colIdx != nullValue) {
      auto& col = _matrix()->get_column(colIdx);
      if ((dim == Master_matrix::template get_null_value<Dimension>() || col.get_dimension() == dim) &&
          (!col.is_paired() || get_position(i) < get_position(_matrix()->get_pivot(col.get_paired_chain_index())))) {
        birthToCycle_[get_position(i)] = c;
        ++c;
      }
    }
  }
 
  representativeCycles_.resize(c);
  tbb::parallel_for(static_cast<Index>(0), nberColumns, [&](Index i) {
    Index colIdx = get_col_index(i);
    // condition order matters
    if (colIdx != nullValue && birthToCycle_[get_position(i)] != nullValue) {
      auto& col = _matrix()->get_column(colIdx);
      representativeCycles_[birthToCycle_[get_position(i)]] =
          Master_matrix::build_cycle_from_range(col.get_non_zero_content_range());
    }
  });
#else
  for (ID_index i = 0; i < nberColumns; i++) {
    Index colIdx = get_col_index(i);
    if (colIdx != nullValue) {
      auto& col = _matrix()->get_column(colIdx);
      if ((dim == Master_matrix::template get_null_value<Dimension>() || col.get_dimension() == dim) &&
          (!col.is_paired() || get_position(i) < get_position(_matrix()->get_pivot(col.get_paired_chain_index())))) {
        representativeCycles_.push_back(Master_matrix::build_cycle_from_range(col.get_non_zero_content_range()));
        birthToCycle_[get_position(i)] = representativeCycles_.size() - 1;
      }
    }
  }
#endif
}
 
template <class Master_matrix>
inline void Chain_representative_cycles<Master_matrix>::update_representative_cycle(const Bar& bar)
{
  auto nberColumns = _matrix()->get_number_of_columns();
  auto get_position = [&](ID_index pivot) {
    if constexpr (Master_matrix::Option_list::has_vine_update) {
      if constexpr (Master_matrix::Option_list::has_map_column_container) {
        return _matrix()->map_.at(pivot);
      } else {
        return _matrix()->map_[pivot];
      }
    } else {
      return pivot;
    }
  };
 
  Index nullValue = Master_matrix::template get_null_value<Index>();
 
  if (birthToCycle_.size() <= bar.birth) {
    birthToCycle_.resize(bar.birth + 1, nullValue);
  }
  if (birthToCycle_[bar.birth] == nullValue) {
    birthToCycle_[bar.birth] = representativeCycles_.size();
    representativeCycles_.resize(representativeCycles_.size() + 1);
  }
 
  ID_index pivot = bar.birth;
 
  if constexpr (Master_matrix::Option_list::has_vine_update) {
    ID_index i = 0;
    while (i < nberColumns && get_position(i) != bar.birth) ++i;
    GUDHI_CHECK(i < nberColumns, std::invalid_argument("Given birth value is not valid."));
    pivot = i;
  }
 
  // if bar.birth does not correspond to a valid birth time, get_column_with_pivot will throw an out of range.
  auto& col = _matrix()->get_column(_matrix()->get_column_with_pivot(pivot));
  representativeCycles_[birthToCycle_[bar.birth]] =
      Master_matrix::build_cycle_from_range(col.get_non_zero_content_range());
}
 
template <class Master_matrix>
inline const std::vector<typename Chain_representative_cycles<Master_matrix>::Cycle>&
Chain_representative_cycles<Master_matrix>::get_all_representative_cycles() const
{
  return representativeCycles_;
}
 
template <class Master_matrix>
inline const typename Chain_representative_cycles<Master_matrix>::Cycle&
Chain_representative_cycles<Master_matrix>::get_representative_cycle(const Bar& bar) const
{
  return representativeCycles_[birthToCycle_[bar.birth]];
}
 
template <class Master_matrix>
inline void Chain_representative_cycles<Master_matrix>::_reset()
{
  representativeCycles_.clear();
  birthToCycle_.clear();
}
 
}  // namespace persistence_matrix
}  // namespace Gudhi
 
#endif  // PM_CHAIN_REP_CYCLES_H