pol2terspecial.hpp

Go to the documentation of this file.

//#include"pol2bisspecial.hpp"
//
//
//
#include<future>
#include<string>
#include<iostream>
#include<iterator>
#include"borderbasix/auxfun3.hpp"
#ifndef LEVELSELECT
#define LEVELSELECT 1e-6
#endif
//#define NBFUTURE 7
#define NBFUTURE 7

template<typename pol_t>
struct resmul
{
  typedef pol_t pol;
  typedef typename pol_t::monom_t monom_t;
 typedef typename pol_t::monom_t mon;
  pol poly;
  mon * exceed;
  int sizeexceed;
};

struct stack_elem {
  typedef pol2ter<mon,COEFF> *typPk;
  COEFF coeff;
  pol2ter<mon,COEFF> *p;
  mon *exceed;
  int sizeexceed;
};

template<typename typPk_t, typename mon>
struct component
{
  typedef typPk_t typPk;
  typedef typename typPk::value_type pol;
  typedef typename pol::coeff_t coeff;

  typPk_t poly;
  vector<mon *> exceed;
  vector<int> sizeexceed;
};
template<typename typPk_t, typename mon>
bool operator==(const component<typPk_t,mon> &a, const component<typPk_t,mon> &b)
{
  bool res=true;
  return (a.poly==b.poly)&&(a.exceed==b.exceed)&&(a.sizeexceed==b.sizeexceed);
}

//
//Routine a remlacer plus tard par une propapagtion des ancestors 
//lors du mult et une definition lors newcrochetchoix
//
template <typename typcomp, typename workspace, typename typserv>
void  update_pol(typcomp &comp_to_update,typcomp &updating_comp, workspace *w,
                 typserv &serv);
#ifndef ASYNC
template<typename typPk, typename Base>
void compute_ancestors(typPk &tmpPk, Base & b)
{
 //
  //
  //cet algo repose sur le fait que les composantes sont placées dans l'ordre d'inclusion
  // on a pas de petite composante avant une plus grosse
  //a changer par la suite
  //
  //
  typedef typename typPk::value_type::monom_t mon;
  for (auto &iter : tmpPk)
    {
     
      mon tmp=iter.ind;
      //Pour etre ancetre il faut diviser 
      //et etre dans la meme composante de B
      mon *candidate_ancestor=NULL;
      tmp.SetCoeff(1);
      for(auto &iterbase : b)
    {
      mon *compancestor=NULL;
      int isincomponent=0;
      for(int i=0;i<iterbase.taille1;i++)
        if(isdivisible(iterbase.accept[i],tmp,b))
          {
        isincomponent=1; 
#ifdef DEBUG
        cout<<"je break ici"<<endl;
#endif
        break;//pas la peine d'aller plus loin on est dans la composante
          }
      if(isincomponent)//sinon ca sert a rien
        {
          for(int i=0;i<iterbase.taille2;i++)
        if(isdivisible(iterbase.refuse[i],tmp,b))
          {
            if(compancestor)//on a deux monomomes qui divisent
              {
            candidate_ancestor=NULL;
#ifdef DEBUG
            cout<<"je break ici"<<endl;
#endif
            break;//pas la peine d'aller plus loin
              }
            else
              {
#ifdef DEBUG
            cout<<"g trouvé ici un gus "<<iterbase.refuse[i]<<endl;
#endif
            compancestor=&iterbase.refuse[i];
              }
          }
        }
      if(compancestor)
        {
#ifdef DEBUG
          cout<<"je passe ici compancestor non null"<<endl;
#endif
          candidate_ancestor=compancestor;
        }
    }
      iter.ancestor=candidate_ancestor;
    }
#ifdef DEBUG
  for(auto &&iter : tmpPk)
    {
      cout<<"ancestor ---- "<<iter.ind<<endl;
      if (iter.ancestor)
    cout<<"\t\t -> "<<*(iter.ancestor)<<endl;
      else
    cout<<"\e[31mcaca\e(B\e[m"<<endl;
    }
#endif
}
#endif
template<typename typgraph, typename typPk, 
     typename typexceed, typename typsizeexceed, typename typmon2face,
     typename Base>
void makegraph(typgraph & G, const typPk tmpPk, const typexceed &exceed,
           const typsizeexceed & sizeexceed, 
           typmon2face &mon2face,
           Base &b)
{
  //on commence par numeroter les faces
  //Et on associe les ind et leur face
  //->>>>>>>>>>>>>on utilise une map
  typedef typename typPk::value_type::monom_t mon;
  //  typedef typename typPk::value_type::order_t ord;
  typedef my_Dlex<mon> ord;
  int numface=0;
  map<typename mon::rep_t,int> root2face;//associee le monome ancetre à sa face
  //map<typename mon::rep_t,int> mon2face;//associe l'ind a sa face
  //le root c'est un monome pres de l'ancetre
  //qui represente l'ind sur la face
  for(auto iterpk=tmpPk.begin();iterpk!=tmpPk.end();iterpk++)
    {
      mon root=*(iterpk->ancestor);
      int count=0;
      for(int i=0;i<b.nbvars;i++)
    if(iterpk->ind.GetDegree(i)>root.GetDegree(i))
      {
        root*=mon(i,1);
        count++;
      }
      //ICI il faut repasser pour s'assurer que le root a bien exactement
      // nvars-1 position augmentées par rapport a l'ancetre
      //en cas de multiappartenance on choisi celle qui est la plus petite
      //dans l'ordre lexico (la multiappartenance c qd l'ind  est sur une 
      //intersefction de deux faces
      for(int i=0;(count<(b.nbvars-1)) && (i<b.nbvars);i++)
    if(iterpk->ancestor->GetDegree(i)==root.GetDegree(i))
      {
        root*=mon(i,1);
        count++;
      }
      if(root2face.count(root.rep)==0) 
    root2face[root.rep]=numface++;
      mon2face[iterpk->ind.rep]=root2face[root.rep];
    }
#ifdef DEBUG
  cout<<"roottoface"<<endl;
  for(auto &i :root2face)
    cout<<i.first<<" --- "<<i.second<<endl;
  cout<<"mon2face"<<endl;
  for(auto &i :mon2face)
    cout<<i.first<<" --- "<<i.second<<endl;
#endif
  //on allour le graph
  G.resize(numface);
  for(auto &i : G)
    {
      i.resize(numface);
      for(auto &j :i)
    j=0;
    }
  //Ensuite on construit le graph
  auto iterexceed=exceed.begin();
  auto itersize=sizeexceed.begin();
  auto iterpk=tmpPk.begin();
  for(;iterpk!=tmpPk.end();iterpk++,iterexceed++,itersize++)
    {
//     cout<<"ici on traite "<<iterpk->ind<<endl;
      for(int i=0;i<*itersize;i++)
    G[mon2face[iterpk->ind.rep]][mon2face[(*iterexceed)[i].rep]]=1;
      // On considere le graph dense
    }
  
}

template<typename typstack, typename typgraph, typename typvec>
void strongconnect(typstack &S, typgraph &G,int i, typvec &vec,
           int &index, typvec &lowlink,int startindex)
{
  S.push_back(i);
  vec[i]=index;
  lowlink[i]=index;
  index++;
//  cout<<"je visite "<<i<<" stack "<<S.size()<<endl;
//   cout<<"i "<<i<<" lowlink[i] "<<lowlink[i]<<endl;
  for(int j=0;j<G[i].size();j++)
    {
//      cout<<"j "<<j<<endl;
      if((i!=j)&& G[i][j])
    {
//    cout<<"dans le if !"<<endl;
      if(vec[j]==-1) 
        {
          strongconnect(S,G,j,vec,index,lowlink,startindex);
          lowlink[i]=min(lowlink[i],lowlink[j]);
        }
      else
        {
          if(vec[j]>=startindex)
            lowlink[i]=min(lowlink[i], vec[j]);
          //      cout<<"lowlink[i] "<<lowlink[i]<<endl;
        }
    }
    }
  //Est ce une racine
  //  cout<<"racine "<<index-1<<" "<<i<<" ::vec[i]:: "<<vec[i]<<" ::lowlink[i]:: "<<
  //  lowlink[i]<<endl;
  if(lowlink[i]==vec[i])
    {
      int tmpindex;
      //Oui alors on créé une composante fortement connexe
    //  cout<<"S.size() "<<S.size()<<endl;
     // cout<<"vec[i] "<<vec[i]<<" i "<<i<<endl;
     // copy(S.begin(),S.end(),ostream_iterator<int>(cout," "));cout<<endl;
      do
    {
      tmpindex=S.back();
//    cout<<"pour la comp "<<tmpindex<<endl;
      S.pop_back();
      vec[tmpindex]=i;
      
    }
      while(tmpindex!=i);//ici i designe un sommet pas un index
    }
}

template<typename typgraph, typename typvec>
void tarjan(typgraph &G, typvec & vec)
/*
algorithm tarjan is
  input: graph G = (V, E)
  output: set of strongly connected components (sets of vertices)

  index := 0
  S := empty
  for each v in V do
    if (v.index is undefined) then
      strongconnect(v)
    end if
  repeat

  function strongconnect(v)
    // Set the depth index for v to the smallest unused index
    v.index := index
    v.lowlink := index
    index := index + 1
    S.push(v)

    // Consider successors of v
    for each (v, w) in E do
      if (w.index is undefined) then
        // Successor w has not yet been visited; recurse on it
        strongconnect(w)
        v.lowlink := min(v.lowlink, w.lowlink)
      else if (w is in S) then
        // Successor w is in stack S and hence in the current SCC
        v.lowlink := min(v.lowlink, w.index)
      end if
    repeat

    // If v is a root node, pop the stack and generate an SCC
    if (v.lowlink = v.index) then
      start a new strongly connected component
      repeat
        w := S.pop()
        add w to current strongly connected component
      until (w = v)
      output the current strongly connected component
    end if
  end function
*/
{
  int index=0;
  vector<int> lowlink(G.size(),-1);
  list<int> S;

  for(int i=0;i<G.size();i++)
    {
      if(vec[i]==-1)
    {
//    cout<<"j'attaque une nvlle composante "<<i<<endl;
      strongconnect(S,G,i,vec,index,lowlink,i);
        
    }
    }
  
  //On remet les numeros comme il faut 0,1,2,3 pour eviter les trous
  vector<int> translator(G.size(),-1);
  int numcomp=0;
  for(int i=0;i<vec.size();i++)
  {
    if(translator[vec[i]]==-1)
      translator[vec[i]]=numcomp++;
     vec[i]=translator[vec[i]];
  }
  copy(vec.begin(),vec.end(),ostream_iterator<int>(cout," "));
  cout<<endl;
  
}

//forward declare 
void reconstruct_pr_pc(int *&pr,int *&pc,int *tmppr,int *tmppc,int sizetmp,int indexinside);
//

template<typename typMat, typename typPk, typename typMk, typename typexceed,
     typename typsizeexceed, typename typP>
void construct_matrix_from_newpol_nofree(typMat &mat, typPk &tmpPk,
                  typMk &tmpMk, typexceed &exceed,
                  typsizeexceed &sizeexceed, typP &P)
 {
   typedef typename typP::value_type::order_t ord;
   typedef typename typPk::value_type pol;
   typedef typename pol::coeff_t coeff;
   typedef typename typP::value_type polyalp;
   int nnz=0,*ncol;
   //memset((void*)ncol,0,(*mat)->ncol*sizeof(int));
   map<mon,int,ord> mapint;
   typename typMk::iterator itermk=tmpMk.begin();
   typename typPk::iterator iterpk=tmpPk.begin();
   //   typename list<mon*>::iterator iterex=exceed.begin();
   typename typexceed::iterator iterex=exceed.begin();
   typename typsizeexceed::iterator itersize=sizeexceed.begin();
   //#ifdef modifnewpk
   typsizeexceed tmpsizeex;
   //#endif
   int line_num=0;
   mat.nrow=tmpPk.size();
   mat.ncol=tmpMk.size();
   ncol=(int*)calloc(mat.ncol,sizeof(int));
   //
   //On crée une map pour indexer les colonnes
   //   
     for(int i=0;itermk!=tmpMk.end();i++,itermk++)
       {
     mapint[*itermk]=i;
       } 
     iterpk=tmpPk.begin();
     nnz=0;
     for(;iterpk!=tmpPk.end();iterpk++,iterex++,itersize++,line_num++)
       {
     //#ifdef modifnewpk
     int tmp=0;
     //#endif
    {
      // coeff tmpcoeff=(iterpk->ind).GetCoeff()-(coeff)1;
     for(int i=0;*itersize>i;i++)
       {
         //cout<<"iterpk icici ->GetCoeff()"<<iterpk->ind.GetCoeff()<<endl;
         if (mapint.count((*iterex)[i]))
           {
         ++ncol[mapint[(*iterex)[i]]];
           }
         //#ifdef modifnewpk
         else
           {
         tmp++;
           }
         //#endif
       }
     //#ifdef modifnewpk
     nnz+=*itersize-tmp;
     tmpsizeex.push_back(tmp);
     //#endif
    }
       }
     //fin rajout c pour prendre en compte que lors de 
     //de la proj le monome de tete se trouve dans exceed
     //fin enlevement
     cout<<"nnz of matrix "<<nnz<<endl;
     //en general la matrice est rectangulaire
     void * toto=MAC_REV_MALLOC<coeff>(nnz*sizeof(typename typMat::coeff_t));
     if (mat.nzval!=NULL) MAC_REV_FREE<typename typMat::coeff_t>
                (mat.nzval,mat.nnz
                 *sizeof(typename typMat::coeff_t));
     mat.nnz=nnz;
     mat.nzval=toto;
     
     if(mat.rowind!=NULL) free(mat.rowind);
     mat.rowind=(int*)malloc(nnz*sizeof(int));
     if(mat.colptr!=NULL) free(mat.colptr);
     mat.colptr=(int*)malloc((mat.ncol+1)*sizeof(int));
     cout<<"matrix of size "<<mat.nrow<<"x"<<mat.ncol<<endl;
   //matrice cree rest plus qu'a ecrire dedans ca
     typename typMat::coeff_t* nzval=(typename typMat::coeff_t*)mat.nzval;
     int *rowind=mat.rowind;
     int *colptr=mat.colptr;
    //mise  en place des colptr
    *colptr=0;
    for(int i=0;i<mat.ncol-1;i++)
      {
    colptr[i+1]=colptr[i]+ncol[i];
      }
    colptr[mat.ncol]=mat.nnz;
    
    //fin mise en place
    //la on fait des manip symboliques pour savoir ou ecrire 
    //les coeffs de chaques colonnes ca evite d'avoir a les
    //recopier 50000 fois pour construire la matrice 
    //c bcp plus rapide
    itermk=tmpMk.begin();
    iterpk=tmpPk.begin();
    iterex=exceed.begin();
    itersize=sizeexceed.begin();
    line_num=0;
    //#ifdef modifnewpk
    auto itertmp=tmpsizeex.begin();
    //#endif
    for(;(iterpk!=tmpPk.end());iterpk++,iterex++,itersize++,line_num++,itertmp++)
      {
    //#ifdef modifnewpk
    int ktmp=0;
    mon *tmpex=(mon*)MAC_REV_MALLOC<mon>(*itertmp*sizeof(mon));
    //#endif
      {
        for(int i=0;*itersize>i;i++)
          {
        //#ifdef modifnewpk
        if(mapint.count((*iterex)[i]))
          {
            //#endif
            const coeff stock=(*iterex)[i].GetCoeff();
            (*iterex)[i].SetCoeff(1);
            const unsigned int tmpplusex=mapint[(*iterex)[i]]+1;
            const unsigned int tmpex=mapint[(*iterex)[i]];
            nzval[colptr[tmpplusex]-ncol[tmpex]]
              =stock;
            rowind[colptr[tmpplusex]-ncol[tmpex]]=line_num; 
            ncol[tmpex]--;
            //#ifdef modifnewpk
              }
        else
          {
            tmpex[ktmp++]=(*iterex)[i];
          }
        //#endif
          }
      }
      MAC_REV_FREE<mon>(*iterex,*itersize*sizeof(mon));
      //#ifdef modifnewpk
      *iterex=tmpex;
//    cout<<"itersize "<<*itersize<<" itertmp "<<*itertmp<<endl;
      *itersize=*itertmp;
//    cout<<"regarde ici "<<iterpk->ind<<endl;
//    for(int i=0;i<*itertmp;i++)
//      cout<<"ccccccccccccc "<<(*iterex)[i]<<endl;
      //#endif
      }
    free(ncol);
    
   for(typename typPk::iterator iterpk=tmpPk.begin()
     ;iterpk!=tmpPk.end();iterpk++)
     iterpk->ind.SetCoeff(1);
 }

#ifndef ASYNC
template<typename typMat,typename typPk,typename typP,
         typename typMk,typename Base,typename typdump,typename typexceed,
         typename typsizeexceed, typename workspace, typename monomial_server>
