//
// $Source: /cvsroot/gambit/gambit/sources/tools/enumpoly/ideal.imp,v $
// $Date: 2006/01/07 06:37:21 $
// $Revision: 1.5 $
//
// DESCRIPTION:
// Implementation of gIdeal
//
// This file is part of Gambit
// Copyright (c) 2002, The Gambit Project
//
// 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 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, write to the Free Software
// Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
//

#include "ideal.h"

//---------------------------------------------------------------
//                      gIdeal
//---------------------------------------------------------------

//---------------------------
// Constructors / Destructors
//---------------------------


template<class T>  gIdeal<T>::gIdeal(const gSpace * ls,
				     const term_order * ordr)
: Space(ls), order(ordr), basis(ls,ordr)
{
}

template<class T>  gIdeal<T>::gIdeal(const gSpace * ls, 
				     const term_order * ordr,
				     const Gambit::List< gPoly<T> *> & polys)
: Space(ls), order(ordr), basis(ls,ordr,polys)
{
}

template<class T>  gIdeal<T>::gIdeal(const term_order * ordr,
				     const gPolyList<T> & bss)
: Space(bss.AmbientSpace()), 
  order(ordr), 
  basis(bss)
{
  basis.ToSortedReducedGrobner(*ordr);
}

template<class T>  gIdeal<T>::gIdeal(const gIdeal<T> & ideal)
: Space(ideal.Space), order(ideal.order), basis(ideal.basis)
{
}

template<class T>  gIdeal<T>::~gIdeal()
{
}

//----------------------------------
//        Operators
//----------------------------------


template<class T> gIdeal<T>& gIdeal<T>::operator=(const gIdeal<T> & rhs)
{
  assert (Space == rhs.Space);

  if (*this != rhs) {
    this->order = rhs.order;
    this->basis = rhs.basis;
  }
  return *this;
}

template<class T>  bool gIdeal<T>::operator==(const gIdeal<T> & rhs) const
{
  if (Space != rhs.Space) { return false; }
  else {
    if (*order == *(rhs.order)) {
      if (basis == rhs.basis) return true;  
      else                    return false;
    }
    if (basis.Length() != rhs.basis.Length()) return false;
    else {
//      gPolyList<T> reordered_rhs(rhs.basis,*order);
//      gPolyList<T> reordered_rhs(rhs.basis);
      const gIdeal<T> reordered_rhs(order,rhs.basis);
      if (basis == reordered_rhs.basis) return true;
      else                              return false;
    }
  }
}

template<class T>  bool gIdeal<T>::operator!=(const gIdeal<T> & rhs) const
{
  return !(*this == rhs);
}

template<class T>  gIdeal<T> gIdeal<T>::operator+ (const gIdeal<T> & rhs) const
{
  Gambit::List<gPoly<T> *> combined;

  for (int i = 1; i <= basis.Length(); i++) {
    gPoly<T>* temp = new gPoly<T>(basis[i]);
    combined.Append(temp);
  }
  for (int j = 1; j <= rhs.basis.Length(); j++) {
    gPoly<T>* temp = new gPoly<T>(rhs.basis[j]);
    combined.Append(temp);
  }

  return gIdeal<T>(Space,order,combined);
}

template<class T>  gIdeal<T> gIdeal<T>::operator* (const gIdeal<T> & rhs) const
{
  Gambit::List<gPoly<T> *> basis_products;

  for (int i = 1; i <= basis.Length(); i++) 
    for (int j = 1; j <= rhs.basis.Length(); j++) {
      gPoly<T>* temp = new gPoly<T>(basis[i] * rhs.basis[j]);
      basis_products.Append(temp);
    }

  return gIdeal<T>(Space,order,basis_products);
}

//----------------------------------
//          Information
//----------------------------------
  
template<class T> gIdeal<T> gIdeal<T>::MonomialIdeal() const
{
  gPolyList<T> mon_bss(Space,order);
  for (int i = 1; i <= basis.Length(); i++) 
    mon_bss += basis[i].LeadingTerm(*order);

  return gIdeal<T>(order,mon_bss);
}
  
template<class T> Gambit::List<exp_vect> gIdeal<T>::MonomialBasis() const
{
  Gambit::List<exp_vect> answer;

  Gambit::Array<int> MinIndices(Dmnsn()), MaxIndices(Dmnsn());
  for (int i = 1; i <= Dmnsn(); i++) MinIndices[i] = 0;

  for (int j = 1; j <= basis.Length(); j++) {
    if (basis[j].LeadingPowerProduct(*order).IsUnivariate())
      MaxIndices[basis[j].LeadingPowerProduct(*order).SoleActiveVariable()] 
	= basis[j].LeadingPowerProduct(*order).TotalDegree() - 1;
  }

  gIndexOdometer odometer(MinIndices,MaxIndices);
  while (odometer.Turn()) {
    exp_vect candidate(Space,odometer.CurrentIndices());
    bool add = true;
    for (int j = 1; j <= basis.Length(); j++) 
      if (basis[j].LeadingPowerProduct(*order) <= candidate) add = false;
    if (add) answer.Append(candidate);
  }

  return answer;
}
  
template<class T> bool gIdeal<T>::IsRoot(const Gambit::Vector<T>& v) const
{
  for (int i = 1; i <= basis.Length(); i++)
    if (basis[i].Evaluate(v) != (T)0) return false;
  return true;
}
  
template<class T> bool gIdeal<T>::ZeroDimensional() const
{
  if (Dmnsn() == 0) return true;

  Gambit::Array<bool> HasLeadingTermThatIsPowerOfVariable(Dmnsn());
  int i;
  for (i = 1; i <= Dmnsn(); i++)
    HasLeadingTermThatIsPowerOfVariable[i] = false;

  for (int j = 1; j <= basis.Length(); j++) {

    if (basis[j].LeadingPowerProduct(*order).IsConstant()) return false;

    if (basis[j].LeadingPowerProduct(*order).IsUnivariate())
      HasLeadingTermThatIsPowerOfVariable[
	 basis[j].LeadingPowerProduct(*order).SoleActiveVariable()
	    ] = true;
  }

  for (i = 1; i <= Dmnsn(); i++)
    if (!HasLeadingTermThatIsPowerOfVariable[i]) return false;
  return true;
}

template<class T> bool gIdeal<T>::Contains(gPoly<T> & f) const
{
  gPoly<T> zero(Space,(T)0,order);
  return (basis.ReductionOf(f,*order) == zero);
}
  
template<class T> bool gIdeal<T>::IsEntireRing() const
{
  gPoly<T> one(Space,(T)1,order);
  return Contains(one);
}