Skeleton_blocker_sub_complex.h

/*    This file is part of the Gudhi Library. The Gudhi library
 *    (Geometric Understanding in Higher Dimensions) is a generic C++
 *    library for computational topology.
 *
 *    Author(s):       David Salinas
 *
 *    Copyright (C) 2014 Inria
 *
 *    This program is free software: you can redistribute it and/or modify
 *    it under the terms of the GNU General Public License as published by
 *    the Free Software Foundation, either version 3 of the License, or
 *    (at your option) any later version.
 *
 *    This program is distributed in the hope that it will be useful,
 *    but WITHOUT ANY WARRANTY; without even the implied warranty of
 *    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *    GNU General Public License for more details.
 *
 *    You should have received a copy of the GNU General Public License
 *    along with this program.  If not, see <http://www.gnu.org/licenses/>.
 */

#ifndef SKELETON_BLOCKER_SKELETON_BLOCKER_SUB_COMPLEX_H_
#define SKELETON_BLOCKER_SKELETON_BLOCKER_SUB_COMPLEX_H_

#include <gudhi/Skeleton_blocker_complex.h>
#include <gudhi/Skeleton_blocker/Skeleton_blocker_simplex.h>
#include <gudhi/Debug_utils.h>

#include <map>
#include <vector>

namespace Gudhi {

namespace skeleton_blocker {

template<typename ComplexType>
class Skeleton_blocker_sub_complex : public ComplexType {
 protected:
  template<class T> friend class Skeleton_blocker_link_complex;

  typedef typename ComplexType::Graph Graph;
  typedef typename ComplexType::Edge_handle Edge_handle;

  typedef typename ComplexType::boost_vertex_handle boost_vertex_handle;

 public:
  using ComplexType::add_vertex;
  using ComplexType::add_edge_without_blockers;
  using ComplexType::add_blocker;

  typedef typename ComplexType::Vertex_handle Vertex_handle;
  typedef typename ComplexType::Root_vertex_handle Root_vertex_handle;
  typedef typename ComplexType::Simplex Simplex;
  typedef typename ComplexType::Root_simplex_handle Root_simplex_handle;

 protected:
  typedef std::map<Root_vertex_handle, Vertex_handle> IdAddressMap;
  typedef typename IdAddressMap::value_type AddressPair;
  typedef typename IdAddressMap::iterator IdAddressMapIterator;
  typedef typename IdAddressMap::const_iterator IdAddressMapConstIterator;
  std::map<Root_vertex_handle, Vertex_handle> adresses;

 public:
  Vertex_handle add_vertex(Root_vertex_handle global) {
    assert(!this->contains_vertex(global));
    Vertex_handle address(boost::add_vertex(this->skeleton));
    this->num_vertices_++;
    (*this)[address].activate();
    (*this)[address].set_id(global);
    adresses.insert(AddressPair(global, address));
    this->degree_.push_back(0);
    return address;
  }

  void add_edge_without_blockers(Root_vertex_handle v1_root, Root_vertex_handle v2_root) {
    auto v1_sub(this->get_address(v1_root));
    auto v2_sub(this->get_address(v2_root));
    assert(v1_sub && v2_sub);
    this->ComplexType::add_edge_without_blockers(*v1_sub, *v2_sub);
  }

  void add_blocker(const Root_simplex_handle& blocker_root) {
    auto blocker_sub = this->get_address(blocker_root);
    assert(blocker_sub);
    this->add_blocker(new Simplex(*blocker_sub));
  }

 public:
  void make_restricted_complex(const ComplexType & parent_complex,
                               const Simplex& simplex) {
    this->clear();
    // add vertices to the sub complex
    for (auto x : simplex) {
      assert(parent_complex.contains_vertex(x));
      auto x_local = this->add_vertex(parent_complex[x].get_id());
      (*this)[x_local] = parent_complex[x];
    }

    // add edges to the sub complex
    for (auto x : simplex) {
      // x_neigh is the neighbor of x intersected with vertices_simplex
      Simplex x_neigh;
      parent_complex.add_neighbours(x, x_neigh, true);
      x_neigh.intersection(simplex);
      for (auto y : x_neigh) {
        this->add_edge_without_blockers(parent_complex[x].get_id(), parent_complex[y].get_id());
      }
    }

    // add blockers to the sub complex
    for (auto blocker : parent_complex.const_blocker_range()) {
      // check if it is the first time we encounter the blocker
      if (simplex.contains(*blocker)) {
        Root_simplex_handle blocker_root(parent_complex.get_id(*(blocker)));
        Simplex blocker_restr(
                              *(this->get_simplex_address(blocker_root)));
        this->add_blocker(new Simplex(blocker_restr));
      }
    }
  }

