/* Copyright (c) 1997-2006
   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_GENERIC_INCIDENCE_MATRIX_H
#define _POLYMAKE_GENERIC_INCIDENCE_MATRIX_H "$Project: polymake $$Id: GenericIncidenceMatrix.h,v 1.21 2006/03/09 16:21:35 gawrilow Exp $"

#ifndef _POLYMAKE_MATRIX_METHODS_H
#  include <matrix_methods.h>
#endif
#ifndef _POLYMAKE_GENERIC_IO_H
#  include <GenericIO.h>
#endif

namespace pm {

/* ------------------------
 *  GenericIncidenceMatrix
 * ------------------------ */

template <bool symmetric=false> class IncidenceMatrix;

template <typename Matrix, typename _sym_discr=typename Matrix::sym_discr>
class GenericIncidenceMatrix
   : public Generic<Matrix>,  public matrix_methods<Matrix, bool>,  public operators::base {
protected:
   GenericIncidenceMatrix() { }
   GenericIncidenceMatrix(const GenericIncidenceMatrix&) { }
public:
   typedef bool element_type;
   typedef _sym_discr sym_discr;
   static const bool is_symmetric=sym_discr::value,
                     prefer_concat=false;
   typedef IncidenceMatrix<is_symmetric> persistent_type;
   typedef GenericIncidenceMatrix generic_type;
   typedef bool value_type;
   typedef bool reference;
   typedef const bool const_reference;
protected:
   template <typename Matrix2>
   void assign(const GenericIncidenceMatrix<Matrix2>& M, True*)
   {
      DIMENSION_CHECK((M.rows() != this->rows() || M.cols() != this->cols()), "operator= - matrix dimension mismatch");
      std::copy(pm::rows(M).begin(), entire(pm::rows(*this)));
   }

   template <typename Matrix2>
   bool trivial_assignment(const GenericIncidenceMatrix<Matrix2,sym_discr>&) const { return false; }

   bool trivial_assignment(const GenericIncidenceMatrix& M) const { return this==&M; }

public:
   Matrix& operator= (const GenericIncidenceMatrix& M)
   {
      if (!trivial_assignment(M))
	 this->top().assign(M, (True*)0);
      return this->top();
   }

   template <typename Matrix2>
   Matrix& operator= (const GenericIncidenceMatrix<Matrix2>& M)
   {
      this->top().assign(M, (bool2type<!is_symmetric || Matrix2::is_symmetric>*)0);
      return this->top();
   }

   template <typename Matrix2>
   void swap(GenericIncidenceMatrix<Matrix2,sym_discr>& M)
   {
      if (trivial_assignment(M)) return;
      DIMENSION_CHECK((this->rows() != M.rows() || this->cols() != M.cols()), "swap - matrix dimension mismatch");
      std::swap_ranges(entire(pm::rows(*this)), pm::rows(M).begin());
   }

   Matrix& transitive_closure();

   template <typename Result>
   struct rebind_generic {
      typedef GenericIncidenceMatrix<Result, sym_discr> type;
   };
};

struct is_incidence_matrix;

template <typename Matrix, typename _sym_discr>
struct spec_object_traits< GenericIncidenceMatrix<Matrix,_sym_discr> >
   : spec_or_model_traits<Matrix,is_container> {
   static const int dimension=2;
   typedef is_incidence_matrix generic_tag;

   typedef Rows<Matrix> serialized;
   static serialized& serialize(Matrix& M) { return rows(M); }
   static const serialized& serialize(const Matrix& M) { return rows(M); }
};

template <typename Matrix, typename _sym_discr>
Matrix& GenericIncidenceMatrix<Matrix,_sym_discr>::transitive_closure()
{
   DIMENSION_CHECK(this->rows() != this->cols(), "transitive_closure - matrix not square");
   bool changed;
   do {
      changed=false;
      for (typename Entire< Rows<Matrix> >::iterator r_i=entire(pm::rows(this->top()));  !r_i.at_end();  ++r_i) {
	 const int ones=r_i->size();
	 for (typename Matrix::row_type::iterator e=r_i->begin(); !e.at_end(); ++e)
	    if (e.index() != r_i.index())
	       (*r_i) += this->row(e.index());
	 if (r_i->size() > ones) changed=true;
      }
   } while (changed);
   return this->top();
}

/* ----------------------------------------
 *  SingleIncidenceRow, SingleIncidenceCol
 * ---------------------------------------- */

template <typename _Set>
class SingleIncidenceRow
   : public GenericIncidenceMatrix< SingleIncidenceRow<_Set>, False> { };

template <typename _Set>
class SingleIncidenceCol
   : public GenericIncidenceMatrix< SingleIncidenceCol<_Set>, False> { };

class SingleElementIncidenceLine
   : public modified_container_impl< SingleElementIncidenceLine,
	       cons< Hidden< single_elem_container<int&, true> >,
		     Accessor< BuildUnaryIt<operations::index2element> > >
            >,
     public GenericSet< SingleElementIncidenceLine, int, operations::cmp<> > { };


template <typename Iterator1Ref, typename Iterator2Ref>
struct SingleElementIncidenceLine_factory {
   typedef Iterator1Ref first_argument_type;
   typedef Iterator2Ref second_argument_type;
   typedef SingleElementIncidenceLine line;
   typedef const line& result_type;

   result_type operator() (Iterator1Ref, Iterator2Ref) const
   {
      return reinterpret_cast<result_type>(*this);
   }
   template <typename Iterator1, typename Iterator2>
   result_type operator() (operations::partial, const Iterator1&, const Iterator2&) const
   {
      return reinterpret_cast<result_type>(line::hidden_type::create_empty());
   }
};

template <typename _Set>
class IncidenceLines_across
   : public modified_container_pair_impl< IncidenceLines_across<_Set>,
               cons< Container1< sequence >,
               cons< Container2< _Set >,
               cons< Coupler< zipping_coupler<operations::cmp<>, set_union_zipper> >,
               cons< Accessor< BuildBinaryIt<SingleElementIncidenceLine_factory> >,
                     Hidden< _Set > > > > >
            > {
   typedef modified_container_pair_impl<IncidenceLines_across> _super;
public:
   sequence get_container1() const { return sequence(0, size()); }
   int size() const { return this->hidden().dim(); }
};

template <typename _Set>
class Rows< SingleIncidenceRow<_Set> > : public single_elem_container<const _Set&> {
protected:
   ~Rows();
};

template <typename _Set>
class Cols< SingleIncidenceRow<_Set> > : public IncidenceLines_across<_Set> {
protected:
   ~Cols();
};

template <typename _Set>
class Rows< SingleIncidenceCol<_Set> > : public IncidenceLines_across<_Set> {
protected:
   ~Rows();
};

template <typename _Set>
class Cols< SingleIncidenceCol<_Set> > : public single_elem_container<const _Set&> {
protected:
   ~Cols();
};

template <typename ContainerRef1, typename ContainerRef2>
class IncidenceLineChain
   : public modified_container_impl< IncidenceLineChain<ContainerRef1,ContainerRef2>,
	       cons< Container< ContainerChain<ContainerRef1,ContainerRef2> >,
	       cons< Accessor< BuildUnaryIt<operations::index2element> >,
		     ExpectedFeatures<indexed> > >
            >,
     public GenericSet< IncidenceLineChain<ContainerRef1,ContainerRef2>, int, operations::cmp<> > {
   typedef modified_container_impl<IncidenceLineChain> _super;
protected:
   typename _super::container chain;
public:
   const typename _super::container& get_container() const { return chain; }

   IncidenceLineChain(typename _super::container::first_arg_type arg1,
		      typename _super::container::second_arg_type arg2)
      : chain(arg1,arg2) { }

   int dim() const
   {
      return chain.get_container1().dim() + chain.get_container2().dim();
   }
};

template <typename ContainerRef1, typename ContainerRef2>
struct check_container_feature<IncidenceLineChain<ContainerRef1, ContainerRef2>, sparse_compatible> : True { };

template <typename ContainerRef1, typename ContainerRef2>
struct spec_object_traits< IncidenceLineChain<ContainerRef1, ContainerRef2> >
   : spec_object_traits<is_container> {
   static const bool is_temporary=true, is_always_const=true;
};

/* ----------------------------
 *  SameElementIncidenceMatrix
 * ---------------------------- */

template <bool elem>
class SameElementIncidenceMatrix
   : public GenericIncidenceMatrix<SameElementIncidenceMatrix<elem>, False> {
public:
   SameElementIncidenceMatrix(int r_arg, int c_arg) : r(r_arg), c(c_arg) { }
protected:
   int r,c;

   friend class Rows<SameElementIncidenceMatrix>;
   friend class Cols<SameElementIncidenceMatrix>;
};

template <bool elem>
class SameElementIncidenceLine
   : public modified_container_impl< SameElementIncidenceLine<elem>,
               cons< Container< sequence >,
	       cons< Accessor< pair<nothing, operations::identity<int> > >,
                     ExpectedFeatures<end_sensitive> > >
            >,
     public GenericSet< SameElementIncidenceLine<elem>, int, operations::cmp<> > {
public:
   int dim() const { return reinterpret_cast<const int&>(*this); }
   int size() const { return elem ? dim() : 0; }
   bool empty() const { return !elem; }

   sequence get_container() const { return sequence(0, size()); }
};

template <bool elem>
struct check_container_feature< SameElementIncidenceLine<elem>, sparse_compatible> : True { };

template <bool elem>
class Rows< SameElementIncidenceMatrix<elem> >
   : public constant_container< Rows< SameElementIncidenceMatrix<elem> >, SameElementIncidenceLine<elem> > {
   typedef constant_container<Rows> _super;
protected:
   ~Rows();
public:
   const SameElementIncidenceMatrix<elem>& hidden() const
   {
      return reinterpret_cast<const SameElementIncidenceMatrix<elem>&>(*this);
   }
   typename _super::const_reference front() const
   {
      return reinterpret_cast<typename _super::const_reference>(hidden().c);
   }
   int size() const { return hidden().r; }
   bool empty() const { return !size(); }
};

template <bool elem>
class Cols< SameElementIncidenceMatrix<elem> >
   : public constant_container< Cols< SameElementIncidenceMatrix<elem> >, SameElementIncidenceLine<elem> > {
   typedef constant_container<Cols> _super;
protected:
   ~Cols();
public:
   const SameElementIncidenceMatrix<elem>& hidden() const
   {
      return reinterpret_cast<const SameElementIncidenceMatrix<elem>&>(*this);
   }
   typename _super::const_reference front() const
   {
      return reinterpret_cast<typename _super::const_reference>(hidden().r);
   }
   int size() const { return hidden().c; }
   bool empty() const { return !size(); }
};

template <typename Matrix1, typename Matrix2> inline
const RowChain< ColChain<const Matrix1&, SameElementIncidenceMatrix<false> >,
                ColChain<SameElementIncidenceMatrix<false>, const Matrix2&> >
diag(const GenericIncidenceMatrix<Matrix1>& m1, const GenericIncidenceMatrix<Matrix2>& m2)
{
   return RowChain< ColChain<const Matrix1&, SameElementIncidenceMatrix<false> >,
                    ColChain<SameElementIncidenceMatrix<false>, const Matrix2&> >
             ( ColChain<const Matrix1&, SameElementIncidenceMatrix<false> >
	          ( m1.top(), SameElementIncidenceMatrix<false>(m1.rows(), m2.cols()) ),
	       ColChain<SameElementIncidenceMatrix<false>, const Matrix2&>
	          ( SameElementIncidenceMatrix<false>(m2.rows(), m1.cols()), m2.top() )
	     );
}

/// the same, but with corner blocks filled with 1
template <typename Matrix1, typename Matrix2> inline
const RowChain< ColChain<const Matrix1&, SameElementIncidenceMatrix<true> >,
                ColChain<SameElementIncidenceMatrix<true>, const Matrix2&> >
diag_1(const GenericIncidenceMatrix<Matrix1>& m1, const GenericIncidenceMatrix<Matrix2>& m2)
{
   return RowChain< ColChain<const Matrix1&, SameElementIncidenceMatrix<true> >,
                    ColChain<SameElementIncidenceMatrix<true>, const Matrix2&> >
             ( ColChain<const Matrix1&, SameElementIncidenceMatrix<true> >
	          ( m1.top(), SameElementIncidenceMatrix<true>(m1.rows(), m2.cols()) ),
	       ColChain<SameElementIncidenceMatrix<true>, const Matrix2&>
	          ( SameElementIncidenceMatrix<true>(m2.rows(), m1.cols()), m2.top() )
	     );
}

/* ------------
 *  complement
 * ------------ */

template <typename Matrix>
class ComplementIncidenceMatrix
   : public GenericIncidenceMatrix< ComplementIncidenceMatrix<Matrix>, typename Matrix::sym_discr > {
protected:
   ComplementIncidenceMatrix();
   ~ComplementIncidenceMatrix();
public:
   typedef Matrix base_type;
   const Matrix& get_base() const { return reinterpret_cast<const Matrix&>(*this); }
};

template <typename Matrix>
class matrix_random_access_methods< ComplementIncidenceMatrix<Matrix> > {
public:
   bool operator() (int i, int j) const
   {
      return ! reinterpret_cast<const Matrix&>(*this)(i,j);
   }
};

template <typename LineRef>
struct ComplementIncidenceLine_factory {
   typedef LineRef argument_type;
   typedef LazySet2<sequence, LineRef, set_difference_zipper> result_type;

   result_type operator() (typename function_argument<LineRef>::const_type l) const
   {
      return result_type(sequence(0,l.dim()), l);
   }
};

template <typename Matrix>
class Rows< ComplementIncidenceMatrix<Matrix> >
   : public modified_container_impl< Rows< ComplementIncidenceMatrix<Matrix> >,
			       cons< Hidden< Rows<Matrix> >,
				     Accessor< BuildUnary<ComplementIncidenceLine_factory> > > > {
protected:
   ~Rows();
};

template <typename Matrix>
class Cols< ComplementIncidenceMatrix<Matrix> >
   : public modified_container_impl< Cols< ComplementIncidenceMatrix<Matrix> >,
			       cons< Hidden< Cols<Matrix> >,
				     Accessor< BuildUnary<ComplementIncidenceLine_factory> > > > {
protected:
   ~Cols();
};

/* ------------------------------------
 *  operations building block matrices
 * ------------------------------------ */

namespace operations {

template <typename OpRef>
struct bitwise_inv_impl<OpRef, is_incidence_matrix>
   : reinterpret_impl<OpRef, const ComplementIncidenceMatrix<typename deref<OpRef>::type>&> { };

template <typename Matrix>
struct bitwise_inv_impl<const ComplementIncidenceMatrix<Matrix>&, is_incidence_matrix> {
   typedef const ComplementIncidenceMatrix<Matrix>& argument_type;
   typedef const Matrix& result_type;

   result_type operator() (argument_type x) const
   {
      return x.base();
   }
};

template <typename LeftRef, typename RightRef, bool check>
struct concat_impl<LeftRef, RightRef, check, cons<is_set, is_set> > {
   typedef LeftRef first_argument_type;
   typedef RightRef second_argument_type;
   typedef typename attrib<LeftRef>::plus_const Left;
   typedef typename attrib<RightRef>::plus_const Right;
   typedef subset_classifier::data_helper<Left, true, is_set> left_helper;
   typedef subset_classifier::data_helper<Right, true, is_set> right_helper;
   typedef IncidenceLineChain<typename left_helper::container_ref, typename right_helper::container_ref> chain;
   typedef const chain result_type;

   result_type operator() (typename function_argument<LeftRef>::const_type l,
			   typename function_argument<RightRef>::const_type r) const
   {
      return chain(left_helper::get(l), right_helper::get(r));
   }
};

template <typename LeftRef, typename RightRef, bool check>
struct div_impl<LeftRef, RightRef, check, cons<is_incidence_matrix, is_incidence_matrix> > {
   typedef LeftRef first_argument_type;
   typedef RightRef second_argument_type;
   typedef typename attrib<LeftRef>::plus_const Left;
   typedef typename attrib<RightRef>::plus_const Right;
   typedef RowChain<Left, Right> chain;
   typedef const chain result_type;

   result_type operator() (typename function_argument<LeftRef>::const_type l,
			   typename function_argument<RightRef>::const_type r) const
   {
      if (check) DIMENSION_CHECK((l.cols() && r.cols() && l.cols() != r.cols()),
				 "concat_rows - column dimensions mismatch");
      return chain(l,r);
   }

   typename function_argument<LeftRef>::lvalue
   assign(typename function_argument<LeftRef>::lvalue l,
	  typename function_argument<RightRef>::const_type r)
   {
      return static_cast<typename function_argument<LeftRef>::lvalue> (l/=r);
   }
};

template <typename LeftRef, typename RightRef, bool check>
struct div_impl<LeftRef, RightRef, check, cons<is_incidence_matrix, is_set> > {
   typedef LeftRef first_argument_type;
   typedef RightRef second_argument_type;
   typedef typename attrib<LeftRef>::plus_const Left;
   typedef typename attrib<RightRef>::plus_const Right;
   typedef Set_with_dim_helper<Right> right_helper;
   typedef RowChain<Left, pm::masquerade<SingleIncidenceRow, typename right_helper::container> > chain;
   typedef const chain result_type;

   result_type operator() (typename function_argument<LeftRef>::const_type l,
			   typename function_argument<RightRef>::const_type r)
   {
      if (check) DIMENSION_CHECK(!right_helper::check(r, l.cols()),
				 "concat_rows - column dimensions mismatch");
      return chain(l, right_helper::create(r, l.cols()));
   }

   typename function_argument<LeftRef>::lvalue
   assign(typename function_argument<LeftRef>::lvalue l,
	  typename function_argument<RightRef>::const_type r)
   {
      return static_cast<typename function_argument<LeftRef>::lvalue> (l/=r);
   }
};

template <typename LeftRef, typename RightRef, bool check>
struct div_impl<LeftRef, RightRef, check, cons<is_set, is_incidence_matrix> > {
   typedef LeftRef first_argument_type;
   typedef RightRef second_argument_type;
   typedef typename attrib<LeftRef>::plus_const Left;
   typedef typename attrib<RightRef>::plus_const Right;
   typedef Set_with_dim_helper<Left> left_helper;
   typedef RowChain<pm::masquerade<SingleIncidenceRow, typename left_helper::container>, Right> chain;
   typedef const chain result_type;

   result_type operator() (typename function_argument<LeftRef>::const_type l,
			   typename function_argument<RightRef>::const_type r)
   {
      if (check) DIMENSION_CHECK(!left_helper::check(l, r.cols()),
				 "concat_rows - column dimensions mismatch");
      return chain(left_helper::create(l, r.cols()), r);
   }
};

template <typename LeftRef, typename RightRef, bool check>
struct bitwise_or_impl<LeftRef, RightRef, check, cons<is_incidence_matrix, is_incidence_matrix> > {
   typedef LeftRef first_argument_type;
   typedef RightRef second_argument_type;
   typedef typename attrib<LeftRef>::plus_const Left;
   typedef typename attrib<RightRef>::plus_const Right;
   typedef ColChain<Left, Right> chain;
   typedef const chain result_type;

   result_type operator() (typename function_argument<LeftRef>::const_type l,
			   typename function_argument<RightRef>::const_type r) const
   {
      if (check) DIMENSION_CHECK((l.rows() && r.rows() && l.rows() != r.rows()),
				 "concat_cols - row dimensions mismatch");
      return chain(l,r);
   }

