Loading...
Searching...
No Matches
cfUnivarGcd.cc
Go to the documentation of this file.
#include "config.h"
#include "debug.h"
#include "cf_algorithm.h"
#include "templates/ftmpl_functions.h"
#include "cf_primes.h"
#include "cfGcdUtil.h"
#include "cfUnivarGcd.h"
#ifdef HAVE_NTL
#include "NTLconvert.h"
#endif
#ifdef HAVE_FLINT
#include "FLINTconvert.h"
#endif
#ifdef HAVE_NTL
#ifndef HAVE_FLINT
CanonicalForm
gcd_univar_ntl0( const CanonicalForm & F, const CanonicalForm & G )
{
ZZX F1=convertFacCF2NTLZZX(F);
ZZX G1=convertFacCF2NTLZZX(G);
ZZX R=GCD(F1,G1);
return convertNTLZZX2CF(R,F.mvar());
}
CanonicalForm
gcd_univar_ntlp( const CanonicalForm & F, const CanonicalForm & G )
{
if (fac_NTL_char!=getCharacteristic())
{
fac_NTL_char=getCharacteristic();
zz_p::init(getCharacteristic());
}
zz_pX F1=convertFacCF2NTLzzpX(F);
zz_pX G1=convertFacCF2NTLzzpX(G);
zz_pX R=GCD(F1,G1);
return convertNTLzzpX2CF(R,F.mvar());
}
#endif
#endif
#ifdef HAVE_FLINT
CanonicalForm
gcd_univar_flintp (const CanonicalForm& F, const CanonicalForm& G)
{
nmod_poly_t F1, G1;
convertFacCF2nmod_poly_t (F1, F);
convertFacCF2nmod_poly_t (G1, G);
nmod_poly_gcd (F1, F1, G1);
CanonicalForm result= convertnmod_poly_t2FacCF (F1, F.mvar());
nmod_poly_clear (F1);
nmod_poly_clear (G1);
return result;
}
CanonicalForm
gcd_univar_flint0( const CanonicalForm & F, const CanonicalForm & G )
{
fmpz_poly_t F1, G1;
convertFacCF2Fmpz_poly_t(F1, F);
convertFacCF2Fmpz_poly_t(G1, G);
fmpz_poly_gcd (F1, F1, G1);
CanonicalForm result= convertFmpz_poly_t2FacCF (F1, F.mvar());
fmpz_poly_clear (F1);
fmpz_poly_clear (G1);
return result;
}
#endif
#ifndef HAVE_NTL
CanonicalForm gcd_poly_univar0( const CanonicalForm & F, const CanonicalForm & G, bool primitive )
{
CanonicalForm f, g, c, cg, cl, BB, B, M, q, Dp, newD, D, newq;
int p, i;
if ( primitive )
{
f = F;
g = G;
c = 1;
}
else
{
CanonicalForm cF = content( F ), cG = content( G );
f = F / cF;
g = G / cG;
c = bgcd( cF, cG );
}
cg = gcd( f.lc(), g.lc() );
cl = ( f.lc() / cg ) * g.lc();
// B = 2 * cg * tmin(
// maxnorm(f)*power(CanonicalForm(2),f.degree())*isqrt(f.degree()+1),
// maxnorm(g)*power(CanonicalForm(2),g.degree())*isqrt(g.degree()+1)
// )+1;
M = tmin( maxNorm(f), maxNorm(g) );
BB = power(CanonicalForm(2),tmin(f.degree(),g.degree()))*M;
q = 0;
i = cf_getNumSmallPrimes() - 1;
while ( true )
{
B = BB;
while ( i >= 0 && q < B )
{
p = cf_getSmallPrime( i );
i--;
while ( i >= 0 && mod( cl, p ) == 0 )
{
p = cf_getSmallPrime( i );
i--;
}
setCharacteristic( p );
Dp = gcd( mapinto( f ), mapinto( g ) );
Dp = ( Dp / Dp.lc() ) * mapinto( cg );
setCharacteristic( 0 );
if ( Dp.degree() == 0 )
return c;
if ( q.isZero() )
{
D = mapinto( Dp );
q = p;
B = power(CanonicalForm(2),D.degree())*M+1;
}
else
{
if ( Dp.degree() == D.degree() )
{
chineseRemainder( D, q, mapinto( Dp ), p, newD, newq );
q = newq;
D = newD;
}
else if ( Dp.degree() < D.degree() )
{
// all previous p's are bad primes
q = p;
D = mapinto( Dp );
B = power(CanonicalForm(2),D.degree())*M+1;
}
// else p is a bad prime
}
}
if ( i >= 0 )
{
// now balance D mod q
D = pp( balance_p( D, q ) );
if ( fdivides( D, f ) && fdivides( D, g ) )
return D * c;
else
q = 0;
}
else
return gcd_poly( F, G );
DEBOUTLN( cerr, "another try ..." );
}
}
#endif
/** CanonicalForm extgcd ( const CanonicalForm & f, const CanonicalForm & g, CanonicalForm & a, CanonicalForm & b )
*
* extgcd() - returns polynomial extended gcd of f and g.
*
* Returns gcd(f, g) and a and b sucht that f*a+g*b=gcd(f, g).
* The gcd is calculated using an extended euclidean polynomial
* remainder sequence, so f and g should be polynomials over an
* euclidean domain. Normalizes result.
*
* Note: be sure that f and g have the same level!
*
**/
CanonicalForm
extgcd ( const CanonicalForm & f, const CanonicalForm & g, CanonicalForm & a, CanonicalForm & b )
{
if (f.isZero())
{
a= 0;
b= 1;
return g;
}
else if (g.isZero())
{
a= 1;
b= 0;
return f;
}
#ifdef HAVE_FLINT
if (( getCharacteristic() > 0 ) && (CFFactory::gettype() != GaloisFieldDomain)
&& (f.level()==g.level()) && isPurePoly(f) && isPurePoly(g))
{
nmod_poly_t F1, G1, A, B, R;
convertFacCF2nmod_poly_t (F1, f);
convertFacCF2nmod_poly_t (G1, g);
nmod_poly_init (R, getCharacteristic());
nmod_poly_init (A, getCharacteristic());
nmod_poly_init (B, getCharacteristic());
nmod_poly_xgcd (R, A, B, F1, G1);
a= convertnmod_poly_t2FacCF (A, f.mvar());
b= convertnmod_poly_t2FacCF (B, f.mvar());
CanonicalForm r= convertnmod_poly_t2FacCF (R, f.mvar());
nmod_poly_clear (F1);
nmod_poly_clear (G1);
nmod_poly_clear (A);
nmod_poly_clear (B);
nmod_poly_clear (R);
return r;
}
#elif defined(HAVE_NTL)
if (( getCharacteristic() > 0 ) && (CFFactory::gettype() != GaloisFieldDomain)
&& (f.level()==g.level()) && isPurePoly(f) && isPurePoly(g))
{
if (fac_NTL_char!=getCharacteristic())
{
fac_NTL_char=getCharacteristic();
zz_p::init(getCharacteristic());
}
zz_pX F1=convertFacCF2NTLzzpX(f);
zz_pX G1=convertFacCF2NTLzzpX(g);
zz_pX R;
zz_pX A,B;
XGCD(R,A,B,F1,G1);
a=convertNTLzzpX2CF(A,f.mvar());
b=convertNTLzzpX2CF(B,f.mvar());
return convertNTLzzpX2CF(R,f.mvar());
}
#endif
#ifdef HAVE_FLINT
if (( getCharacteristic() ==0) && (f.level()==g.level())
&& isPurePoly(f) && isPurePoly(g))
{
fmpq_poly_t F1, G1;
convertFacCF2Fmpq_poly_t (F1, f);
convertFacCF2Fmpq_poly_t (G1, g);
fmpq_poly_t R, A, B;
fmpq_poly_init (R);
fmpq_poly_init (A);
fmpq_poly_init (B);
fmpq_poly_xgcd (R, A, B, F1, G1);
a= convertFmpq_poly_t2FacCF (A, f.mvar());
b= convertFmpq_poly_t2FacCF (B, f.mvar());
CanonicalForm r= convertFmpq_poly_t2FacCF (R, f.mvar());
fmpq_poly_clear (F1);
fmpq_poly_clear (G1);
fmpq_poly_clear (A);
fmpq_poly_clear (B);
fmpq_poly_clear (R);
return r;
}
#elif defined(HAVE_NTL)
if (( getCharacteristic() ==0)
&& (f.level()==g.level()) && isPurePoly(f) && isPurePoly(g))
{
CanonicalForm fc=bCommonDen(f);
CanonicalForm gc=bCommonDen(g);
ZZX F1=convertFacCF2NTLZZX(f*fc);
ZZX G1=convertFacCF2NTLZZX(g*gc);
ZZX R=GCD(F1,G1);
CanonicalForm r=convertNTLZZX2CF(R,f.mvar());
ZZ RR;
ZZX A,B;
if (r.inCoeffDomain())
{
XGCD(RR,A,B,F1,G1,1);
CanonicalForm rr=convertZZ2CF(RR);
if(!rr.isZero())
{
a=convertNTLZZX2CF(A,f.mvar())*fc/rr;
b=convertNTLZZX2CF(B,f.mvar())*gc/rr;
return CanonicalForm(1);
}
else
{
F1 /= R;
G1 /= R;
XGCD (RR, A,B,F1,G1,1);
rr=convertZZ2CF(RR);
a=convertNTLZZX2CF(A,f.mvar())*(fc/rr);
b=convertNTLZZX2CF(B,f.mvar())*(gc/rr);
}
}
else
{
XGCD(RR,A,B,F1,G1,1);
CanonicalForm rr=convertZZ2CF(RR);
if (!rr.isZero())
{
a=convertNTLZZX2CF(A,f.mvar())*fc;
b=convertNTLZZX2CF(B,f.mvar())*gc;
}
else
{
F1 /= R;
G1 /= R;
XGCD (RR, A,B,F1,G1,1);
rr=convertZZ2CF(RR);
a=convertNTLZZX2CF(A,f.mvar())*(fc/rr);
b=convertNTLZZX2CF(B,f.mvar())*(gc/rr);
}
return r;
}
}
#endif
// may contain bug in the co-factors, see track 107
CanonicalForm contf = content( f );
CanonicalForm contg = content( g );
CanonicalForm p0 = f / contf, p1 = g / contg;
CanonicalForm f0 = 1, f1 = 0, g0 = 0, g1 = 1, q, r;
while ( ! p1.isZero() )
{
divrem( p0, p1, q, r );
p0 = p1; p1 = r;
r = g0 - g1 * q;
g0 = g1; g1 = r;
r = f0 - f1 * q;
f0 = f1; f1 = r;
}
CanonicalForm contp0 = content( p0 );
a = f0 / ( contf * contp0 );
b = g0 / ( contg * contp0 );
p0 /= contp0;
if ( p0.sign() < 0 )
{
p0 = -p0;
a = -a;
b = -b;
}
return p0;
}
convertFmpq_poly_t2FacCF
CanonicalForm convertFmpq_poly_t2FacCF(const fmpq_poly_t p, const Variable &x)
conversion of a FLINT poly over Q to CanonicalForm
Definition FLINTconvert.cc:294
convertnmod_poly_t2FacCF
CanonicalForm convertnmod_poly_t2FacCF(const nmod_poly_t poly, const Variable &x)
conversion of a FLINT poly over Z/p to CanonicalForm
Definition FLINTconvert.cc:208
convertFmpz_poly_t2FacCF
CanonicalForm convertFmpz_poly_t2FacCF(const fmpz_poly_t poly, const Variable &x)
conversion of a FLINT poly over Z to CanonicalForm
Definition FLINTconvert.cc:172
convertFacCF2Fmpq_poly_t
void convertFacCF2Fmpq_poly_t(fmpq_poly_t result, const CanonicalForm &f)
conversion of a factory univariate polynomials over Q to fmpq_poly_t
Definition FLINTconvert.cc:320
convertFacCF2Fmpz_poly_t
void convertFacCF2Fmpz_poly_t(fmpz_poly_t result, const CanonicalForm &f)
conversion of a factory univariate polynomial over Z to a fmpz_poly_t
Definition FLINTconvert.cc:146
FLINTconvert.h
This file defines functions for conversion to FLINT (www.flintlib.org) and back.
convertZZ2CF
CanonicalForm convertZZ2CF(const ZZ &a)
NAME: convertZZ2CF.
Definition NTLconvert.cc:495
convertFacCF2NTLZZX
ZZX convertFacCF2NTLZZX(const CanonicalForm &f)
Definition NTLconvert.cc:691
convertNTLzzpX2CF
CanonicalForm convertNTLzzpX2CF(const zz_pX &poly, const Variable &x)
Definition NTLconvert.cc:255
convertNTLZZX2CF
CanonicalForm convertNTLZZX2CF(const ZZX &polynom, const Variable &x)
Definition NTLconvert.cc:285
convertFacCF2NTLzzpX
zz_pX convertFacCF2NTLzzpX(const CanonicalForm &f)
Definition NTLconvert.cc:105
fac_NTL_char
VAR long fac_NTL_char
Definition NTLconvert.cc:46
NTLconvert.h
Conversion to and from NTL.
bgcd
CanonicalForm bgcd(const CanonicalForm &f, const CanonicalForm &g)
CanonicalForm bgcd ( const CanonicalForm & f, const CanonicalForm & g ).
Definition canonicalform.cc:1648
divrem
void divrem(const CanonicalForm &f, const CanonicalForm &g, CanonicalForm &q, CanonicalForm &r)
Definition canonicalform.cc:1032
power
CanonicalForm power(const CanonicalForm &f, int n)
exponentiation
Definition canonicalform.cc:1896
content
CanonicalForm FACTORY_PUBLIC content(const CanonicalForm &)
CanonicalForm content ( const CanonicalForm & f ).
Definition cf_gcd.cc:603
mapinto
CanonicalForm mapinto(const CanonicalForm &f)
Definition canonicalform.h:355
mod
CF_NO_INLINE FACTORY_PUBLIC CanonicalForm mod(const CanonicalForm &, const CanonicalForm &)
pp
CanonicalForm FACTORY_PUBLIC pp(const CanonicalForm &)
CanonicalForm pp ( const CanonicalForm & f ).
Definition cf_gcd.cc:676
setCharacteristic
void FACTORY_PUBLIC setCharacteristic(int c)
Definition cf_char.cc:28
gcd_poly
CanonicalForm FACTORY_PUBLIC gcd_poly(const CanonicalForm &f, const CanonicalForm &g)
CanonicalForm gcd_poly ( const CanonicalForm & f, const CanonicalForm & g ).
Definition cf_gcd.cc:492
getCharacteristic
int FACTORY_PUBLIC getCharacteristic()
Definition cf_char.cc:70
i
int i
Definition cfEzgcd.cc:132
cfGcdUtil.h
coprimality check and change of representation mod n
p
int p
Definition cfModGcd.cc:4086
cl
cl
Definition cfModGcd.cc:4108
g
g
Definition cfModGcd.cc:4098
b
CanonicalForm b
Definition cfModGcd.cc:4111
cg
CanonicalForm cg
Definition cfModGcd.cc:4091
gcd_univar_flint0
CanonicalForm gcd_univar_flint0(const CanonicalForm &F, const CanonicalForm &G)
Definition cfUnivarGcd.cc:61
extgcd
CanonicalForm extgcd(const CanonicalForm &f, const CanonicalForm &g, CanonicalForm &a, CanonicalForm &b)
CanonicalForm extgcd ( const CanonicalForm & f, const CanonicalForm & g, CanonicalForm & a,...
Definition cfUnivarGcd.cc:174
gcd_univar_flintp
CanonicalForm gcd_univar_flintp(const CanonicalForm &F, const CanonicalForm &G)
Definition cfUnivarGcd.cc:48
cfUnivarGcd.h
univariate Gcd over finite fields and Z, extended GCD over finite fields and Q
bCommonDen
CanonicalForm bCommonDen(const CanonicalForm &f)
void psqr ( const CanonicalForm & f, const CanonicalForm & g, CanonicalForm & q, CanonicalForm & r,...
Definition cf_algorithm.cc:293
maxNorm
CanonicalForm maxNorm(const CanonicalForm &f)
bool fdivides ( const CanonicalForm & f, const CanonicalForm & g )
Definition cf_algorithm.cc:536
cf_algorithm.h
declarations of higher level algorithms.
isPurePoly
bool isPurePoly(const CanonicalForm &f)
Definition cf_factor.cc:248
fdivides
bool fdivides(const CanonicalForm &f, const CanonicalForm &g)
chineseRemainder
void FACTORY_PUBLIC chineseRemainder(const CanonicalForm &x1, const CanonicalForm &q1, const CanonicalForm &x2, const CanonicalForm &q2, CanonicalForm &xnew, CanonicalForm &qnew)
void chineseRemainder ( const CanonicalForm & x1, const CanonicalForm & q1, const CanonicalForm & x2,...
Definition cf_chinese.cc:57
GaloisFieldDomain
#define GaloisFieldDomain
Definition cf_defs.h:18
gcd_poly_univar0
static CanonicalForm gcd_poly_univar0(const CanonicalForm &F, const CanonicalForm &G, bool primitive)
Definition cf_gcd.cc:178
balance_p
static CanonicalForm balance_p(const CanonicalForm &f, const CanonicalForm &q, const CanonicalForm &qh)
Definition cf_gcd.cc:149
cf_getNumSmallPrimes
int cf_getNumSmallPrimes()
Definition cf_primes.cc:34
cf_getSmallPrime
int cf_getSmallPrime(int i)
Definition cf_primes.cc:28
cf_primes.h
access to prime tables
f
FILE * f
Definition checklibs.c:9
CFFactory::gettype
static int gettype()
Definition cf_factory.h:28
CanonicalForm
factory's main class
Definition canonicalform.h:86
CanonicalForm::isZero
CF_NO_INLINE bool isZero() const
CanonicalForm::degree
int degree() const
Returns -1 for the zero polynomial and 0 if CO is in a base domain.
Definition canonicalform.cc:414
CanonicalForm::sign
int sign() const
int CanonicalForm::sign () const
Definition canonicalform.cc:1360
CanonicalForm::mvar
Variable mvar() const
mvar() returns the main variable of CO or Variable() if CO is in a base domain.
Definition canonicalform.cc:593
CanonicalForm::inCoeffDomain
bool inCoeffDomain() const
Definition canonicalform.cc:122
CanonicalForm::lc
CanonicalForm lc() const
CanonicalForm CanonicalForm::lc (), Lc (), LC (), LC ( v ) const.
Definition canonicalform.cc:337
debug.h
functions to print debug output
DEBOUTLN
#define DEBOUTLN(stream, objects)
Definition debug.h:49
result
return result
Definition facAbsBiFact.cc:76
B
b *CanonicalForm B
Definition facBivar.cc:52
nmod_poly_init
nmod_poly_init(FLINTmipo, getCharacteristic())
convertFacCF2nmod_poly_t
convertFacCF2nmod_poly_t(FLINTmipo, M)
nmod_poly_clear
nmod_poly_clear(FLINTmipo)
ftmpl_functions.h
some useful template functions.
tmin
template CanonicalForm tmin(const CanonicalForm &, const CanonicalForm &)
D
#define D(A)
Definition gentable.cc:128
G
STATIC_VAR TreeM * G
Definition janet.cc:31
R
#define R
Definition sirandom.c:27
A
#define A
Definition sirandom.c:24
M
#define M
Definition sirandom.c:25
gcd
int gcd(int a, int b)
Definition walkSupport.cc:836