/* Copyright (c) 1997-2007
   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.
*/

#ifndef _POLYMAKE_PERMUTATIONS_H
#define _POLYMAKE_PERMUTATIONS_H "$Project: polymake $$Id: permutations.h 7556 2007-01-12 17:36:36Z gawrilow $"

#include <GenericMatrix.h>
#include <Integer.h>
#include <Map.h>
#include <Set.h>
#include <Array.h>
#include <Bitset.h>

#include <algorithm>
#include <vector>
#include <sstream>

namespace pm {

template <typename Permutation>
int permutation_sign (const Permutation& v)
{
   const int l=v.size();
   if (l<=1) return 1;

   std::vector<int> w(v.size());
   std::copy(v.begin(), entire(w));

   int sign=1;
   for (int i=0; i<l; ) {
      if (w[i]!=i) {
	 int k=w[i];
	 w[i]=w[k];
	 w[k]=k;
	 sign = -sign;
      } else {
	 ++i;
      }
   }
   return sign;
}

template <typename Permutation, typename InvPermutation> inline
void inverse_permutation (const Permutation& perm, InvPermutation& inv_perm)
{
   inv_perm.resize(perm.size());
   int pos=0;
   for (typename Entire<Permutation>::const_iterator i=entire(perm);  !i.at_end();  ++i, ++pos)
      inv_perm[*i]=pos;
}

template <typename Iterator>
class permutation_cycle_iterator {
public:
   typedef forward_iterator_tag iterator_category;
   typedef int value_type;
   typedef const value_type& reference;
   typedef const value_type* pointer;
   typedef ptrdiff_t difference_type;
   typedef permutation_cycle_iterator iterator;
   typedef iterator const_iterator;

   permutation_cycle_iterator() { }

   permutation_cycle_iterator(Iterator pos_arg, value_type cur_arg, Bitset *visited_arg)
      : pos(pos_arg), cur(cur_arg), visited(visited_arg) { }

   permutation_cycle_iterator(typename alt_constructor<Iterator>::arg_type pos_arg, value_type cur_arg, Bitset *visited_arg)
      : pos(pos_arg), cur(cur_arg), visited(visited_arg) { }

   reference operator* () const { return cur; }
   pointer operator-> () const { return &cur; }

   iterator& operator++ ()
   {
      *visited+=cur;
      value_type prev=cur;
      cur=*pos;
      std::advance(pos, cur-prev);
      if (visited->contains(cur)) visited=0;
      return *this;
   }
   const iterator operator++ (int) { iterator copy=*this; operator++(); return copy; }

   bool operator== (const permutation_cycle_iterator& it) const { return pos==it.pos; }
   bool operator!= (const permutation_cycle_iterator& it) const { return !operator==(it); }

   bool at_end() const { return visited==0; }

protected:
   Iterator pos;
   value_type cur;
   Bitset *visited;
};

template <typename Permutation>
class permutation_cycles_iterator;

template <typename Permutation>
class PermutationCycle {
public:
   typedef typename Permutation::const_iterator perm_iterator;
   typedef permutation_cycle_iterator<perm_iterator> iterator;
   typedef iterator const_iterator;
   typedef int value_type;
   typedef typename iterator::reference reference;
   typedef reference const_reference;

   PermutationCycle(perm_iterator start_arg, value_type start_val_arg, int size_arg)
      : start(start_arg), start_val(start_val_arg), limit(size_arg), visited(size_arg)
   {
      valid_position();
   }

   iterator begin() const { return iterator(start,start_val,&visited); }
   iterator end() const { return iterator(start,start_val,0); }

   bool empty() const { return false; }
   int size() const { return count_it(begin()); }
   reference front() const { return start_val; }

protected:
   perm_iterator start;
   value_type start_val, limit;
   mutable Bitset visited;

   void incr() { ++start;  ++start_val; }

   void valid_position()
   {
      while (start_val<limit && (start_val==*start || visited.contains(start_val)))
	 incr();
   }

   friend class permutation_cycles_iterator<Permutation>;
};

template <typename Permutation>
class permutation_cycles_iterator {
public:
   typedef forward_iterator_tag iterator_category;
   typedef PermutationCycle<Permutation> value_type;
   typedef const value_type& reference;
   typedef const value_type* pointer;
   typedef ptrdiff_t difference_type;
   typedef permutation_cycles_iterator iterator;
   typedef iterator const_iterator;

   permutation_cycles_iterator() { }
   permutation_cycles_iterator(const Permutation& perm)
      : value(perm.begin(), 0, perm.size()) { }
   permutation_cycles_iterator(const Permutation& perm, bool _end)
      : value(perm.end(), perm.size(), perm.size()) { }

   reference operator* () const { return value; }
   pointer operator-> () const { return &value; }

   iterator& operator++ ()
   {
      value.incr();
      value.valid_position();
      return *this;
   }
   const iterator operator++ (int) { iterator copy=*this; operator++(); return copy; }

   bool operator== (const permutation_cycles_iterator& it) const { return value.start==it.value.start; }
   bool operator!= (const permutation_cycles_iterator& it) const { return !operator==(it); }

   bool at_end() const { return value.start_val>=value.limit; }
protected:
   value_type value;
};

template <typename Permutation>
class PermutationCycles {
protected:
   PermutationCycles();
   ~PermutationCycles();

   const Permutation& hidden() const { return reinterpret_cast<const Permutation&>(*this); }
public:
   typedef permutation_cycles_iterator<Permutation> iterator;
   typedef iterator const_iterator;
   typedef typename iterator::value_type value_type;
   typedef typename iterator::reference reference;
   typedef reference const_reference;

   iterator begin() const { return hidden(); }
   iterator end() const { return iterator(hidden(), true); }
   value_type front() const { return value_type(hidden().begin(), 0, hidden().size()); }

   bool empty() const { return begin().at_end(); }
   int size() const { return count_it(begin()); }
};

template <typename Permutation> inline
const PermutationCycles<Permutation>& permutation_cycles(const Permutation& p)
{
   return reinterpret_cast<const PermutationCycles<Permutation>&>(p);
}

template <typename Permutation>
struct check_iterator_feature<permutation_cycle_iterator<Permutation>, end_sensitive> : True { };

template <typename Permutation>
struct check_iterator_feature<permutation_cycles_iterator<Permutation>, end_sensitive> : True { };

template <typename Permutation>
struct spec_object_traits< PermutationCycle<Permutation> >
   : spec_object_traits<is_container> {
   static const bool is_temporary=true, is_always_const=true;
   static const IO_separator_kind IO_separator=IO_sep_inherit;
};

template <typename Permutation>
struct spec_object_traits< PermutationCycles<Permutation> >
   : spec_object_traits<is_container> {
   static const bool is_always_const=true;
   static const IO_separator_kind IO_separator=IO_sep_inherit;
};

template <typename Input1, typename Input2, typename Output, typename Comparator>
void find_permutation (Input1 src1, Input2 src2, Output dst, const Comparator& comparator)
{
   typedef Map<typename iterator_traits<Input1>::value_type, int, Comparator> map_type;
   map_type index_map(comparator);

   for (int i=0; !src1.at_end(); ++src1, ++i)
      index_map[*src1]=i;

   for (; !src2.at_end(); ++src2, ++dst) {
      typename map_type::iterator where=index_map.find(*src2);
      if (where.at_end()) {
	 std::string reason;
	 if (index_map.empty()) {
	    reason="not a permutation: first sequence is shorter";
	 } else {
	    std::ostringstream os;
	    wrap(os) << "not a permutation: no match for <" << *src2 << ">";
	    reason=os.str();
	 }
	 throw no_match(reason);
      }
      *dst=where->second;
      index_map.erase(where);
   }

   if (!index_map.empty())
      throw no_match("not a permutation: second sequence is shorter");
}

template <typename Container1, typename Container2, typename Comparator> inline
Array<int> find_permutation (const Container1& c1, const Container2& c2, const Comparator& comparator)
{
   Array<int> perm(c1.size());
   find_permutation(entire(c1), entire(c2), perm.begin(), comparator);
   return perm;
}

template <typename Container1, typename Container2> inline
Array<int> find_permutation (const Container1& c1, const Container2& c2)
{
   return find_permutation(c1,c2,polymake::operations::cmp());
}


enum permutation_sequence {
   permutations_heap,
   permutations_adjacent,
   permutations_lex
};

template <permutation_sequence kind=permutations_heap>
class AllPermutations;

template <permutation_sequence kind>
class permutation_iterator;

class permutation_iterator_base {
public:
   typedef forward_iterator_tag iterator_category;
   typedef std::vector<int> value_type;
   typedef const value_type* pointer;
   typedef const value_type& reference;
   typedef ptrdiff_t difference_type;
protected:
   value_type perm;

