1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115
|
// polys.cc : implements interface to NTL polynomials
//////////////////////////////////////////////////////////////////////////
//
// Copyright 1990-2012 John Cremona
//
// This file is part of the eclib package.
//
// eclib is free software; you can redistribute it and/or modify it
// under the terms of the GNU General Public License as published by the
// Free Software Foundation; either version 2 of the License, or (at your
// option) any later version.
//
// eclib is distributed in the hope that it will be useful, but WITHOUT
// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
// for more details.
//
// You should have received a copy of the GNU General Public License
// along with eclib; if not, write to the Free Software Foundation,
// Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA
//
//////////////////////////////////////////////////////////////////////////
#include <eclib/marith.h>
#include <eclib/polys.h>
FqPoly reduce(const ZPoly& f, const galois_field& Fq)
{
NewFqPoly(Fq,fmodq);
SetDegree(fmodq,Degree(f));
for(int i=0; i<=Degree(f); i++)
SetCoeff(fmodq,i,ZtoGF(Fq,PolyCoeff(f,i)));
return fmodq;
}
vector<gf_element> roots(const FqPoly& f)
{
// make f monic:
FqPoly f1=f;
MakeMonic(f1);
// reduce to distinct roots case:
ZZ_pX X; SetX(X);
ZZ_pX g = PowerXMod(ZZ_p::modulus(),f1)-X;
vec_ZZ_p r; FindRoots(r,GCD(f1,g));
vector<gf_element>ans;
for(int i=0; i<r.length(); i++) ans.push_back(r[i]);
return ans;
}
vector<bigint> rootsmod(const vector<bigint>& coeffs, bigint q)
{
galois_field Fq(q);
NewFqPoly(Fq,f);
unsigned long i, deg = coeffs.size()-1;
SetDegree(f,deg);
for (i=0; i<=deg; i++) SetCoeff(f,i,ZtoGF(Fq,coeffs[i]));
vector<gf_element> r = roots(f);
vector<bigint>ans;
for(i=0; i<r.size(); i++) ans.push_back(LiftGF(r[i]));
sort(ans.begin(),ans.end());
return ans;
}
//#define TRACE_ROOTS
vector<bigrational> roots(const vector<bigint>& coeffs)
{
#ifdef TRACE_ROOTS
cout<<"Finding rational roots of polynomial f with coefficients "<<coeffs<<endl;
#endif
ZZX f,g; ZZ c;
vector<bigrational> ans;
int i, d = coeffs.size()-1; // degree
if(d<1) return ans;
for(i=0; i<=d; i++)
SetCoeff(f,d-i,coeffs[i]);
#ifdef TRACE_ROOTS
cout<<"f = "<<f<<endl;
#endif
vec_pair_ZZX_long factors;
factor(c,factors,f);
#ifdef TRACE_ROOTS
cout<<"f has "<<factors.length()<<" factors"<<endl;
#endif
for(i=0; i<factors.length(); i++)
{
g = factors[i].a;
#ifdef TRACE_ROOTS
cout<<"factor "<<g<<" has degree "<<deg(g)<<endl;
#endif
if(deg(g)==1)
{
bigrational root = bigrational(-coeff(g,0),coeff(g,1));
#ifdef TRACE_ROOTS
cout<<"root "<<root<<endl;
#endif
ans.push_back(root);
}
}
sort(ans.begin(), ans.end());
return ans;
}
// find the number of roots of X^3 + bX^2 + cX + d = 0 (mod p)
int nrootscubic(const bigint& b, const bigint& c, const bigint& d, const bigint& p)
{
vector<bigint> coeffs;
coeffs.push_back(d);
coeffs.push_back(c);
coeffs.push_back(b);
coeffs.push_back(BIGINT(1));
return rootsmod(coeffs,p).size();
}
|