/* Copyright (c) 1997-2004
   Ewgenij Gawrilow, Michael Joswig (Technische Universitaet Berlin, Germany)
   http://www.math.tu-berlin.de/polymake,  mailto:polymake@math.tu-berlin.de

   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 2, or (at your option) any
   later version: http://www.gnu.org/licenses/gpl.txt.

   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.
*/

#ident "$Project: polymake $$Id: complex_tools.tcc 7477 2006-11-28 22:27:27Z witte $"

#include <ext/hash_map>

#include <FaceMap.h>
#include <Bitset.h>

namespace polymake { namespace topaz {

template <typename VertexSet, typename Complex>
HasseDiagram<VertexSet> hasse_diagram(const Complex& C, const int d, const int end_dim)
{
   // sort facets according to theire dimension
   Array< std::list<VertexSet> > facets_of_dim(d+1);
   int dim=-1;

   for (typename Entire<Complex>::const_iterator c_it=entire(C); !c_it.at_end(); ++c_it) {
      typename Complex::value_type f=*c_it;
      while (f.size() > facets_of_dim.size())
	 facets_of_dim.resize(2*facets_of_dim.size());
      if (f.size()-1 > dim)
	 dim = f.size()-1;

      facets_of_dim[f.size()-1].push_back(f);
   } 

   // fill the hasse diagram
   HasseDiagram<VertexSet> HD;
   typename HasseDiagram<VertexSet>::_filler HD_filler(HD);
   HD_filler.add_node(sequence(0,0));	        // node #0 is the top node representing the whole complex
   HD_filler.increase_dim();
   int face_end_this_dim=1, face_index=1;
  
   if (C.empty())                                // if c is empty, the HD is a single node only
      return HD;
  
   FaceMap<> FM;
   const int last_level= end_dim<0 ? dim + end_dim : end_dim;
   for (int d=dim; d>=last_level; --d) {

      // add facets of dimension d
      int i=face_end_this_dim;
      for (typename Entire< std::list<VertexSet> >::const_iterator f_it=entire(facets_of_dim[d]);
	   !f_it.at_end(); ++f_it, ++i) {
	 HD_filler.add_node(*f_it);
	 HD_filler.add_edge(i,0);
      }

      // add faces of dimension d which are not a facet
      if (d < dim) {
	 while (face_index < face_end_this_dim) {
	    const Subsets_less_1<VertexSet> faces_1_below(HD.node(face_index));
	    // add faces one below
	    for (typename Entire< Subsets_less_1<VertexSet> >::const_iterator fbi=entire(faces_1_below);
		 !fbi.at_end(); ++fbi) {
	       int &face = FM[*fbi];
	       if (face==-1) {
		  face=HD_filler.add_node(*fbi);
	       }
	       HD_filler.add_edge(face, face_index);
	    }
	    ++face_index;
	 }
      }

      face_end_this_dim=HD.nodes();
      HD_filler.increase_dim();
   }
  
   HD_filler.add_node(sequence(0,0));  // the bottom node (empty face)
   while (face_index < face_end_this_dim) {
      HD_filler.add_edge(face_end_this_dim, face_index);
      ++face_index;
   }
  
   return HD;
}

template <typename VertexSet, typename Complex>
HasseDiagram<VertexSet> pure_hasse_diagram(const Complex& C, const int end_dim)
{
  // works for pure complexes only
  HasseDiagram<VertexSet> HD;
  typename HasseDiagram<VertexSet>::_filler HD_filler(HD);
  HD_filler.add_node(sequence(0,0));	        // node #0 is the top node representing the whole complex
  HD_filler.increase_dim();
  
  if (C.empty())                                // if c is empty, the HD is a single node only
    return HD;
  
  const int dim=C.front().size()-1;
  
  // add facets of C
  HD_filler.add_nodes(C.size(), C.begin());
  HD_filler.increase_dim();
  int face_index=1, face_end_this_dim=C.size()+1;
  for (int i=1; i<face_end_this_dim; ++i)
    HD_filler.add_edge(i,0);
  
  // add faces of dimension d= dim-1,..0
  FaceMap<> FM;

  const int last_level= end_dim<0 ? dim + end_dim : end_dim;
  for (int d=dim-1; d>=last_level; --d) { 
    while (face_index < face_end_this_dim) {
      const Subsets_less_1<VertexSet> faces_1_below(HD.node(face_index));
      // add faces one below
      for (typename Entire< Subsets_less_1<VertexSet> >::const_iterator fbi=entire(faces_1_below);
	   !fbi.at_end(); ++fbi) {
	int &face = FM[*fbi];
	if (face==-1) {
	  face=HD_filler.add_node(*fbi);
	}
	HD_filler.add_edge(face, face_index);
      }
      ++face_index;
    }
    face_end_this_dim=HD.nodes();
    HD_filler.increase_dim();
  }
  
  HD_filler.add_node(sequence(0,0));  // the bottom node (empty face)
  while (face_index < face_end_this_dim) {
    HD_filler.add_edge(face_end_this_dim, face_index);
    ++face_index;
  }

  return HD;
}

template <typename VertexSet, typename _Set>
void remove_facet (HasseDiagram<VertexSet>& HD, const GenericSet< _Set,int >& F)
{
   // find the specified facet in the hasse diagram 
   int start_node = -1;
   for (Entire<sequence>::iterator it=entire(HD.nodes_of_dim(-1)); !it.at_end(); ++it)
      if (incl(HD.node(*it),F) == 0)
	 start_node = *it;
   

#ifdef POLY_DEBUG
   if (start_node == -1)
      throw std::runtime_error("remove_facet: Given set does not specify a facet.\ninvalid input");
#endif
   
   // remove all nodes of the HD that can be reached from start_node only
   Bitset visited(HD.nodes());
   std::list<int> node_queue;
   Set<int> nodes_2B_removed;
   
   node_queue.push_back(start_node);
   visited+=start_node;
   HD.out_edges(start_node).clear();

   while (!node_queue.empty()) {
      const int n = node_queue.front();    // remove first node from the queue
      node_queue.pop_front();
      
      if (HD.out_degree(n) == 0) {          // node is contained in F only
	 nodes_2B_removed += n;
	 for (typename Entire<typename HasseDiagram<VertexSet>::in_edge_list>::const_iterator e=entire(HD.in_edges(n));
	      !e.at_end(); ++e) {
	    const int nn=e.from_node();
	    if (!visited.contains(nn) && nn!=HD.bottom_node() ) {    // nn not been visited yet and nn != bottom node
	       visited += nn;                                        // -> add nn to the queue
	       node_queue.push_back(nn);
	    }
	 }
	 HD.in_edges(n).clear();
      }
   }
   
   // remove nodes to be removed
   HD.delete_nodes(nodes_2B_removed);
}
	
template <typename Complex, typename Set>
bool adj_numbering (Complex& C, const Set& V)
{
   if (V.empty())
      return false;

   const bool renumber= V.front()!=0 || V.back()+1!=V.size();

   if (renumber) {
      std_ext::hash_map<int, int> vertex_map(V.size());
      int count=0;
      for (typename Entire< Set >::const_iterator s_it=entire(V); !s_it.at_end(); ++s_it, ++count)
	 vertex_map[*s_it]=count;
      
      for (typename Entire<Complex>::iterator c_it=entire(C); !c_it.at_end(); ++c_it) {
	 typedef typename Complex::value_type Facet;
	 Facet f; 
	 for (typename Entire<Facet>::iterator s_it=entire(*c_it); !s_it.at_end(); ++s_it)
	    f += vertex_map[*s_it];
	 *c_it = f;
      }
   }

   return renumber;
}

template <typename VertexSet, typename OutputIterator>
bool is_pseudo_manifold (const HasseDiagram<VertexSet>& HD, OutputIterator boundary_consumer, const bool verbose) {
  if (HD.in_edges(HD.top_node()).empty())
    return true;

  int test_dim = -1;
  for (typename Entire<typename HasseDiagram<VertexSet>::in_edge_list>::const_iterator it=entire(HD.in_edges(HD.top_node()));
       !it.at_end(); ++it) {
      const int n = it.from_node();
    
      if (test_dim == -1)  // frist facet
	test_dim = HD.node(n).size()-1;

      if ( HD.node(n).size()-1 != test_dim ) {  // complex is not pure   
	if (verbose)
	  cout << "is_pseudo_manifold: complex is not PURE." << endl;
	return false;
      }
   }
  
   bool is_pmf=true;
   for (typename Entire<sequence>::const_iterator it=entire(HD.nodes_of_dim(-2)); !it.at_end(); ++it) {
      const int d = HD.out_degree(*it);
      if ( d > 2 ) {
	 if (!verbose)
	    return false;
	 cout << "is_pseudo_manifold: " << HD.node(*it) << " is contained in " << d << ">2 facets." << endl;
	 is_pmf=false;
      }
      if (!pm::derived_from_instance<OutputIterator, pm::black_hole>::value && d == 1 )
	 *boundary_consumer++ = HD.node(*it);
   }
   
   return is_pmf;
}

} }


// Local Variables:
// mode:C++
// c-basic-offset:3
// End: