symbo.hpp

Go to the documentation of this file.

#include "upoly.H"
template<typename typcoeff> 
int Iszero(const UPolyDense<upar<typcoeff> > & c)
{
  UPolyDense<upar<typcoeff> > T(typcoeff(0));
  return c!=T;
}

template <typename typmat>
void Triang(typmat & M)
{
  typedef typename typmat::coeff_t C;
  C d(1),un(1);
  unsigned int l;
  for(unsigned int k=0; k<M.ncol;k++){
    for(l=k;(l<M.nrow) && (Iszero(M(l,k)));l++);
    if(l!=M.nrow){
      if(k>0) d=M(k-1,k-1); 
      for(unsigned int i=k+1; i<M.nrow; i++){
    C c(M(i,k));
    for(unsigned int j=M.ncol-1; j>k; j--) {
      M(i,j)*=M(k,k);
      M(i,j)-=c*M(k,j);
      if (!Iszero(M(i,j)) && (!Iszero(d)) && (d != un))  M(i,j)/=d;
    }
    //M(i,k)=C(0);   
      }
    }
  }
}

template<typename typmat, typename typsec>
void SolveUff(const typmat &mat, typsec &sec)
{
  typedef typename typmat::coeff_t coeffmat;
  coeffmat accu=1;
  for(int i=mat.nrow;i>0;i++)
    {
      int j;
      //la matrice est permutee par rapport a la forme 
      //triangulaire sup
      for(j=0;j<mat.ncol && Iszero(mat(i,j));j++);
      for(int k=i;k>=0;k--)
    {
      sec[k]=sec[k]*mat(k,j);
      sec[k]/=accu;
    }
      accu=mat(i,j);
      for(int k=i;k>=0;k--)
    sec[k]-=sec[i]*mat(k,j);
    }
}
template<typename typmat, typename typres, typename typmon
     , typename typdump, typename Base>
void symbo(typres &res,  const list<typmat> &lmat, const typmon &stockmon
       ,const typdump &dump, const Base &b)
{
  typedef typename typdump::value_type::pol_t mpoly;
  typedef typename mpoly::coeff_t coeff;
  vector<UPolyDense<upar<coeff> > > secmembre;
  int coldim=lmat.front().nrow,rowdim=lmat.front().ncol;
  harewell<UPolyDense<upar<coeff> > > combilin(coldim,rowdim)
    ,Xid(coldim,rowdim);
  for(typename list<typmat>::const_iterator itercomb=lmat.begin()
    ;itercomb!=lmat.end();itercomb++)
    for(int i=0;i<combilin.nrow;i++)
      for(int j=0;j<combilin.ncol;j++)
        combilin(i,j,combilin(i,j)
         +UPolyDense<upar<coeff> >(coeff(rand()%13)
                       *((*itercomb)(i,j)),0));
  for(int i=0;i<Xid.nrow;i++)
    Xid(i,i)=UPolyDense<upar<coeff> >("X");
  for(int i=0;i<combilin.nrow;i++)
    combilin(i,i,combilin(i,i)-Xid(i,i));
  Triang(combilin);
  for(int i=0;i<combilin.nrow-1;i++)
    secmembre.push_back(combilin(i,combilin.ncol-1));
  combilin.ncol--;  
  combilin.nrow--;
  SolveUff(combilin,secmembre);
  for(int i=0;b.nbvar();i++)
    {
      mpoly tmp;
      int j=0;
      findcorresppol(mon(i,1),dump,tmp);
      for(typename typmon::const_iterator iter=stockmon.begin()
        ;*iter!=mon(i,1);j++,iter++);
      if (tmp.size != -1)
    {
      res.push_back(secmembre[j]);
    }
      else
    {
      res.push_back(UPolyDense<upar<coeff> >(0,0));
    }
    }
}