#ifndef NTL_ZZX__H #define NTL_ZZX__H #include #include #include NTL_OPEN_NNS class ZZX { public: vec_ZZ rep; /*************************************************************** Constructors, Destructors, and Assignment ****************************************************************/ ZZX() // initial value 0 { } ZZX(INIT_SIZE_TYPE, long n) // initial value 0, but space is pre-allocated for n coefficients { rep.SetMaxLength(n); } ZZX(const ZZX& a) : rep(a.rep) { } // initial value is a ZZX& operator=(const ZZX& a) { rep = a.rep; return *this; } ~ZZX() { } void normalize(); // strip leading zeros void SetMaxLength(long n) // pre-allocate space for n coefficients. // Value is unchanged { rep.SetMaxLength(n); } void kill() // free space held by this polynomial. Value becomes 0. { rep.kill(); } static const ZZX& zero(); inline ZZX(long i, const ZZ& c); inline ZZX(long i, long c); inline ZZX& operator=(long a); inline ZZX& operator=(const ZZ& a); ZZX(ZZX& x, INIT_TRANS_TYPE) : rep(x.rep, INIT_TRANS) { } }; /******************************************************************** input and output I/O format: [a_0 a_1 ... a_n], represents the polynomial a_0 + a_1*X + ... + a_n*X^n. On output, all coefficients will be integers between 0 and p-1, amd a_n not zero (the zero polynomial is [ ]). Leading zeroes are stripped. *********************************************************************/ NTL_SNS istream& operator>>(NTL_SNS istream& s, ZZX& x); NTL_SNS ostream& operator<<(NTL_SNS ostream& s, const ZZX& a); /********************************************************** Some utility routines ***********************************************************/ inline long deg(const ZZX& a) { return a.rep.length() - 1; } // degree of a polynomial. // note that the zero polynomial has degree -1. const ZZ& coeff(const ZZX& a, long i); // zero if i not in range void GetCoeff(ZZ& x, const ZZX& a, long i); // x = a[i], or zero if i not in range const ZZ& LeadCoeff(const ZZX& a); // zero if a == 0 const ZZ& ConstTerm(const ZZX& a); // zero if a == 0 void SetCoeff(ZZX& x, long i, const ZZ& a); // x[i] = a, error is raised if i < 0 void SetCoeff(ZZX& x, long i, long a); inline ZZX::ZZX(long i, const ZZ& a) { SetCoeff(*this, i, a); } inline ZZX::ZZX(long i, long a) { SetCoeff(*this, i, a); } void SetCoeff(ZZX& x, long i); // x[i] = 1, error is raised if i < 0 void SetX(ZZX& x); // x is set to the monomial X long IsX(const ZZX& a); // test if x = X inline void clear(ZZX& x) // x = 0 { x.rep.SetLength(0); } inline void set(ZZX& x) // x = 1 { x.rep.SetLength(1); set(x.rep[0]); } inline void swap(ZZX& x, ZZX& y) // swap x & y (only pointers are swapped) { swap(x.rep, y.rep); } void trunc(ZZX& x, const ZZX& a, long m); // x = a % X^m inline ZZX trunc(const ZZX& a, long m) { ZZX x; trunc(x, a, m); NTL_OPT_RETURN(ZZX, x); } void RightShift(ZZX& x, const ZZX& a, long n); // x = a/X^n inline ZZX RightShift(const ZZX& a, long n) { ZZX x; RightShift(x, a, n); NTL_OPT_RETURN(ZZX, x); } void LeftShift(ZZX& x, const ZZX& a, long n); // x = a*X^n inline ZZX LeftShift(const ZZX& a, long n) { ZZX x; LeftShift(x, a, n); NTL_OPT_RETURN(ZZX, x); } #ifndef NTL_TRANSITION inline ZZX operator>>(const ZZX& a, long n) { ZZX x; RightShift(x, a, n); NTL_OPT_RETURN(ZZX, x); } inline ZZX operator<<(const ZZX& a, long n) { ZZX x; LeftShift(x, a, n); NTL_OPT_RETURN(ZZX, x); } inline ZZX& operator<<=(ZZX& x, long n) { LeftShift(x, x, n); return x; } inline ZZX& operator>>=(ZZX& x, long n) { RightShift(x, x, n); return x; } #endif void diff(ZZX& x, const ZZX& a); // x = derivative of a inline ZZX diff(const ZZX& a) { ZZX x; diff(x, a); NTL_OPT_RETURN(ZZX, x); } void InvTrunc(ZZX& x, const ZZX& a, long m); // computes x = a^{-1} % X^m // constant term must be non-zero inline ZZX InvTrunc(const ZZX& a, long m) { ZZX x; InvTrunc(x, a, m); NTL_OPT_RETURN(ZZX, x); } void MulTrunc(ZZX& x, const ZZX& a, const ZZX& b, long n); // x = a * b % X^n inline ZZX MulTrunc(const ZZX& a, const ZZX& b, long n) { ZZX x; MulTrunc(x, a, b, n); NTL_OPT_RETURN(ZZX, x); } void SqrTrunc(ZZX& x, const ZZX& a, long n); // x = a^2 % X^n inline ZZX SqrTrunc(const ZZX& a, long n) { ZZX x; SqrTrunc(x, a, n); NTL_OPT_RETURN(ZZX, x); } void reverse(ZZX& c, const ZZX& a, long hi); inline ZZX reverse(const ZZX& a, long hi) { ZZX x; reverse(x, a, hi); NTL_OPT_RETURN(ZZX, x); } inline void reverse(ZZX& c, const ZZX& a) { reverse(c, a, deg(a)); } inline ZZX reverse(const ZZX& a) { ZZX x; reverse(x, a); NTL_OPT_RETURN(ZZX, x); } inline void VectorCopy(vec_ZZ& x, const ZZX& a, long n) { VectorCopy(x, a.rep, n); } inline vec_ZZ VectorCopy(const ZZX& a, long n) { return VectorCopy(a.rep, n); } /******************************************************************* conversion routines ********************************************************************/ void conv(ZZX& x, long a); inline ZZX to_ZZX(long a) { ZZX x; conv(x, a); NTL_OPT_RETURN(ZZX, x); } inline ZZX& ZZX::operator=(long a) { conv(*this, a); return *this; } void conv(ZZX& x, const ZZ& a); inline ZZX to_ZZX(const ZZ& a) { ZZX x; conv(x, a); NTL_OPT_RETURN(ZZX, x); } inline ZZX& ZZX::operator=(const ZZ& a) { conv(*this, a); return *this; } void conv(ZZX& x, const vec_ZZ& a); inline ZZX to_ZZX(const vec_ZZ& a) { ZZX x; conv(x, a); NTL_OPT_RETURN(ZZX, x); } void conv(zz_pX& x, const ZZX& a); inline zz_pX to_zz_pX(const ZZX& a) { zz_pX x; conv(x, a); NTL_OPT_RETURN(zz_pX, x); } void conv(ZZ_pX& x, const ZZX& a); inline ZZ_pX to_ZZ_pX(const ZZX& a) { ZZ_pX x; conv(x, a); NTL_OPT_RETURN(ZZ_pX, x); } void conv(ZZX& x, const ZZ_pX& a); inline ZZX to_ZZX(const ZZ_pX& a) { ZZX x; conv(x, a); NTL_OPT_RETURN(ZZX, x); } void conv(ZZX& x, const zz_pX& a); inline ZZX to_ZZX(const zz_pX& a) { ZZX x; conv(x, a); NTL_OPT_RETURN(ZZX, x); } /************************************************************* Comparison **************************************************************/ long IsZero(const ZZX& a); long IsOne(const ZZX& a); long operator==(const ZZX& a, const ZZX& b); inline long operator!=(const ZZX& a, const ZZX& b) { return !(a == b); } long operator==(const ZZX& a, const ZZ& b); long operator==(const ZZX& a, long b); inline long operator==(const ZZ& a, const ZZX& b) { return b == a; } inline long operator==(long a, const ZZX& b) { return b == a; } inline long operator!=(const ZZX& a, const ZZ& b) { return !(a == b); } inline long operator!=(const ZZX& a, long b) { return !(a == b); } inline long operator!=(const ZZ& a, const ZZX& b) { return !(a == b); } inline long operator!=(long a, const ZZX& b) { return !(a == b); } /*************************************************************** Addition ****************************************************************/ void add(ZZX& x, const ZZX& a, const ZZX& b); // x = a + b void sub(ZZX& x, const ZZX& a, const ZZX& b); // x = a - b void negate(ZZX& x, const ZZX& a); // x = -a // scalar versions void add(ZZX & x, const ZZX& a, const ZZ& b); // x = a + b void add(ZZX& x, const ZZX& a, long b); inline void add(ZZX& x, const ZZ& a, const ZZX& b) { add(x, b, a); } inline void add(ZZX& x, long a, const ZZX& b) { add(x, b, a); } void sub(ZZX & x, const ZZX& a, const ZZ& b); // x = a - b void sub(ZZX& x, const ZZX& a, long b); void sub(ZZX& x, const ZZ& a, const ZZX& b); void sub(ZZX& x, long a, const ZZX& b); inline ZZX operator+(const ZZX& a, const ZZX& b) { ZZX x; add(x, a, b); NTL_OPT_RETURN(ZZX, x); } inline ZZX operator+(const ZZX& a, const ZZ& b) { ZZX x; add(x, a, b); NTL_OPT_RETURN(ZZX, x); } inline ZZX operator+(const ZZX& a, long b) { ZZX x; add(x, a, b); NTL_OPT_RETURN(ZZX, x); } inline ZZX operator+(const ZZ& a, const ZZX& b) { ZZX x; add(x, a, b); NTL_OPT_RETURN(ZZX, x); } inline ZZX operator+(long a, const ZZX& b) { ZZX x; add(x, a, b); NTL_OPT_RETURN(ZZX, x); } inline ZZX operator-(const ZZX& a, const ZZX& b) { ZZX x; sub(x, a, b); NTL_OPT_RETURN(ZZX, x); } inline ZZX operator-(const ZZX& a, const ZZ& b) { ZZX x; sub(x, a, b); NTL_OPT_RETURN(ZZX, x); } inline ZZX operator-(const ZZX& a, long b) { ZZX x; sub(x, a, b); NTL_OPT_RETURN(ZZX, x); } inline ZZX operator-(const ZZ& a, const ZZX& b) { ZZX x; sub(x, a, b); NTL_OPT_RETURN(ZZX, x); } inline ZZX operator-(long a, const ZZX& b) { ZZX x; sub(x, a, b); NTL_OPT_RETURN(ZZX, x); } inline ZZX& operator+=(ZZX& x, const ZZX& b) { add(x, x, b); return x; } inline ZZX& operator+=(ZZX& x, const ZZ& b) { add(x, x, b); return x; } inline ZZX& operator+=(ZZX& x, long b) { add(x, x, b); return x; } inline ZZX& operator-=(ZZX& x, const ZZX& b) { sub(x, x, b); return x; } inline ZZX& operator-=(ZZX& x, const ZZ& b) { sub(x, x, b); return x; } inline ZZX& operator-=(ZZX& x, long b) { sub(x, x, b); return x; } inline ZZX operator-(const ZZX& a) { ZZX x; negate(x, a); NTL_OPT_RETURN(ZZX, x); } inline ZZX& operator++(ZZX& x) { add(x, x, 1); return x; } inline void operator++(ZZX& x, int) { add(x, x, 1); } inline ZZX& operator--(ZZX& x) { sub(x, x, 1); return x; } inline void operator--(ZZX& x, int) { sub(x, x, 1); } /***************************************************************** Multiplication ******************************************************************/ void mul(ZZX& x, const ZZX& a, const ZZX& b); // x = a * b void sqr(ZZX& x, const ZZX& a); inline ZZX sqr(const ZZX& a) { ZZX x; sqr(x, a); NTL_OPT_RETURN(ZZX, x); } // x = a^2 void PlainMul(ZZX& x, const ZZX& a, const ZZX& b); void PlainSqr(ZZX& x, const ZZX& a); void KarMul(ZZX& x, const ZZX& a, const ZZX& b); void KarSqr(ZZX& x, const ZZX& a); void HomMul(ZZX& x, const ZZX& a, const ZZX& b); void HomSqr(ZZX& x, const ZZX& a); void SSMul(ZZX& x, const ZZX& a, const ZZX& b); void SSSqr(ZZX& x, const ZZX& a); double SSRatio(long na, long maxa, long nb, long maxb); void mul(ZZX & x, const ZZX& a, const ZZ& b); void mul(ZZX& x, const ZZX& a, long b); inline void mul(ZZX& x, const ZZ& a, const ZZX& b) { mul(x, b, a); } inline void mul(ZZX& x, long a, const ZZX& b) { mul(x, b, a); } inline ZZX operator*(const ZZX& a, const ZZX& b) { ZZX x; mul(x, a, b); NTL_OPT_RETURN(ZZX, x); } inline ZZX operator*(const ZZX& a, const ZZ& b) { ZZX x; mul(x, a, b); NTL_OPT_RETURN(ZZX, x); } inline ZZX operator*(const ZZX& a, long b) { ZZX x; mul(x, a, b); NTL_OPT_RETURN(ZZX, x); } inline ZZX operator*(const ZZ& a, const ZZX& b) { ZZX x; mul(x, a, b); NTL_OPT_RETURN(ZZX, x); } inline ZZX operator*(long a, const ZZX& b) { ZZX x; mul(x, a, b); NTL_OPT_RETURN(ZZX, x); } inline ZZX& operator*=(ZZX& x, const ZZX& b) { mul(x, x, b); return x; } inline ZZX& operator*=(ZZX& x, const ZZ& b) { mul(x, x, b); return x; } inline ZZX& operator*=(ZZX& x, long b) { mul(x, x, b); return x; } /************************************************************* Division **************************************************************/ // "plain" versions void PlainPseudoDivRem(ZZX& q, ZZX& r, const ZZX& a, const ZZX& b); void PlainPseudoDiv(ZZX& q, const ZZX& a, const ZZX& b); void PlainPseudoRem(ZZX& r, const ZZX& a, const ZZX& b); // "homomorphic imaging" versions void HomPseudoDivRem(ZZX& q, ZZX& r, const ZZX& a, const ZZX& b); void HomPseudoDiv(ZZX& q, const ZZX& a, const ZZX& b); void HomPseudoRem(ZZX& r, const ZZX& a, const ZZX& b); inline void PseudoDivRem(ZZX& q, ZZX& r, const ZZX& a, const ZZX& b) // performs pseudo-division: computes q and r // with deg(r) < deg(b), and LeadCoeff(b)^(deg(a)-deg(b)+1) a = b q + r. // current implementation always defaults to "plain" { PlainPseudoDivRem(q, r, a, b); } inline void PseudoDiv(ZZX& q, const ZZX& a, const ZZX& b) { PlainPseudoDiv(q, a, b); } inline void PseudoRem(ZZX& r, const ZZX& a, const ZZX& b) { PlainPseudoRem(r, a, b); } inline ZZX PseudoDiv(const ZZX& a, const ZZX& b) { ZZX x; PseudoDiv(x, a, b); NTL_OPT_RETURN(ZZX, x); } inline ZZX PseudoRem(const ZZX& a, const ZZX& b) { ZZX x; PseudoRem(x, a, b); NTL_OPT_RETURN(ZZX, x); } #ifndef NTL_TRANSITION void DivRem(ZZX& q, ZZX& r, const ZZX& a, const ZZX& b); void div(ZZX& q, const ZZX& a, const ZZX& b); void div(ZZX& q, const ZZX& a, const ZZ& b); void div(ZZX& q, const ZZX& a, long b); void rem(ZZX& r, const ZZX& a, const ZZX& b); inline ZZX operator/(const ZZX& a, const ZZX& b) { ZZX x; div(x, a, b); NTL_OPT_RETURN(ZZX, x); } inline ZZX operator/(const ZZX& a, const ZZ& b) { ZZX x; div(x, a, b); NTL_OPT_RETURN(ZZX, x); } inline ZZX operator/(const ZZX& a, long b) { ZZX x; div(x, a, b); NTL_OPT_RETURN(ZZX, x); } inline ZZX& operator/=(ZZX& x, const ZZ& b) { div(x, x, b); return x; } inline ZZX& operator/=(ZZX& x, long b) { div(x, x, b); return x; } inline ZZX& operator/=(ZZX& x, const ZZX& b) { div(x, x, b); return x; } inline ZZX operator%(const ZZX& a, const ZZX& b) { ZZX x; rem(x, a, b); NTL_OPT_RETURN(ZZX, x); } inline ZZX& operator%=(ZZX& x, const ZZX& b) { rem(x, x, b); return x; } #endif // Modular arithemtic---f must be monic, and other args // must have degree less than that of f void MulMod(ZZX& x, const ZZX& a, const ZZX& b, const ZZX& f); inline ZZX MulMod(const ZZX& a, const ZZX& b, const ZZX& f) { ZZX x; MulMod(x, a, b, f); NTL_OPT_RETURN(ZZX, x); } void SqrMod(ZZX& x, const ZZX& a, const ZZX& f); inline ZZX SqrMod(const ZZX& a, const ZZX& f) { ZZX x; SqrMod(x, a, f); NTL_OPT_RETURN(ZZX, x); } void MulByXMod(ZZX& x, const ZZX& a, const ZZX& f); inline ZZX MulByXMod(const ZZX& a, const ZZX& f) { ZZX x; MulByXMod(x, a, f); NTL_OPT_RETURN(ZZX, x); } // these always use "plain" division long PlainDivide(ZZX& q, const ZZX& a, const ZZX& b); long PlainDivide(const ZZX& a, const ZZX& b); // these always use "homomorphic imaging" long HomDivide(ZZX& q, const ZZX& a, const ZZX& b); long HomDivide(const ZZX& a, const ZZX& b); long divide(ZZX& q, const ZZX& a, const ZZX& b); // if b | a, sets q = a/b and returns 1; otherwise returns 0 long divide(const ZZX& a, const ZZX& b); long divide(ZZX& q, const ZZX& a, const ZZ& b); // if b | a, sets q = a/b and returns 1; otherwise returns 0 long divide(const ZZX& a, const ZZ& b); // if b | a, returns 1; otherwise returns 0 //single-precision versions long divide(ZZX& q, const ZZX& a, long b); long divide(const ZZX& a, long b); void content(ZZ& d, const ZZX& f); // d = content of f, sign(d) == sign(LeadCoeff(f)) inline ZZ content(const ZZX& f) { ZZ x; content(x, f); NTL_OPT_RETURN(ZZ, x); } void PrimitivePart(ZZX& pp, const ZZX& f); // pp = primitive part of f, LeadCoeff(pp) >= 0 inline ZZX PrimitivePart(const ZZX& f) { ZZX x; PrimitivePart(x, f); NTL_OPT_RETURN(ZZX, x); } void GCD(ZZX& d, const ZZX& a, const ZZX& b); // d = gcd(a, b), LeadCoeff(d) >= 0 inline ZZX GCD(const ZZX& a, const ZZX& b) { ZZX x; GCD(x, a, b); NTL_OPT_RETURN(ZZX, x); } long MaxBits(const ZZX& f); // returns max NumBits of coefficients of f long CharPolyBound(const ZZX& a, const ZZX& f); /*************************************************************** traces, norms, resultants ****************************************************************/ void TraceVec(vec_ZZ& S, const ZZX& f); // S[i] = Trace(X^i mod f), for i = 0..deg(f)-1. // f must be a monic polynomial. inline vec_ZZ TraceVec(const ZZX& f) { vec_ZZ x; TraceVec(x, f); NTL_OPT_RETURN(vec_ZZ, x); } void TraceMod(ZZ& res, const ZZX& a, const ZZX& f); inline ZZ TraceMod(const ZZX& a, const ZZX& f) { ZZ x; TraceMod(x, a, f); NTL_OPT_RETURN(ZZ, x); } // res = trace of (a mod f) // f must be monic void resultant(ZZ& res, const ZZX& a, const ZZX& b, long deterministic=0); inline ZZ resultant(const ZZX& a, const ZZX& b, long deterministic=0) { ZZ x; resultant(x, a, b, deterministic); NTL_OPT_RETURN(ZZ, x); } // res = resultant of a and b // if !deterministic, then it may use a randomized strategy // that errs with probability no more than 2^{-80}. void NormMod(ZZ& res, const ZZX& a, const ZZX& f, long deterministic=0); inline ZZ NormMod(const ZZX& a, const ZZX& f, long deterministic=0) { ZZ x; NormMod(x, a, f, deterministic); NTL_OPT_RETURN(ZZ, x); } // res = norm of (a mod f) // f must be monic // if !deterministic, then it may use a randomized strategy // that errs with probability no more than 2^{-80}. void discriminant(ZZ& d, const ZZX& a, long deterministic=0); inline ZZ discriminant(const ZZX& a, long deterministic=0) { ZZ x; discriminant(x, a, deterministic); NTL_OPT_RETURN(ZZ, x); } // d = discriminant of a // = (-1)^{m(m-1)/2} resultant(a, a')/lc(a), // where m = deg(a) // if !deterministic, then it may use a randomized strategy // that errs with probability no more than 2^{-80}. void CharPolyMod(ZZX& g, const ZZX& a, const ZZX& f, long deterministic=0); inline ZZX CharPolyMod(const ZZX& a, const ZZX& f, long deterministic=0) { ZZX x; CharPolyMod(x, a, f, deterministic); NTL_OPT_RETURN(ZZX, x); } // g = char poly of (a mod f) // f must be monic // if !deterministic, then it may use a randomized strategy // that errs with probability no more than 2^{-80}. void MinPolyMod(ZZX& g, const ZZX& a, const ZZX& f); inline ZZX MinPolyMod(const ZZX& a, const ZZX& f) { ZZX x; MinPolyMod(x, a, f); NTL_OPT_RETURN(ZZX, x); } // g = min poly of (a mod f) // f must be monic // may use a probabilistic strategy that errs with // probability no more than 2^{-80} void XGCD(ZZ& r, ZZX& s, ZZX& t, const ZZX& a, const ZZX& b, long deterministic=0); // r = resultant of a and b; // if r != 0, then computes s and t such that: // a*s + b*t = r; // otherwise s and t not affected. // if !deterministic, then resultant computation may use a randomized strategy // that errs with probability no more than 2^{-80}. /****************************************************** Incremental Chinese Remaindering *******************************************************/ long CRT(ZZX& a, ZZ& prod, const zz_pX& A); long CRT(ZZX& a, ZZ& prod, const ZZ_pX& A); // If p is the current modulus with (p, prod) = 1; // Computes b such that b = a mod prod and b = A mod p, // with coefficients in the interval (-p*prod/2, p*prod/2]; // Sets a = b, prod = p*prod, and returns 1 if a's value changed. // vectors NTL_vector_decl(ZZX,vec_ZZX) NTL_eq_vector_decl(ZZX,vec_ZZX) NTL_io_vector_decl(ZZX,vec_ZZX) NTL_CLOSE_NNS #endif