   typename function_argument<LeftRef>::lvalue
   assign(typename function_argument<LeftRef>::lvalue l,
	  typename function_argument<RightRef>::const_type r)
   {
      return static_cast<typename function_argument<LeftRef>::lvalue> (l|=r);
   }
};

template <typename LeftRef, typename RightRef, bool check>
struct bitwise_or_impl<LeftRef, RightRef, check, cons<is_incidence_matrix, is_set> > {
   typedef LeftRef first_argument_type;
   typedef RightRef second_argument_type;
   typedef typename attrib<LeftRef>::plus_const Left;
   typedef typename attrib<RightRef>::plus_const Right;
   typedef Set_with_dim_helper<Right> right_helper;
   typedef ColChain<Left, pm::masquerade<SingleIncidenceCol, typename right_helper::container> > chain;
   typedef const chain result_type;

   result_type operator() (typename function_argument<LeftRef>::const_type l,
			   typename function_argument<RightRef>::const_type r)
   {
      if (check) DIMENSION_CHECK(!right_helper::check(r, l.rows()),
				 "concat_cols - row dimensions mismatch");
      return chain(l, right_helper::create(r, l.rows()));
   }

   typename function_argument<LeftRef>::lvalue
   assign(typename function_argument<LeftRef>::lvalue l,
	  typename function_argument<RightRef>::const_type r)
   {
      return static_cast<typename function_argument<LeftRef>::lvalue> (l|=r);
   }
};

template <typename LeftRef, typename RightRef, bool check>
struct bitwise_or_impl<LeftRef, RightRef, check, cons<is_set, is_incidence_matrix> > {
   typedef LeftRef first_argument_type;
   typedef RightRef second_argument_type;
   typedef typename attrib<LeftRef>::plus_const Left;
   typedef typename attrib<RightRef>::plus_const Right;
   typedef Set_with_dim_helper<Left> left_helper;
   typedef ColChain<pm::masquerade<SingleIncidenceCol, typename left_helper::container>, Right> chain;
   typedef const chain result_type;

   result_type operator() (typename function_argument<LeftRef>::const_type l,
			   typename function_argument<RightRef>::const_type r)
   {
      if (check) DIMENSION_CHECK(!Set_with_dim_helper<LeftRef>::check(l, r.rows()),
				 "concat_cols - row dimensions mismatch");
      return chain(Set_with_dim_helper<LeftRef>::create(l, r.rows()), r);
   }
};

} // end namespace operations
namespace operators {

template <typename Matrix1, typename Matrix2> inline
bool operator== (const GenericIncidenceMatrix<Matrix1>& l, const GenericIncidenceMatrix<Matrix2>& r)
{
   DIMENSION_CHECK((l.rows() != r.rows() ||
		    l.cols() != r.cols()),  "operator== - matrix dimension mismatch");
   operations::eq<const Matrix1&, const Matrix2&> op;
   return op(l.top(), r.top());
}

template <typename Matrix1, typename Matrix2> inline
bool operator< (const GenericIncidenceMatrix<Matrix1>& l, const GenericIncidenceMatrix<Matrix2>& r)
{
   DIMENSION_CHECK((l.rows() != r.rows() ||
		    l.cols() != r.cols()),  "operator< - matrix dimension mismatch");
   operations::lt<const Matrix1&, const Matrix2&> op;
   return op(l.top(), r.top());
}

} // end namespace operators
} // end namespace pm

namespace polymake {
   using pm::GenericIncidenceMatrix;
}

namespace std {
   template <typename Matrix1, typename Matrix2, typename _sym_discr> inline
   void swap(pm::GenericIncidenceMatrix<Matrix1,_sym_discr>& M1, pm::GenericIncidenceMatrix<Matrix2,_sym_discr>& M2)
   {
      M1.top().swap(M2.top());
   }

   template <typename Matrix1, typename Matrix2, typename _sym_discr> inline
   void swap(const pm::GenericIncidenceMatrix<Matrix1,_sym_discr>& M1, const pm::GenericIncidenceMatrix<Matrix2,_sym_discr>& M2)
   {
      pm::non_const(M1.top()).swap(pm::non_const(M2.top()));
   }
}

#ifdef _POLYMAKE_GENERIC_MATRIX_H
#  include <Incidence_and_SparseMatrix.h>
#endif
#endif // _POLYMAKE_GENERIC_INCIDENCE_MATRIX_H

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