   permutation_iterator_base() { }
   permutation_iterator_base(int n) : perm(n) { reset(n); }
   void reset(int n) { std::copy(entire(sequence(0,n)), perm.begin()); }

public:
   reference operator* () const { return perm; }
   pointer operator-> () const { return &perm; }
};

/// Implementation of the Heap's algorithm by R. Sedgewick
template <>
class permutation_iterator<permutations_heap>
   : public permutation_iterator_base {
   friend class AllPermutations<permutations_heap>;
public:
   typedef permutation_iterator iterator;
   typedef iterator const_iterator;

protected:
   std::vector<int> cnt;
   int n, pos;

   permutation_iterator(int n_arg, bool)
      : permutation_iterator_base(n_arg),
	cnt(n_arg,0), n(n_arg), pos(n_arg) { }

public:
   permutation_iterator()
      : n(0), pos(0) { }

   explicit permutation_iterator(int n_arg)
      : permutation_iterator_base(n_arg),
	cnt(size_t(n_arg),0),
	n(n_arg), pos(n_arg>1) { }

   iterator& operator++ ()
   {
      do {
	 if (cnt[pos]<pos) {
	    std::swap(perm[pos], perm[pos%2 * cnt[pos]]);
	    ++cnt[pos];
	    pos=1;
	    break;
	 }
	 cnt[pos]=0;
      } while (++pos < n);
      return *this;
   }
   const iterator operator++ (int) { iterator copy(*this); operator++(); return copy; }

   bool at_end() const { return pos>=n; }

   bool operator== (const iterator& it) const
   {
      return n==it.n && pos==it.pos && std::equal(entire(cnt), it.cnt.begin());
   }
   bool operator!= (const iterator& it) const
   {
      return !operator==(it);
   }

   void rewind()
   {
      reset(n);
      pos=1;
      std::fill(entire(cnt), 0);
   }
};

template <>
class permutation_iterator<permutations_adjacent>
   : public permutation_iterator_base {
   friend class AllPermutations<permutations_adjacent>;
public:
   typedef permutation_iterator iterator;
   typedef iterator const_iterator;

protected:
   std::vector<int> move;
   int n, next;

   void incr()
   {
      if (next==0) {
	 next=-1; return;
      }
      int lim=n, offset=0;
      std::vector<int>::iterator pos=move.begin();
      for (;;) {
	 int m=*pos;
	 if (m>0) {
	    if (m<lim) {
	       next=m+offset;
	       ++*pos;
	       break;
	    }
	 } else {
	    if (m<-1) {
	       next=offset-m-1;
	       ++*pos;
	       break;
	    }
	    ++offset;
	 }
	 if (--lim > 1) {
	    *pos=-m;
	    ++pos;
	 } else {
	    next=0;
	    break;
	 }
      }
   }

public:
   permutation_iterator() : n(0) { }

   explicit permutation_iterator(int n_arg)
      : permutation_iterator_base(n_arg),
	move(size_t(n_arg),1), n(n_arg), next(n-1)
   {
      if (next>0) {
	 next=1; move[0]=2;
      }
   }

   iterator& operator++ ()
   {
      std::swap(perm[next-1], perm[next]);
      incr();
      return *this;
   }
   const iterator operator++ (int) { iterator copy(*this); operator++(); return copy; }

   int next_swap() const { return next; }

   bool at_end() const { return next<0; }

   bool operator== (const iterator& it) const
   {
      return next==it.next && n==it.n && (!next || std::equal(entire(move), it.move.begin()));
   }
   bool operator!= (const iterator& it) const
   {
      return !operator==(it);
   }

   void rewind()
   {
      reset(n);
      std::fill(entire(move),1);
      next=n-1;
      if (next>1) {
	 next=1; move[0]=2;
      }
   }
};

// permutations in lexicographic order via STL
template <>
class permutation_iterator<permutations_lex>
   : public permutation_iterator_base {
   friend class AllPermutations<permutations_lex>;
public:
   typedef permutation_iterator iterator;
   typedef iterator const_iterator;

protected:
   bool _at_end;

public:
   permutation_iterator() { }

   explicit permutation_iterator(int n, bool end_arg=false)
      : permutation_iterator_base(n), _at_end(end_arg || n==0) { }

   iterator& operator++ ()
   {
      _at_end=!std::next_permutation(perm.begin(), perm.end());
      return *this;
   }
   const iterator operator++ (int) { iterator copy(*this); operator++(); return copy; }

   bool at_end() const { return _at_end; }

   bool operator== (const iterator& it) const
   {
      return _at_end ? it._at_end : !it._at_end && perm.size()==it.perm.size() && std::equal(entire(perm), it.perm.begin());
   }
   bool operator!= (const iterator& it) const
   {
      return !operator==(it);
   }

   void rewind()
   {
      reset(perm.size());
      _at_end=perm.empty();
   }
};

template <permutation_sequence kind>
struct check_iterator_feature<permutation_iterator<kind>, end_sensitive> : True { };
template <permutation_sequence kind>
struct check_iterator_feature<permutation_iterator<kind>, rewindable> : True { };

template <permutation_sequence kind>
class AllPermutations {
protected:
   int n;
public:
   explicit AllPermutations(int n_arg=0) : n(n_arg) { }

   typedef permutation_iterator<kind> iterator;
   typedef iterator const_iterator;
   typedef typename iterator::value_type value_type;
   typedef typename iterator::reference reference;
   typedef reference const_reference;

