/* Copyright (c) 1997-2010
   Ewgenij Gawrilow, Michael Joswig (Technische Universitaet Darmstadt, Germany)
   http://www.polymake.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.
--------------------------------------------------------------------------------
   $Project: polymake $$Id: intersection.cc 9590 2010-03-22 01:46:42Z gawrilow $
*/

#include <polymake.h>
#include <ListMatrix.h>
#include <Matrix.h>
#include <Vector.h>
#include <Array.h>

namespace polymake { namespace polytope {

   template <typename Scalar>
   perl::Object intersection(const Array<perl::Object>& pp_in)
   {
      Entire< Array<perl::Object> >::const_iterator p_in=entire(pp_in);
      ListMatrix< Vector<Scalar> > Inequalities, Equations;
      const int dim=p_in->give("AMBIENT_DIM");
      try {
	 p_in->give("FACETS | INEQUALITIES") >> Inequalities;
	 p_in->lookup("AFFINE_HULL | EQUATIONS") >> Equations;
      } catch (perl::undefined) {
	 p_in->give("AFFINE_HULL | EQUATIONS") >> Equations;
      }

      perl::Object p_out(p_in->type());
      std::string descr_names=p_in->name();

      while (! (++p_in).at_end()) {
	 const int dim2=p_in->give("AMBIENT_DIM");
	 if (dim != dim2) throw std::runtime_error("dimension mismatch");
	 try {
	    const Matrix<Scalar> F=p_in->give("FACETS | INEQUALITIES"),
	                        AH=p_in->lookup("AFFINE_HULL | EQUATIONS");
	    Inequalities /= F;
	    Equations /= AH;
	 } catch (perl::undefined) {
	    const Matrix<Scalar> AH=p_in->give("AFFINE_HULL | EQUATIONS");
	    Equations /= AH;
	 }
	 descr_names+=", ";
	 descr_names+=p_in->name();
      }

      p_out.set_description() << "Intersection of polytopes " << descr_names << endl;
      p_out.take("AMBIENT_DIM") << dim;
      p_out.take("INEQUALITIES") << Inequalities;
      if (Equations.rows())
	 p_out.take("EQUATIONS") << Equations;

      return p_out;
   }

   UserFunctionTemplate4perl("# @category Producing a new polyhedron from others"
			     "# Construct a new polyhedron as the intersection of given polyhedra."
			     "# @param Polytope P ... polyhedra to be intersected"
			     "# @return Polytope",
			     "intersection<Scalar>(Polytope<Scalar> +)");
} }
