/* Functions to manipulate integer-sized mod-2 polynomials  */

typedef unsigned int IntMod2Poly; /* accommodates degree 31 or less */


/* *****************************************************
 * Convert an Integer-Size Mod-2 Polynomial to an Integer
 *
 * mod2polyint( p[], dp ):
 * [0] Initialize ip = 0
 * [1] For d=dp down to 0, let ip = 2*ip + p[d]
 * [2] Return ip
 *
 */
IntMod2Poly
mod2polyint( int p[], int dp )
{
  IntMod2Poly ip = 0;
  int d;
  for(d=dp; d>=0; d--)
    ip = 2*ip+p[d];
  return ip;
}

/* *****************************************************
 * Integer-Size Mod-2 Polynomial Sums
 *
 * intmod2polysum( p, q ):
 * [0] sum = p^q	  
 * [1] Return sum	       
 *
 */
IntMod2Poly
intmod2polysum(IntMod2Poly p, IntMod2Poly q )
{
  return p^q;
}

/* *****************************************************
 * Integer-Size Mod-2 Polynomial Products
 *
 * intmod2polyproduct( p, q ):
 * [0] Initialize prod = 0
 * [1] While p>0, do [2] to [4]
 * [2]    If p is odd, then replace prod ^= q
 * [3]    Replace q <<= 1 by its bitshift left
 * [4]    Replace p >>= 1 by its bitshift right
 * [5] Return prod
 *
 */
IntMod2Poly
intmod2polyproduct(IntMod2Poly p, IntMod2Poly q )
{
  IntMod2Poly prod=0;
  while(p>0){
    if(p&1) prod ^= q;
    q<<=1;
    p>>=1;
  }
  return prod;
}

/* *****************************************************
 * Integer-Size Mod-2 Polynomial Degree
 *
 * intmod2polydegree( p ):
 * [0] Initialize degree = -1
 * [1] While p>0, do [2] to [3]
 * [2]    Bitshift right p >>= 1
 * [3]    Increment degree += 1
 * [4] Return degree
 *
 */
int
intmod2polydegree( IntMod2Poly p )
{
  int degree=-1;
  while(p>0) {
    p >>= 1;
    ++degree;
  }
  return degree;
}

/* *****************************************************
 * Integer-Size Mod-2 Polynomial Quotient and Remainder
 *
 * intmod2polydivision( p1, p2 ):
 * [0] Initialize q=0 and  r=p1
 * [1] Let sh = intmod2polydegree(r)-intmod2polydegree(p2)
 * [2] If sh >= 0, then do [3] to [4]
 * [3]    Replace r ^= (p2<<sh) and replace q ^= (1<<sh)
 * [4]    Go to [1]
 * [5] Return q, r
 *
 * Return q,r such that p1 = p2*q+r as mod-2 polynomials.
 *
 * In this implementation, the function's return value is the remainder r,
 * while quotient q is returned by reference.
 */
IntMod2Poly
intmod2polydivision( IntMod2Poly *q, IntMod2Poly p1, IntMod2Poly p2 )
{
  int sh;
  *q = 0;
  while((sh = intmod2polydegree(p1)-intmod2polydegree(p2)) >= 0) {
    p1 ^= (p2<<sh);
    *q ^= (1<<sh);
  }
  return p1;
}

/* *****************************************************
 * 
 * Euclid's Algorithm for Integer-Size Mod-2 Polynomials
 * 
 * intmod2polygcd( x, y ):
 * [0] Let z = x
 * [1] Call intmod2polydivision(y,x), replace x = remainder
 * [2] Let y = z
 * [3] If x>0, then go to [0]
 * [4] Return y
 *
 *
 */
IntMod2Poly
intmod2polygcd( IntMod2Poly x, IntMod2Poly y )
{
  IntMod2Poly z, q;
  while(x > 0) {
    z = x;
    x = intmod2polydivision(&q, y, x);
    y = z;
  }

  return y;
}