int newpkmkMatrixofSolve(typMat &mat,typMat & L,typPk &Pk, typP &P,typMk &Mk,
                    int &k,Base &b, typdump &dump,typexceed &exceed_ret, 
             typsizeexceed & sizeexceed_ret, workspace &w, monomial_server &serv)
{
  
typedef typename typMk::value_type mon;
  typedef typename typP::value_type::order_t ord;
  typedef typename typPk::value_type pol;
  typedef typename pol::coeff_t coeff;
  typedef typename typP::value_type polyalp; 
  typedef component<vector<pol>,mon> compo;
  typPk recall;//to be erased as soon as treat_crochet et compute_... 
  //have been patched
  map<typename mon::rep_t,int> mon2face, mon2col; 
  vector<vector<int> > G;
  vector<int> vec;
  typPk tmpPk,workvect,convcrochet;
  typMk tmpMk;
  //list<mon *> exceed,exceedcrochet;
  vector<mon *> exceed,exceedcrochet;
  list<int> sizeexceed,sizeexceedcrochet;
  int nbproj=0,flagdoplus=1;
  int nbcompo=0,*pr,*pc;
  list<compo>listcomp;
  compo crochet;//ceci est le regroupement des crochet
  //  list<list<mon *> > listexceed;
  // list<list<int> > listsizeexceed;
  int indexinside=0;
  //
  mat.destroystore();mat.ncol=0;mat.nrow=0;mat.nnz=0;
  exceed_ret.erase(exceed_ret.begin(),exceed_ret.end());
  sizeexceed_ret.erase(sizeexceed_ret.begin(),sizeexceed_ret.end());
 //prise en compte des pol pas encore introduits
  cerr<<"k on entering newpkmkmatrixof "<<k<<endl;
  typP ttmpP=Crochet(P,k+1);//
    if (!ttmpP.empty())
    {
//flagdo plus est positinne dans cett fonction au cas ou il ait chute de degre  
    treat_noempty_crochet(P, ttmpP, tmpPk, tmpMk, Pk,
                            Mk, dump, b, flagdoplus, k, recall, exceed, sizeexceed, serv);
    cout<<"en sortie tmpPk.size "<<tmpPk.size()<<endl;
     }
    if(flagdoplus)
      {
    cout<<"newpol"<<endl;
    compute_newpol(Pk,tmpPk,tmpMk,dump,exceed,sizeexceed, b, k, recall, w,serv);
    cout<<"fin newpol"<<endl;
      }      
 //Ces polynome la n'ont pas d'ind à proprement parler 
 //Il ne faut pas les prendre en compte pour faire le graph 
 //
    //
    //Il faut effacer tmpMk et le reconstruire apres coup en
    //fonction des composantes fortement connexes decouvertes
    //
   
  conv_merge(convcrochet,ttmpP,exceedcrochet,sizeexceedcrochet,b,serv);
  crochet.poly.resize(convcrochet.size());
  crochet.sizeexceed.resize(convcrochet.size());
  copy(convcrochet.begin(),convcrochet.end(),crochet.poly.begin());
  crochet.exceed=exceedcrochet;
  copy(sizeexceedcrochet.begin(),sizeexceedcrochet.end(),crochet.sizeexceed.begin());
#ifdef DEBUG1 
 cout<<"crochet "<<endl;
  for(int i=0;i<crochet.poly.size();i++)
    {
      cout<<"ind "<<crochet.poly[i].ind<<endl;
      copy(crochet.exceed[i],crochet.exceed[i]+crochet.sizeexceed[i],ostream_iterator<mon>(cout," + "));
      cout<<invconv<Poly>(crochet.poly[i],serv)<<endl;
    }
#endif
 // 
 // cout<<"icccccccccccccccccccccccccccccccccccccccccccccccccccccci"<<endl;
  auto itersizeexceed=sizeexceed.begin();
 // for(auto &p : exceed)
  //  copy(p,p+*(itersizeexceed++),ostream_iterator<mon>(cout," "));
  //crochet.sizeexceed=sizeexceedcrochet;
  //relire ce suit pour verifier les matrices pas carres

  //On constuit le graphe
  //Revoir l'alloc de G 
     
compute_ancestors(tmpPk,b);
  makegraph(G, tmpPk, exceed, sizeexceed, mon2face,b);
  //On construit le grpahe et
  //on conserve l'appairage entre les ind et le faces
  vec.resize(G.size(),-1);
#ifdef DEBUG1
  cout<<"le graph"<<endl;
  for(auto &v : G)
    {
      for(auto i: v)
    cout<<i<<", ";
      cout<<endl;
    }
#endif
  cout<<"je rentre ds tarjan "<<G.size()<<" "<<tmpMk.size()<<endl;
 tarjan( G, vec);
  pr=(int *)malloc(tmpPk.size()*sizeof(int));
  pc=(int *)malloc(tmpMk.size()*sizeof(int));
   tmpMk.erase(tmpMk.begin(),tmpMk.end());
 //
  //On construit les listes de pol a partir des Comp Connexes
  for_each(vec.begin(),vec.end(),[&nbcompo](int i){nbcompo=(nbcompo<i?i:nbcompo);});
  //nbcompo est l'index de la derniere composante
  for(int i=0;i<=nbcompo;i++)
   {
    typPk tmp;
    list<mon*> tmpexceed;
    list<int> tmpsizeexceed;
    compo tmpcomp;
    int j=0;
    auto itersizeexceed=sizeexceed.begin();
    auto iterexceed=exceed.begin();

    for_each(tmpPk.begin(),tmpPk.end(),
        [&]
        (typename typPk::value_type &iter)
    {
#ifdef DEBUG
      cout<<"mon2face "<<iter.ind.rep<<" "<<mon2face[iter.ind.rep]<<endl;
      cout<<"vec[mon2face[iter.ind.rep]] "<<vec[mon2face[iter.ind.rep]]<<
        " "<<i<<endl;
#endif
      if(vec[mon2face[iter.ind.rep]]==i) 
        {
          //tmp.push_back(iter);
          tmpcomp.poly.push_back(iter);
          //tmpsizeexceed.push_back(*itersizeexceed++);
#ifdef DEBUG
          for(int i=0;i<*itersizeexceed;i++)
        cout<<" "<<(*iterexceed)[i];
          cout <<endl;
#endif
          tmpcomp.sizeexceed.push_back(*itersizeexceed);
          //        tmpexceed.push_back($iterexceed++);

          tmpcomp.exceed.push_back(*iterexceed);
        }
            itersizeexceed++;
        iterexceed++;
    });//fin de la lambda et du foreach
    //  listcomp.push_back(tmp);
    //  listexceed.push_back(tmpexceed);
    //  listsizeexceed.push_back(tmpsizeexceed);
    listcomp.push_back(tmpcomp);
   }  
#ifdef DEBUG
  for(auto iter=listcomp.begin();iter!=listcomp.end();iter++)
  {
    for(auto &m : iter->poly)
      cout<<"ind compo "<<m.ind<<endl;
   cout<<"----------------------------------------------------------"<<endl;
  }
#endif
//Pour chaque liste matrice + solve
//  itersizeexceed=sizeexceed.begin();
//  iterexceed=exceed!.begin();
  int prout=0;
  int numcompo=0;
  cout<<"nb total de compo dans la loiste "<<listcomp.size()<<endl;
  for(auto iter=listcomp.begin();iter!=listcomp.end();iter++,numcompo++)
   {
    //il faut construit tmpMK pour la composante
     int *tmppr,*tmppc;
     typMk verytmpMk;
     typMat tmpmat;
#ifdef DEBUG1
     cout<<"pour la composante "<<prout++<<endl;
     {
     auto toto=iter->sizeexceed.begin();
     typename vector<pol>::iterator tutu=iter->poly.begin();
    for(auto m : iter->exceed)
      {
    cout<<"sizeexceed ici "<<*toto<<endl;
    copy(m,m+*(toto++),ostream_iterator<mon>(cout," + ")),cout<<endl;;
    cout<<tutu->ind<<endl;
    if(tutu->size)
      copy(tutu->nf,tutu->nf+tutu->size,ostream_iterator<coeff>(cout," ")),cout<<endl;
    cout<<invconv<Poly>(*(tutu++),serv)<<endl;
    cout<<"------------------------------------ "<<endl;
      }
     } 
#endif
     for(auto &&iterind : iter->poly)
     {
       verytmpMk.push_back(iterind.ind);
       //   cout<<"je push "<<iterind.ind<<endl;
     }
     construct_matrix_from_newpol_nofree(tmpmat, iter->poly, verytmpMk, 
                     iter->exceed, iter->sizeexceed, P);

#ifdef DEBUG
     cout<<"matrix of comp"<<endl;
     affdebug(tmpmat);
     cout<<"\e[32mpour la composante \e[37m"<<prout<<endl;
     {
       auto toto=iter->sizeexceed.begin();
       typename vector<pol>::iterator tutu=iter->poly.begin();
       for(auto m : iter->exceed)
     {
       copy(m,m+*(toto++),ostream_iterator<mon>(cout," - ")),cout<<endl;;
       cout<<tutu->ind<<endl;
       cout<<"size tutu "<<tutu->size<<endl;
        if(tutu->size)
         copy(tutu->nf,tutu->nf+tutu->size,ostream_iterator<coeff>(cout," "));
       cout<<invconv<Poly>(*(tutu++),serv)<<endl;
       cout<<"------------------------------------ "<<endl;
     }
     }
#endif 
#ifdef DEBUG
     cout<<"la petite matrice"<<endl;
     affdebug(tmpmat);
#endif
        /*
     modifier constructmatrix pour ressortir les exceed a cote de la composante
     //FAIT
    */


    //P est la pour forwareder les types ord et pol 
    //On resoud 
    //            pr et pc sont alloués avec malloc
    //            dans la fonction solve
    //            il faut les desallouer
    // 
    //Ceci est une fonction solve a ecrire pour reporter
    //les operations sur les poly aussi sur les exceed
    //
     Solve(tmpmat,L,*iter,&tmppr,&tmppc,verytmpMk,w,serv);
#ifdef DEBUG1
    cout<<"\e[31mApres SOLVE \e[37m"<<prout<<endl;
    {
      auto toto=iter->sizeexceed.begin();
      typename vector<pol>::iterator tutu=iter->poly.begin();
      for(auto m : iter->exceed)
        {
          copy(m,m+*(toto++),ostream_iterator<mon>(cout," - ")),cout<<endl;;
          cout<<tutu->ind<<endl;
          cout<<invconv<Poly>(*(tutu++),serv)<<endl;
        }
    } 
#endif

//  copy(tmppc,tmppc+verytmpMk.size(),ostream_iterator<int>(cout," "));
//  cout<<endl;
    // inverse(tmppc,verytmpMk.size());
    // Dopermu(verytmpMk,tmppc);
     for(auto &m :verytmpMk)
       {
         m.SetCoeff(1);
         tmpMk.push_back(m);
       } 
#if 0
     // ces operation sont faites dans reorder ind entre solveL et solveU
    {
    //ici on remet les ind a la bonne valeur 
    //a inclure dans une petite fonction inline apres
    int j=0;
    typename typMk::iterator itermon=verytmpMk.begin();

    for(auto iter_in_poly=iter->poly.begin();
        itermon!=verytmpMk.end();iter_in_poly++,itermon++,j++)
     {
       cout<<" dans verytmp "<<*itermon<<
         invconv<typename typP::value_type>(*iter_in_poly,serv)<<endl;
     
       if(!Iszero(tmpmat(j,j)))
         {
           //cout<<"ligne "<<j<<" je mets le monome "<<*itermon<<" "<<iter->ind.GetCoeff()<<endl:
           iter_in_poly->ind=*itermon;//*(iter_in_poly->ind.GetCoeff());
           iter_in_poly->ind.SetCoeff(1);
         }
       else
         {
           //On rajoute aux exceed des autre polyn\çomes de la composante les
           //les monomes correspondant a cette colonne 
           for(int i=tmpmat.colptr[j];i<tmpmat.colptr[j+1];i++)
         if(!(((coeff*)tmpmat.nzval)[i]==0))
           {
             //il y a eu permiutation des colonnes il
             // faut donc se fier a itermon plutot qu'aux ind
             mon tmpmon=*itermon;//iter_in_poly->ind;
             coeff tmpcoeff=((coeff*)tmpmat.nzval)[i];
             //Bon maintenant on regarde si le coefficient diagonal de U de la 
             //ligne tmpmat.rowind[i] est non nul et si c'est le cas
             //il fautr diviser par lui 
             if(!Iszero(tmpmat(tmpmat.rowind[i],tmpmat.rowind[i])))
               tmpcoeff/=tmpmat(tmpmat.rowind[i],tmpmat.rowind[i]);
             //Apres on peut affecter
             tmpmon.SetCoeff(tmpcoeff);
             iter->exceed[tmpmat.rowind[i]]=(mon*)
               MAC_REV_REALLOC<mon>(iter->exceed[tmpmat.rowind[i]],
                        iter->sizeexceed[tmpmat.rowind[i]]*sizeof(mon),
                        (iter->sizeexceed[tmpmat.rowind[i]]+1)*sizeof(mon));
             
             iter->exceed[tmpmat.rowind[i]][iter->sizeexceed[tmpmat.rowind[i]]]=tmpmon;
             iter->sizeexceed[tmpmat.rowind[i]]++;
           }
           iter_in_poly->ind=mon(0);
           
         }
     }
    }
#endif
#ifdef DEBUG11
    cout<<"apres solve les poly de la compo"<<endl;
    for(int i=0;i<iter->poly.size();i++)
      {
        cout<<"ind "<<iter->poly[i].ind<<endl;
        copy(iter->exceed[i],iter->exceed[i]+iter->sizeexceed[i],
         ostream_iterator<mon>(cout," + "));cout<<endl;
        cout<<invconv<typename typP::value_type>(iter->poly[i],serv)<<endl;
        cout<<"--------------------********************"<<endl;
      }
    for(auto i : iter->poly)
      cout<<"ind "<<i.ind<<endl;
    copy(iter->sizeexceed.begin(),iter->sizeexceed.end(),ostream_iterator<int>(cout," "));
    cout<<"iciiiiiiiiiiiiiiiiiiiiiiiiiiiiiii"<<endl;


#endif
    cout<<"update "<<endl;
    //mise a jour des autres polynomes et recnstruction de la matrice complete
#if 1
    cout<<"la map"<<endl;
    map<typename mon::rep_t,resmul<typename compo::typPk::value_type> > mon2pol;
    {
      auto iterexceed=iter->exceed.begin();
          auto itersizeexceed=iter->sizeexceed.begin();
      auto iterpoly=iter->poly.begin();
      int i=0;
      for(auto iterpoly=iter->poly.begin();iterpoly!=iter->poly.end();
          iterpoly++,iterexceed++,itersizeexceed++,i++)
        {
          resmul<typename compo::typPk::value_type> tmp;
          tmp.poly=*iterpoly;
          tmp.exceed=*iterexceed;
          tmp.sizeexceed=*itersizeexceed;
          if(!(tmpmat(i,i)==0)) mon2pol[iterpoly->ind.rep]=tmp;
          //on augmente le workspace pour update apres
          for(;w.sizenf<iterpoly->sizenf;increase_space(&w));
          }
    }
    cout<<"le graphe"<<endl;
#endif
    //on regarde pour chaque sommet de la composante les arcs sortant
    //ceux qui aboutissent a des sommets qui ne sont pas dans la compo
    //indiquent un update a faire
    vector<char> toupdate(nbcompo+1,0);
    int indexcompo=0;
#if 1
    for(int i_vec=0;i_vec<G.size();i_vec++)
    {
        if(vec[i_vec]==numcompo)
        {
          //on a trouvé un sommet dans la composante
          //en cours de traitement
          for(int j=0;j<G.size();j++)
          {
            //          cout<<"G["<<j<<"]["<<i_vec<<"] "<<G[j][i_vec]<<endl;
            // cout<<"<<<<<<"<<vec[j]<<" "<<numcompo<<endl;
//          if(G[i_vec][j] && vec[j]!=numcompo)
            if(G[j][i_vec] && vec[j]!=numcompo)
              {
            toupdate[vec[j]]=1;
            //On rajoute les composantes qui vont m'updater
            //a celle trouvee
            for(int j_graph=0;j_graph<G.size();j_graph++)
              G[j][j_graph]|=G[i_vec][j_graph];
            //on a plus besoin re parser cet update
            G[j][i_vec]=0;

              }
          }
        }
    }
#endif
    cout<<"calcul poly"<<endl;
    for(auto itercomp=listcomp.begin();itercomp!=listcomp.end();
        itercomp++,indexcompo++)
      //le test est fait dans update pol
       if(toupdate[indexcompo])
         {
#ifdef DEBUG    
           cout<<"la compo a updater avant"<<endl;
           for(int i=0;i<itercomp->poly.size();i++)
         {
           cout<<"ind "<<itercomp->poly[i].ind<<endl;
           copy(itercomp->exceed[i],itercomp->exceed[i]+itercomp->sizeexceed[i],
            ostream_iterator<mon>(cout," + "));cout<<endl;
           cout<<invconv<typename typP::value_type>(itercomp->poly[i],serv)<<endl;
           cout<<"--------------------********************"<<endl;
         }  
#endif 
           update_pol(*itercomp,*iter,mon2pol,w,serv);

#ifdef DEBUG    
           cout<<"la compo apres update"<<endl;
           for(int i=0;i<itercomp->poly.size();i++)
         {
           cout<<"ind "<<itercomp->poly[i].ind<<endl;
           copy(itercomp->exceed[i],itercomp->exceed[i]+itercomp->sizeexceed[i],
            ostream_iterator<mon>(cout," + "));cout<<endl;
           cout<<invconv<typename typP::value_type>(itercomp->poly[i],serv)<<endl;
           cout<<"--------------------********************"<<endl;
         }
#endif
         }

      //update_pol(*itercomp,*iter,w,serv);
    //   cout<<"fin update "<<endl;
     update_pol(crochet,*iter,mon2pol,w,serv);//on fait uassi agir sur les polynomes du crochet
     cout<<"fin des updates "<<endl;
#ifdef DEBUG
     cout<<"apres update du crochet "<<endl;
     for(int i=0;i<crochet.poly.size();i++)
       {
         cout<<"ind "<<crochet.poly[i].ind<<endl;
         copy(crochet.exceed[i],crochet.exceed[i]+crochet.sizeexceed[i],ostream_iterator<mon>(cout," + "));
         cout<<invconv<Poly>(crochet.poly[i],serv)<<endl;
       }
#endif
     reconstruct_mat(mat,tmpmat);
#ifdef DEBUG1
    affdebug(mat);
    //il faut reconstruire les pr et pc globaux
    cout<<"indexinside "<<indexinside<<endl;
#endif
    reconstruct_pr_pc(pr,pc,tmppr,tmppc,tmpmat.ncol,indexinside);
        indexinside+=verytmpMk.size();
    free(tmppr);
    free(tmppc);
   }
  //Reconstruction des tmpPk tmpMk classiques
  // Pour la matrice
  // On gere aussi ici le cas des coefficient des colonnes de pivot nul 
  //comme ca on connait les indices de ligne et colonne
  //
  typPk verytmpPk;
#ifdef DEBUG11  
  cout<<"les pol apres process newpkmk"<<endl;
  for(auto &c : listcomp)
    for(int i=0;i<c.poly.size();i++)
      {
    cout<<"ind "<<c.poly[i].ind<<endl;
    copy(c.exceed[i],c.exceed[i]+c.sizeexceed[i],
         ostream_iterator<mon>(cout," + "));cout<<endl;
    cout<<invconv<typename typP::value_type>(c.poly[i],serv)<<endl;
    cout<<"--------------------********************------------"<<endl;
      }
#endif

  for(auto c : listcomp)
    {
      for(auto &p: c.poly)
    verytmpPk.push_back(p);
      for(auto &p: c.exceed)
    exceed_ret.push_back(p);
      for(auto &p: c.sizeexceed)
    {
      //  cout<<"sizeexceed "<<p<<endl;
      sizeexceed_ret.push_back(p);
    }
     
    }
  cout<<"apres reconstruction des Pk "<<endl;
#ifdef DEBUG
  copy(sizeexceed_ret.begin(),sizeexceed_ret.end(),ostream_iterator<int>(cout," **--** "));cout<<endl;
#endif
  //On rajoute les poly de crochet dans Pk
  //ET (c'est important pour addB on incrmente nrow de mat en fonction
   for(auto &p : crochet.poly)
     {
       p.ind=mon(0);
       verytmpPk.push_back(p);
     }
   for(auto &p : crochet.exceed)
     exceed_ret.push_back(p);
   for(auto &p : crochet.sizeexceed)
     sizeexceed_ret.push_back(p);
   
#if 1
  int indexmat=0;
  for(auto m : tmpMk)
    {
      //      cout<<"mon2col "<<m.rep<<" "<<indexmat<<endl;
      mon2col[m.rep]=indexmat++;
    }
  indexmat=0;
  cout<<"avant updatematrix "<<endl;
#ifdef DEBUG 
  copy(sizeexceed_ret.begin(),sizeexceed_ret.end(),ostream_iterator<int>(cout," **--** "));cout<<endl;
#endif
#if 1
  //#ifdef GLOBALUPDATE
        global_update_matrix(mat,listcomp,mon2col);
#else
  for(auto c : listcomp)
    {
      update_matrix(mat,c,mon2col,indexmat);
      indexmat+=c.poly.size();
    }
#endif
  cout<<"apres updatematrix"<<endl;
  //affdebug(mat);
  //copy(sizeexceed_ret.begin(),sizeexceed_ret.end(),ostream_iterator<int>(cout," **--** "));cout<<endl;
#endif 


   mat.nrow+=crochet.poly.size();
#ifdef DEBUG
 cout<<"verytmp"<<endl;
  for(auto &p :verytmpPk)
    cout<<invconv<typename typP::value_type>(p,serv)<<endl;
#endif
  Pk=verytmpPk;
  Mk=tmpMk;
  return flagdoplus;
}
#endif
//Fonction qui met a jour la matrice globale lorsque des sous matrices on ete de rang 
//deficient
template<typename typmat, typename typlistcomp, typename typmap>
void global_update_matrix(typmat &m, typlistcomp &listcomp, typmap &mon2col)
{
  typedef typename typmat::coeff_t coeff;
  int *shiftcol=new int[m.ncol], *addcol=new int[m.ncol];
  int compnnz=0;
  int  tmpligne=0;
#ifdef DEBUG
  cout<<"en entree"<<endl;
  for(int i=0;i<m.nnz;i++)
    cout<<((coeff*)m.nzval)[i]<<"/"<<i<<" ";
  cout<<endl;
  for(int i=0;i<m.nnz;i++)
    cout<<m.rowind[i]<<"/"<<i<<" ";
  cout<<endl;
  cout<<"----------------------------------------"<<endl;
#endif
  for(int i=0;i<m.ncol;i++)
    shiftcol[i]=0;
  //  cout<<"boucle1"<<endl;
  for(auto c : listcomp)
    {
      auto iterexceed=c.exceed.begin();
      auto itersizeexceed=c.sizeexceed.begin();
      auto iterpoly=c.poly.begin();
      for(;iterexceed!=c.exceed.end();iterexceed++,itersizeexceed++, iterpoly++)
    {
      //on compte colonne par colonne les ajouts a faire
      for(int i=0;i<*itersizeexceed;i++)
        if (!Iszero((*iterexceed)[i].GetCoeff()))
          {
        shiftcol[mon2col[(*iterexceed)[i].rep]]++;
        //on compte le nombre total de nonzeor rajoutés
        compnnz++;
          }
    }
    }

  //On realloue la matrice
  m.nzval=MAC_REV_REALLOC<coeff>(m.nzval,m.nnz*sizeof(coeff),
                 (m.nnz+compnnz)*sizeof(coeff));
  m.rowind=(int*)MAC_REV_REALLOC<int>((int*)m.rowind,m.nnz*sizeof(int),
                      (m.nnz+compnnz)*sizeof(int));
  // if(compnnz==0) return;
  //on calcul les shift absolu de chaque colonne
  //  cout<<"boucle2"<<endl;
  for(int i=0;i<m.ncol;i++)
    addcol[i]=shiftcol[i];
  for(int i=1;i<m.ncol;i++)
      shiftcol[i]+=shiftcol[i-1];
  //On decale les colonnes
  for(int i=m.ncol-1;i>0;i--)
    {
      for(int j=m.colptr[i+1]-1;j>=m.colptr[i];j--)
    {
      // cout<<j<<" "<<  ((coeff*)m.nzval)[j]<<endl;
      ((coeff*)m.nzval)[j+shiftcol[i-1]]=((coeff*)m.nzval)[j];
      m.rowind[j+shiftcol[i-1]]=m.rowind[j];
    }
    }
  //  cout<<"boucle3"<<endl;
  //on reparcour tout pour faire les mises a jour
  for(auto c : listcomp)
    {
      auto iterexceed=c.exceed.begin();
      auto itersizeexceed=c.sizeexceed.begin();
      auto iterpoly=c.poly.begin();
      for(;iterexceed!=c.exceed.end();iterexceed++,itersizeexceed++, iterpoly++,tmpligne++)
    {
      for(int i=0;i<*itersizeexceed;i++)
        {
          if(!Iszero((*iterexceed)[i].GetCoeff()))
        {
          int numcol=mon2col[(*iterexceed)[i].rep];
          m.rowind[m.colptr[numcol+1]+shiftcol[numcol]-1]=tmpligne;
          ((coeff*)m.nzval)[m.colptr[numcol+1]+shiftcol[numcol]-- -1]=(*iterexceed)[i].GetCoeff();
          //  cout<<"G ecrit "<<(*iterexceed)[i].GetCoeff()<<endl;
          //cout<<"en colonne "<<numcol<<" ligne "<<m.colptr[numcol+1]+shiftcol[numcol]<<" "<<tmpligne<<endl;
          
        }
        }
    }
    }
  m.nnz+=compnnz;
#ifdef DEBUG
  cout<<"compnnz "<<compnnz<<endl;
  for(int i=0;i<m.nnz;i++)
    cout<<((coeff*)m.nzval)[i]<<"/"<<i<<" ";
  cout<<endl;
  for(int i=0;i<m.nnz;i++)
    cout<<m.rowind[i]<<"/"<<i<<" ";
  cout<<endl;
#endif
  for(int i=1;i<m.ncol;i++)
    addcol[i]+=addcol[i-1];
  for(int i=1;i<=m.ncol;i++)
    m.colptr[i]+=addcol[i-1];
  delete[] shiftcol;
  delete[] addcol;
}

