polynomialt.cc.txt

压缩文件中是Error Correction Coding - Mathematical Methods and Algorithms(Wiley 2005)作者:(Todd K. Moon )的配

// polynomialT.C -- templatized polynomial arithmetic
// Todd Moon, with modifications due to Stewart Weber

// Copyright 2004 by Todd K. Moon
// Permission is granted to use this program/data
// for educational/research only

#include <assert.h>
#include "polynomialT.h"
extern "C" {
#include <string.h>  // for strlen
}

#ifndef MAX
#define MAX(a,b) (a>b?(a):(b))
#endif

#ifndef MIN
#define MIN(a,b) (a<b?(a):(b))
#endif


template <class T> int polynomialT<T>::decreasingprint = 0;
template <class T> char polynomialT<T>::defvarname[6] = "x";
template <class T> char polynomialT<T>::beforeafterprint[3][6] = {{0},{0},{0}};

// default constructor
// same as constant assignment constructor
template <class T>
polynomialT<T>:: polynomialT(const T coeff0, const char *invarname) 
   : degree(0) 
{
   coeff = new T[1];
   coeff[0] = coeff0;
   
   varname = NULL;
   if(invarname) {
	  varname = new char[strlen(invarname)+1];
	  strcpy(varname,invarname);
   }
}


// normal constructor
template <class T>
polynomialT<T>:: polynomialT(int indegree, const T *incoeff, 
							 const char *invarname) 
{
   // check indegree
   for( ; (indegree>0) && (incoeff[indegree]==0); indegree--); 

   degree = indegree;
   coeff = new T[degree+1];
   for( ; indegree>=0; indegree--) // copy over coefficients
		coeff[indegree] = incoeff[indegree];
   polytemp<T>::resettemplist();

   varname = NULL;
   if(invarname) {
	  varname = new char[strlen(invarname)+1];
	  strcpy(varname,invarname);
   }
}


// copy constructor
template <class T>
polynomialT<T>:: polynomialT(const polynomialT<T> &inpoly, 
							 const char *invarname) 
{
   int i;

	i = degree = inpoly.degree;
	coeff = new T[degree+1];
	for( ; i>=0; i--)
		coeff[i] = inpoly.coeff[i];
	polytemp<T>::resettemplist();

	varname = NULL;
	if(inpoly.varname) {
	   varname = new char[strlen(inpoly.varname)+1];
	   strcpy(varname,inpoly.varname);
	}
}

template <class T>
polynomialT<T>::~polynomialT()
{ 
   delete [] coeff;
   if(varname) {
	  delete[] varname;
   }
}

// Protected member functions

// resize the polynomial (no copy)
template <class T> void
polynomialT<T>:: resize(const int newdegree)
{
   if(degree != newdegree) {
	  degree = newdegree;
	  delete [] coeff;
	  coeff = new T[newdegree+1];
   }
}

// resize and copy over the old 
// If the polynomial is made bigger, insert zeros in the 
// newly allocated spaces.  This makes things easier for 
// the functions that call resizecopy.

template <class T> void
polynomialT<T>:: resizecopy(const int newdegree)
{
   T *newcoeff;
   int i;
   T z(0);						// build a zero
   if(degree != newdegree) {
	  newcoeff = new T[newdegree+1];
	  for(i=newdegree; i>degree; i--)
		 newcoeff[i] = z;
	  for( ; i>=0; i--)
		 newcoeff[i] = coeff[i];
	  degree = newdegree;
	  delete [] coeff;
	  coeff = newcoeff;
   }
}


// Public member functions

// set a coefficient (can be used to change degree)
template <class T> polynomialT<T>&
polynomialT<T>:: setc(int idx, const T &cval) 
{
	int i;
	
	assert(idx>=0);
	
	if(idx <= degree) {
	   coeff[idx] = cval;
	   // fix poly if leading zeros
	   for(i=degree; (i>0) && (coeff[i]==0); i--);
	   if(i != degree) resizecopy(i);
	} else if(cval != 0) {
	   // note that idx>degree here
	   resizecopy(idx);
	   coeff[idx] = cval;
	}
	return *this;
}


