#include <NTL/ZZX.h>#include <NTL/new.h>NTL_START_IMPLconst ZZX& ZZX::zero(){   static ZZX z;   return z;}void conv(ZZ_pX& x, const ZZX& a){   conv(x.rep, a.rep);   x.normalize();}void conv(ZZX& x, const ZZ_pX& a){   conv(x.rep, a.rep);   x.normalize();}istream& operator>>(istream& s, ZZX& x){   s >> x.rep;   x.normalize();   return s;}ostream& operator<<(ostream& s, const ZZX& a){   return s << a.rep;}void ZZX::normalize(){   long n;   const ZZ* p;   n = rep.length();   if (n == 0) return;   p = rep.elts() + (n-1);   while (n > 0 && IsZero(*p)) {      p--;       n--;   }   rep.SetLength(n);}long IsZero(const ZZX& a){   return a.rep.length() == 0;}long IsOne(const ZZX& a){    return a.rep.length() == 1 && IsOne(a.rep[0]);}long operator==(const ZZX& a, const ZZX& b){   long i, n;   const ZZ *ap, *bp;   n = a.rep.length();   if (n != b.rep.length()) return 0;   ap = a.rep.elts();   bp = b.rep.elts();   for (i = 0; i < n; i++)      if (ap[i] != bp[i]) return 0;   return 1;}long operator==(const ZZX& a, long b){   if (b == 0)      return IsZero(a);   if (deg(a) != 0)      return 0;   return a.rep[0] == b;}long operator==(const ZZX& a, const ZZ& b){   if (IsZero(b))      return IsZero(a);   if (deg(a) != 0)      return 0;   return a.rep[0] == b;}void GetCoeff(ZZ& x, const ZZX& a, long i){   if (i < 0 || i > deg(a))      clear(x);   else      x = a.rep[i];}void SetCoeff(ZZX& x, long i, const ZZ& a){   long j, m;   if (i < 0)       Error("SetCoeff: negative index");   if (i >= (1L << (NTL_BITS_PER_LONG-4)))      Error("overflow in SetCoeff");   m = deg(x);   if (i > m) {      long pos = x.rep.position(a);      x.rep.SetLength(i+1);      if (pos != -1)         x.rep[i] = x.rep.RawGet(pos);      else         x.rep[i] = a;      for (j = m+1; j < i; j++)         clear(x.rep[j]);   }   else      x.rep[i] = a;   x.normalize();}void SetCoeff(ZZX& x, long i){   long j, m;   if (i < 0)       Error("coefficient index out of range");   if (i >= (1L << (NTL_BITS_PER_LONG-4)))      Error("overflow in SetCoeff");   m = deg(x);   if (i > m) {      x.rep.SetLength(i+1);      for (j = m+1; j < i; j++)         clear(x.rep[j]);   }   set(x.rep[i]);   x.normalize();}void SetX(ZZX& x){   clear(x);   SetCoeff(x, 1);}long IsX(const ZZX& a){   return deg(a) == 1 && IsOne(LeadCoeff(a)) && IsZero(ConstTerm(a));}            const ZZ& coeff(const ZZX& a, long i){   if (i < 0 || i > deg(a))      return ZZ::zero();   else      return a.rep[i];}const ZZ& LeadCoeff(const ZZX& a){   if (IsZero(a))      return ZZ::zero();   else      return a.rep[deg(a)];}const ZZ& ConstTerm(const ZZX& a){   if (IsZero(a))      return ZZ::zero();   else      return a.rep[0];}void conv(ZZX& x, const ZZ& a){   if (IsZero(a))      x.rep.SetLength(0);   else {      x.rep.SetLength(1);      x.rep[0] = a;   }}void conv(ZZX& x, long a){   if (a == 0)       x.rep.SetLength(0);   else {      x.rep.SetLength(1);      conv(x.rep[0], a);   }}void conv(ZZX& x, const vec_ZZ& a){   x.rep = a;   x.normalize();}void add(ZZX& x, const ZZX& a, const ZZX& b){   long da = deg(a);   long db = deg(b);   long minab = min(da, db);   long maxab = max(da, db);   x.rep.SetLength(maxab+1);   long i;   const ZZ *ap, *bp;    ZZ* xp;   for (i = minab+1, ap = a.rep.elts(), bp = b.rep.elts(), xp = x.rep.elts();        i; i--, ap++, bp++, xp++)      add(*xp, (*ap), (*bp));   if (da > minab && &x != &a)      for (i = da-minab; i; i--, xp++, ap++)         *xp = *ap;   else if (db > minab && &x != &b)      for (i = db-minab; i; i--, xp++, bp++)         *xp = *bp;   else      x.normalize();}void add(ZZX& x, const ZZX& a, const ZZ& b){   long n = a.rep.length();   if (n == 0) {      conv(x, b);   }   else if (&x == &a) {      add(x.rep[0], a.rep[0], b);      x.normalize();   }   else if (x.rep.MaxLength() == 0) {      x = a;      add(x.rep[0], a.rep[0], b);      x.normalize();   }   else {      // ugly...b could alias a coeff of x      ZZ *xp = x.rep.elts();      add(xp[0], a.rep[0], b);      x.rep.SetLength(n);      xp = x.rep.elts();      const ZZ *ap = a.rep.elts();      long i;      for (i = 1; i < n; i++)         xp[i] = ap[i];      x.normalize();   }}void add(ZZX& x, const ZZX& a, long b){   if (a.rep.length() == 0) {      conv(x, b);   }   else {      if (&x != &a) x = a;      add(x.rep[0], x.rep[0], b);      x.normalize();   }}void sub(ZZX& x, const ZZX& a, const ZZX& b){   long da = deg(a);   long db = deg(b);   long minab = min(da, db);   long maxab = max(da, db);   x.rep.SetLength(maxab+1);   long i;   const ZZ *ap, *bp;    ZZ* xp;   for (i = minab+1, ap = a.rep.elts(), bp = b.rep.elts(), xp = x.rep.elts();        i; i--, ap++, bp++, xp++)      sub(*xp, (*ap), (*bp));   if (da > minab && &x != &a)      for (i = da-minab; i; i--, xp++, ap++)         *xp = *ap;   else if (db > minab)      for (i = db-minab; i; i--, xp++, bp++)         negate(*xp, *bp);   else      x.normalize();}void sub(ZZX& x, const ZZX& a, const ZZ& b){   long n = a.rep.length();   if (n == 0) {      conv(x, b);      negate(x, x);   }   else if (&x == &a) {      sub(x.rep[0], a.rep[0], b);      x.normalize();   }   else if (x.rep.MaxLength() == 0) {      x = a;      sub(x.rep[0], a.rep[0], b);      x.normalize();   }   else {      // ugly...b could alias a coeff of x      ZZ *xp = x.rep.elts();      sub(xp[0], a.rep[0], b);      x.rep.SetLength(n);      xp = x.rep.elts();      const ZZ *ap = a.rep.elts();      long i;      for (i = 1; i < n; i++)         xp[i] = ap[i];      x.normalize();   }}void sub(ZZX& x, const ZZX& a, long b){   if (b == 0) {      x = a;      return;   }   if (a.rep.length() == 0) {      x.rep.SetLength(1);      conv(x.rep[0], b);      negate(x.rep[0], x.rep[0]);   }   else {      if (&x != &a) x = a;      sub(x.rep[0], x.rep[0], b);   }   x.normalize();}void sub(ZZX& x, long a, const ZZX& b){   negate(x, b);   add(x, x, a);}void sub(ZZX& x, const ZZ& b, const ZZX& a){   long n = a.rep.length();   if (n == 0) {      conv(x, b);   }   else if (x.rep.MaxLength() == 0) {      negate(x, a);      add(x.rep[0], a.rep[0], b);      x.normalize();   }   else {      // ugly...b could alias a coeff of x      ZZ *xp = x.rep.elts();      sub(xp[0], b, a.rep[0]);      x.rep.SetLength(n);      xp = x.rep.elts();      const ZZ *ap = a.rep.elts();      long i;      for (i = 1; i < n; i++)         negate(xp[i], ap[i]);      x.normalize();   }}void negate(ZZX& x, const ZZX& a){   long n = a.rep.length();   x.rep.SetLength(n);   const ZZ* ap = a.rep.elts();   ZZ* xp = x.rep.elts();   long i;   for (i = n; i; i--, ap++, xp++)      negate((*xp), (*ap));}long MaxBits(const ZZX& f){   long i, m;   m = 0;   for (i = 0; i <= deg(f); i++) {      m = max(m, NumBits(f.rep[i]));   }   return m;}void PlainMul(ZZX& x, const ZZX& a, const ZZX& b){   if (&a == &b) {      PlainSqr(x, a);      return;   }   long da = deg(a);   long db = deg(b);   if (da < 0 || db < 0) {      clear(x);      return;   }   long d = da+db;   const ZZ *ap, *bp;   ZZ *xp;      ZZX la, lb;   if (&x == &a) {      la = a;      ap = la.rep.elts();   }   else      ap = a.rep.elts();