template<typename typmat, typename typcomp, typename typmap>
void update_matrix(typmat &m, typcomp &c, typmap &mon2col, int indexmat)
{
  //On parcourt la composante
  //s'il reste des monomes dans les exceeds c'est des monomes de colonnes 
  //de pivot nul
  typedef typename typmat::coeff_t coeff;
  auto iterexceed=c.exceed.begin();
  auto itersizeexceed=c.sizeexceed.begin();
  auto iterpoly=c.poly.begin();
  int numligne=0;
  for(;iterexceed!=c.exceed.end();iterexceed++,itersizeexceed++, iterpoly++,numligne++)
    {
      //  cout<<"itersizeexceed "<<*itersizeexceed<<endl;
      m.nzval=MAC_REV_REALLOC<coeff>(m.nzval,m.nnz*sizeof(coeff),
                    (m.nnz+*itersizeexceed)*sizeof(coeff));
      m.rowind=(int*)MAC_REV_REALLOC<int>((void*)m.rowind,m.nnz*sizeof(int),
                    (m.nnz+*itersizeexceed)*sizeof(int));
      for(int i=0;i<*itersizeexceed;i++)
    {
      //  cout<<(*iterexceed)[i]<<endl;
      if((*iterexceed)[i].rep!=iterpoly->ind.rep)//Ce test a Virer 
        //Qd debug termineé     {
          for(int j=m.nnz;j>=m.colptr[mon2col[(*iterexceed)[i].rep]+1];j--)
        {
          ((coeff*)m.nzval)[j]=((coeff*)m.nzval)[j-1];
          m.rowind[j]=m.rowind[j-1];
        }  
          //   cout<<"mon2col "<<(*iterexceed)[i].rep<<" "<<m.colptr[mon2col[(*iterexceed)[i].rep]+1]<<endl;
          m.nnz++;
          m.rowind[m.colptr[mon2col[(*iterexceed)[i].rep]+1]]=indexmat+numligne;
          //ICI MODIF DU 11 JANVIER
          //OUILLE OUILLE OUILLE
//        m.rowind[m.colptr[mon2col[(*iterexceed)[i].rep]+1]]=indexmat+i;
          ((coeff*)m.nzval)[m.colptr[mon2col[(*iterexceed)[i].rep]+1]]=(*iterexceed)[i].GetCoeff();
#ifdef DEBUG
          cout<<"G ecrit "<<(*iterexceed)[i].GetCoeff()<<endl;
          cout<<"en colonne "<<mon2col[(*iterexceed)[i].rep]<<" ligne "<<m.colptr[mon2col[(*iterexceed)[i].rep]+1]<<" "<<indexmat+numligne<<endl;
#endif
          for(int j=mon2col[(*iterexceed)[i].rep]+1;j<=m.ncol;j++)
        m.colptr[j]++;
          //les colonnes d'apres commencent maintenant un cran plus ooin
          
        }
    }
      //MAC_REV_FREE<mon>(*iterexceed,*itersizeexceed*sizeof(mon));
      //ICI c dans ADDB qu'on fait le FREE
      //
    }

}
//
// Les polyn\^omes des autres composantes fortement connexes sont mis à jour en 
// fonction du resultat de l'inversion de la composante en cours
//
template <typename typcomp, typename typserv>
void  update_pol(typcomp &comp_to_update,typcomp &updating_comp, typserv &serv)
{
  typedef typename typcomp::typPk typPk;
  typedef typename typserv::monom_t mon;
  typedef typename typcomp::pol pol;
  typedef typename pol::coeff_t coeff_t;
//  cout<<"rentre dans update"<<endl;
  if(comp_to_update==updating_comp)
    return;
//  cout<<"rentre dans compo pas les memes"<<endl;
  auto iterexceed=updating_comp.exceed.begin();
  auto itersizeexceed=updating_comp.sizeexceed.begin();
  auto iterpoly=updating_comp.poly.begin();
  //
  //
  
  auto myitersize=updating_comp.sizeexceed.begin();
  #ifdef DEBUG
  for(auto &i : updating_comp.poly)
    cout<<"ind qui update "<<i.ind<<" "<<*myitersize++<<endl;
#endif
  myitersize=comp_to_update.sizeexceed.begin();
 #ifdef DEBUG
  for(auto &i : comp_to_update.poly)
    cout<<"ind a up "<<i.ind<<" "<<*myitersize++<<endl;
#endif
  //on créé l'association enter les ind et les structures correspondantes
  map<typename mon::rep_t,resmul<typename typPk::value_type> > mon2pol;
  for(;iterpoly!=updating_comp.poly.end();iterpoly++,iterexceed++,itersizeexceed++)
    {
      resmul<typename typPk::value_type> tmp;
      tmp.poly=*iterpoly;
      tmp.exceed=*iterexceed;
      tmp.sizeexceed=*itersizeexceed;
      mon2pol[iterpoly->ind.rep]=tmp;
      //on augmente le workspace pour update apres
      //for(;w->sizenf<iterpoly->sizenf;increase_space(w));
    }
  iterexceed=comp_to_update.exceed.begin();
  itersizeexceed=comp_to_update.sizeexceed.begin();
  for(iterpoly=comp_to_update.poly.begin();
    iterpoly!=comp_to_update.poly.end();
    iterpoly++,iterexceed++,itersizeexceed++)
    //pour chaque polynome dans la composante a mettre a jour
    {
      //On construit une map pour l'exceed correspondant
      int k=0;
      map<typename mon::rep_t,COEFF> exceedmap;
      map<typename mon::rep_t,COEFF> nfmap;//il faudrait utiliser le workspace ici au lieu de cette map
 //     cout<<"itersizexceed ici "<<*itersizeexceed<<endl;
      int pos_nf=0;
      int newsizenfind=iterpoly->sizenf;
      //augmente ici le workspace pour que la mise a jour
      //puisse se faire
      //for(;w->sizenf<iterpoly->sizenf;increase_space(w));
      for(int i=0;i<iterpoly->sizenf;i++)
    if(iterpoly->nfind)
    {
      for(int j=0;j<8;j++)
        if(iterpoly->nfind[i]>>j&1)
        {
          mon tmp;
          serv.int2mon(8*i+j,tmp);
          nfmap[tmp.rep]=iterpoly->nf[pos_nf++];
        }
    }
      for(int i=0;i<*itersizeexceed;i++)
    {
      //On ne rajoute pas dans exceedmap 
      //les monomes de mon2pol!
      if(!mon2pol.count((*iterexceed)[i].rep))
       exceedmap[(*iterexceed)[i].rep]= 
         (*iterexceed)[i].GetCoeff();
    }
      for(int i=0;i<*itersizeexceed;i++)
    {
      if(mon2pol.count((*iterexceed)[i].rep))
        //On a trouvé un monome de la composante 
        //qu'on utilise pour mettre a jour
        {
          const resmul<typename typPk::value_type> &tmpresmul=mon2pol[(*iterexceed)[i].rep];
//        cout<<"je declenche une mise a jour par "<<tmpresmul.poly.ind<<endl;
//        cout<<"   de "<<iterpoly->ind<<endl;
          if(newsizenfind<tmpresmul.poly.sizenf)
        newsizenfind=tmpresmul.poly.sizenf;
          //
          //A VERIFIER QUE
          //les exceedd  sont modifié comme il faut pendant le solve
          //fait
//        cout<<"le coeff maj "<<(*iterexceed)[i].GetCoeff()<<endl;
//        cout<<"exceed avant "<<*itersizeexceed<<" tmpresmul.sizeexceed "<<tmpresmul.sizeexceed<<endl;
//        for(auto &i : exceedmap)
//      cout<<i.first<<" "<<i.second<<" ";
//        cout<<endl;
          for(int j=0;j<tmpresmul.sizeexceed;j++)
        {
          exceedmap[tmpresmul.exceed[j].rep]-=
            tmpresmul.exceed[j].GetCoeff()*(*iterexceed)[i].GetCoeff();
        }
//        cout<<"nf avant"<<endl;
//        for(auto &i : nfmap)
//      cout<<i.first<<" "<<i.second<<" ";
//        cout<<endl;
        int pos_nfiter=0;
        for(int ii=0;ii<tmpresmul.poly.sizenf;ii++)
            if(tmpresmul.poly.nfind)
            {
              for(int j=0;j<8;j++)
              if(tmpresmul.poly.nfind[ii]>>j&1)
                {
                  //ICI PASSER PAR LE WORKSPACE!!
                  //ca evite mon2int et la map derriere!
                  mon tmp;
                  serv.int2mon(8*ii+j,tmp);
                 #if 0
                  cout<<"mon mis a jour "<<tmp<<" "<<tmpresmul.poly.nf[pos_nfiter]<<endl;
                  cout<<"nfmap "<<nfmap[tmp.rep]<<endl;
                  cout<<tmpresmul.poly.nf[pos_nfiter]<<endl;
                  cout<<(*iterexceed)[i].GetCoeff()<<endl;
                  cout<<coeff_t(0)-tmpresmul.poly.nf[pos_nfiter]
                                     *(*iterexceed)[i].GetCoeff()<<endl;
                  cout<<nfmap[tmp.rep]-tmpresmul.poly.nf[pos_nfiter]
                                     *(*iterexceed)[i].GetCoeff()<<endl;
                  #endif
                  nfmap[tmp.rep]-=tmpresmul.poly.nf[pos_nfiter]
                *(*iterexceed)[i].GetCoeff();
              //    cout<<"apres operation "<<nfmap[tmp.rep]<<endl;
                pos_nfiter++;

                }
            }
#ifdef DEBUG
          cout<<"exceed apres"<<endl;
          for(auto &i : exceedmap)
        cout<<i.first<<" "<<i.second<<" ";
          cout<<endl;
          cout<<"nf apres"<<endl;
          for(auto &i : nfmap)
        cout<<i.first<<" "<<i.second<<" ";
          cout<<endl;
#endif       
        }
    }   
#ifdef DEBUG
      cout<<"nfmap.size() "<<nfmap.size()<<endl;
      cout<<"newsizenfind "<<newsizenfind<<" "<<(void*)iterpoly->nfind<<endl;
#endif
      *iterexceed=(mon*)MAC_REV_REALLOC<mon>(*iterexceed,*itersizeexceed*sizeof(mon)
                         ,exceedmap.size()*sizeof(mon));
      iterpoly->nfind=(unsigned char *)realloc(iterpoly->nfind,newsizenfind);
      if(iterpoly->size != nfmap.size())
       iterpoly->nf=(coeff_t*)MAC_REV_REALLOC<coeff_t>(iterpoly->nf,
                        iterpoly->size*sizeof(coeff_t),
                    nfmap.size()*sizeof(coeff_t));
      for(int k=0;k<newsizenfind;k++)
    iterpoly->nfind[k]=0;
      for(auto iter : nfmap)
      {
    mon m(1);
    m.rep=iter.first;
    int place=serv.mon2int(m);
    iterpoly->nfind[place/8]|=(1<<(place%8));
      }
      int posnf=0;
      for(int k=0;k<newsizenfind;k++)
      {
    if(iterpoly->nfind[k])
    {
      for(int kk=0;kk<8;kk++)
        if((iterpoly->nfind[k]>>kk)&1)
         {
           mon m;
           serv.int2mon(8*k+kk,m);
           iterpoly->nf[posnf++]=nfmap[m.rep];
         }
    }
      }
      iterpoly->size=nfmap.size();
      iterpoly->sizenf=newsizenfind;
      for(auto iter : exceedmap)
    {
      mon tmpmon;
      tmpmon.rep=iter.first;
      tmpmon.coeff=iter.second;
      (*iterexceed)[k++]=tmpmon;
    }
      *itersizeexceed=exceedmap.size();
    }
}