// set a whole array of coefficients, setting degree based on the array
template <class T> polynomialT<T>&
polynomialT<T>:: setc(const T *cvals, int indegree) 
{
   int i;
   int d;
   assert(indegree>=0);
   // determine the actual degree
   for(d = indegree; (cvals[d] == 0) && (d > 0); d--) ;
   indegree = d;
   
   // make just enough room for the array
   resize(indegree);
   
   // copy over
   for(i = 0; i <= degree; i++) {
	  coeff[i] = cvals[i];
   }
   return *this;
}

// set a coefficient (can be used to change degree)
template <class T> polynomialT<T>&
polynomialT<T>:: buildspace(int indegree)
{
	int i;
	
	assert(indegree>=0);

	if(indegree > degree) {
	   resize(indegree);
	}
	return *this;
}


// Overloaded operators

// polynomial addition
template <class T> const polynomialT<T>& 
polynomialT<T>:: operator+(const polynomialT<T> &r) const 
{
   polytemp<T> *out;
   int i;
   
   // There are 3 cases for deciding the size of
   // out (and 3 methods will be used to add), 
   // if l>r, then out is the size of 'l' 
   // (because we can be confident that 'l' doesn't
   // have any leading zeros).  If l<r, then out is
   // the size of 'r'.  If l==r, then we need to find
   // the highest nonzero coefficient of l+r to determine
	// the size.
   
   if(degree>r.degree) {
	  out = polytemp<T>::gettemppoly(degree);
	  for(i=degree; i>r.degree; i--)
		 out->coeff[i] = coeff[i];
	  for( ; i>=0; i--)
		 out->coeff[i] = coeff[i] + r.coeff[i];
   } else if(degree<r.degree) {
	  out = polytemp<T>::gettemppoly(r.degree);
	  for(i=r.degree; i>degree; i--)
		 out->coeff[i] = r.coeff[i];
	  for( ; i>=0; i--)
		 out->coeff[i] = coeff[i] + r.coeff[i];
   } else {
	  // degree==r.degree
	  // look for first nonzero coefficient
	  for(i=degree; (i>0) && (coeff[i]+r.coeff[i] == 0); i--);
	  out = polytemp<T>::gettemppoly(i);
	  for( ; i>=0; i--)
		 out->coeff[i] = coeff[i] + r.coeff[i];
   }
   
   return *(out);
} // overloaded +


// polynomial subtraction
// refer to operator+ for explanation of code
template <class T> const polynomialT<T>& 
polynomialT<T>:: operator-(const polynomialT<T> &r) const 
{
   polytemp<T> *out;
   int i;

   if(degree>r.degree) {
	  out = polytemp<T>::gettemppoly(degree);
	  for(i=degree; i>r.degree; i--)
		 out->coeff[i] = coeff[i];
	  for( ; i>=0; i--)
		 out->coeff[i] = coeff[i] - r.coeff[i];
	  
   } else if(degree<r.degree) {
	  out = polytemp<T>::gettemppoly(r.degree);
	  
	  for(i=r.degree; i>degree; i--)
		 out->coeff[i] = -(r.coeff[i]);
	  for( ; i>=0; i--)
		 out->coeff[i] = coeff[i] - r.coeff[i];
	  
   } else {
	  // degree==r.degree
	  // look for first nonzero coefficient
	  for(i=degree; (i>0) && (coeff[i]-r.coeff[i] == 0); i--);
	  
	  out = polytemp<T>::gettemppoly(i);
	  
	  for( ; i>=0; i--)
		 out->coeff[i] = coeff[i] - r.coeff[i];
   }
   
   return *(out);
} // overloaded -