  void clear() {
    adresses.clear();
    ComplexType::clear();
  }

  boost::optional<Vertex_handle> get_address(Root_vertex_handle global) const {
    boost::optional < Vertex_handle > res;
    IdAddressMapConstIterator it = adresses.find(global);
    if (it == adresses.end())
      res.reset();
    else
      res = (*it).second;
    return res;
  }

  // /**
  //  * Allocates a simplex in L corresponding to the simplex s in K
  //  * with its local adresses and returns an AddressSimplex.
  //  */
  // boost::optional<Simplex> get_address(const Root_simplex_handle & s) const;

  // private:
  // @remark should be private but problem with VS

  std::vector<boost::optional<Vertex_handle> > get_addresses(
                                                             const Root_simplex_handle & s) const {
    std::vector < boost::optional<Vertex_handle> > res;
    for (auto i : s) {
      res.push_back(get_address(i));
    }
    return res;
  }
};

template<typename ComplexType>
bool proper_face_in_union(
                          Skeleton_blocker_sub_complex<ComplexType> & link,
                          std::vector<boost::optional<typename ComplexType::Vertex_handle> > & addresses_sigma_in_link,
                          std::size_t vertex_to_be_ignored) {
  // we test that all vertices of 'addresses_sigma_in_link' but 'vertex_to_be_ignored'
  // are in link1 if it is the case we construct the corresponding simplex
  bool vertices_sigma_are_in_link = true;
  typename ComplexType::Simplex sigma_in_link;
  for (std::size_t i = 0; i < addresses_sigma_in_link.size(); ++i) {
    if (i != vertex_to_be_ignored) {
      if (!addresses_sigma_in_link[i]) {
        vertices_sigma_are_in_link = false;
        break;
      } else {
        sigma_in_link.add_vertex(*addresses_sigma_in_link[i]);
      }
    }
  }
  // If one of vertices of the simplex is not in the complex then it returns false
  // Otherwise, it tests if the simplex is in the complex
  return vertices_sigma_are_in_link && link.contains(sigma_in_link);
}

// Remark: this function should be friend in order to leave get_adresses private
// however doing so seemes currently not possible due to a visual studio bug c2668
// "the compiler does not support partial ordering of template functions as specified in the C++ Standard"
// http://www.serkey.com/error-c2668-ambiguous-call-to-overloaded-function-bb45ft.html

template<typename ComplexType>
bool proper_faces_in_union(
                           Skeleton_blocker_simplex<typename ComplexType::Root_vertex_handle> & sigma,
                           Skeleton_blocker_sub_complex<ComplexType> & link1,
                           Skeleton_blocker_sub_complex<ComplexType> & link2) {
  typedef typename ComplexType::Vertex_handle Vertex_handle;
  std::vector < boost::optional<Vertex_handle> > addresses_sigma_in_link1 =
      link1.get_addresses(sigma);
  std::vector < boost::optional<Vertex_handle> > addresses_sigma_in_link2 =
      link2.get_addresses(sigma);

  for (std::size_t current_index = 0; current_index < addresses_sigma_in_link1.size();
       ++current_index) {
    if (!proper_face_in_union(link1, addresses_sigma_in_link1, current_index)
        && !proper_face_in_union(link2, addresses_sigma_in_link2,
                                 current_index)) {
      return false;
    }
  }
  return true;
}

}  // namespace skeleton_blocker

namespace skbl = skeleton_blocker;

}  // namespace Gudhi

#endif  // SKELETON_BLOCKER_SKELETON_BLOCKER_SUB_COMPLEX_H_