template <typename typcomp, typename typserv>
void  wrap_update_pol(typcomp &comp_to_update,typcomp &updating_comp, typserv &serv)
{
  workspace<typename typcomp::pol::coeff_t> w;
  update_pol(comp_to_update,updating_comp,w,serv);
  destroy_space(&w);
}
template <typename typcomp, typename typwork, typename typserv>
void  update_pol(typcomp &comp_to_update,typcomp &updating_comp, typwork &w,
         typserv &serv)
{
  typedef typename typcomp::typPk typPk;
  typedef typename typserv::monom_t mon;
  typedef typename typcomp::pol pol;
  typedef typename pol::coeff_t coeff_t;
//  cout<<"rentre dans update"<<endl;
  if(comp_to_update==updating_comp)
    return;
//  cout<<"rentre dans compo pas les memes"<<endl;
  auto iterexceed=updating_comp.exceed.begin();
  auto itersizeexceed=updating_comp.sizeexceed.begin();
  auto iterpoly=updating_comp.poly.begin();
  //
  //
  
  auto myitersize=updating_comp.sizeexceed.begin();
  #ifdef DEBUG
  for(auto &i : updating_comp.poly)
    cout<<"ind qui update "<<i.ind<<" "<<*myitersize++<<endl;
#endif
  myitersize=comp_to_update.sizeexceed.begin();
 #ifdef DEBUG
  for(auto &i : comp_to_update.poly)
    cout<<"ind a up "<<i.ind<<" "<<*myitersize++<<endl;
#endif
  //on créé l'association enter les ind et les structures correspondantes
  map<typename mon::rep_t,resmul<typename typPk::value_type> > mon2pol;
  for(;iterpoly!=updating_comp.poly.end();iterpoly++,iterexceed++,itersizeexceed++)
    {
      resmul<typename typPk::value_type> tmp;
      tmp.poly=*iterpoly;
      tmp.exceed=*iterexceed;
      tmp.sizeexceed=*itersizeexceed;
      mon2pol[iterpoly->ind.rep]=tmp;
      //on augmente le workspace pour update apres
      for(;w.sizenf<iterpoly->sizenf;increase_space(&w));
    }
  iterexceed=comp_to_update.exceed.begin();
  itersizeexceed=comp_to_update.sizeexceed.begin();
  for(iterpoly=comp_to_update.poly.begin();
    iterpoly!=comp_to_update.poly.end();
    iterpoly++,iterexceed++,itersizeexceed++)
    //pour chaque polynome dans la composante a mettre a jour
    {
      //On construit une map pour l'exceed correspondant
      int k=0;
      map<typename mon::rep_t,COEFF> exceedmap;
      map<typename mon::rep_t,COEFF> nfmap;//il faudrait utiliser le workspace ici au lieu de cette map
 //     cout<<"itersizexceed ici "<<*itersizeexceed<<endl;
      int pos_nf=0;
      int newsizenfind=iterpoly->sizenf;
      //augmente ici le workspace pour que la mise a jour
      //puisse se faire
      for(;w.sizenf<iterpoly->sizenf;increase_space(&w));
      for(int i=0;i<iterpoly->sizenf;i++)
    if(iterpoly->nfind)
    {
      w.tabnf[i]=iterpoly->nfind[i];
      tabfun1[w.tabnf[i]](w.tabcoeff+(8*i),iterpoly->nf+pos_nf);
      pos_nf+=nbbits[iterpoly->nfind[i]];
    }
      for(int i=0;i<*itersizeexceed;i++)
    {
      //On ne rajoute pas dans exceedmap 
      //les monomes de mon2pol!
      if(!mon2pol.count((*iterexceed)[i].rep))
       exceedmap[(*iterexceed)[i].rep]= 
         (*iterexceed)[i].GetCoeff();
    }
      for(int i=0;i<*itersizeexceed;i++)
    {
      if(mon2pol.count((*iterexceed)[i].rep))
        //On a trouvé un monome de la composante 
        //qu'on utilise pour mettre a jour
        {
          resmul<typename typPk::value_type> &tmpresmul=mon2pol[(*iterexceed)[i].rep];
          if(newsizenfind<tmpresmul.poly.sizenf)
        newsizenfind=tmpresmul.poly.sizenf;
          for(int j=0;j<tmpresmul.sizeexceed;j++)
        {
          exceedmap[tmpresmul.exceed[j].rep]-=
            tmpresmul.exceed[j].GetCoeff()*(*iterexceed)[i].GetCoeff();
        }
        int pos_nfiter=0;
        for(int ii=0;ii<tmpresmul.poly.sizenf;ii++)
            if(tmpresmul.poly.nfind)
            {
            coeff_t tmp=(*iterexceed)[i].GetCoeff();    
            w.tabnf[ii]=tmpresmul.poly.nfind[ii];
            tabfun2[w.tabnf[ii]](w.tabcoeff+(8*ii),
                &tmp,
                tmpresmul.poly.nf+pos_nfiter);
            pos_nfiter+=nbbits[tmpresmul.poly.nfind[ii]];

            }
#ifdef DEBUG
          cout<<"exceed apres"<<endl;
          for(auto &i : exceedmap)
        cout<<i.first<<" "<<i.second<<" ";
          cout<<endl;
          cout<<"nf apres"<<endl;
          for(auto &i : nfmap)
        cout<<i.first<<" "<<i.second<<" ";
          cout<<endl;
#endif       
        }
    }   
#ifdef DEBUG
      cout<<"newsizenfind "<<newsizenfind<<" "<<(void*)iterpoly->nfind<<endl;
      cout<<"w.sizenf "<<w.sizenf<<endl;
#endif
      *iterexceed=(mon*)MAC_REV_REALLOC<mon>(*iterexceed,*itersizeexceed*sizeof(mon)
                         ,exceedmap.size()*sizeof(mon));
      iterpoly->nfind=(unsigned char *)realloc(iterpoly->nfind,newsizenfind);
      int newsize=0;
      for(int i=0;i<w.sizenf;newsize+=nbbits[w.tabnf[i++]]);
       iterpoly->nf=(coeff_t*)MAC_REV_REALLOC<coeff_t>(iterpoly->nf,
                        iterpoly->size*sizeof(coeff_t),
                    newsize*sizeof(coeff_t));
      for(int k=0;k<newsizenfind;k++)
    {
      iterpoly->nfind[k]=w.tabnf[k];
      w.tabnf[k]=0;
    }
      for(int k=0;k<newsize;k++)
      {
    iterpoly->nf[k]=w.tabcoeff[k];
    w.tabcoeff[k]=0;
      }
      iterpoly->size=newsize;
      iterpoly->sizenf=newsizenfind;
      for(auto iter : exceedmap)
    {
      mon tmpmon;
      tmpmon.rep=iter.first;
      tmpmon.coeff=iter.second;
      (*iterexceed)[k++]=tmpmon;
    }
      *itersizeexceed=exceedmap.size();
    }
}
template <typename typcomp, typename typmon2pol, typename typserv>
void  wrap_update_pol(typcomp &comp_to_update,const typcomp &updating_comp, 
    const typmon2pol &mon2pol, const typserv &serv)
{
  workspace<typename typcomp::pol::coeff_t> w;
  w.sizenf=0;
  w.size=0;
  w.tabnf=NULL;
  w.tabcoeff=NULL;
  for(auto &p : updating_comp.poly)
    for(;w.sizenf<=p.sizenf;increase_space(&w));
  update_pol(comp_to_update,updating_comp,mon2pol,w,serv);
  destroy_space(&w);
}
template <typename typcomp, typename typmon2pol, typename typwork, typename typserv>
void  update_pol(typcomp &comp_to_update,const typcomp &updating_comp, 
    const typmon2pol &mon2pol, typwork &w, const typserv &serv)
{
  typedef typename typcomp::typPk typPk;
  typedef typename typserv::monom_t mon;
  typedef typename typcomp::pol pol;
  typedef typename pol::coeff_t coeff_t;
//  cout<<"rentre dans update"<<endl;
  if(comp_to_update==updating_comp)
    return;
//  cout<<"rentre dans compo pas les memes"<<endl;
  auto iterexceed=comp_to_update.exceed.begin();
  auto itersizeexceed=comp_to_update.sizeexceed.begin();
  auto iterpoly=comp_to_update.poly.begin();
  //
  //
  
  auto myitersize=updating_comp.sizeexceed.begin();
#ifdef DEBUG_UPDATE
  for(auto &i : updating_comp.poly)
    cout<<"ind qui update "<<i.ind<<" "<<*myitersize++<<endl;
#endif
  myitersize=comp_to_update.sizeexceed.begin();
#ifdef DEBUG_UPDATE
  for(auto &i : comp_to_update.poly)
    cout<<"ind a up "<<i.ind<<" "<<*myitersize++<<endl;
#endif
  //on créé l'association enter les ind et les structures correspondantes
  iterexceed=comp_to_update.exceed.begin();
  itersizeexceed=comp_to_update.sizeexceed.begin();
  for(iterpoly=comp_to_update.poly.begin();
    iterpoly!=comp_to_update.poly.end();
    iterpoly++,iterexceed++,itersizeexceed++)
    //pour chaque polynome dans la composante a mettre a jour
    {
       //On construit une map pour l'exceed correspondant
      int k=0;
      map<typename mon::rep_t,COEFF> exceedmap;
 //     cout<<"itersizexceed ici "<<*itersizeexceed<<endl;
      int pos_nf=0;
      int newsizenfind=iterpoly->sizenf;
      //augmente ici le workspace pour que la mise a jour
      //puisse se faire
      for(;w.sizenf<=iterpoly->sizenf;increase_space(&w));
      for(int i=0;i<iterpoly->sizenf;i++)
    if(iterpoly->nfind)
    {
      w.tabnf[i]=iterpoly->nfind[i];
      tabfun1[w.tabnf[i]](w.tabcoeff+(8*i),iterpoly->nf+pos_nf);
      pos_nf+=nbbits[iterpoly->nfind[i]];
    }
      for(int i=0;i<*itersizeexceed;i++)
    {
      //On ne rajoute pas dans exceedmap 
      //les monomes de mon2pol!
      if(!mon2pol.count((*iterexceed)[i].rep))
       exceedmap[(*iterexceed)[i].rep]= 
         (*iterexceed)[i].GetCoeff();
    }
      for(int i=0;i<*itersizeexceed;i++)
    {
      if(mon2pol.count((*iterexceed)[i].rep))
        //On a trouvé un monome de la composante 
        //qu'on utilise pour mettre a jour
        {
          const resmul<typename typPk::value_type> &tmpresmul=
                             mon2pol.find((*iterexceed)[i].rep)->second;
          if(newsizenfind<tmpresmul.poly.sizenf)
        newsizenfind=tmpresmul.poly.sizenf;
          for(int j=0;j<tmpresmul.sizeexceed;j++)
        {
          exceedmap[tmpresmul.exceed[j].rep]-=
            tmpresmul.exceed[j].GetCoeff()*(*iterexceed)[i].GetCoeff();
        }
        int pos_nfiter=0;
        for(int ii=0;ii<tmpresmul.poly.sizenf;ii++)
            if(tmpresmul.poly.nfind[ii])
            {
            coeff_t tmp=(*iterexceed)[i].GetCoeff();    
            w.tabnf[ii]|=tmpresmul.poly.nfind[ii];
            //          tabfun2[w.tabnf[ii]](w.tabcoeff+(8*ii),
            tabfun2[tmpresmul.poly.nfind[ii]](w.tabcoeff+(8*ii),
                     &tmp,
                     tmpresmul.poly.nf+pos_nfiter);
            pos_nfiter+=nbbits[tmpresmul.poly.nfind[ii]];

            }
        }
    }   
#ifdef DEBUG_UPDATE
      cout<<"newsizenfind "<<newsizenfind<<" "<<(void*)iterpoly->nfind<<endl;
      cout<<"w.sizenf "<<w.sizenf<<endl;
#endif
      *iterexceed=(mon*)MAC_REV_REALLOC<mon>(*iterexceed,*itersizeexceed*sizeof(mon)
                         ,exceedmap.size()*sizeof(mon));
      iterpoly->nfind=(unsigned char *)realloc(iterpoly->nfind,newsizenfind);
      int newsize=0;
      for(int i=0;i<w.sizenf;newsize+=nbbits[w.tabnf[i++]]);
       iterpoly->nf=(coeff_t*)MAC_REV_REALLOC<coeff_t>(iterpoly->nf,
                        iterpoly->size*sizeof(coeff_t),
                    newsize*sizeof(coeff_t));
#ifdef DEBUG_UPDATE
       cout<<"newsize "<<newsize<<endl;
#endif
       for(int k=0;k<newsizenfind;k++)
    {
      iterpoly->nfind[k]=w.tabnf[k];
      w.tabnf[k]=0;
    }
       //      for(int k=0;k<newsize;k++)
       // {
       //   iterpoly->nf[k]=w.tabcoeff[k];
       //   w.tabcoeff[k]=0;
       //}
       int tmpcompt=0;
       for(int k=0;k<newsizenfind;k++)
     {
       if(iterpoly->nfind[k])
         {
           for(int j=0;j<8;j++)
         if((iterpoly->nfind[k]>>j)&1)
           {
             iterpoly->nf[tmpcompt++]=w.tabcoeff[8*k+j];
             w.tabcoeff[8*k+j]=0;
           }
         }
     }
      iterpoly->size=newsize;
      iterpoly->sizenf=newsizenfind;
#ifdef DEBUG_UPDATE
      cout<<"exceed apres "<<exceedmap.size()<<endl;
      for(auto &i : exceedmap)
    cout<<i.first<<" "<<i.second<<" ";
      cout<<endl;
#endif      
      for(auto iter : exceedmap)
    {
      mon tmpmon;
      tmpmon.rep=iter.first;
      tmpmon.coeff=iter.second;
      (*iterexceed)[k++]=tmpmon;
    }
      *itersizeexceed=exceedmap.size();
      // cout<<"---------------------------"<<endl;
    }
}
//est sense reconstruire la matrice mat a partir 
//des matrices correspondant aux composantes fortmeùent connexes

template<typename typMat>
void reconstruct_mat(typMat &mat, const typMat &tmpmat)
{
  typedef typename typMat::coeff_t coeff_t;
#ifdef DEBUG
  cout<<"la petitee "<<tmpmat.nnz<<endl;
  affdebug(tmpmat);
  cout<<"la grosse "<<endl;
  if(mat.ncol)
    affdebug(mat);
#endif
  mat.colptr=(int *)realloc(mat.colptr,(mat.ncol+tmpmat.ncol+1)*sizeof(int));
  mat.nzval=MAC_REV_REALLOC<coeff_t>(mat.nzval,mat.nnz*sizeof(coeff_t),
    (mat.nnz+tmpmat.nnz)*sizeof(coeff_t));
  mat.rowind=(int *)realloc(mat.rowind,(mat.nnz+tmpmat.nnz)*sizeof(int));
  //
  for(int i=0;i<=tmpmat.ncol;i++)
    mat.colptr[mat.ncol+i]=mat.nnz+tmpmat.colptr[i];
  //
  for(int i=0;i<tmpmat.nnz;i++)
    mat.rowind[mat.nnz+i]=mat.nrow+tmpmat.rowind[i];
  //
  for(int i=0;i<tmpmat.nnz;i++)
    ((coeff_t *)mat.nzval)[mat.nnz+i]=((coeff_t *)tmpmat.nzval)[i];
  mat.ncol+=tmpmat.ncol;
  mat.nrow+=tmpmat.nrow;
  mat.nnz+=tmpmat.nnz;
#ifdef DEBUG
   cout<<"dans reconstruct mat.ncol en sortie "<<mat.ncol<<endl;
   cout<<"dans reconstruct mat.nrow en sortie "<<mat.nrow<<endl;
#endif
}


void reconstruct_pr_pc(int *&pr,int *&pc,int *tmppr,int *tmppc,int sizetmp,int indexinside)
{
  for(int i=0;i<sizetmp;i++)
    pr[i+indexinside]=tmppr[i]+indexinside;
  
  for(int i=0;i<sizetmp;i++)
    pc[i+indexinside]=tmppc[i]+indexinside;
}

template<typename typMat,typename typP,typename typPk,typename typMk, typename Base,typename typdump, class monomial_server, typename workspace>
void init_algo(typP &P,int &nouveaumon, int &maxdeg, int & k,
           typPk &Pk, typMk & Mk, Base & b, typdump & dump, workspace &w,
           monomial_server &serv);

template<typename typP,typename typPk,typename typdump,typename Base
     , template<typename  T> class tmptypmat ,
     template<typename mon, typename pol> class monomial_server, typename mon>
void algo(typP P,typdump &dump,Base &b, 
      monomial_server<mon,pol2ter<mon,typename mon::coeff_t> > &serv)
{
  typedef typename typP::value_type polyalp;
  typedef list<typename typPk::value_type::monom_t> typMk;
  //  typedef typename typMk::value_type mon;
  typedef typename typPk::value_type pol;

  typMk Mk;
  typMk Spoldejavu;
  typPk Pk,S;
  typPk recall;
  workspace<typename typP::value_type::coeff_t> w;
  typedef tmptypmat<typename typP::value_type::coeff_t> typMat;
  typMat Mat,L;
  int nouveaumon,allpolused=0,maxdeg,k,*pr,*pc;
  vector<typename typMk::value_type*> exceed;
  list<int> sizeexceed;
  //cout<<"avant init"<<endl;
  
  init_algo<typMat>(P,nouveaumon,maxdeg,k,Pk,Mk,b,dump,w,serv);
  //cout<<"apres init"<<endl;
#ifdef GOTZMANN
  int nextdim=-1;
  initial_mon(k,b,dump,serv);
  while((!critere_got(b,dump,serv))||(k<=maxdeg))
#else  
    while(nouveaumon||(k<=maxdeg))
#endif
    {
    int killspol;
    //    cout<<"Spolordsup"<<endl;
    S=Spolordsup(Pk,b,Spoldejavu);
#ifdef CACABOUDIN    
 for(typename typPk::iterator iter=S.begin();iter!=S.end();iter++)
      {
    //   cout<<"poly de taille "<<iter->size<<endl;
    cout<<invconv<typename typP::value_type>(*iter,serv)<<endl;       
      }
#endif
#if 0//ICICICICICICICICICICCII
cout<<"poly a la fin de la boulcle de l'algo "<<Pk.size()<<endl;
for(typename typPk::iterator iter=Pk.begin();iter!=Pk.end();iter++)
 {
   //   cout<<"poly de taille "<<iter->size<<endl;
   cout<<invconv<typename typP::value_type>(*iter,serv)<<endl;       
}
    cout<<endl<<"newpkmkmatrixof"<<endl;
#endif
    killspol=!newpkmkMatrixofSolve(Mat,L,Pk,P,Mk,k,b,dump,exceed,sizeexceed
                   ,w,serv);
#if 0//ICICICICICICICICCI
    cout<<"poly apres solve (attention exceed)"<<endl;
    copy(sizeexceed.begin(),sizeexceed.end(),ostream_iterator<int>(cout," *** "));cout<<endl;
    for(typename typPk::iterator iter=Pk.begin();iter!=Pk.end();iter++)
      {
    
    cout<<invconv<typename typP::value_type>(*iter,serv)<<endl;       
      }
#endif
    if(killspol) S.erase(S.begin(),S.end());
    if (k+1>=maxdeg) allpolused=1;
    if(Mk.size()==0) nouveaumon=0;
    //    cout<<"le rang "<<my_rank(Mat)<<endl;
    // cout<<"Mk.size() "<<Mk.size()<<endl;
   // copy(Mk.begin(),Mk.end(),ostream_iterator<mon>(cout," "));
   // cout<<endl;
    if ((unsigned)my_rank(Mat)==Mk.size())
      //matrice de rang plein
      //donc tous les exceed sont vides
      {
    typMk tmpSpoldejavu;
    typPk redspol,secd=secondmembre(Mat,Pk);;
    //secondmembre2(Mat,Pk,secd);
    //typename typPk::iterator iter=Pk.begin();
    cout<<"matrice de rang plein "<<endl;
    Dump(Pk,dump);
#ifdef GOTZMANN
    nextmon(k,b,dump,serv);
#endif
    cout<<"secd/size() "<<secd.size()<<endl;
    ReduceSpol(S,dump,redspol,b,w,serv);//tmpSpoldejavu);
    my_merge(redspol,secd);
    cout<<"redspol.size "<<redspol.size()<<endl;
    if(redspol.size()==0)
      {
        k++;
        my_merge(Spoldejavu,tmpSpoldejavu);
        L.destroystore();
        continue;
        //la matrice a rang plein 
        //ET les Spol\secondmembre de ce niveau se ,reduisent a 0
      }
    else
      {
        //on a des spol\secondmembre qui se reduisent pas a 0
        //cad un defaut de commuatation
        int tmpdeg=mindeg(redspol,serv);
        max_degppcm=k+2;
        nouveaumon=1;
        cout<<"redspol.size() "<<redspol.size()<<endl;

        if (tmpdeg==(k+1))
          my_merge(Spoldejavu,tmpSpoldejavu);
        k=tmpdeg;

        recover(Pk,Mk,dump,k);

        forget(b,dump,k,recall);
        cout<<"debut newcrochetchoix"<<endl;
        NewCrochetchoix(P,Pk,Mk,redspol,k,dump,b,recall,w,serv);
        cout<<"fin newcrochetchoix"<<endl;

        NewDestroyredspol(redspol,Pk);
        serv.compress(Pk,dump,b,k);
        Dump(Pk,dump);
#ifdef GOTZMANN
        nextmon(k,b,dump,serv);
#endif
        L.destroystore();
        //misa jour de B dans crochet choix!!
        continue;
      }
      }
    else
      //la matrice n'est pas de rang plein
      {
    cout<<"matrice de rang pas plein "<<endl;
    typPk secd;
    cout<<"addB"<<endl;
#ifdef DEBUG11  
    cout<<"avant addB "<<Pk.size()<<endl;
    auto iterexceed=exceed.begin();
    auto itersizeexceed=sizeexceed.begin();
    for(typename typPk::iterator iter=Pk.begin();iter!=Pk.end();iter++,iterexceed++,itersizeexceed++)
      {
        cout<<"ind "<<iter->ind<<" est sikzeexceed "<<*itersizeexceed<<endl;
        copy(*iterexceed,*iterexceed+*itersizeexceed,ostream_iterator<mon>(cout," + "));
        cout<<endl;
        cout<<invconv<typename typP::value_type>(*iter,serv)<<endl;       
      }
#endif
    AddB(Mat,Pk,Mk,exceed,sizeexceed,dump,b,serv);
#if 0
    cout<<"apres addB "<<Pk.size()<<endl;
    for(typename typPk::iterator iter=Pk.begin();iter!=Pk.end();iter++)
      {
        cout<<invconv<typename typP::value_type>(*iter,serv)<<endl;       
      }
#endif
    //
    //secondmembre supprime des trucs dans Pk!!
    //
    //  cout<<"second membre"<<endl;
    secd=secondmembre(Mat,Pk);
    //  cout<<"voisin"<<endl;
    rajoute_voisins(Pk,Mk,dump,b,w,serv);
    //  cout<<"fin voisin"<<endl;
    //ne selectionne que les non nuls
       if(secd.size()==0)
     {
       typPk redspol;
       Dump(Pk,dump);
#ifdef GOTZMANN
       nextmon(k,b,dump,serv);
#endif
      
       ReduceSpol(S,dump,redspol,b,w,serv);//,Spoldejavu);
#ifdef NO_TWICE    
for(typename typPk::iterator iter=redspol.begin();iter!=redspol.end();iter++)
      P.push_back(invconv<polyalp>(*iter));
#endif
       if(redspol.size()==0)
         {
           k++;
           L.destroystore();
           continue;
           //la matrice a pas rang plein second membre =0 
           //ET les Spol de ce niveau se ,reduisent a 0
         }
       else
         {
           //on a des spol qui se reduisent pas a 0
           //cad un defaut de commutation 
           k=mindeg(redspol,serv);
           recover(Pk,Mk,dump,k);
           forget(b,dump,k,recall);
           cout<<"il ya un passage par mat rang pas"
               <<"plein et pol non zero"<<endl;
           NewCrochetchoix(P,Pk,Mk,redspol,k,dump,b,recall,w,serv);
           NewDestroyredspol(redspol,Pk);
           serv.compress(Pk,dump,b,k);
           Dump(Pk,dump);
#ifdef GOTZMANN
           nextmon(k,b,dump,serv);
#endif
 
           L.destroystore();
           continue;
         }
     }
       else
     {
       //on a des pol qui se reduisent pas a 0
       //cad un defaut de commuatation dans les coins
       typPk redspol;
       Dump(Pk,dump);
#ifdef GOTZMANN
       nextmon(k,b,dump,serv);
#endif
       ReduceSpol(S,dump,redspol,b,w,serv);//,Spoldejavu);
       my_merge(secd,redspol);
#ifdef NO_TWICE    
for(typename typPk::iterator iter=redspol.begin();iter!=redspol.end();iter++)
         P.push_back(invconv<polyalp>(*iter));
#endif
           k=mindeg(secd,serv);
       recover(Pk,Mk,dump,k);
       forget(b,dump,k,recall);
#if 0
        cout<<"poly avant newcrochetchoix (attention exceed)"<<endl;
        copy(sizeexceed.begin(),sizeexceed.end(),ostream_iterator<int>(cout," *** "));cout<<endl;
        for(typename typPk::iterator iter=Pk.begin();iter!=Pk.end();iter++)
          {
        
        cout<<invconv<typename typP::value_type>(*iter,serv)<<endl;       
          }
#endif

       NewCrochetchoix(P,Pk,Mk,secd,k,dump,b,recall,w,serv);
#if 0
        cout<<"poly apres newcrochetchoix (attention exceed)"<<endl;
        copy(sizeexceed.begin(),sizeexceed.end(),ostream_iterator<int>(cout," *** "));cout<<endl;
        for(typename typPk::iterator iter=Pk.begin();iter!=Pk.end();iter++)
          {
        
        cout<<invconv<typename typP::value_type>(*iter,serv)<<endl;       
          }
#endif

       NewDestroyredspol(secd,Pk);
       serv.compress(Pk,dump,b,k);
       Dump(Pk,dump);
#ifdef GOTZMANN
       nextmon(k,b,dump,serv);
#endif

       L.destroystore();
       continue;
     }
      
      }
  }
  
  Mat.destroystore();
  //while nouveau mon etc...
  destroy_space(&w);
  cout<<"sortie d'algo flash et dimension "<<endl;
  return ;//dump;
  //comme ca on conserve le tri par degre des formes normales et 
  //donc la construction des operateurs de multiplication se 
  //fait plus vite.. :-)

}


template<typename mon, typename T,typename Base ,typename monomial_server>
mon choice_naif(const pol2ter<mon,T> & p, const Base &b, monomial_server &serv)
{
  list<mon> tmp;
  mon m; 
  int maxdeg=-1,k=0;
  //int maxtmpdeg=0;
  for(int i=0;i<p.sizenf;i++)
    {
      if (p.nfind[i])
    for(int j=0;j<8;j++)
      if(p.nfind[i]>>j&1)
      {
        mon tmpmon;
        serv.int2mon(8*i+j,tmpmon);
        if(tmpmon.GetDegree()>maxdeg)
          {
        tmp.erase(tmp.begin(),tmp.end());
        maxdeg=tmpmon.GetDegree();
          }
        if(tmpmon.GetDegree()==maxdeg)
          {
        tmpmon.SetCoeff(p.nf[k]);
        tmp.push_back(tmpmon);
          }
        k++;
      }
    }
  m=tmp.front();
  //  cout<<"monom choisi "<<m<<endl;
  return m;
}

template<typename mon, typename T,typename Base, typename monomial_server>
mon choice_num(const pol2ter<mon,T> & p, const Base &b, monomial_server &serv)
{
  list<mon> tmp;
  mon m; 
  int maxdeg=-1,k=0;
  //int maxtmpdeg=0;
  for(int i=0;i<p.sizenf;i++)
    {
      if (p.nfind[i])
    for(int j=0;j<8;j++)
      if(p.nfind[i]>>j&1)
      {
        mon tmpmon;
        serv.int2mon(8*i+j,tmpmon);
        if(tmpmon.GetDegree()>maxdeg)
          {
        tmp.erase(tmp.begin(),tmp.end());
        maxdeg=tmpmon.GetDegree();
          }
        if(tmpmon.GetDegree()==maxdeg)
          {
        tmpmon.SetCoeff(p.nf[k]);
        tmp.push_back(tmpmon);
          }
        k++;
      }
    }
  m=tmp.front();
  for(typename list<mon>::iterator itertmp=tmp.begin();itertmp!=tmp.end();itertmp++)
    {
      if(m.GetCoeff()<itertmp->GetCoeff()) {m=*itertmp;};
    }
  //cout<<"--------------------------------------------------------------------------------"<<endl;
  // cout<<"monom choisi "<<m<<endl;
  //cout<<"pol"<<endl<<invconv<Poly>(p,serv)<<endl;
  //cout<<"--------------------------------------------------------------------------------"<<endl;
  return m;
}

template<typename mon, typename T,typename Base, typename monomial_server>
mon choice_grevlex(const pol2ter<mon,T> & p,const Base &b, 
           monomial_server &serv)
{
  list<mon> tmp;
  mon m; 
  int maxdeg=-1,k=0;
  //int maxtmpdeg=0;
  for(int i=0;i<p.sizenf;i++)
    {
      if (p.nfind[i])
    for(int j=0;j<8;j++)
      if(p.nfind[i]>>j&1)
        
    {
      mon tmpmon;
      serv.int2mon(8*i+j,tmpmon);
      if((tmpmon.GetDegree()>maxdeg) && (p.nf[k]!=0) )
        {
          tmp.erase(tmp.begin(),tmp.end());
          maxdeg=tmpmon.GetDegree();
        }
      if(tmpmon.GetDegree()==maxdeg)
        {
          tmpmon.SetCoeff(p.nf[k]);
          tmp.push_back(tmpmon);
        }
      k++;
    }
    }
  for(;tmp.front().GetCoeff()==0;tmp.pop_front());
  m=tmp.front();
  for(typename list<mon>::iterator itertmp=tmp.begin();itertmp!=tmp.end();itertmp++)
    {
      if((itertmp->GetCoeff()!=0) && (my_ord(*itertmp,m)<0)) {m=*itertmp;};
    }
  // cout<<"monom choisi "<<m<<endl;
  return m;
}

template<typename mon, typename T,typename Base, typename monomial_server>
mon choice_dlex(const pol2ter<mon,T> & p,const Base &b, monomial_server &serv)
{
  list<mon> tmp;
  mon m; 
  int maxdeg=-1,k=0;
  //int maxtmpdeg=0;
  for(int i=0;i<p.sizenf;i++)
    {
      if (p.nfind[i])
    for(int j=0;j<8;j++)
      if(p.nfind[i]>>j&1)
        {
          mon tmpmon;
          serv.int2mon(8*i+j,tmpmon);
          if(tmpmon.GetDegree()>maxdeg)
        {
          tmp.erase(tmp.begin(),tmp.end());
          maxdeg=tmpmon.GetDegree();
        }
          if(tmpmon.GetDegree()==maxdeg)
        {
          tmpmon.SetCoeff(p.nf[k]);
          tmp.push_back(tmpmon);
        }
          k++;
        }
    }
  m=tmp.front();
  for(typename list<mon>::iterator itertmp=tmp.begin();itertmp!=tmp.end();itertmp++)
    {
      if(my_ord(*itertmp,m)>0) {m=*itertmp;};
    }
  //  cout<<"monom choisi "<<m<<endl;
  return m;
}


template<typename mon, typename T,typename Base, typename monomial_server>
mon choice_mac(const pol2ter<mon,T> & p,const Base &b, monomial_server &serv)
{
  list<mon> tmp;
  mon m; 
  int maxdeg=-1,k=0;
  int maxtmpdeg=0;
  for(int i=0;i<p.sizenf;i++)
    {
      if (p.nfind[i])
    for(int j=0;j<8;j++)
      if(p.nfind[i]>>j&1)
        {
          mon tmpmon;
          serv.int2mon(8*i+j,tmpmon);
          if(tmpmon.GetDegree()>maxdeg && !Iszero(p.nf[k]))
        {
          tmp.erase(tmp.begin(),tmp.end());
          maxdeg=tmpmon.GetDegree();
        }
          if(tmpmon.GetDegree()==maxdeg && !Iszero(p.nf[k]))
        {
          tmpmon.SetCoeff(p.nf[k]);
          tmp.push_back(tmpmon);
        }
          k++;
        }
    }
  m=tmp.front();
  for(typename list<mon>::iterator itertmp=tmp.begin();itertmp!=tmp.end();itertmp++)
    {
      //cout<<"monomes regarde "<<*itertmp<<" Lavr  "<<lvar(itertmp->rep)<<endl;
      for(int i=0;i<=lvar(itertmp->rep);i++)
    {
      if(itertmp->GetDegree(i)>=maxtmpdeg)
        {
          maxtmpdeg=itertmp->GetDegree(i);
          m=*itertmp;
        }
    }
    }
  return m;
}


template<typename mon, typename T,typename Base, typename monomial_server>
mon choice_mpq(const pol2ter<mon,T> & p,const Base &b, monomial_server &serv)
{
  int i,degp=Degree(p),k=0;
  int mini=INT_MAX;
  mon tmpmon;
  mon stockmon;
  for(i=0;i<p.sizenf;i++)
    if(p.nfind[i])
      for(int j=0;j<8;j++)
    if(p.nfind[i]>>j&1)
      {
        serv.int2mon(8*i+j,tmpmon);
        if(tmpmon.GetDegree()==degp)
          {
        __mpz_struct  *numref,*denref;
        numref=mpq_numref(&(p.nf[k].rep()));
        denref=mpq_denref(&(p.nf[k].rep()));
        if(mpz_size(numref)+
           mpz_size(denref)<mini) 
          {
            mini=mpz_size(mpq_numref(&(p.nf[k].rep())))+
              mpz_size(mpq_denref(&(p.nf[k].rep())));
            //cout<<"newmax "<<max<<endl;
            //cout<<" les monomes "<<*iter<<endl;
            stockmon=tmpmon*p.nf[k];
          }
          }
        k++;
      }
  //cout<<"je retourne "<<stockmon<<endl;
  return stockmon;
}

template<typename mon, typename T, typename Base, typename monomial_server>
mon choice_fall(const pol2ter<mon,T> & p, const Base &b, monomial_server &serv)
{
  typedef typename mon::coeff_t coeff_t;
  int i,k=0;
  int maxi=-1;
  mon tmpmon;
  mon stockmon;
  mon res=choice_mac(p,b,serv);
  typename Base::const_iterator iterb;
  list<mon> tmplist;
  int degp=Degree(p,serv);
#ifdef DEB 
  cout<<"Degree(p) "<<degp<<endl;
#endif
  for(i=0;i<p.sizenf;i++)
    if(p.nfind[i])
      for(int j=0;j<8;j++)
    if((p.nfind[i]>>j)&1) 
      {
        serv.int2mon(8*i+j,tmpmon);
        if ((tmpmon.GetDegree()==degp) && (!Iszero(p.nf[k])))
          tmplist.push_back(tmpmon*p.nf[k]);
        k++;    
        //    cout<<"k "<<k<<" "<<tmpmon<<endl;
      }
  for(typename list<mon>::iterator iterlist=tmplist.begin()
    ;iterlist!=tmplist.end();iterlist++)
    {
      for(iterb=b.begin();iterb!=b.end();iterb++)
    {
      int isgood=0;
      for(int i=0;i<iterb->taille1;i++)
        if(isdivisible(iterb->accept[i],*iterlist,b))
          {isgood=1;break;}
      for(int i=0;i<iterb->taille2;i++)
        if(isdivisible(iterb->refuse[i],*iterlist,b))
          {isgood=0;break;}
      if(isgood)
        {

          for(int i=0;i<iterb->taille2;i++)
        {
          mon m=lcm(*iterlist,iterb->refuse[i]);
          int dist=0;
          for(int j=0;j<iterb->taille2;j++)
            if ((j!=i) && (isdivisible(iterb->refuse[j],m,b)))
              isgood=0;
          if(isgood)
            {
              for(int j=0;j<b.nbvar();j++)
            if(m.GetDegree(j)>iterb->refuse[i].GetDegree(j))
              dist+=m.GetDegree(j)-iterb->refuse[i].GetDegree(j);
              if(dist>maxi)
            {
              res=*iterlist;
              maxi=dist;
            }
            }
        }
        }
      
    }
    }
#if 1
  cout<<"le pol "<<invconv<Poly>(p,serv)<<endl<<"le mon choisi "<<res<<endl;
#endif
  return res;
}
#ifndef ASYNC
template<typename typdump,typename Base,typename mon, 
     typename workspace, typename monomial_server> 
pol<mon,COEFF> 
inline mult(int i,const pol<mon, COEFF> &p,const typdump &dump,
        const Base &b,mon **exce,int &exsize,workspace &w,
        monomial_server & serv)
{
  static const mon monzero=mon(0);
  static COEFF zero=0;
  static int flag=1;
  pol<mon,COEFF> ret;
  (*exce)=(mon*)MAC_REV_MALLOC<mon>((p.size+1)*sizeof(mon));
  int *stockindicetab=(int*)malloc(p.size*sizeof(int));
  int noquo=0,maxind=0,pp=0,q=0;
  COEFF *tmpspace;
  if(flag)
    {
      init_pointer<COEFF>();
      flag=0;cout<<"pouet pouet"<<endl;
    }
  //#ifdef DEBUG
#ifdef DEBUG
  cout<<"ds mult je traite "<<p.ind<<" ---- "<<i<<endl;
  cout<<invconv<Poly>(p,serv)<<endl;
  cout<<"p.size "<<p.size<<endl;
#endif
  for(int j=0;j<p.sizenf;++j)
    if (p.nfind[j])
     {
       for(int k=0;k<8;k++)
     if((p.nfind[j]>>k)&1)
       {
         int stockindice2=serv.mulind(8*j+k,i,b);
         if(stockindice2<0)
           stockindice2=-1*(8*j+k)-1;
         stockindicetab[pp++]=stockindice2;
         if (stockindice2>=maxind)
           maxind=stockindice2;
         cout<<"stockindice2 "<<8*j+k<<" "<<i<<stockindice2<<endl;
       }
     }
  for(;maxind/8+1>=w.sizenf;increase_space(&w));
  for(int j=0;j<pp;j++)
    if(stockindicetab[j]>=0)
      {
    w.tabcoeff[stockindicetab[j]]+=p.nf[j];
    if(!Iszero(w.tabcoeff[stockindicetab[j]]))
      w.tabnf[stockindicetab[j]/8]|=1<<(stockindicetab[j]%8);
      }
    else  
      {
    pol<mon,COEFF>  tmp;
    tmp=serv.nf(i,-1*(stockindicetab[j]+1),dump,b);
          if (tmp.size!=-1)
        {
          int compteurtmp=0;
          //On fait un plus dense dans tmpspace
          int k,l;
          unsigned long *tmpptr=(unsigned long*)tmp.nfind;
          unsigned long mask=1;
          for(;tmp.sizenf>=w.sizenf;increase_space(&w));

          for(k=0;k<tmp.sizenf/sizeof(unsigned long);k++)
            if(tmpptr[k])
              {
            mask=tmpptr[k];
            ((unsigned long*)w.tabnf)[k]|=mask;
            for(l=0;l<8*sizeof(unsigned long);l+=8,mask>>=8)
              {
                if(mask&255)
                  {
                tabfun2[mask&255](
                          w.tabcoeff+(8*sizeof(unsigned long)*k+l)
                          ,&(p.nf[j]),tmp.nf+compteurtmp);
                compteurtmp+=nbbits[mask&255];
                  }
#if 0
                for(int m=0;m<8;m++)
                  if(Iszero(w.tabcoeff[8*k*sizeof(unsigned long)+l+m]))
                ((unsigned long*)w.tabnf)[k]&=~(1<<(l+m));
#endif
              }
              }

          for(k=tmp.sizenf-(tmp.sizenf%sizeof(unsigned long))
            ;k<tmp.sizenf;k++)
            if(tmp.nfind[k])
              {
            w.tabnf[k]|=tmp.nfind[k];
            mask=1;
            for(l=0;l<8;l++,mask<<=1)
              {
                if(tmp.nfind[k]&mask)
                  w.tabcoeff[8*k+l]-=p.nf[j]*tmp.nf[compteurtmp++];
                if(Iszero(w.tabcoeff[8*k+l]))
                  w.tabnf[k]&=~mask;

              }
              }
        }
          else  
        {
          //On est tombé sur un monome du degre courant
          //qui est sur le bord
          mon tmpmon;
          (*exce)[noquo]=mon(i,1);
          serv.int2mon(-1*(stockindicetab[j]+1),tmpmon);
          (*exce)[noquo]*=tmpmon;//int2mon<mon>(j);
          (*exce)[noquo]*=p.nf[j];
          ++noquo;
        };
      }
#ifdef ANCESTOR
  //On met juste a jour l'ancetre pour la construction du grpahe apres
  ret.ancestor=p.ancestor;
#endif
  ret.ind=p.ind;
  ret.ind*=mon(i,1);
  (*exce)[noquo++]=ret.ind;
  exsize=noquo;
  (*exce)=(mon*)MAC_REV_REALLOC<mon>((*exce),(p.size+1)*sizeof(mon)
                   ,(noquo)*sizeof(mon));
  ret.ind.SetCoeff(noquo+1);
  ret.size=0;
  pp=0;
  for(int j=0;j<w.sizenf;j++)
    if(w.tabnf[j])
      {
    for(int m=0;m<8;m++)
      if(Iszero(w.tabcoeff[8*j+m]))
        w.tabnf[j]&=~(1<<m);
    if(w.tabnf[j])
      {
        ret.size+=nbbits[w.tabnf[j]];
        pp=j;
      }
      }
  ret.nf=(COEFF*)MAC_REV_MALLOC<COEFF>(ret.size*sizeof(COEFF));
  ret.sizenf=pp+1;
  ret.nfind=(unsigned char *)malloc(ret.sizenf);
  pp=0;
  for(int j=0;j<ret.sizenf;j++)
    {
      if(w.tabnf[j])
    for(int k=0;k<8;k++)
      if((w.tabnf[j]>>k)&1)
        {
          ret.nf[pp++]=w.tabcoeff[8*j+k];
          w.tabcoeff[8*j+k]=zero;
        }
      ret.nfind[j]=w.tabnf[j];
      w.tabnf[j]=0;
    }
  free(stockindicetab);
#ifdef DEBUG
  cout<<"en sortie"<<endl;
  cout<<invconv<Poly>(ret,serv)<<endl;
  copy(*exce,*exce+exsize,ostream_iterator<mon>(cout," + "));cout<<endl;
  cout<<"----------"<<endl;
#endif
  return ret;
}
#endif

#if 1
template<typename typmat,typename typPk, typename typMk,typename typdump,
     typename typexceed, typename typsizeexceed, typename typB,typename monomial_server>
     void AddB(const typmat & mat, typPk & Pk,const typMk &Mk,const typexceed & exceed,
           const typsizeexceed & sizeexceed, const typdump & dump,
      typB &b, monomial_server &serv)
{
  typedef typename typPk::value_type pol;
  typedef typename pol::monom_t mon;
  typedef typename pol::coeff_t coeff_t;
  typename typPk::iterator iterpk=Pk.begin();;
  auto iterexceed=exceed.begin();
  auto itersizeexceed=sizeexceed.begin();
  predicat<mon> tmp;
  typename typMk::const_iterator iter=Mk.begin();
  int r=my_rank(mat),comprefuse=0,compaccept=0;
  unsigned int i;
  //cout<<"DANS ADDB"<<endl;
  //cout<<"rang de la matrice m "<<r<<" Mk.size "<<Mk.size()<<endl;
  //cout<<"taille de la matrice "<<mat.nrow<<"x"<<mat.ncol<<endl;
  tmp.refuse=(mon*)MAC_REV_MALLOC<mon>(r*sizeof(mon));
  tmp.taille2=r;
  tmp.accept=(mon*)MAC_REV_MALLOC<mon>((Mk.size()-r)*sizeof(mon));
  tmp.taille1=Mk.size()-r;
  //for(i=0;i<(unsigned)r;i++,iter++)
  //  tmp.refuse[i]=*iter;
  iter=Mk.begin();
  for(i=0;i<Mk.size();i++,iter++)
    {
      if(Iszero(mat(i,i))) 
    {
      tmp.accept[compaccept++]=*iter;
      //cout<<"monome rajoute au quotient "<<*iter<<endl;
    }
      else
    {
      tmp.refuse[comprefuse++]=*iter;
    }
    }
  b.def.push_back(tmp);
  //bon maintenant que l'on vient de rajouter des 
  //monomes au quotient il faut remttre les polynomes comme il faut
  iter=Mk.begin();
//cout<<"report"<<endl;
//On fait une passe pour determiner les polynomes a modifier et on effectue 
//les modifications en lisant les exceed (au passage on FREE)

  
  for(unsigned int i=0;i<Pk.size();i++,iterpk++,iter++,iterexceed++,itersizeexceed++)
    {
      vector<sortstruct<coeff_t> > decal;
      int maxindex=-1;
      int expend=0;
      for(int j=0;j<*itersizeexceed;j++)
    if(!Iszero((*iterexceed)[j].GetCoeff()))
      //si des monomes nuls sortent comme coeff de la matrice
      //et bien on les prends pas
      {
        sortstruct<coeff_t> tmp;
        expend++;
        maxindex=max(maxindex,serv.mon2int((*iterexceed)[j]));
        tmp.index=serv.mon2int((*iterexceed)[j]);
        tmp.coeff=(*iterexceed)[j].GetCoeff();
        decal.push_back(tmp);
      }
      // Si on a rien a rajouter on arrete
      // cout<<"decal.empty "<<decal.empty()<<endl;
      if(decal.empty()) continue;
      //on augmente la taille du polynome  considere
      iterpk->nf=(typename mon::coeff_t*)
    MAC_REV_REALLOC<typename mon::coeff_t>(iterpk->nf
           ,iterpk->size*sizeof(typename mon::coeff_t)
           ,(iterpk->size+expend)*sizeof(typename mon::coeff_t));
      if((maxindex/8)>=iterpk->sizenf)
    {
      iterpk->nfind=(unsigned char*)
        MAC_REV_REALLOC<unsigned char>(iterpk->nfind
                       ,iterpk->sizenf,maxindex/8+1);
      for(int i=iterpk->sizenf;i<maxindex/8+1;i++)
        iterpk->nfind[i]=0;
      iterpk->sizenf=maxindex/8+1;
      
    }
// la place est faite il reste plus qu'a decaler les monomes de iterpk->nf
//pour intercaler eventuellment les nouveaux monomes
       stable_sort(decal.begin(),decal.end(),sortstructorder<coeff_t>());
       int nbcoeff=iterpk->size-1;
      // cout<<"iterpk->sizenf "<<iterpk->sizenf<<endl;
       for(int j=iterpk->sizenf-1;j>=0;j--)
        {
          //cout<<"decal.index "<<decal.rbegin()->index<<endl;
          if(8*j>decal.rbegin()->index)
            for (int jj=0;jj<nbbits[iterpk->nfind[j]];jj++)
              {
                iterpk->nf[nbcoeff+decal.size()]=
                  iterpk->nf[nbcoeff];
                nbcoeff--;
            }
          else
            {
              int indexnf=decal.back().index;
              typename typMk::value_type tmpmon;
              //cout<<"8*j"<<8*j<<" indexnf "<<indexnf<<" valeur "<<decal.back().index<<endl;
              int jj=7;
              while(8*j<=decal.back().index)
              {
              indexnf=decal.back().index % 8;
              for (;jj>(indexnf);jj--)
              if(iterpk->nfind[j]&(1<<jj))
                {
                 iterpk->nf[nbcoeff+decal.size()]=
                 iterpk->nf[nbcoeff];
                 nbcoeff--;
                }
                      //cout<<"je rajoute "<<decal.back().coeff<<" a la place "<<nbcoeff+decal.size()<<endl;
//          if(nbcoeff+decal.size()+1<iterpk->size+expend)
            //  cout<<"coeff apres "<<iterpk->nf[nbcoeff+decal.size()+1]<<endl;
              iterpk->nf[nbcoeff+decal.size()]=decal.back().coeff;
              iterpk->nfind[j]|=(1<<(indexnf));
          //   serv.int2mon(8*j+indexnf,tmpmon);
        //  cout<<"monome a la place "<<decal.back().index<<" "<<tmpmon<<endl;    
          //   serv.int2mon(decal.back().index,tmpmon);
        //  cout<<"monome a la place "<<decal.back().index<<" "<<tmpmon<<endl;    
              decal.pop_back();
             jj--; 
            if(decal.empty()) break;
            }   
             for (;jj>=0;jj--)
              if(iterpk->nfind[j]&(1<<jj))
               {
                 iterpk->nf[nbcoeff+decal.size()]=
                        iterpk->nf[nbcoeff];
                  nbcoeff--;
                 }
            if(decal.empty()) break;
            }
                
        }
       iterpk->size+=expend;
       MAC_REV_FREE<mon>(*iterexceed,*itersizeexceed*sizeof(mon));
    }   
  return;
}
#endif
#ifdef ASYNC
#include "pol2terspecial-async.hpp"
#endif

        {
          for(int j=m.nnz;j>=m.colptr[mon2col[(*iterexceed)[i].rep]+1];j--)
        {
          ((coeff*)m.nzval)[j]=((coeff*)m.nzval)[j-1];
          m.rowind[j]=m.rowind[j-1];
        }  
          //   cout<<"mon2col "<<(*iterexceed)[i].rep<<" "<<m.colptr[mon2col[(*iterexceed)[i].rep]+1]<<endl;
          m.nnz++;
          m.rowind[m.colptr[mon2col[(*iterexceed)[i].rep]+1]]=indexmat+numligne;
          //ICI MODIF DU 11 JANVIER
          //OUILLE OUILLE OUILLE
//        m.rowind[m.colptr[mon2col[(*iterexceed)[i].rep]+1]]=indexmat+i;
          ((coeff*)m.nzval)[m.colptr[mon2col[(*iterexceed)[i].rep]+1]]=(*iterexceed)[i].GetCoeff();
#ifdef DEBUG
          cout<<"G ecrit "<<(*iterexceed)[i].GetCoeff()<<endl;
          cout<<"en colonne "<<mon2col[(*iterexceed)[i].rep]<<" ligne "<<m.colptr[mon2col[(*iterexceed)[i].rep]+1]<<" "<<indexmat+numligne<<endl;
#endif
          for(int j=mon2col[(*iterexceed)[i].rep]+1;j<=m.ncol;j++)
        m.colptr[j]++;
          //les colonnes d'apres commencent maintenant un cran plus ooin
          
        }
    }
      //MAC_REV_FREE<mon>(*iterexceed,*itersizeexceed*sizeof(mon));
      //ICI c dans ADDB qu'on fait le FREE
      //
    }

}
//
// Les polyn\^omes des autres composantes fortement connexes sont mis à jour en 
// fonction du resultat de l'inversion de la composante en cours
//
template <typename typcomp, typename typserv>
void  update_pol(typcomp &comp_to_update,typcomp &updating_comp, typserv &serv)
{
  typedef typename typcomp::typPk typPk;
  typedef typename typserv::monom_t mon;
  typedef typename typcomp::pol pol;
  typedef typename pol::coeff_t coeff_t;
//  cout<<"rentre dans update"<<endl;
  if(comp_to_update==updating_comp)
    return;
//  cout<<"rentre dans compo pas les memes"<<endl;
  auto iterexceed=updating_comp.exceed.begin();
  auto itersizeexceed=updating_comp.sizeexceed.begin();
  auto iterpoly=updating_comp.poly.begin();
  //
  //
  
  auto myitersize=updating_comp.sizeexceed.begin();
  #ifdef DEBUG
  for(auto &i : updating_comp.poly)
    cout<<"ind qui update "<<i.ind<<" "<<*myitersize++<<endl;
#endif
  myitersize=comp_to_update.sizeexceed.begin();
 #ifdef DEBUG
  for(auto &i : comp_to_update.poly)
    cout<<"ind a up "<<i.ind<<" "<<*myitersize++<<endl;
#endif
  //on créé l'association enter les ind et les structures correspondantes
  map<typename mon::rep_t,resmul<typename typPk::value_type> > mon2pol;
  for(;iterpoly!=updating_comp.poly.end();iterpoly++,iterexceed++,itersizeexceed++)
    {
      resmul<typename typPk::value_type> tmp;
      tmp.poly=*iterpoly;
      tmp.exceed=*iterexceed;
      tmp.sizeexceed=*itersizeexceed;
      mon2pol[iterpoly->ind.rep]=tmp;
      //on augmente le workspace pour update apres
      //for(;w->sizenf<iterpoly->sizenf;increase_space(w));
    }
  iterexceed=comp_to_update.exceed.begin();
  itersizeexceed=comp_to_update.sizeexceed.begin();
  for(iterpoly=comp_to_update.poly.begin();
    iterpoly!=comp_to_update.poly.end();
    iterpoly++,iterexceed++,itersizeexceed++)
    //pour chaque polynome dans la composante a mettre a jour
    {
      //On construit une map pour l'exceed correspondant
      int k=0;
      map<typename mon::rep_t,COEFF> exceedmap;
      map<typename mon::rep_t,COEFF> nfmap;//il faudrait utiliser le workspace ici au lieu de cette map
 //     cout<<"itersizexceed ici "<<*itersizeexceed<<endl;
      int pos_nf=0;
      int newsizenfind=iterpoly->sizenf;
      //augmente ici le workspace pour que la mise a jour
      //puisse se faire
      //for(;w->sizenf<iterpoly->sizenf;increase_space(w));
      for(int i=0;i<iterpoly->sizenf;i++)
    if(iterpoly->nfind)
    {
      for(int j=0;j<8;j++)
        if(iterpoly->nfind[i]>>j&1)
        {
          mon tmp;
          serv.int2mon(8*i+j,tmp);
          nfmap[tmp.rep]=iterpoly->nf[pos_nf++];
        }
    }
      for(int i=0;i<*itersizeexceed;i++)
    {
      //On ne rajoute pas dans exceedmap 
      //les monomes de mon2pol!
      if(!mon2pol.count((*iterexceed)[i].rep))
       exceedmap[(*iterexceed)[i].rep]= 
         (*iterexceed)[i].GetCoeff();
    }
      for(int i=0;i<*itersizeexceed;i++)
    {
      if(mon2pol.count((*iterexceed)[i].rep))
        //On a trouvé un monome de la composante 
        //qu'on utilise pour mettre a jour
        {
          const resmul<typename typPk::value_type> &tmpresmul=mon2pol[(*iterexceed)[i].rep];
//        cout<<"je declenche une mise a jour par "<<tmpresmul.poly.ind<<endl;
//        cout<<"   de "<<iterpoly->ind<<endl;
          if(newsizenfind<tmpresmul.poly.sizenf)
        newsizenfind=tmpresmul.poly.sizenf;
          //
          //A VERIFIER QUE
          //les exceedd  sont modifié comme il faut pendant le solve
          //fait
//        cout<<"le coeff maj "<<(*iterexceed)[i].GetCoeff()<<endl;
//        cout<<"exceed avant "<<*itersizeexceed<<" tmpresmul.sizeexceed "<<tmpresmul.sizeexceed<<endl;
//        for(auto &i : exceedmap)
//      cout<<i.first<<" "<<i.second<<" ";
//        cout<<endl;
          for(int j=0;j<tmpresmul.sizeexceed;j++)
        {
          exceedmap[tmpresmul.exceed[j].rep]-=
            tmpresmul.exceed[j].GetCoeff()*(*iterexceed)[i].GetCoeff();
        }
//        cout<<"nf avant"<<endl;
//        for(auto &i : nfmap)
//      cout<<i.first<<" "<<i.second<<" ";
//        cout<<endl;
        int pos_nfiter=0;
        for(int ii=0;ii<tmpresmul.poly.sizenf;ii++)
            if(tmpresmul.poly.nfind)
            {
              for(int j=0;j<8;j++)
              if(tmpresmul.poly.nfind[ii]>>j&1)
                {
                  //ICI PASSER PAR LE WORKSPACE!!
                  //ca evite mon2int et la map derriere!
                  mon tmp;
                  serv.int2mon(8*ii+j,tmp);
                 #if 0
                  cout<<"mon mis a jour "<<tmp<<" "<<tmpresmul.poly.nf[pos_nfiter]<<endl;
                  cout<<"nfmap "<<nfmap[tmp.rep]<<endl;
                  cout<<tmpresmul.poly.nf[pos_nfiter]<<endl;
                  cout<<(*iterexceed)[i].GetCoeff()<<endl;
                  cout<<coeff_t(0)-tmpresmul.poly.nf[pos_nfiter]
                                     *(*iterexceed)[i].GetCoeff()<<endl;
                  cout<<nfmap[tmp.rep]-tmpresmul.poly.nf[pos_nfiter]
                                     *(*iterexceed)[i].GetCoeff()<<endl;
                  #endif
                  nfmap[tmp.rep]-=tmpresmul.poly.nf[pos_nfiter]
                *(*iterexceed)[i].GetCoeff();
              //    cout<<"apres operation "<<nfmap[tmp.rep]<<endl;
                pos_nfiter++;

                }
            }
#ifdef DEBUG
          cout<<"exceed apres"<<endl;
          for(auto &i : exceedmap)
        cout<<i.first<<" "<<i.second<<" ";
          cout<<endl;
          cout<<"nf apres"<<endl;
          for(auto &i : nfmap)
        cout<<i.first<<" "<<i.second<<" ";
          cout<<endl;
#endif       
        }
    }   
#ifdef DEBUG
      cout<<"nfmap.size() "<<nfmap.size()<<endl;
      cout<<"newsizenfind "<<newsizenfind<<" "<<(void*)iterpoly->nfind<<endl;
#endif
      *iterexceed=(mon*)MAC_REV_REALLOC<mon>(*iterexceed,*itersizeexceed*sizeof(mon)
                         ,exceedmap.size()*sizeof(mon));
      iterpoly->nfind=(unsigned char *)realloc(iterpoly->nfind,newsizenfind);
      if(iterpoly->size != nfmap.size())
       iterpoly->nf=(coeff_t*)MAC_REV_REALLOC<coeff_t>(iterpoly->nf,
                        iterpoly->size*sizeof(coeff_t),
                    nfmap.size()*sizeof(coeff_t));
      for(int k=0;k<newsizenfind;k++)
    iterpoly->nfind[k]=0;
      for(auto iter : nfmap)
      {
    mon m(1);
    m.rep=iter.first;
    int place=serv.mon2int(m);
    iterpoly->nfind[place/8]|=(1<<(place%8));
      }
      int posnf=0;
      for(int k=0;k<newsizenfind;k++)
      {
    if(iterpoly->nfind[k])
    {
      for(int kk=0;kk<8;kk++)
        if((iterpoly->nfind[k]>>kk)&1)
         {
           mon m;
           serv.int2mon(8*k+kk,m);
           iterpoly->nf[posnf++]=nfmap[m.rep];
         }
    }
      }
      iterpoly->size=nfmap.size();
      iterpoly->sizenf=newsizenfind;
      for(auto iter : exceedmap)
    {
      mon tmpmon;
      tmpmon.rep=iter.first;
      tmpmon.coeff=iter.second;
      (*iterexceed)[k++]=tmpmon;
    }
      *itersizeexceed=exceedmap.size();
    }
}

template <typename typcomp, typename typserv>
void  wrap_update_pol(typcomp &comp_to_update,typcomp &updating_comp, typserv &serv)
{
  workspace<typename typcomp::pol::coeff_t> w;
  update_pol(comp_to_update,updating_comp,w,serv);
  destroy_space(&w);
}
template <typename typcomp, typename typwork, typename typserv>
void  update_pol(typcomp &comp_to_update,typcomp &updating_comp, typwork &w,
         typserv &serv)
{
  typedef typename typcomp::typPk typPk;
  typedef typename typserv::monom_t mon;
  typedef typename typcomp::pol pol;
  typedef typename pol::coeff_t coeff_t;
//  cout<<"rentre dans update"<<endl;
  if(comp_to_update==updating_comp)
    return;
//  cout<<"rentre dans compo pas les memes"<<endl;
  auto iterexceed=updating_comp.exceed.begin();
  auto itersizeexceed=updating_comp.sizeexceed.begin();
  auto iterpoly=updating_comp.poly.begin();
  //
  //
  
  auto myitersize=updating_comp.sizeexceed.begin();
  #ifdef DEBUG
  for(auto &i : updating_comp.poly)
    cout<<"ind qui update "<<i.ind<<" "<<*myitersize++<<endl;
#endif
  myitersize=comp_to_update.sizeexceed.begin();
 #ifdef DEBUG
  for(auto &i : comp_to_update.poly)
    cout<<"ind a up "<<i.ind<<" "<<*myitersize++<<endl;
#endif
  //on créé l'association enter les ind et les structures correspondantes
  map<typename mon::rep_t,resmul<typename typPk::value_type> > mon2pol;
  for(;iterpoly!=updating_comp.poly.end();iterpoly++,iterexceed++,itersizeexceed++)
    {
      resmul<typename typPk::value_type> tmp;
      tmp.poly=*iterpoly;
      tmp.exceed=*iterexceed;
      tmp.sizeexceed=*itersizeexceed;
      mon2pol[iterpoly->ind.rep]=tmp;
      //on augmente le workspace pour update apres
      for(;w.sizenf<iterpoly->sizenf;increase_space(&w));
    }
  iterexceed=comp_to_update.exceed.begin();
  itersizeexceed=comp_to_update.sizeexceed.begin();
  for(iterpoly=comp_to_update.poly.begin();
    iterpoly!=comp_to_update.poly.end();
    iterpoly++,iterexceed++,itersizeexceed++)
    //pour chaque polynome dans la composante a mettre a jour
    {
      //On construit une map pour l'exceed correspondant
      int k=0;
      map<typename mon::rep_t,COEFF> exceedmap;
      map<typename mon::rep_t,COEFF> nfmap;//il faudrait utiliser le workspace ici au lieu de cette map
 //     cout<<"itersizexceed ici "<<*itersizeexceed<<endl;
      int pos_nf=0;
      int newsizenfind=iterpoly->sizenf;
      //augmente ici le workspace pour que la mise a jour
      //puisse se faire
      for(;w.sizenf<iterpoly->sizenf;increase_space(&w));
      for(int i=0;i<iterpoly->sizenf;i++)
    if(iterpoly->nfind)
    {
      w.tabnf[i]=iterpoly->nfind[i];
      tabfun1[w.tabnf[i]](w.tabcoeff+(8*i),iterpoly->nf+pos_nf);
      pos_nf+=nbbits[iterpoly->nfind[i]];
    }
      for(int i=0;i<*itersizeexceed;i++)
    {
      //On ne rajoute pas dans exceedmap 
      //les monomes de mon2pol!
      if(!mon2pol.count((*iterexceed)[i].rep))
       exceedmap[(*iterexceed)[i].rep]= 
         (*iterexceed)[i].GetCoeff();
    }
      for(int i=0;i<*itersizeexceed;i++)
    {
      if(mon2pol.count((*iterexceed)[i].rep))
        //On a trouvé un monome de la composante 
        //qu'on utilise pour mettre a jour
        {
          resmul<typename typPk::value_type> &tmpresmul=mon2pol[(*iterexceed)[i].rep];
          if(newsizenfind<tmpresmul.poly.sizenf)
        newsizenfind=tmpresmul.poly.sizenf;
          for(int j=0;j<tmpresmul.sizeexceed;j++)
        {
          exceedmap[tmpresmul.exceed[j].rep]-=
            tmpresmul.exceed[j].GetCoeff()*(*iterexceed)[i].GetCoeff();
        }
        int pos_nfiter=0;
        for(int ii=0;ii<tmpresmul.poly.sizenf;ii++)
            if(tmpresmul.poly.nfind)
            {
            coeff_t tmp=(*iterexceed)[i].GetCoeff();    
            w.tabnf[ii]=tmpresmul.poly.nfind[ii];
            tabfun2[w.tabnf[ii]](w.tabcoeff+(8*ii),
                &tmp,
                tmpresmul.poly.nf+pos_nfiter);
            pos_nfiter+=nbbits[tmpresmul.poly.nfind[ii]];

            }
#ifdef DEBUG
          cout<<"exceed apres"<<endl;
          for(auto &i : exceedmap)
        cout<<i.first<<" "<<i.second<<" ";
          cout<<endl;
          cout<<"nf apres"<<endl;
          for(auto &i : nfmap)
        cout<<i.first<<" "<<i.second<<" ";
          cout<<endl;
#endif       
        }
    }   
#ifdef DEBUG
      cout<<"newsizenfind "<<newsizenfind<<" "<<(void*)iterpoly->nfind<<endl;
      cout<<"w.sizenf "<<w.sizenf<<endl;
#endif
      *iterexceed=(mon*)MAC_REV_REALLOC<mon>(*iterexceed,*itersizeexceed*sizeof(mon)
                         ,exceedmap.size()*sizeof(mon));
      iterpoly->nfind=(unsigned char *)realloc(iterpoly->nfind,newsizenfind);
      int newsize=0;
      for(int i=0;i<w.sizenf;newsize+=nbbits[w.tabnf[i++]]);
       iterpoly->nf=(coeff_t*)MAC_REV_REALLOC<coeff_t>(iterpoly->nf,
                        iterpoly->size*sizeof(coeff_t),
                    newsize*sizeof(coeff_t));
      for(int k=0;k<newsizenfind;k++)
    {
      iterpoly->nfind[k]=w.tabnf[k];
      w.tabnf[k]=0;
    }
      for(int k=0;k<newsize;k++)
      {
    iterpoly->nf[k]=w.tabcoeff[k];
    w.tabcoeff[k]=0;
      }
      iterpoly->size=newsize;
      iterpoly->sizenf=newsizenfind;
      for(auto iter : exceedmap)
    {
      mon tmpmon;
      tmpmon.rep=iter.first;
      tmpmon.coeff=iter.second;
      (*iterexceed)[k++]=tmpmon;
    }
      *itersizeexceed=exceedmap.size();
    }
}
template <typename typcomp, typename typmon2pol, typename typserv>
void  wrap_update_pol(typcomp &comp_to_update,const typcomp &updating_comp, 
    const typmon2pol &mon2pol, const typserv &serv)
{
  workspace<typename typcomp::pol::coeff_t> w;
  w.sizenf=0;
  w.size=0;
  w.tabnf=NULL;
  w.tabcoeff=NULL;
  for(auto &p : updating_comp.poly)
    for(;w.sizenf<=p.sizenf;increase_space(&w));
  update_pol(comp_to_update,updating_comp,mon2pol,w,serv);
  destroy_space(&w);
}
template <typename typcomp, typename typmon2pol, typename typwork, typename typserv>
void  update_pol(typcomp &comp_to_update,const typcomp &updating_comp, 
    const typmon2pol &mon2pol, typwork &w, const typserv &serv)
{
  typedef typename typcomp::typPk typPk;
  typedef typename typserv::monom_t mon;
  typedef typename typcomp::pol pol;
  typedef typename pol::coeff_t coeff_t;
//  cout<<"rentre dans update"<<endl;
  if(comp_to_update==updating_comp)
    return;
//  cout<<"rentre dans compo pas les memes"<<endl;
  auto iterexceed=comp_to_update.exceed.begin();
  auto itersizeexceed=comp_to_update.sizeexceed.begin();
  auto iterpoly=comp_to_update.poly.begin();
  //
  //
  
  auto myitersize=updating_comp.sizeexceed.begin();
#ifdef DEBUG_UPDATE
  for(auto &i : updating_comp.poly)
    cout<<"ind qui update "<<i.ind<<" "<<*myitersize++<<endl;
#endif
  myitersize=comp_to_update.sizeexceed.begin();
#ifdef DEBUG_UPDATE
  for(auto &i : comp_to_update.poly)
    cout<<"ind a up "<<i.ind<<" "<<*myitersize++<<endl;
#endif
  //on créé l'association enter les ind et les structures correspondantes
  iterexceed=comp_to_update.exceed.begin();
  itersizeexceed=comp_to_update.sizeexceed.begin();
  for(iterpoly=comp_to_update.poly.begin();
    iterpoly!=comp_to_update.poly.end();
    iterpoly++,iterexceed++,itersizeexceed++)
    //pour chaque polynome dans la composante a mettre a jour
    {
       //On construit une map pour l'exceed correspondant
      int k=0;
      map<typename mon::rep_t,COEFF> exceedmap;
 //     cout<<"itersizexceed ici "<<*itersizeexceed<<endl;
      int pos_nf=0;
      int newsizenfind=iterpoly->sizenf;
      //augmente ici le workspace pour que la mise a jour
      //puisse se faire
      for(;w.sizenf<=iterpoly->sizenf;increase_space(&w));
      for(int i=0;i<iterpoly->sizenf;i++)
    if(iterpoly->nfind)
    {
      w.tabnf[i]=iterpoly->nfind[i];
      tabfun1[w.tabnf[i]](w.tabcoeff+(8*i),iterpoly->nf+pos_nf);
      pos_nf+=nbbits[iterpoly->nfind[i]];
    }
      for(int i=0;i<*itersizeexceed;i++)
    {
      //On ne rajoute pas dans exceedmap 
      //les monomes de mon2pol!
      if(!mon2pol.count((*iterexceed)[i].rep))
       exceedmap[(*iterexceed)[i].rep]= 
         (*iterexceed)[i].GetCoeff();
    }
      for(int i=0;i<*itersizeexceed;i++)
    {
      if(mon2pol.count((*iterexceed)[i].rep))
        //On a trouvé un monome de la composante 
        //qu'on utilise pour mettre a jour
        {
          const resmul<typename typPk::value_type> &tmpresmul=
                             mon2pol.find((*iterexceed)[i].rep)->second;
          if(newsizenfind<tmpresmul.poly.sizenf)
        newsizenfind=tmpresmul.poly.sizenf;
          for(int j=0;j<tmpresmul.sizeexceed;j++)
        {
          exceedmap[tmpresmul.exceed[j].rep]-=
            tmpresmul.exceed[j].GetCoeff()*(*iterexceed)[i].GetCoeff();
        }
        int pos_nfiter=0;
        for(int ii=0;ii<tmpresmul.poly.sizenf;ii++)
            if(tmpresmul.poly.nfind[ii])
            {
            coeff_t tmp=(*iterexceed)[i].GetCoeff();    
            w.tabnf[ii]|=tmpresmul.poly.nfind[ii];
            //          tabfun2[w.tabnf[ii]](w.tabcoeff+(8*ii),
            tabfun2[tmpresmul.poly.nfind[ii]](w.tabcoeff+(8*ii),
                     &tmp,
                     tmpresmul.poly.nf+pos_nfiter);
            pos_nfiter+=nbbits[tmpresmul.poly.nfind[ii]];

            }
        }
    }   
#ifdef DEBUG_UPDATE
      cout<<"newsizenfind "<<newsizenfind<<" "<<(void*)iterpoly->nfind<<endl;
      cout<<"w.sizenf "<<w.sizenf<<endl;
#endif
      *iterexceed=(mon*)MAC_REV_REALLOC<mon>(*iterexceed,*itersizeexceed*sizeof(mon)
                         ,exceedmap.size()*sizeof(mon));
      iterpoly->nfind=(unsigned char *)realloc(iterpoly->nfind,newsizenfind);
      int newsize=0;
      for(int i=0;i<w.sizenf;newsize+=nbbits[w.tabnf[i++]]);
       iterpoly->nf=(coeff_t*)MAC_REV_REALLOC<coeff_t>(iterpoly->nf,
                        iterpoly->size*sizeof(coeff_t),
                    newsize*sizeof(coeff_t));
#ifdef DEBUG_UPDATE
       cout<<"newsize "<<newsize<<endl;
#endif
       for(int k=0;k<newsizenfind;k++)
    {
      iterpoly->nfind[k]=w.tabnf[k];
      w.tabnf[k]=0;
    }
       //      for(int k=0;k<newsize;k++)
       // {
       //   iterpoly->nf[k]=w.tabcoeff[k];
       //   w.tabcoeff[k]=0;
       //}
       int tmpcompt=0;
       for(int k=0;k<newsizenfind;k++)
     {
       if(iterpoly->nfind[k])
         {
           for(int j=0;j<8;j++)
         if((iterpoly->nfind[k]>>j)&1)
           {
             iterpoly->nf[tmpcompt++]=w.tabcoeff[8*k+j];
             w.tabcoeff[8*k+j]=0;
           }
         }
     }
      iterpoly->size=newsize;
      iterpoly->sizenf=newsizenfind;
#ifdef DEBUG_UPDATE
      cout<<"exceed apres "<<exceedmap.size()<<endl;
      for(auto &i : exceedmap)
    cout<<i.first<<" "<<i.second<<" ";
      cout<<endl;
#endif      
      for(auto iter : exceedmap)
    {
      mon tmpmon;
      tmpmon.rep=iter.first;
      tmpmon.coeff=iter.second;
      (*iterexceed)[k++]=tmpmon;
    }
      *itersizeexceed=exceedmap.size();
      // cout<<"---------------------------"<<endl;
    }
}
//est sense reconstruire la matrice mat a partir 
//des matrices correspondant aux composantes fortmeùent connexes

template<typename typMat>
void reconstruct_mat(typMat &mat, const typMat &tmpmat)
{
  typedef typename typMat::coeff_t coeff_t;
#ifdef DEBUG
  cout<<"la petitee "<<tmpmat.nnz<<endl;
  affdebug(tmpmat);
  cout<<"la grosse "<<endl;
  if(mat.ncol)
    affdebug(mat);
#endif
  mat.colptr=(int *)realloc(mat.colptr,(mat.ncol+tmpmat.ncol+1)*sizeof(int));
  mat.nzval=MAC_REV_REALLOC<coeff_t>(mat.nzval,mat.nnz*sizeof(coeff_t),
    (mat.nnz+tmpmat.nnz)*sizeof(coeff_t));
  mat.rowind=(int *)realloc(mat.rowind,(mat.nnz+tmpmat.nnz)*sizeof(int));
  //
  for(int i=0;i<=tmpmat.ncol;i++)
    mat.colptr[mat.ncol+i]=mat.nnz+tmpmat.colptr[i];
  //
  for(int i=0;i<tmpmat.nnz;i++)
    mat.rowind[mat.nnz+i]=mat.nrow+tmpmat.rowind[i];
  //
  for(int i=0;i<tmpmat.nnz;i++)
    ((coeff_t *)mat.nzval)[mat.nnz+i]=((coeff_t *)tmpmat.nzval)[i];
  mat.ncol+=tmpmat.ncol;
  mat.nrow+=tmpmat.nrow;
  mat.nnz+=tmpmat.nnz;
#ifdef DEBUG
   cout<<"dans reconstruct mat.ncol en sortie "<<mat.ncol<<endl;
   cout<<"dans reconstruct mat.nrow en sortie "<<mat.nrow<<endl;
#endif
}


void reconstruct_pr_pc(int *&pr,int *&pc,int *tmppr,int *tmppc,int sizetmp,int indexinside)
{
  for(int i=0;i<sizetmp;i++)
    pr[i+indexinside]=tmppr[i]+indexinside;
  
  for(int i=0;i<sizetmp;i++)
    pc[i+indexinside]=tmppc[i]+indexinside;
}

template<typename typMat,typename typP,typename typPk,typename typMk, typename Base,typename typdump, class monomial_server, typename workspace>
void init_algo(typP &P,int &nouveaumon, int &maxdeg, int & k,
           typPk &Pk, typMk & Mk, Base & b, typdump & dump, workspace &w,
           monomial_server &serv);

template<typename typP,typename typPk,typename typdump,typename Base
     , template<typename  T> class tmptypmat ,
     template<typename mon, typename pol> class monomial_server, typename mon>
void algo(typP P,typdump &dump,Base &b, 
      monomial_server<mon,pol2ter<mon,typename mon::coeff_t> > &serv)
{
  typedef typename typP::value_type polyalp;
  typedef list<typename typPk::value_type::monom_t> typMk;
  //  typedef typename typMk::value_type mon;
  typedef typename typPk::value_type pol;

  typMk Mk;
  typMk Spoldejavu;
  typPk Pk,S;
  typPk recall;
  workspace<typename typP::value_type::coeff_t> w;
  typedef tmptypmat<typename typP::value_type::coeff_t> typMat;
  typMat Mat,L;
  int nouveaumon,allpolused=0,maxdeg,k,*pr,*pc;
  vector<typename typMk::value_type*> exceed;
  list<int> sizeexceed;
  //cout<<"avant init"<<endl;
  
  init_algo<typMat>(P,nouveaumon,maxdeg,k,Pk,Mk,b,dump,w,serv);
  //cout<<"apres init"<<endl;
#ifdef GOTZMANN
  int nextdim=-1;
  initial_mon(k,b,dump,serv);
  while((!critere_got(b,dump,serv))||(k<=maxdeg))
#else  
    while(nouveaumon||(k<=maxdeg))
#endif
    {
    int killspol;
    //    cout<<"Spolordsup"<<endl;
    S=Spolordsup(Pk,b,Spoldejavu);
#ifdef CACABOUDIN    
 for(typename typPk::iterator iter=S.begin();iter!=S.end();iter++)
      {
    //   cout<<"poly de taille "<<iter->size<<endl;
    cout<<invconv<typename typP::value_type>(*iter,serv)<<endl;       
      }
#endif
#if 0//ICICICICICICICICICICCII
cout<<"poly a la fin de la boulcle de l'algo "<<Pk.size()<<endl;
for(typename typPk::iterator iter=Pk.begin();iter!=Pk.end();iter++)
 {
   //   cout<<"poly de taille "<<iter->size<<endl;
   cout<<invconv<typename typP::value_type>(*iter,serv)<<endl;       
}
    cout<<endl<<"newpkmkmatrixof"<<endl;
#endif
    killspol=!newpkmkMatrixofSolve(Mat,L,Pk,P,Mk,k,b,dump,exceed,sizeexceed
                   ,w,serv);
#if 0//ICICICICICICICICCI
    cout<<"poly apres solve (attention exceed)"<<endl;
    copy(sizeexceed.begin(),sizeexceed.end(),ostream_iterator<int>(cout," *** "));cout<<endl;
    for(typename typPk::iterator iter=Pk.begin();iter!=Pk.end();iter++)
      {
    
    cout<<invconv<typename typP::value_type>(*iter,serv)<<endl;       
      }
#endif
    if(killspol) S.erase(S.begin(),S.end());
    if (k+1>=maxdeg) allpolused=1;
    if(Mk.size()==0) nouveaumon=0;
    //    cout<<"le rang "<<my_rank(Mat)<<endl;
    // cout<<"Mk.size() "<<Mk.size()<<endl;
   // copy(Mk.begin(),Mk.end(),ostream_iterator<mon>(cout," "));
   // cout<<endl;
    if ((unsigned)my_rank(Mat)==Mk.size())
      //matrice de rang plein
      //donc tous les exceed sont vides
      {
    typMk tmpSpoldejavu;
    typPk redspol,secd=secondmembre(Mat,Pk);;
    //secondmembre2(Mat,Pk,secd);
    //typename typPk::iterator iter=Pk.begin();
    cout<<"matrice de rang plein "<<endl;
    Dump(Pk,dump);
#ifdef GOTZMANN
    nextmon(k,b,dump,serv);
#endif
    cout<<"secd/size() "<<secd.size()<<endl;
    ReduceSpol(S,dump,redspol,b,w,serv);//tmpSpoldejavu);
    my_merge(redspol,secd);
    cout<<"redspol.size "<<redspol.size()<<endl;
    if(redspol.size()==0)
      {
        k++;
        my_merge(Spoldejavu,tmpSpoldejavu);
        L.destroystore();
        continue;
        //la matrice a rang plein 
        //ET les Spol\secondmembre de ce niveau se ,reduisent a 0
      }
    else
      {
        //on a des spol\secondmembre qui se reduisent pas a 0
        //cad un defaut de commuatation
        int tmpdeg=mindeg(redspol,serv);
        max_degppcm=k+2;
        nouveaumon=1;
        cout<<"redspol.size() "<<redspol.size()<<endl;

        if (tmpdeg==(k+1))
          my_merge(Spoldejavu,tmpSpoldejavu);
        k=tmpdeg;

        recover(Pk,Mk,dump,k);

        forget(b,dump,k,recall);
        cout<<"debut newcrochetchoix"<<endl;
        NewCrochetchoix(P,Pk,Mk,redspol,k,dump,b,recall,w,serv);
        cout<<"fin newcrochetchoix"<<endl;

        NewDestroyredspol(redspol,Pk);
        serv.compress(Pk,dump,b,k);
        Dump(Pk,dump);
#ifdef GOTZMANN
        nextmon(k,b,dump,serv);
#endif
        L.destroystore();
        //misa jour de B dans crochet choix!!
        continue;
      }
      }
    else
      //la matrice n'est pas de rang plein
      {
    cout<<"matrice de rang pas plein "<<endl;
    typPk secd;
    cout<<"addB"<<endl;
#ifdef DEBUG11  
    cout<<"avant addB "<<Pk.size()<<endl;
    auto iterexceed=exceed.begin();
    auto itersizeexceed=sizeexceed.begin();
    for(typename typPk::iterator iter=Pk.begin();iter!=Pk.end();iter++,iterexceed++,itersizeexceed++)
      {
        cout<<"ind "<<iter->ind<<" est sikzeexceed "<<*itersizeexceed<<endl;
        copy(*iterexceed,*iterexceed+*itersizeexceed,ostream_iterator<mon>(cout," + "));
        cout<<endl;
        cout<<invconv<typename typP::value_type>(*iter,serv)<<endl;       
      }
#endif
    AddB(Mat,Pk,Mk,exceed,sizeexceed,dump,b,serv);
#if 0
    cout<<"apres addB "<<Pk.size()<<endl;
    for(typename typPk::iterator iter=Pk.begin();iter!=Pk.end();iter++)
      {
        cout<<invconv<typename typP::value_type>(*iter,serv)<<endl;       
      }
#endif
    //
    //secondmembre supprime des trucs dans Pk!!
    //
    //  cout<<"second membre"<<endl;
    secd=secondmembre(Mat,Pk);
    //  cout<<"voisin"<<endl;
    rajoute_voisins(Pk,Mk,dump,b,w,serv);
    //  cout<<"fin voisin"<<endl;
    //ne selectionne que les non nuls
       if(secd.size()==0)
     {
       typPk redspol;
       Dump(Pk,dump);
#ifdef GOTZMANN
       nextmon(k,b,dump,serv);
#endif
      
       ReduceSpol(S,dump,redspol,b,w,serv);//,Spoldejavu);
#ifdef NO_TWICE    
for(typename typPk::iterator iter=redspol.begin();iter!=redspol.end();iter++)
      P.push_back(invconv<polyalp>(*iter));
#endif
       if(redspol.size()==0)
         {
           k++;
           L.destroystore();
           continue;
           //la matrice a pas rang plein second membre =0 
           //ET les Spol de ce niveau se ,reduisent a 0
         }
       else
         {
           //on a des spol qui se reduisent pas a 0
           //cad un defaut de commutation 
           k=mindeg(redspol,serv);
           recover(Pk,Mk,dump,k);
           forget(b,dump,k,recall);
           cout<<"il ya un passage par mat rang pas"
               <<"plein et pol non zero"<<endl;
           NewCrochetchoix(P,Pk,Mk,redspol,k,dump,b,recall,w,serv);
           NewDestroyredspol(redspol,Pk);
           serv.compress(Pk,dump,b,k);
           Dump(Pk,dump);
#ifdef GOTZMANN
           nextmon(k,b,dump,serv);
#endif
 
           L.destroystore();
           continue;
         }
     }
       else
     {
       //on a des pol qui se reduisent pas a 0
       //cad un defaut de commuatation dans les coins
       typPk redspol;
       Dump(Pk,dump);
#ifdef GOTZMANN
       nextmon(k,b,dump,serv);
#endif
       ReduceSpol(S,dump,redspol,b,w,serv);//,Spoldejavu);
       my_merge(secd,redspol);
#ifdef NO_TWICE    
for(typename typPk::iterator iter=redspol.begin();iter!=redspol.end();iter++)
         P.push_back(invconv<polyalp>(*iter));
#endif
           k=mindeg(secd,serv);
       recover(Pk,Mk,dump,k);
       forget(b,dump,k,recall);
#if 0
        cout<<"poly avant newcrochetchoix (attention exceed)"<<endl;
        copy(sizeexceed.begin(),sizeexceed.end(),ostream_iterator<int>(cout," *** "));cout<<endl;
        for(typename typPk::iterator iter=Pk.begin();iter!=Pk.end();iter++)
          {
        
        cout<<invconv<typename typP::value_type>(*iter,serv)<<endl;       
          }
#endif

       NewCrochetchoix(P,Pk,Mk,secd,k,dump,b,recall,w,serv);
#if 0
        cout<<"poly apres newcrochetchoix (attention exceed)"<<endl;
        copy(sizeexceed.begin(),sizeexceed.end(),ostream_iterator<int>(cout," *** "));cout<<endl;
        for(typename typPk::iterator iter=Pk.begin();iter!=Pk.end();iter++)
          {
        
        cout<<invconv<typename typP::value_type>(*iter,serv)<<endl;       
          }
#endif

       NewDestroyredspol(secd,Pk);
       serv.compress(Pk,dump,b,k);
       Dump(Pk,dump);
#ifdef GOTZMANN
       nextmon(k,b,dump,serv);
#endif

       L.destroystore();
       continue;
     }
      
      }
  }
  
  Mat.destroystore();
  //while nouveau mon etc...
  destroy_space(&w);
  cout<<"sortie d'algo flash et dimension "<<endl;
  return ;//dump;
  //comme ca on conserve le tri par degre des formes normales et 
  //donc la construction des operateurs de multiplication se 
  //fait plus vite.. :-)

}


template<typename mon, typename T,typename Base ,typename monomial_server>
mon choice_naif(const pol2ter<mon,T> & p, const Base &b, monomial_server &serv)
{
  list<mon> tmp;
  mon m; 
  int maxdeg=-1,k=0;
  //int maxtmpdeg=0;
  for(int i=0;i<p.sizenf;i++)
    {
      if (p.nfind[i])
    for(int j=0;j<8;j++)
      if(p.nfind[i]>>j&1)
      {
        mon tmpmon;
        serv.int2mon(8*i+j,tmpmon);
        if(tmpmon.GetDegree()>maxdeg)
          {
        tmp.erase(tmp.begin(),tmp.end());
        maxdeg=tmpmon.GetDegree();
          }
        if(tmpmon.GetDegree()==maxdeg)
          {
        tmpmon.SetCoeff(p.nf[k]);
        tmp.push_back(tmpmon);
          }
        k++;
      }
    }
  m=tmp.front();
  //  cout<<"monom choisi "<<m<<endl;
  return m;
}

template<typename mon, typename T,typename Base, typename monomial_server>
mon choice_num(const pol2ter<mon,T> & p, const Base &b, monomial_server &serv)
{
  list<mon> tmp;
  mon m; 
  int maxdeg=-1,k=0;
  //int maxtmpdeg=0;
  for(int i=0;i<p.sizenf;i++)
    {
      if (p.nfind[i])
    for(int j=0;j<8;j++)
      if(p.nfind[i]>>j&1)
      {
        mon tmpmon;
        serv.int2mon(8*i+j,tmpmon);
        if(tmpmon.GetDegree()>maxdeg)
          {
        tmp.erase(tmp.begin(),tmp.end());
        maxdeg=tmpmon.GetDegree();
          }
        if(tmpmon.GetDegree()==maxdeg)
          {
        tmpmon.SetCoeff(p.nf[k]);
        tmp.push_back(tmpmon);
          }
        k++;
      }
    }
  m=tmp.front();
  for(typename list<mon>::iterator itertmp=tmp.begin();itertmp!=tmp.end();itertmp++)
    {
      if(m.GetCoeff()<itertmp->GetCoeff()) {m=*itertmp;};
    }
  //cout<<"--------------------------------------------------------------------------------"<<endl;
  // cout<<"monom choisi "<<m<<endl;
  //cout<<"pol"<<endl<<invconv<Poly>(p,serv)<<endl;
  //cout<<"--------------------------------------------------------------------------------"<<endl;
  return m;
}