// polynomial multiplication
template <class T> const polynomialT<T>& 
polynomialT<T>:: operator*(const polynomialT<T> &r) const
{
   polytemp<T> *out;
   T z(0);						// build a zero
   T endvalue = z;
   int i = degree + r.degree;
   int j, jdone;
   
   // Multiplication is just convolution.
   // Figure out degree of output by doing the
   // convolution of the endpoint first, then work
   // backwords.  The degree is the first nonzero coefficient
   
   // convolution without storage
   for( ; (i>=0) && (endvalue==0); i--) {
	  jdone = MIN(i, degree);
	  for(j=MAX(0, i-r.degree); j<=jdone; j++)
		 endvalue += (coeff[j]*r.coeff[i-j]);
   }
   out = polytemp<T>::gettemppoly(i+1);
   out->coeff[i+1] = endvalue;
   
   // finish convolution with storage
   for( ; i>=0; i--) {
	  out->coeff[i] = z;
	  
	  jdone = MIN(i, degree);
	  for(j=MAX(0, i-r.degree); j<=jdone; j++)
		 out->coeff[i] += (coeff[j]*r.coeff[i-j]);
   }
   return *out;
} // overloaded *

// polynomial division

// Since the operation of division produces both a quotient
// and a remainder, and sometimes both are wanted, both the quotient
// and remainder are saved in temporary variables.  They can then be
// retrieved using the functions getlastquotient and getlastremainder.

template <class T> const void 
polynomialT<T>:: divide(const polynomialT<T> &b) const
// compute this/b, with quotient and remainder
{
   int i,j,k;
   T z(0);						// build a zero
   if(degree < b.degree) {  // improper division
	  polytemp<T> *quot = polytemp<T>::gettempquot(0);
	  polytemp<T> *rem = polytemp<T>::gettemprem(degree);
	  quot->coeff[0] = z;
	  for(i = 0; i <= degree; i++) {
		 rem->coeff[i] = coeff[i];
	  }
   }
   else if(b.degree==0) {
	  T gninv = T(1)/b.coeff[0];
	  polytemp<T> *quot = polytemp<T>::gettempquot(degree);
	  polytemp<T> *rem = polytemp<T>::gettemprem(0);
	  for(i = 0; i <= degree; i++) {
		 (*quot)[i] = coeff[i] * gninv;
	  }
	  (*rem)[0] = z;
   }
   else {
	  int quotdeg = degree - b.degree;
	  int nr1 = b.degree-1;
	  T q, gninv;
	  int i,j,k;
	  polytemp<T> *quot = polytemp<T>::gettempquot(quotdeg);
	  polytemp<T> *rem = polytemp<T>::gettemprem(nr1);
	  // set up initial
	  gninv = T(1)/b.coeff[b.degree];
	  for(j = degree, i=nr1; i>=0; j--, i--) {
		 rem->coeff[i] = coeff[j];
	  }
	  for(k=0; k <= quotdeg; j--,k++) {  // work the rest of them
		 q = rem->coeff[nr1]*gninv;
		 quot->coeff[quotdeg-k] = q;
		 for(i = nr1; i>0; i--) {
			rem->coeff[i] = rem->coeff[i-1] - q*b.coeff[i];
		 }
		 rem->coeff[0] = coeff[j]-q*b.coeff[0];
	  }
	  // determine the actual degree of the remainder
	  rem->degree = 0;
	  for(i = nr1; i >=0; i--) {
		 if(rem->coeff[i] != 0) {
			rem->degree = i;
			break;
		 }
	  }
   }
}


template <class T> const polynomialT<T>& 
polynomialT<T>:: operator/(const polynomialT<T> &r) const
{
   divide(r);
   return *polytemp<T>::tempquotient;
}


template <class T> const polynomialT<T>& 
polynomialT<T>:: operator%(const polynomialT<T> &r) const
{  
   divide(r);
   return *polytemp<T>::temprem;
} // overloaded %


// polynomial div (/=)
template <class T> const polynomialT<T>& 
polynomialT<T>:: operator/=(const polynomialT<T> &r) 
{
   (*this)=(*this)/r;
   polytemp<T>::resettemplist();
   return *this;
} // overloaded /=

// polynomial div (%=)
template <class T> const polynomialT<T>& 
polynomialT<T>:: operator%=(const polynomialT<T> &r) 
{
   (*this)=(*this)%r;
   polytemp<T>::resettemplist();
   return *this;
} // overloaded %=

template <class T> const polynomialT<T>&
polynomialT<T>:: getlastquotient(void) {
   return *polytemp<T>::tempquotient;
}

template <class T> const polynomialT<T>&
polynomialT<T>:: getlastremainder(void) {
   return *polytemp<T>::temprem;
}


// unary -
template <class T> const polynomialT<T>& 
polynomialT<T>:: operator-(void) const {
   polytemp<T>* out;
   int i;

   out = polytemp<T>::gettemppoly(degree);
   for(i=degree ; i>=0; i--)
	  out->coeff[i] = -(coeff[i]);
   return *(out);
} // overloaded unary -


// assignment
template <class T> const polynomialT<T>&
polynomialT<T>:: operator=(const polynomialT<T> &r) 
{
   int i = r.degree;

   if(this != &r) {
	  resize(r.degree);
	  
	  for( ; i>=0; i--)
		 coeff[i] = r.coeff[i];
	  
   }
   polytemp<T>::resettemplist();
   varname = NULL;
   if(r.varname) {
	  varname = new char[strlen(r.varname)+1];
	  strcpy(varname,r.varname);
   }
   return *this;
} // overloaded =


// polynomial increment (+=)
template <class T> const polynomialT<T>& 
polynomialT<T>:: operator+=(const polynomialT<T> &r) 
{
   int i = r.degree;

	if(degree==r.degree) {
	   // look for first nonzero coefficient
	   for( ; (i>0) && (coeff[i]+r.coeff[i]==0); i--);
	   resizecopy(i);
	   
	} else if(degree<r.degree) resizecopy(r.degree);
	
	for( ; i>=0; i--) coeff[i] += r.coeff[i];
	
	polytemp<T>::resettemplist();
	return *this;
} // overloaded +=


// polynomial decrement (-=)
template <class T> const polynomialT<T>& 
polynomialT<T>:: operator-=(const polynomialT<T> &r) 
{
   int i = r.degree;

   if(degree==r.degree) {
	  // look for first nonzero coefficient
	  for( ; (i>0) && (coeff[i]-r.coeff[i]==0); i--);
	  resizecopy(i);
	  
   } else if(degree<r.degree) resizecopy(r.degree);
   
   for( ; i>=0; i--) coeff[i] -= r.coeff[i];
   
   polytemp<T>::resettemplist();
   return *this;
} // overloaded -=


// polynomial scale (*=)
template <class T> const polynomialT<T>& 
polynomialT<T>:: operator*=(const polynomialT<T> &r) 
{
   (*this)=(*this)*r;
   polytemp<T>::resettemplist();
   return *this;
} // overloaded *=


// divide by x -- reduce the degree (non-cyclic shift)
template <class T> const polynomialT<T>& 
polynomialT<T>:: operator>>(const int nshift) const 
{
   polytemp<T> *out;
   T save;
   int i;
   T z(0);						// build a zero

   if(nshift == 0) {
	  return *this;
   }
   else if(nshift < 0) {		// shifting the other way
	  return (*this) << -nshift;
   }
   else {
	  if(nshift > degree) {
		 out = polytemp<T>::gettemppoly(0);
		 out->coeff[0] = z;
	  }
	  else {
		 out = polytemp<T>::gettemppoly(degree-nshift);
		 for(i=degree; i>nshift-1; i--) {
			out->coeff[i-nshift] = coeff[i];
		 }
	  }
   }
   return *out;
} // overloaded >>

// multiply by x -- increase the degree (non-cyclic shift)
template <class T> const polynomialT<T>& 
polynomialT<T>:: operator<<(const int nshift) const 
{
   polytemp<T> *out;
   int i;
   T save;
   T z(0);						// build a zero