   iterator begin() const { return iterator(n); }
   iterator end() const { return iterator(n,true); }

   value_type front() const
   {
      return value_type(sequence(0,n).begin(), sequence(0,n).end());
   }

   Integer size() const { return n ? Integer::fac(n) : Integer(0); }
   bool empty() const { return n==0; }
};

template <permutation_sequence kind>
struct spec_object_traits< AllPermutations<kind> >
   : spec_object_traits<is_container> {
   static const bool is_always_const=true;
};

inline
AllPermutations<> all_permutations(int n)
{
   return AllPermutations<>(n);
}

template <permutation_sequence kind>
AllPermutations<kind> all_permutations(int n)
{
   return AllPermutations<kind>(n);
}


template <typename PermutationRef, typename Element=int>
class PermutationMatrix
: public GenericMatrix< PermutationMatrix<PermutationRef,Element>, Element > {
protected:
   alias<PermutationRef> perm;
   mutable std::vector<int> inv_perm;
public:
   typedef Element value_type;
   typedef const Element& reference;
   typedef reference const_reference;

   typedef typename alias<PermutationRef>::arg_type arg_type;
   PermutationMatrix(arg_type perm_arg) : perm(perm_arg) { }

   typename alias<PermutationRef>::const_reference get_perm() const { return perm; }

   const std::vector<int>& get_inv_perm() const
   {
      if (inv_perm.empty() && !get_perm().empty())
	 inverse_permutation(get_perm(), inv_perm);
      return inv_perm;
   }

   const Element& one() const
   {
      static const Element One(1);
      return One;
   }
};

template <typename PermutationRef, typename Element>
struct spec_object_traits< PermutationMatrix<PermutationRef,Element> >
   : spec_object_traits<is_container> {
   static const bool is_temporary=true, is_always_const=true;
};

template <typename PermutationRef, typename Element>
struct check_container_feature< PermutationMatrix<PermutationRef,Element>, pure_sparse> : True { };

template <typename PermutationRef, typename Element>
class matrix_random_access_methods< PermutationMatrix<PermutationRef,Element> > {
   typedef PermutationMatrix<PermutationRef,Element> master;
public:
   const Element& operator() (int i, int j) const
   {
      const master& me=*static_cast<const master*>(this);
      INDEX_CHECK((i<0 || i>=me.rows()), "operator() - row index out of range");
      INDEX_CHECK((j<0 || j>=me.rows()), "operator() - col index out of range");
      operations::clear<Element> zero;
      return me.get_perm()[i]==j ? me.one() : zero();
   }
};

template <typename PermutationRef, typename Element>
class Rows< PermutationMatrix<PermutationRef,Element> >
   : public modified_container_pair_impl< Rows< PermutationMatrix<PermutationRef,Element> >,
					  list( Container1< PermutationRef >,
						Container2< unlimited_constant_container<const Element&> >,
						Operation< UnitVector_factory<> >,
						MasqueradedTop ) > {
   typedef modified_container_pair_impl<Rows> _super;
protected:
   ~Rows();
public:
   const typename _super::container1& get_container1() const
   {
      return this->hidden().get_perm();
   }
   const typename _super::container2& get_container2() const
   {
      return constant(this->hidden().one());
   }
   typename _super::operation get_operation() const
   {
      return this->size();
   }
};

template <typename PermutationRef, typename Element>
class Cols< PermutationMatrix<PermutationRef,Element> >
   : public modified_container_pair_impl< Cols< PermutationMatrix<PermutationRef,Element> >,
					  list( Container1< std::vector<int> >,
						Container2< unlimited_constant_container<const Element&> >,
						Operation< UnitVector_factory<> >,
						MasqueradedTop ) > {
   typedef modified_container_pair_impl<Cols> _super;
protected:
   ~Cols();
public:
   const typename _super::container1& get_container1() const
   {
      return this->hidden().get_inv_perm();
   }
   const typename _super::container2& get_container2() const
   {
      return constant(this->hidden().one());
   }
   typename _super::operation get_operation() const
   {
      return this->size();
   }
};

template <typename Element, typename Permutation> inline
PermutationMatrix<const Permutation&, Element>
permutation_matrix(const Permutation& perm)
{
   return perm;
}

} // end namespace pm

namespace polymake {
   using pm::permutation_sign;
   using pm::inverse_permutation;
   using pm::permutation_cycles;
   using pm::find_permutation;
   using pm::permutation_sequence;
   using pm::permutations_heap;
   using pm::permutations_adjacent;
   using pm::permutations_lex;
   using pm::AllPermutations;
   using pm::all_permutations;
   using pm::PermutationMatrix;
   using pm::permutation_matrix;
}

#endif // _POLYMAKE_PERMUTATIONS_H

// Local Variables:
// mode:C++
// End:
