powerprod.cpp

#include "powerprod.H"
using namespace std;


// methods for the PowProd class

AutoAlphabet PowProd::default_alphabet;
Alphabet* PowProd::alphabet = &PowProd::default_alphabet;
int PowProd::order = LEX;

PowProd::PowProd(){}

PowProd::PowProd(long var_num){
  long i;

  if(var_num > 0){
    powers.resize(var_num);
    for(i=0; i < var_num -1; i++)
      powers[i]=0;
    powers[var_num -1]=1;
  }
}


PowProd::PowProd(long var_num, long thepower){
  long i;

  if(var_num > 0 && thepower > 0){
    powers.resize(var_num);
    for(i=0; i < var_num -1; i++)
      powers[i]=0;
    powers[var_num -1]= thepower;
  }
}



PowProd::PowProd(const PowProd& that){
  powers=that.powers;
}

long PowProd::nvars() const {
  return powers.size();
}

long PowProd::degree() const {
  long n=0;

  for(long i=0; i < powers.size(); i++)
    n += powers[i];

  return n;
}

long PowProd::power_of(long i) const {
  if( (i > powers.size()) || (i<1) )
    return 0;
  else
    return powers[i-1];
}

PowProd& PowProd::operator=(const PowProd& that){
  powers=that.powers;
  return *this;
}

PowProd& PowProd::operator*=(const PowProd& that){
  long n, i;

  if( that.powers.size() <= powers.size() ) {
    for(i=0; i < that.powers.size(); i++)
      powers[i] += that.powers[i];
  }
  else {
    n= powers.size();
    powers.resize( that.powers.size() );
    for(i=0; i < n; i++)
      powers[i] += that.powers[i];
    for( ; i < powers.size(); i++)
      powers[i]=that.powers[i];
  }
  return *this;
}

PowProd PowProd::operator*( const PowProd& that) const {
  PowProd ans;
  // we try to avoid resizing ans too many times:
  if(powers.size() >= that.powers.size()){
    ans= (*this);
    ans *= that;
  }
  else{
    ans= that;
    ans *= (*this);
  }
  return ans;
}

bool PowProd::lex_lessthan(const PowProd& that) const {
  long n,m;

  n= powers.size();
  m= that.powers.size();

  if( n != m )
    return n < m;
  else{
    do {n--;}
    while( (n >= 0) && (powers[n]==that.powers[n]) );
    if( n < 0 )
      return false;
    else
      return (powers[n] < that.powers[n]);
  }
}

bool PowProd::lexsv_lessthan(const PowProd& that, long sv) const {
  long i,j;
  PowProd a,b;

  i= power_of(sv);
  j= that.power_of(sv);
  if( i != j )
    return i<j;
  //else...
  a= *this;
  b= that;
  if( sv <= powers.size() )
    a.powers[sv-1]= 0;
  if( sv <= that.powers.size() )
    b.powers[sv-1]= 0;
  return a.lex_lessthan(b);
}



bool PowProd::deglex_lessthan(const PowProd& that) const {
  long n,m;

  n= degree();
  m= that.degree();

  if( n != m )
    return n < m;
  else
    return this->lex_lessthan(that);
}

bool PowProd::deglexsv_lessthan(const PowProd& that, long sv) const {
  long i,j;
  PowProd a,b;

  i= power_of(sv);
  j= that.power_of(sv);

  if( i != j )
    return i<j;
  //else...
  a= *this;
  b= that;
  if( sv <= powers.size() )
    a.powers[sv-1]= 0;
  if( sv <= that.powers.size() )
    b.powers[sv-1]= 0;
  return a.deglex_lessthan(b);
}

bool PowProd::degrevlex_lessthan(const PowProd& that) const {
  long n,m,i;

  n= degree();
  m= that.degree();

  if( n != m )
    return n < m;
  else{
    n=powers.size();
    m=that.powers.size();
    n= (n < m) ? n : m;
    for(i=0; i < n; i++){
      if(powers[i]!=that.powers[i])
        break;
    }
    if(i==n)
      return false; // this==that in this case
    else
      return powers[i] > that.powers[i];
  }
}

bool PowProd::degrevlexsv_lessthan(const PowProd& that, long sv) const {
  long i,j;
  PowProd a,b;

  i= power_of(sv);
  j= that.power_of(sv);
  if( i != j )
    return i<j;
  //else...
  a= *this;
  b= that;
  if( sv <= powers.size() )
    a.powers[sv-1]= 0;
  if( sv <= that.powers.size() )
    b.powers[sv-1]= 0;
  return a.degrevlex_lessthan(b);
}




bool PowProd::operator<(const PowProd& that) const {
  switch(order){
  case LEX:
    return this->lex_lessthan(that);
  case DEGLEX:
    return this->deglex_lessthan(that);
  case DEGREVLEX:
    return this->degrevlex_lessthan(that);
  }
  return false; // in case order has a funny value
  // also keeps the compiler happy.
}

void PowProd::use_order(int theorder){
  order=theorder;
}

int PowProd::current_order(){
  return order;
}

bool PowProd::operator==(const PowProd& that) const {
  return powers == that.powers;
}


bool PowProd::divides(const PowProd& that) const {
  bool ans;

  if(powers.size() > that.powers.size())
    return false;
  else{
    ans=true;
    for(long i=0; i < powers.size(); i++)
      if( powers[i] > that.powers[i]){
        ans=false;
        break;
      }
  }
  return ans;
}

PowProd& PowProd::operator/=(const PowProd& that){
  Tuple<long> temp;
  long i,j;
  // set i=nb of substractions
  if(powers.size() > that.powers.size()){
    i=that.powers.size();
    temp= powers;
  }
  else{
    for(i=powers.size() -1; i>=0;i--){
      if(powers[i]!=that.powers[i])
        break;
    }
    i++;
    temp.resize(i);
  }
  //compute the answer
  for(j=0; j<i; j++)
    temp[j]=powers[j] - that.powers[j];
  //return
  powers= temp;
  return *this;
}

PowProd PowProd::operator/(const PowProd& that) const {
  PowProd ans;

  ans= *this;
  ans /= that;
  return ans;
}

bool PowProd::prime_to(const PowProd& that) const {
  long i,min;

  min= powers.size() > that.powers.size() ? that.powers.size() : powers.size();

  for(i=0; i<min; i++)
    if( (powers[i] != 0 ) && (that.powers[i] != 0 ) )
      return false;

  return true;

}



void PowProd::swap(long i, long j){
  long max, min, k, dummy, ii, jj;

  ii=i-1;
  jj=j-1;
  max= i > j ? ii : jj;
  min= i > j ? jj : ii;
  // easy case: everything in scope, size will not change
  if( max < powers.size() ){
    if( (max < powers.size()-1) || powers[min] !=0 ){
      dummy= powers[ii];
      powers[ii]= powers[jj];
      powers[jj]= dummy;
      return;
    }// if not, we need to reduce the size
    // in this case max= powersize -1, powers[max] is 
    // nonzero, powers[min] is zero, so min < max
    powers[min]= powers[max];
    for(k= max-1; powers[k]==0; k--);
    powers.resize(k+1);
    return;
  }
  // if none of this happened, we have max >= powersize
  // check if eveything is out of scope
    if( min >= powers.size() || powers[min]==0 )
      return; // in this case the swap involves variables
              // which are not there in the powprod: we return.
    // final case: we resize
    k= powers.size();
    powers.resize(max+1);
    for(; k<= max; k++)
      powers[k]=0;
    
    dummy= powers[ii];
    powers[ii]= powers[jj];
    powers[jj]= dummy;
    return; //phew !
}

void PowProd::shift(long n) {
  Tuple<long> temp;
  long i;

  if(powers.size() == 0) // empty?
    return;
  //else...
  temp.resize( powers.size() + n );
  for(i=0; i < n; i++)
    temp[i]=0;
  for(i=0; i < powers.size(); i++)
    temp[n + i]= powers[i];
  powers= temp;
}



Alphabet* PowProd::get_current_alphabet() {
  return PowProd::alphabet;
}

Alphabet* PowProd::get_default_alphabet() {
  return &PowProd::default_alphabet;
}


void PowProd::set_alphabet(Alphabet *thealphabet){
  PowProd::alphabet= thealphabet;
}




PowProd lcm(const PowProd& a, const PowProd& b){
  long i,n,m;
  PowProd ans;

  n= a.powers.size();
  m= b.powers.size();
  if(n>m){
    ans.powers= a.powers;
  }
  else{
    m= n;
    ans.powers= b.powers;
  }// now m=min
  //  ans.powers.resize(n); ????? what was this ????
  for(i=0; i<m; i++)
    ans.powers[i]= a.powers[i] > b.powers[i] ? a.powers[i] : b.powers[i];
  return ans;
}

bool powprod_compare(const PowProd& one, const PowProd& two, int theorder, long priority){
  switch(theorder){
  case LEX:
    return one.lex_lessthan(two);
  case LEXSV:
    return one.lexsv_lessthan(two, priority);
  case DEGLEX:
    return one.deglex_lessthan(two);
  case DEGLEXSV:
    return one.deglexsv_lessthan(two, priority);
  case DEGREVLEX:
    return one.degrevlex_lessthan(two);
  case DEGREVLEXSV:
    return one.degrevlexsv_lessthan(two, priority);
  }
  return false; // in case order has a funny value


}




ostream& operator<<(ostream& os, const PowProd& P){
  long i;

  if(P.nvars() == 0){
    os << "I";
    return os;
  }

  for(i=0; i < P.nvars(); i++){
    if(P.powers[i] != 0){
      os << PowProd::alphabet->nameof(i+1);
      if(P.powers[i] > 1) 
        os << "^" << P.powers[i];
    }
  }
  return os;
}


ofstream& operator<<(ofstream& file, PowProd& x){

  file << x.powers;
  return file;
}

void operator>>(ifstream& file, PowProd& x){

  file >> x.powers;
}

---