template<typename mon, typename T,typename Base, typename monomial_server>
mon choice_grevlex(const pol2ter<mon,T> & p,const Base &b, 
           monomial_server &serv)
{
  list<mon> tmp;
  mon m; 
  int maxdeg=-1,k=0;
  //int maxtmpdeg=0;
  for(int i=0;i<p.sizenf;i++)
    {
      if (p.nfind[i])
    for(int j=0;j<8;j++)
      if(p.nfind[i]>>j&1)
        
    {
      mon tmpmon;
      serv.int2mon(8*i+j,tmpmon);
      if((tmpmon.GetDegree()>maxdeg) && (p.nf[k]!=0) )
        {
          tmp.erase(tmp.begin(),tmp.end());
          maxdeg=tmpmon.GetDegree();
        }
      if(tmpmon.GetDegree()==maxdeg)
        {
          tmpmon.SetCoeff(p.nf[k]);
          tmp.push_back(tmpmon);
        }
      k++;
    }
    }
  for(;tmp.front().GetCoeff()==0;tmp.pop_front());
  m=tmp.front();
  for(typename list<mon>::iterator itertmp=tmp.begin();itertmp!=tmp.end();itertmp++)
    {
      if((itertmp->GetCoeff()!=0) && (my_ord(*itertmp,m)<0)) {m=*itertmp;};
    }
  // cout<<"monom choisi "<<m<<endl;
  return m;
}

template<typename mon, typename T,typename Base, typename monomial_server>
mon choice_dlex(const pol2ter<mon,T> & p,const Base &b, monomial_server &serv)
{
  list<mon> tmp;
  mon m; 
  int maxdeg=-1,k=0;
  //int maxtmpdeg=0;
  for(int i=0;i<p.sizenf;i++)
    {
      if (p.nfind[i])
    for(int j=0;j<8;j++)
      if(p.nfind[i]>>j&1)
        {
          mon tmpmon;
          serv.int2mon(8*i+j,tmpmon);
          if(tmpmon.GetDegree()>maxdeg)
        {
          tmp.erase(tmp.begin(),tmp.end());
          maxdeg=tmpmon.GetDegree();
        }
          if(tmpmon.GetDegree()==maxdeg)
        {
          tmpmon.SetCoeff(p.nf[k]);
          tmp.push_back(tmpmon);
        }
          k++;
        }
    }
  m=tmp.front();
  for(typename list<mon>::iterator itertmp=tmp.begin();itertmp!=tmp.end();itertmp++)
    {
      if(my_ord(*itertmp,m)>0) {m=*itertmp;};
    }
  //  cout<<"monom choisi "<<m<<endl;
  return m;
}


template<typename mon, typename T,typename Base, typename monomial_server>
mon choice_mac(const pol2ter<mon,T> & p,const Base &b, monomial_server &serv)
{
  list<mon> tmp;
  mon m; 
  int maxdeg=-1,k=0;
  int maxtmpdeg=0;
  for(int i=0;i<p.sizenf;i++)
    {
      if (p.nfind[i])
    for(int j=0;j<8;j++)
      if(p.nfind[i]>>j&1)
        {
          mon tmpmon;
          serv.int2mon(8*i+j,tmpmon);
          if(tmpmon.GetDegree()>maxdeg && !Iszero(p.nf[k]))
        {
          tmp.erase(tmp.begin(),tmp.end());
          maxdeg=tmpmon.GetDegree();
        }
          if(tmpmon.GetDegree()==maxdeg && !Iszero(p.nf[k]))
        {
          tmpmon.SetCoeff(p.nf[k]);
          tmp.push_back(tmpmon);
        }
          k++;
        }
    }
  m=tmp.front();
  for(typename list<mon>::iterator itertmp=tmp.begin();itertmp!=tmp.end();itertmp++)
    {
      //cout<<"monomes regarde "<<*itertmp<<" Lavr  "<<lvar(itertmp->rep)<<endl;
      for(int i=0;i<=lvar(itertmp->rep);i++)
    {
      if(itertmp->GetDegree(i)>=maxtmpdeg)
        {
          maxtmpdeg=itertmp->GetDegree(i);
          m=*itertmp;
        }
    }
    }
  return m;
}


template<typename mon, typename T,typename Base, typename monomial_server>
mon choice_mpq(const pol2ter<mon,T> & p,const Base &b, monomial_server &serv)
{
  int i,degp=Degree(p),k=0;
  int mini=INT_MAX;
  mon tmpmon;
  mon stockmon;
  for(i=0;i<p.sizenf;i++)
    if(p.nfind[i])
      for(int j=0;j<8;j++)
    if(p.nfind[i]>>j&1)
      {
        serv.int2mon(8*i+j,tmpmon);
        if(tmpmon.GetDegree()==degp)
          {
        __mpz_struct  *numref,*denref;
        numref=mpq_numref(&(p.nf[k].rep()));
        denref=mpq_denref(&(p.nf[k].rep()));
        if(mpz_size(numref)+
           mpz_size(denref)<mini) 
          {
            mini=mpz_size(mpq_numref(&(p.nf[k].rep())))+
              mpz_size(mpq_denref(&(p.nf[k].rep())));
            //cout<<"newmax "<<max<<endl;
            //cout<<" les monomes "<<*iter<<endl;
            stockmon=tmpmon*p.nf[k];
          }
          }
        k++;
      }
  //cout<<"je retourne "<<stockmon<<endl;
  return stockmon;
}

template<typename mon, typename T, typename Base, typename monomial_server>
mon choice_fall(const pol2ter<mon,T> & p, const Base &b, monomial_server &serv)
{
  typedef typename mon::coeff_t coeff_t;
  int i,k=0;
  int maxi=-1;
  mon tmpmon;
  mon stockmon;
  mon res=choice_mac(p,b,serv);
  typename Base::const_iterator iterb;
  list<mon> tmplist;
  int degp=Degree(p,serv);
#ifdef DEB 
  cout<<"Degree(p) "<<degp<<endl;
#endif
  for(i=0;i<p.sizenf;i++)
    if(p.nfind[i])
      for(int j=0;j<8;j++)
    if((p.nfind[i]>>j)&1) 
      {
        serv.int2mon(8*i+j,tmpmon);
        if ((tmpmon.GetDegree()==degp) && (!Iszero(p.nf[k])))
          tmplist.push_back(tmpmon*p.nf[k]);
        k++;    
        //    cout<<"k "<<k<<" "<<tmpmon<<endl;
      }
  for(typename list<mon>::iterator iterlist=tmplist.begin()
    ;iterlist!=tmplist.end();iterlist++)
    {
      for(iterb=b.begin();iterb!=b.end();iterb++)
    {
      int isgood=0;
      for(int i=0;i<iterb->taille1;i++)
        if(isdivisible(iterb->accept[i],*iterlist,b))
          {isgood=1;break;}
      for(int i=0;i<iterb->taille2;i++)
        if(isdivisible(iterb->refuse[i],*iterlist,b))
          {isgood=0;break;}
      if(isgood)
        {

          for(int i=0;i<iterb->taille2;i++)
        {
          mon m=lcm(*iterlist,iterb->refuse[i]);
          int dist=0;
          for(int j=0;j<iterb->taille2;j++)
            if ((j!=i) && (isdivisible(iterb->refuse[j],m,b)))
              isgood=0;
          if(isgood)
            {
              for(int j=0;j<b.nbvar();j++)
            if(m.GetDegree(j)>iterb->refuse[i].GetDegree(j))
              dist+=m.GetDegree(j)-iterb->refuse[i].GetDegree(j);
              if(dist>maxi)
            {
              res=*iterlist;
              maxi=dist;
            }
            }
        }
        }
      
    }
    }
#if 1
  cout<<"le pol "<<invconv<Poly>(p,serv)<<endl<<"le mon choisi "<<res<<endl;
#endif
  return res;
}
#ifndef ASYNC
template<typename typdump,typename Base,typename mon, 
     typename workspace, typename monomial_server> 
pol<mon,COEFF> 
inline mult(int i,const pol<mon, COEFF> &p,const typdump &dump,
        const Base &b,mon **exce,int &exsize,workspace &w,
        monomial_server & serv)
{
  static const mon monzero=mon(0);
  static COEFF zero=0;
  static int flag=1;
  pol<mon,COEFF> ret;
  (*exce)=(mon*)MAC_REV_MALLOC<mon>((p.size+1)*sizeof(mon));
  int *stockindicetab=(int*)malloc(p.size*sizeof(int));
  int noquo=0,maxind=0,pp=0,q=0;
  COEFF *tmpspace;
  if(flag)
    {
      init_pointer<COEFF>();
      flag=0;cout<<"pouet pouet"<<endl;
    }
  //#ifdef DEBUG
#ifdef DEBUG
  cout<<"ds mult je traite "<<p.ind<<" ---- "<<i<<endl;
  cout<<invconv<Poly>(p,serv)<<endl;
  cout<<"p.size "<<p.size<<endl;
#endif
  for(int j=0;j<p.sizenf;++j)
    if (p.nfind[j])
     {
       for(int k=0;k<8;k++)
     if((p.nfind[j]>>k)&1)
       {
         int stockindice2=serv.mulind(8*j+k,i,b);
         if(stockindice2<0)
           stockindice2=-1*(8*j+k)-1;
         stockindicetab[pp++]=stockindice2;
         if (stockindice2>=maxind)
           maxind=stockindice2;
         cout<<"stockindice2 "<<8*j+k<<" "<<i<<stockindice2<<endl;
       }
     }
  for(;maxind/8+1>=w.sizenf;increase_space(&w));
  for(int j=0;j<pp;j++)
    if(stockindicetab[j]>=0)
      {
    w.tabcoeff[stockindicetab[j]]+=p.nf[j];
    if(!Iszero(w.tabcoeff[stockindicetab[j]]))
      w.tabnf[stockindicetab[j]/8]|=1<<(stockindicetab[j]%8);
      }
    else  
      {
    pol<mon,COEFF>  tmp;
    tmp=serv.nf(i,-1*(stockindicetab[j]+1),dump,b);
          if (tmp.size!=-1)
        {
          int compteurtmp=0;
          //On fait un plus dense dans tmpspace
          int k,l;
          unsigned long *tmpptr=(unsigned long*)tmp.nfind;
          unsigned long mask=1;
          for(;tmp.sizenf>=w.sizenf;increase_space(&w));

          for(k=0;k<tmp.sizenf/sizeof(unsigned long);k++)
            if(tmpptr[k])
              {
            mask=tmpptr[k];
            ((unsigned long*)w.tabnf)[k]|=mask;
            for(l=0;l<8*sizeof(unsigned long);l+=8,mask>>=8)
              {
                if(mask&255)
                  {
                tabfun2[mask&255](
                          w.tabcoeff+(8*sizeof(unsigned long)*k+l)
                          ,&(p.nf[j]),tmp.nf+compteurtmp);
                compteurtmp+=nbbits[mask&255];
                  }
#if 0
                for(int m=0;m<8;m++)
                  if(Iszero(w.tabcoeff[8*k*sizeof(unsigned long)+l+m]))
                ((unsigned long*)w.tabnf)[k]&=~(1<<(l+m));
#endif
              }
              }

          for(k=tmp.sizenf-(tmp.sizenf%sizeof(unsigned long))
            ;k<tmp.sizenf;k++)
            if(tmp.nfind[k])
              {
            w.tabnf[k]|=tmp.nfind[k];
            mask=1;
            for(l=0;l<8;l++,mask<<=1)
              {
                if(tmp.nfind[k]&mask)
                  w.tabcoeff[8*k+l]-=p.nf[j]*tmp.nf[compteurtmp++];
                if(Iszero(w.tabcoeff[8*k+l]))
                  w.tabnf[k]&=~mask;

              }
              }
        }
          else  
        {
          //On est tombé sur un monome du degre courant
          //qui est sur le bord
          mon tmpmon;
          (*exce)[noquo]=mon(i,1);
          serv.int2mon(-1*(stockindicetab[j]+1),tmpmon);
          (*exce)[noquo]*=tmpmon;//int2mon<mon>(j);
          (*exce)[noquo]*=p.nf[j];
          ++noquo;
        };
      }
#ifdef ANCESTOR
  //On met juste a jour l'ancetre pour la construction du grpahe apres
  ret.ancestor=p.ancestor;
#endif
  ret.ind=p.ind;
  ret.ind*=mon(i,1);
  (*exce)[noquo++]=ret.ind;
  exsize=noquo;
  (*exce)=(mon*)MAC_REV_REALLOC<mon>((*exce),(p.size+1)*sizeof(mon)
                   ,(noquo)*sizeof(mon));
  ret.ind.SetCoeff(noquo+1);
  ret.size=0;
  pp=0;
  for(int j=0;j<w.sizenf;j++)
    if(w.tabnf[j])
      {
    for(int m=0;m<8;m++)
      if(Iszero(w.tabcoeff[8*j+m]))
        w.tabnf[j]&=~(1<<m);
    if(w.tabnf[j])
      {
        ret.size+=nbbits[w.tabnf[j]];
        pp=j;
      }
      }
  ret.nf=(COEFF*)MAC_REV_MALLOC<COEFF>(ret.size*sizeof(COEFF));
  ret.sizenf=pp+1;
  ret.nfind=(unsigned char *)malloc(ret.sizenf);
  pp=0;
  for(int j=0;j<ret.sizenf;j++)
    {
      if(w.tabnf[j])
    for(int k=0;k<8;k++)
      if((w.tabnf[j]>>k)&1)
        {
          ret.nf[pp++]=w.tabcoeff[8*j+k];
          w.tabcoeff[8*j+k]=zero;
        }
      ret.nfind[j]=w.tabnf[j];
      w.tabnf[j]=0;
    }
  free(stockindicetab);
#ifdef DEBUG
  cout<<"en sortie"<<endl;
  cout<<invconv<Poly>(ret,serv)<<endl;
  copy(*exce,*exce+exsize,ostream_iterator<mon>(cout," + "));cout<<endl;
  cout<<"----------"<<endl;
#endif
  return ret;
}
#endif

#if 1
template<typename typmat,typename typPk, typename typMk,typename typdump,
     typename typexceed, typename typsizeexceed, typename typB,typename monomial_server>
     void AddB(const typmat & mat, typPk & Pk,const typMk &Mk,const typexceed & exceed,
           const typsizeexceed & sizeexceed, const typdump & dump,
      typB &b, monomial_server &serv)
{
  typedef typename typPk::value_type pol;
  typedef typename pol::monom_t mon;
  typedef typename pol::coeff_t coeff_t;
  typename typPk::iterator iterpk=Pk.begin();;
  auto iterexceed=exceed.begin();
  auto itersizeexceed=sizeexceed.begin();
  predicat<mon> tmp;
  typename typMk::const_iterator iter=Mk.begin();
  int r=my_rank(mat),comprefuse=0,compaccept=0;
  unsigned int i;
  //cout<<"DANS ADDB"<<endl;
  //cout<<"rang de la matrice m "<<r<<" Mk.size "<<Mk.size()<<endl;
  //cout<<"taille de la matrice "<<mat.nrow<<"x"<<mat.ncol<<endl;
  tmp.refuse=(mon*)MAC_REV_MALLOC<mon>(r*sizeof(mon));
  tmp.taille2=r;
  tmp.accept=(mon*)MAC_REV_MALLOC<mon>((Mk.size()-r)*sizeof(mon));
  tmp.taille1=Mk.size()-r;
  //for(i=0;i<(unsigned)r;i++,iter++)
  //  tmp.refuse[i]=*iter;
  iter=Mk.begin();
  for(i=0;i<Mk.size();i++,iter++)
    {
      if(Iszero(mat(i,i))) 
    {
      tmp.accept[compaccept++]=*iter;
      //cout<<"monome rajoute au quotient "<<*iter<<endl;
    }
      else
    {
      tmp.refuse[comprefuse++]=*iter;
    }
    }
  b.def.push_back(tmp);
  //bon maintenant que l'on vient de rajouter des 
  //monomes au quotient il faut remttre les polynomes comme il faut
  iter=Mk.begin();
//cout<<"report"<<endl;
//On fait une passe pour determiner les polynomes a modifier et on effectue 
//les modifications en lisant les exceed (au passage on FREE)

  
  for(unsigned int i=0;i<Pk.size();i++,iterpk++,iter++,iterexceed++,itersizeexceed++)
    {
      vector<sortstruct<coeff_t> > decal;
      int maxindex=-1;
      int expend=0;
      for(int j=0;j<*itersizeexceed;j++)
    if(!Iszero((*iterexceed)[j].GetCoeff()))
      //si des monomes nuls sortent comme coeff de la matrice
      //et bien on les prends pas
      {
        sortstruct<coeff_t> tmp;
        expend++;
        maxindex=max(maxindex,serv.mon2int((*iterexceed)[j]));
        tmp.index=serv.mon2int((*iterexceed)[j]);
        tmp.coeff=(*iterexceed)[j].GetCoeff();
        decal.push_back(tmp);
      }
      // Si on a rien a rajouter on arrete
      // cout<<"decal.empty "<<decal.empty()<<endl;
      if(decal.empty()) continue;
      //on augmente la taille du polynome  considere
      iterpk->nf=(typename mon::coeff_t*)
    MAC_REV_REALLOC<typename mon::coeff_t>(iterpk->nf
           ,iterpk->size*sizeof(typename mon::coeff_t)
           ,(iterpk->size+expend)*sizeof(typename mon::coeff_t));
      if((maxindex/8)>=iterpk->sizenf)
    {
      iterpk->nfind=(unsigned char*)
        MAC_REV_REALLOC<unsigned char>(iterpk->nfind
                       ,iterpk->sizenf,maxindex/8+1);
      for(int i=iterpk->sizenf;i<maxindex/8+1;i++)
        iterpk->nfind[i]=0;
      iterpk->sizenf=maxindex/8+1;
      
    }
// la place est faite il reste plus qu'a decaler les monomes de iterpk->nf
//pour intercaler eventuellment les nouveaux monomes
       stable_sort(decal.begin(),decal.end(),sortstructorder<coeff_t>());
       int nbcoeff=iterpk->size-1;
      // cout<<"iterpk->sizenf "<<iterpk->sizenf<<endl;
       for(int j=iterpk->sizenf-1;j>=0;j--)
        {
          //cout<<"decal.index "<<decal.rbegin()->index<<endl;
          if(8*j>decal.rbegin()->index)
            for (int jj=0;jj<nbbits[iterpk->nfind[j]];jj++)
              {
                iterpk->nf[nbcoeff+decal.size()]=
                  iterpk->nf[nbcoeff];
                nbcoeff--;
            }
          else
            {
              int indexnf=decal.back().index;
              typename typMk::value_type tmpmon;
              //cout<<"8*j"<<8*j<<" indexnf "<<indexnf<<" valeur "<<decal.back().index<<endl;
              int jj=7;
              while(8*j<=decal.back().index)
              {
              indexnf=decal.back().index % 8;
              for (;jj>(indexnf);jj--)
              if(iterpk->nfind[j]&(1<<jj))
                {
                 iterpk->nf[nbcoeff+decal.size()]=
                 iterpk->nf[nbcoeff];
                 nbcoeff--;
                }
                      //cout<<"je rajoute "<<decal.back().coeff<<" a la place "<<nbcoeff+decal.size()<<endl;
//          if(nbcoeff+decal.size()+1<iterpk->size+expend)
            //  cout<<"coeff apres "<<iterpk->nf[nbcoeff+decal.size()+1]<<endl;
              iterpk->nf[nbcoeff+decal.size()]=decal.back().coeff;
              iterpk->nfind[j]|=(1<<(indexnf));
          //   serv.int2mon(8*j+indexnf,tmpmon);
        //  cout<<"monome a la place "<<decal.back().index<<" "<<tmpmon<<endl;    
          //   serv.int2mon(decal.back().index,tmpmon);
        //  cout<<"monome a la place "<<decal.back().index<<" "<<tmpmon<<endl;    
              decal.pop_back();
             jj--; 
            if(decal.empty()) break;
            }   
             for (;jj>=0;jj--)
              if(iterpk->nfind[j]&(1<<jj))
               {
                 iterpk->nf[nbcoeff+decal.size()]=
                        iterpk->nf[nbcoeff];
                  nbcoeff--;
                 }
            if(decal.empty()) break;
            }
                
        }
       iterpk->size+=expend;
       MAC_REV_FREE<mon>(*iterexceed,*itersizeexceed*sizeof(mon));
    }   
  return;
}
#endif
#ifdef ASYNC
#include "pol2terspecial-async.hpp"
#endif