Loading...
Searching...
No Matches
cfGcdAlgExt.cc File Reference
`#include "config.h"`
`#include <vector>`
`#include <stdio.h>`
`#include <iostream>`
`#include "cf_assert.h"`
`#include "timing.h"`
`#include "templates/ftmpl_functions.h"`
`#include "cf_defs.h"`
`#include "canonicalform.h"`
`#include "cf_iter.h"`
`#include "cf_primes.h"`
`#include "cf_algorithm.h"`
`#include "cfGcdAlgExt.h"`
`#include "cfUnivarGcd.h"`
`#include "cf_map.h"`
`#include "cf_generator.h"`
`#include "facMul.h"`
`#include "cfNTLzzpEXGCD.h"`
`#include "NTLconvert.h"`
`#include "FLINTconvert.h"`

Go to the source code of this file.

## Functions

TIMING_DEFINE_PRINT (alg_content_p) TIMING_DEFINE_PRINT(alg_content) TIMING_DEFINE_PRINT(alg_compress) TIMING_DEFINE_PRINT(alg_termination) TIMING_DEFINE_PRINT(alg_termination_p) TIMING_DEFINE_PRINT(alg_reconstruction) TIMING_DEFINE_PRINT(alg_newton_p) TIMING_DEFINE_PRINT(alg_recursion_p) TIMING_DEFINE_PRINT(alg_gcd_p) TIMING_DEFINE_PRINT(alg_euclid_p) static int myCompress(const CanonicalForm &F

| compressing two polynomials F and G, M is used for compressing, N to reverse the compression
| for (int i=0;i<=n;i++) degsf[i]
| if (topLevel)
| DELETE_ARRAY (degsg)
void | tryInvert (const CanonicalForm &F, const CanonicalForm &M, CanonicalForm &inv, bool &fail)
static CanonicalForm | trycontent (const CanonicalForm &f, const Variable &x, const CanonicalForm &M, bool &fail)
static CanonicalForm | tryvcontent (const CanonicalForm &f, const Variable &x, const CanonicalForm &M, bool &fail)
static CanonicalForm | trycf_content (const CanonicalForm &f, const CanonicalForm &g, const CanonicalForm &M, bool &fail)
static CanonicalForm | tryNewtonInterp (const CanonicalForm &alpha, const CanonicalForm &u, const CanonicalForm &newtonPoly, const CanonicalForm &oldInterPoly, const Variable &x, const CanonicalForm &M, bool &fail)
static void | leadDeg (const CanonicalForm &f, int degs[])
void | tryBrownGCD (const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &M, CanonicalForm &result, bool &fail, bool topLevel)
| modular gcd over F_p[x]/(M) for not necessarily irreducible M. If a zero divisor is encountered fail is set to true.
static CanonicalForm | myicontent (const CanonicalForm &f, const CanonicalForm &c)
static CanonicalForm | myicontent (const CanonicalForm &f)
CanonicalForm | QGCD (const CanonicalForm &F, const CanonicalForm &G)
| gcd over Q(a)
bool | isLess (int *a, int *b, int lower, int upper)
bool | isEqual (int *a, int *b, int lower, int upper)
CanonicalForm | firstLC (const CanonicalForm &f)

## Variables

const CanonicalForm & | G
const CanonicalForm CFMap & | M
const CanonicalForm CFMap CFMap & | N
const CanonicalForm CFMap CFMap bool | topLevel
int * | degsf = NEW_ARRAY(int,n + 1)
int * | degsg = NEW_ARRAY(int,n + 1)
int | both_non_zero = 0
int | f_zero = 0
int | g_zero = 0
int | both_zero = 0
int | Flevel =F.level()
int | Glevel =G.level()
| else
| return

## Function Documentation

## DELETE_ARRAY()

DELETE_ARRAY | ( | degsg |  | ) |

## firstLC()

CanonicalForm firstLC | ( | const CanonicalForm & | f | ) |

Definition at line 956 of file cfGcdAlgExt.cc.

{ // returns the leading coefficient (LC) of level <= 1
CanonicalForm ret = f;
while(ret.level() > 1)
ret = LC(ret);
return ret;
}
LC
CanonicalForm LC(const CanonicalForm &f)
Definition canonicalform.h:310
f
FILE * f
Definition checklibs.c:9
CanonicalForm
factory's main class
Definition canonicalform.h:86
CanonicalForm::level
int level() const
level() returns the level of CO.
Definition canonicalform.cc:576

## for()

for | ( | int | i = 0;i<=n;i++ | ) |

Definition at line 72 of file cfEzgcd.cc.

{
if (degsf[i] != 0 && degsg[i] != 0)
{
both_non_zero++;
continue;
}
if (degsf[i] == 0 && degsg[i] != 0 && i <= G.level())
{
f_zero++;
continue;
}
if (degsg[i] == 0 && degsf[i] && i <= F.level())
{
g_zero++;
continue;
}
}
degsf
int * degsf
Definition cfEzgcd.cc:59
f_zero
int f_zero
Definition cfEzgcd.cc:69
both_non_zero
const CanonicalForm CFMap CFMap int & both_non_zero
Definition cfEzgcd.cc:57
i
int i
Definition cfEzgcd.cc:132
g_zero
int g_zero
Definition cfEzgcd.cc:70
degsg
int * degsg
Definition cfEzgcd.cc:60
G
STATIC_VAR TreeM * G
Definition janet.cc:31

## if()

if | ( | topLevel |  | ) |

Definition at line 76 of file cfGcdAlgExt.cc.

{
for (int i= 1; i <= n; i++)
{
if (degsf[i] != 0 && degsg[i] != 0)
{
both_non_zero++;
continue;
}
if (degsf[i] == 0 && degsg[i] != 0 && i <= Glevel)
{
f_zero++;
continue;
}
if (degsg[i] == 0 && degsf[i] && i <= Flevel)
{
g_zero++;
continue;
}
}
if (both_non_zero == 0)
{
DELETE_ARRAY(degsf);
DELETE_ARRAY(degsg);
return 0;
}
// map Variables which do not occur in both polynomials to higher levels
int k= 1;
int l= 1;
for (int i= 1; i <= n; i++)
{
if (degsf[i] != 0 && degsg[i] == 0 && i <= Flevel)
{
if (k + both_non_zero != i)
{
M.newpair (Variable (i), Variable (k + both_non_zero));
N.newpair (Variable (k + both_non_zero), Variable (i));
}
k++;
}
if (degsf[i] == 0 && degsg[i] != 0 && i <= Glevel)
{
if (l + g_zero + both_non_zero != i)
{
M.newpair (Variable (i), Variable (l + g_zero + both_non_zero));
N.newpair (Variable (l + g_zero + both_non_zero), Variable (i));
}
l++;
}
}
// sort Variables x_{i} in increasing order of
// min(deg_{x_{i}}(f),deg_{x_{i}}(g))
int m= tmax (Flevel, Glevel);
int min_max_deg;
k= both_non_zero;
l= 0;
int i= 1;
while (k > 0)
{
if (degsf [i] != 0 && degsg [i] != 0)
min_max_deg= tmax (degsf[i], degsg[i]);
else
min_max_deg= 0;
while (min_max_deg == 0)
{
i++;
min_max_deg= tmax (degsf[i], degsg[i]);
if (degsf [i] != 0 && degsg [i] != 0)
min_max_deg= tmax (degsf[i], degsg[i]);
else
min_max_deg= 0;
}
for (int j= i + 1; j <= m; j++)
{
if (tmax (degsf[j],degsg[j]) <= min_max_deg && degsf[j] != 0 && degsg [j] != 0)
{
min_max_deg= tmax (degsf[j], degsg[j]);
l= j;
}
}
if (l != 0)
{
if (l != k)
{
M.newpair (Variable (l), Variable(k));
N.newpair (Variable (k), Variable(l));
degsf[l]= 0;
degsg[l]= 0;
l= 0;
}
else
{
degsf[l]= 0;
degsg[l]= 0;
l= 0;
}
}
else if (l == 0)
{
if (i != k)
{
M.newpair (Variable (i), Variable (k));
N.newpair (Variable (k), Variable (i));
degsf[i]= 0;
degsg[i]= 0;
}
else
{
degsf[i]= 0;
degsg[i]= 0;
}
i++;
}
k--;
}
}
N
const CanonicalForm CFMap CFMap & N
Definition cfEzgcd.cc:56
l
int l
Definition cfEzgcd.cc:100
m
int m
Definition cfEzgcd.cc:128
Flevel
int Flevel
Definition cfEzgcd.cc:101
Glevel
int Glevel
Definition cfEzgcd.cc:102
k
int k
Definition cfEzgcd.cc:99
DELETE_ARRAY
#define DELETE_ARRAY(P)
Definition cf_defs.h:65
Variable
factory's class for variables
Definition factory.h:127
j
int j
Definition facHensel.cc:110
tmax
template CanonicalForm tmax(const CanonicalForm &, const CanonicalForm &)
M
#define M
Definition sirandom.c:25

## isEqual()

bool isEqual | ( | int * | a,
|  | int * | b,
|  | int | lower,
|  | int | upper )

Definition at line 947 of file cfGcdAlgExt.cc.

{ // compares the degree vectors a,b on the specified part. Note: x(i+1) > x(i)
for(int i=lower; i<=upper; i++)
if(a[i] != b[i])
return false;
return true;
}
b
CanonicalForm b
Definition cfModGcd.cc:4111

## isLess()

bool isLess | ( | int * | a,
|  | int * | b,
|  | int | lower,
|  | int | upper )

Definition at line 936 of file cfGcdAlgExt.cc.

{ // compares the degree vectors a,b on the specified part. Note: x(i+1) > x(i)
for(int i=upper; i>=lower; i--)
if(a[i] == b[i])
continue;
else
return a[i] < b[i];
return true;
}

## leadDeg()

void leadDeg | ( | const CanonicalForm & | f,
|  | int | degs[] )

static

Definition at line 373 of file cfGcdAlgExt.cc.

{ // leading degree vector w.r.t. lex. monomial order x(i+1) > x(i)
// 'a' should point to an array of sufficient size level(f)+1
if(f.inCoeffDomain())
return;
CanonicalForm tmp = f;
do
{
degs[tmp.level()] = tmp.degree();
tmp = LC(tmp);
}
while(!tmp.inCoeffDomain());
}
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::inCoeffDomain
bool inCoeffDomain() const
Definition canonicalform.cc:122

## myicontent() [1/2]

CanonicalForm myicontent | ( | const CanonicalForm & | f | ) |

static

Definition at line 722 of file cfGcdAlgExt.cc.

{
#if defined(HAVE_NTL) || defined(HAVE_FLINT)
return myicontent( f, 0 );
#else
return 1;
#endif
}
myicontent
static CanonicalForm myicontent(const CanonicalForm &f, const CanonicalForm &c)
Definition cfGcdAlgExt.cc:672

## myicontent() [2/2]

CanonicalForm myicontent | ( | const CanonicalForm & | f,
|  | const CanonicalForm & | c )

static

Definition at line 672 of file cfGcdAlgExt.cc.

{
#if defined(HAVE_NTL) || defined(HAVE_FLINT)
if (f.isOne() || c.isOne())
return 1;
if ( f.inBaseDomain() && c.inBaseDomain())
{
if (c.isZero()) return abs(f);
return bgcd( f, c );
}
else if ( (f.inCoeffDomain() && c.inCoeffDomain()) ||
(f.inCoeffDomain() && c.inBaseDomain()) ||
(f.inBaseDomain() && c.inCoeffDomain()))
{
if (c.isZero()) return abs (f);
#ifdef HAVE_FLINT
fmpz_poly_t FLINTf, FLINTc;
convertFacCF2Fmpz_poly_t (FLINTf, f);
convertFacCF2Fmpz_poly_t (FLINTc, c);
fmpz_poly_gcd (FLINTc, FLINTc, FLINTf);
CanonicalForm result;
if (f.inCoeffDomain())
result= convertFmpz_poly_t2FacCF (FLINTc, f.mvar());
else
result= convertFmpz_poly_t2FacCF (FLINTc, c.mvar());
fmpz_poly_clear (FLINTc);
fmpz_poly_clear (FLINTf);
return result;
#else
ZZX NTLf= convertFacCF2NTLZZX (f);
ZZX NTLc= convertFacCF2NTLZZX (c);
NTLc= GCD (NTLc, NTLf);
if (f.inCoeffDomain())
return convertNTLZZX2CF(NTLc,f.mvar());
else
return convertNTLZZX2CF(NTLc,c.mvar());
#endif
}
else
{
CanonicalForm g = c;
for ( CFIterator i = f; i.hasTerms() && ! g.isOne(); i++ )
g = myicontent( i.coeff(), g );
return g;
}
#else
return 1;
#endif
}
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
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
abs
Rational abs(const Rational &a)
Definition GMPrat.cc:436
convertFacCF2NTLZZX
ZZX convertFacCF2NTLZZX(const CanonicalForm &f)
Definition NTLconvert.cc:691
convertNTLZZX2CF
CanonicalForm convertNTLZZX2CF(const ZZX &polynom, const Variable &x)
Definition NTLconvert.cc:285
bgcd
CanonicalForm bgcd(const CanonicalForm &f, const CanonicalForm &g)
CanonicalForm bgcd ( const CanonicalForm & f, const CanonicalForm & g ).
Definition canonicalform.cc:1648
g
g
Definition cfModGcd.cc:4098
CFIterator
class to iterate through CanonicalForm's
Definition cf_iter.h:44
CanonicalForm::isZero
CF_NO_INLINE bool isZero() const
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::isOne
CF_NO_INLINE bool isOne() const
CanonicalForm::inBaseDomain
bool inBaseDomain() const
Definition canonicalform.cc:104
result
return result
Definition facAbsBiFact.cc:76

## QGCD()

CanonicalForm QGCD | ( | const CanonicalForm & | F,
|  | const CanonicalForm & | G )

gcd over Q(a)

Definition at line 731 of file cfGcdAlgExt.cc.

{ // f,g in Q(a)[x1,...,xn]
if(F.isZero())
{
if(G.isZero())
return G; // G is zero
if(G.inCoeffDomain())
return CanonicalForm(1);
CanonicalForm lcinv= 1/Lc (G);
return G*lcinv; // return monic G
}
if(G.isZero()) // F is non-zero
{
if(F.inCoeffDomain())
return CanonicalForm(1);
CanonicalForm lcinv= 1/Lc (F);
return F*lcinv; // return monic F
}
if(F.inCoeffDomain() || G.inCoeffDomain())
return CanonicalForm(1);
// here: both NOT inCoeffDomain
CanonicalForm f, g, tmp, M, q, D, Dp, cl, newq, mipo;
int p, i;
int *bound, *other; // degree vectors
bool fail;
bool off_rational=!isOn(SW_RATIONAL);
On( SW_RATIONAL ); // needed by bCommonDen
f = F * bCommonDen(F);
g = G * bCommonDen(G);
TIMING_START (alg_content)
CanonicalForm contf= myicontent (f);
CanonicalForm contg= myicontent (g);
f /= contf;
g /= contg;
CanonicalForm gcdcfcg= myicontent (contf, contg);
TIMING_END_AND_PRINT (alg_content, "time for content in alg gcd: ")
Variable a, b;
if(hasFirstAlgVar(f,a))
{
if(hasFirstAlgVar(g,b))
{
if(b.level() > a.level())
a = b;
}
}
else
{
if(!hasFirstAlgVar(g,a))// both not in extension
{
Off( SW_RATIONAL );
Off( SW_USE_QGCD );
tmp = gcdcfcg*gcd( f, g );
On( SW_USE_QGCD );
if (off_rational) Off(SW_RATIONAL);
return tmp;
}
}
// here: a is the biggest alg. var in f and g AND some of f,g is in extension
setReduce(a,false); // do not reduce expressions modulo mipo
tmp = getMipo(a);
M = tmp * bCommonDen(tmp);
// here: f, g in Z[a][x1,...,xn], M in Z[a] not necessarily monic
Off( SW_RATIONAL ); // needed by mod
// calculate upper bound for degree vector of gcd
int mv = f.level(); i = g.level();
if(i > mv)
mv = i;
// here: mv is level of the largest variable in f, g
bound = new int[mv+1]; // 'bound' could be indexed from 0 to mv, but we will only use from 1 to mv
other = new int[mv+1];
for(int i=1; i<=mv; i++) // initialize 'bound', 'other' with zeros
bound[i] = other[i] = 0;
leadDeg(f,bound); // 'bound' is set the leading degree vector of f
leadDeg(g,other);
for(int i=1; i<=mv; i++) // entry at i=0 not visited
if(other[i] < bound[i])
bound[i] = other[i];
// now 'bound' is the smaller vector
cl = lc(M) * lc(f) * lc(g);
q = 1;
D = 0;
CanonicalForm test= 0;
bool equal= false;
for( i=cf_getNumBigPrimes()-1; i>-1; i-- )
{
p = cf_getBigPrime(i);
if( mod( cl, p ).isZero() ) // skip lc-bad primes
continue;
fail = false;
setCharacteristic(p);
mipo = mapinto(M);
mipo /= mipo.lc();
// here: mipo is monic
TIMING_START (alg_gcd_p)
tryBrownGCD( mapinto(f), mapinto(g), mipo, Dp, fail );
TIMING_END_AND_PRINT (alg_gcd_p, "time for alg gcd mod p: ")
if( fail ) // mipo splits in char p
{
setCharacteristic(0);
continue;
}
if( Dp.inCoeffDomain() ) // early termination
{
tryInvert(Dp,mipo,tmp,fail); // check if zero divisor
setCharacteristic(0);
if(fail)
continue;
setReduce(a,true);
if (off_rational) Off(SW_RATIONAL); else On(SW_RATIONAL);
delete[] other;
delete[] bound;
return gcdcfcg;
}
setCharacteristic(0);
// here: Dp NOT inCoeffDomain
for(int i=1; i<=mv; i++)
other[i] = 0; // reset (this is necessary, because some entries may not be updated by call to leadDeg)
leadDeg(Dp,other);
if(isEqual(bound, other, 1, mv)) // equal
{
chineseRemainder( D, q, mapinto(Dp), p, tmp, newq );
// tmp = Dp mod p
// tmp = D mod q
// newq = p*q
q = newq;
if( D != tmp )
D = tmp;
On( SW_RATIONAL );
TIMING_START (alg_reconstruction)
tmp = Farey( D, q ); // Farey
tmp *= bCommonDen (tmp);
TIMING_END_AND_PRINT (alg_reconstruction,
"time for rational reconstruction in alg gcd: ")
setReduce(a,true); // reduce expressions modulo mipo
On( SW_RATIONAL ); // needed by fdivides
if (test != tmp)
test= tmp;
else
equal= true; // modular image did not add any new information
TIMING_START (alg_termination)
#ifdef HAVE_NTL
#ifdef HAVE_FLINT
if (equal && tmp.isUnivariate() && f.isUnivariate() && g.isUnivariate()
&& f.level() == tmp.level() && tmp.level() == g.level())
{
CanonicalForm Q, R;
newtonDivrem (f, tmp, Q, R);
if (R.isZero())
{
newtonDivrem (g, tmp, Q, R);
if (R.isZero())
{
Off (SW_RATIONAL);
setReduce (a,true);
if (off_rational) Off(SW_RATIONAL); else On(SW_RATIONAL);
TIMING_END_AND_PRINT (alg_termination,
"time for successful termination test in alg gcd: ")
delete[] other;
delete[] bound;
return tmp*gcdcfcg;
}
}
}
else
#endif
#endif
if(equal && fdivides( tmp, f ) && fdivides( tmp, g )) // trial division
{
Off( SW_RATIONAL );
setReduce(a,true);
if (off_rational) Off(SW_RATIONAL); else On(SW_RATIONAL);
TIMING_END_AND_PRINT (alg_termination,
"time for successful termination test in alg gcd: ")
delete[] other;
delete[] bound;
return tmp*gcdcfcg;
}
TIMING_END_AND_PRINT (alg_termination,
"time for unsuccessful termination test in alg gcd: ")
Off( SW_RATIONAL );
setReduce(a,false); // do not reduce expressions modulo mipo
continue;
}
if( isLess(bound, other, 1, mv) ) // current prime unlucky
continue;
// here: isLess(other, bound, 1, mv) ) ==> all previous primes unlucky
q = p;
D = mapinto(Dp); // shortcut CRA // shortcut CRA
for(int i=1; i<=mv; i++) // tighten bound
bound[i] = other[i];
}
// hopefully, we never reach this point
setReduce(a,true);
delete[] other;
delete[] bound;
Off( SW_USE_QGCD );
D = gcdcfcg*gcd( f, g );
On( SW_USE_QGCD );
if (off_rational) Off(SW_RATIONAL); else On(SW_RATIONAL);
return D;
}
isOn
bool isOn(int sw)
switches
Definition canonicalform.cc:1971
On
void On(int sw)
switches
Definition canonicalform.cc:1957
Off
void Off(int sw)
switches
Definition canonicalform.cc:1964
mapinto
CanonicalForm mapinto(const CanonicalForm &f)
Definition canonicalform.h:355
lc
CanonicalForm lc(const CanonicalForm &f)
Definition canonicalform.h:304
mod
CF_NO_INLINE FACTORY_PUBLIC CanonicalForm mod(const CanonicalForm &, const CanonicalForm &)
setCharacteristic
void FACTORY_PUBLIC setCharacteristic(int c)
Definition cf_char.cc:28
hasFirstAlgVar
bool hasFirstAlgVar(const CanonicalForm &f, Variable &a)
check if poly f contains an algebraic variable a
Definition cf_ops.cc:679
Lc
CanonicalForm Lc(const CanonicalForm &f)
Definition canonicalform.h:307
for
for(int i=0;i<=n;i++) degsf[i]
Definition cfEzgcd.cc:72
else
else
Definition cfGcdAlgExt.cc:196
isLess
bool isLess(int *a, int *b, int lower, int upper)
Definition cfGcdAlgExt.cc:936
leadDeg
static void leadDeg(const CanonicalForm &f, int degs[])
Definition cfGcdAlgExt.cc:373
isEqual
bool isEqual(int *a, int *b, int lower, int upper)
Definition cfGcdAlgExt.cc:947
return
return
Definition cfGcdAlgExt.cc:219
tryInvert
void tryInvert(const CanonicalForm &F, const CanonicalForm &M, CanonicalForm &inv, bool &fail)
Definition cfGcdAlgExt.cc:222
p
int p
Definition cfModGcd.cc:4086
false
return false
Definition cfModGcd.cc:85
cl
cl
Definition cfModGcd.cc:4108
gcdcfcg
CanonicalForm gcdcfcg
Definition cfModGcd.cc:4109
equal
bool equal
Definition cfModGcd.cc:4134
test
CanonicalForm test
Definition cfModGcd.cc:4104
bCommonDen
CanonicalForm bCommonDen(const CanonicalForm &f)
void psqr ( const CanonicalForm & f, const CanonicalForm & g, CanonicalForm & q, CanonicalForm & r,...
Definition cf_algorithm.cc:293
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
SW_USE_QGCD
static const int SW_USE_QGCD
set to 1 to use Encarnacion GCD over Q(a)
Definition cf_defs.h:43
SW_RATIONAL
static const int SW_RATIONAL
set to 1 for computations over Q
Definition cf_defs.h:31
bound
static CanonicalForm bound(const CFMatrix &M)
Definition cf_linsys.cc:460
cf_getBigPrime
int cf_getBigPrime(int i)
Definition cf_primes.cc:39
cf_getNumBigPrimes
int cf_getNumBigPrimes()
Definition cf_primes.cc:45
Variable::level
int level() const
Definition factory.h:143
mipo
CanonicalForm mipo
Definition facAlgExt.cc:57
alg_content
CanonicalForm alg_content(const CanonicalForm &f, const CFList &as)
Definition facAlgFunc.cc:42
isZero
bool isZero(const CFArray &A)
checks if entries of A are zero
Definition facSparseHensel.h:468
setReduce
void setReduce(const Variable &alpha, bool reduce)
Definition variable.cc:238
getMipo
CanonicalForm getMipo(const Variable &alpha, const Variable &x)
Definition variable.cc:207
if
if(!FE_OPT_NO_SHELL_FLAG)
Definition fehelp.cc:1000
Farey
static number Farey(number, number, const coeffs)
Definition flintcf_Q.cc:495
D
#define D(A)
Definition gentable.cc:128
TIMING_START
#define TIMING_START(t)
Definition timing.h:92
TIMING_END_AND_PRINT
#define TIMING_END_AND_PRINT(t, msg)
Definition timing.h:94
gcd
int gcd(int a, int b)
Definition walkSupport.cc:836

## TIMING_DEFINE_PRINT()

TIMING_DEFINE_PRINT | ( | alg_content_p |  | ) | const &

compressing two polynomials F and G, M is used for compressing, N to reverse the compression

## tryBrownGCD()

void tryBrownGCD | ( | const CanonicalForm & | F,
|  | const CanonicalForm & | G,
|  | const CanonicalForm & | M,
|  | CanonicalForm & | result,
|  | bool & | fail,
|  | bool | topLevel )

modular gcd over F_p[x]/(M) for not necessarily irreducible M. If a zero divisor is encountered fail is set to true.

Definition at line 387 of file cfGcdAlgExt.cc.

{ // assume F,G are multivariate polys over Z/p(a) for big prime p, M "univariate" polynomial in an algebraic variable
// M is assumed to be monic
if(F.isZero())
{
if(G.isZero())
{
result = G; // G is zero
return;
}
if(G.inCoeffDomain())
{
tryInvert(G,M,result,fail);
if(fail)
return;
result = 1;
return;
}
// try to make G monic modulo M
CanonicalForm inv;
tryInvert(Lc(G),M,inv,fail);
if(fail)
return;
result = inv*G;
result= reduce (result, M);
return;
}
if(G.isZero()) // F is non-zero
{
if(F.inCoeffDomain())
{
tryInvert(F,M,result,fail);
if(fail)
return;
result = 1;
return;
}
// try to make F monic modulo M
CanonicalForm inv;
tryInvert(Lc(F),M,inv,fail);
if(fail)
return;
result = inv*F;
result= reduce (result, M);
return;
}
// here: F,G both nonzero
if(F.inCoeffDomain())
{
tryInvert(F,M,result,fail);
if(fail)
return;
result = 1;
return;
}
if(G.inCoeffDomain())
{
tryInvert(G,M,result,fail);
if(fail)
return;
result = 1;
return;
}
TIMING_START (alg_compress)
CFMap MM,NN;
int lev= myCompress (F, G, MM, NN, topLevel);
if (lev == 0)
{
result= 1;
return;
}
CanonicalForm f=MM(F);
CanonicalForm g=MM(G);
TIMING_END_AND_PRINT (alg_compress, "time to compress in alg gcd: ")
// here: f,g are compressed
// compute largest variable in f or g (least one is Variable(1))
int mv = f.level();
if(g.level() > mv)
mv = g.level();
// here: mv is level of the largest variable in f, g
Variable v1= Variable (1);
#ifdef HAVE_NTL
Variable v= M.mvar();
int ch=getCharacteristic();
if (fac_NTL_char != ch)
{
fac_NTL_char= ch;
zz_p::init (ch);
}
zz_pX NTLMipo= convertFacCF2NTLzzpX (M);
zz_pE::init (NTLMipo);
zz_pEX NTLResult;
zz_pEX NTLF;
zz_pEX NTLG;
#endif
if(mv == 1) // f,g univariate
{
TIMING_START (alg_euclid_p)
#ifdef HAVE_NTL
NTLF= convertFacCF2NTLzz_pEX (f, NTLMipo);
NTLG= convertFacCF2NTLzz_pEX (g, NTLMipo);
tryNTLGCD (NTLResult, NTLF, NTLG, fail);
if (fail)
return;
result= convertNTLzz_pEX2CF (NTLResult, f.mvar(), v);
#else
tryEuclid(f,g,M,result,fail);
if(fail)
return;
#endif
TIMING_END_AND_PRINT (alg_euclid_p, "time for euclidean alg mod p: ")
result= NN (reduce (result, M)); // do not forget to map back
return;
}
TIMING_START (alg_content_p)
// here: mv > 1
CanonicalForm cf = tryvcontent(f, Variable(2), M, fail); // cf is univariate poly in var(1)
if(fail)
return;
CanonicalForm cg = tryvcontent(g, Variable(2), M, fail);
if(fail)
return;
CanonicalForm c;
#ifdef HAVE_NTL
NTLF= convertFacCF2NTLzz_pEX (cf, NTLMipo);
NTLG= convertFacCF2NTLzz_pEX (cg, NTLMipo);
tryNTLGCD (NTLResult, NTLF, NTLG, fail);
if (fail)
return;
c= convertNTLzz_pEX2CF (NTLResult, v1, v);
#else
tryEuclid(cf,cg,M,c,fail);
if(fail)
return;
#endif
// f /= cf
f.tryDiv (cf, M, fail);
if(fail)
return;
// g /= cg
g.tryDiv (cg, M, fail);
if(fail)
return;
TIMING_END_AND_PRINT (alg_content_p, "time for content in alg gcd mod p: ")
if(f.inCoeffDomain())
{
tryInvert(f,M,result,fail);
if(fail)
return;
result = NN(c);
return;
}
if(g.inCoeffDomain())
{
tryInvert(g,M,result,fail);
if(fail)
return;
result = NN(c);
return;
}
std::vector<int> L(mv + 1, 0);
std::vector<int> N(mv + 1, 0);
for(int i=2; i<=mv; i++)
L[i] = N[i] = 0;
leadDeg(f, L.data());
leadDeg(g, N.data());
CanonicalForm gamma;
TIMING_START (alg_euclid_p)
#ifdef HAVE_NTL
NTLF= convertFacCF2NTLzz_pEX (firstLC (f), NTLMipo);
NTLG= convertFacCF2NTLzz_pEX (firstLC (g), NTLMipo);
tryNTLGCD (NTLResult, NTLF, NTLG, fail);
if (fail)
return;
gamma= convertNTLzz_pEX2CF (NTLResult, v1, v);
#else
tryEuclid( firstLC(f), firstLC(g), M, gamma, fail );
if(fail)
return;
#endif
TIMING_END_AND_PRINT (alg_euclid_p, "time for gcd of lcs in alg mod p: ")
for(int i=2; i<=mv; i++) // entries at i=0,1 not visited
if(N[i] < L[i])
L[i] = N[i];
// L is now upper bound for degrees of gcd
std::vector<int> dg_im(mv+1); // for the degree vector of the image we don't need any entry at i=1
for(int i=2; i<=mv; i++)
dg_im[i] = 0; // initialize
CanonicalForm gamma_image, m=1;
CanonicalForm gm=0;
CanonicalForm g_image, alpha, gnew;
FFGenerator gen = FFGenerator();
Variable x= Variable (1);
bool divides= true;
for(FFGenerator gen = FFGenerator(); gen.hasItems(); gen.next())
{
alpha = gen.item();
gamma_image = reduce(gamma(alpha, x),M); // plug in alpha for var(1)
if(gamma_image.isZero()) // skip lc-bad points var(1)-alpha
continue;
TIMING_START (alg_recursion_p)
tryBrownGCD( f(alpha, x), g(alpha, x), M, g_image, fail, false ); // recursive call with one var less
TIMING_END_AND_PRINT (alg_recursion_p,
"time for recursive calls in alg gcd mod p: ")
if(fail)
return;
g_image = reduce(g_image, M);
if(g_image.inCoeffDomain()) // early termination
{
tryInvert(g_image,M,result,fail);
if(fail)
return;
result = NN(c);
return;
}
for(int i=2; i<=mv; i++)
dg_im[i] = 0; // reset (this is necessary, because some entries may not be updated by call to leadDeg)
leadDeg(g_image, dg_im.data());
if(isEqual(dg_im.data(), L.data(), 2, mv))
{
CanonicalForm inv;
tryInvert (firstLC (g_image), M, inv, fail);
if (fail)
return;
g_image *= inv;
g_image *= gamma_image; // multiply by multiple of image lc(gcd)
g_image= reduce (g_image, M);
TIMING_START (alg_newton_p)
gnew= tryNewtonInterp (alpha, g_image, m, gm, x, M, fail);
TIMING_END_AND_PRINT (alg_newton_p,
"time for Newton interpolation in alg gcd mod p: ")
// gnew = gm mod m
// gnew = g_image mod var(1)-alpha
// mnew = m * (var(1)-alpha)
if(fail)
return;
m *= (x - alpha);
if((firstLC(gnew) == gamma) || (gnew == gm)) // gnew did not change
{
TIMING_START (alg_termination_p)
cf = tryvcontent(gnew, Variable(2), M, fail);
if(fail)
return;
divides = true;
g_image= gnew;
g_image.tryDiv (cf, M, fail);
if(fail)
return;
divides= tryFdivides (g_image,f, M, fail); // trial division (f)
if(fail)
return;
if(divides)
{
bool divides2= tryFdivides (g_image,g, M, fail); // trial division (g)
if(fail)
return;
if(divides2)
{
result = NN(reduce (c*g_image, M));
TIMING_END_AND_PRINT (alg_termination_p,
"time for successful termination test in alg gcd mod p: ")
return;
}
}
TIMING_END_AND_PRINT (alg_termination_p,
"time for unsuccessful termination test in alg gcd mod p: ")
}
gm = gnew;
continue;
}
if(isLess(L.data(), dg_im.data(), 2, mv)) // dg_im > L --> current point unlucky
continue;
// here: isLess(dg_im, L, 2, mv) --> all previous points were unlucky
m = CanonicalForm(1); // reset
gm = 0; // reset
for(int i=2; i<=mv; i++) // tighten bound
L[i] = dg_im[i];
}
// we are out of evaluation points
fail = true;
}
convertFacCF2NTLzz_pEX
zz_pEX convertFacCF2NTLzz_pEX(const CanonicalForm &f, const zz_pX &mipo)
Definition NTLconvert.cc:1064
convertNTLzz_pEX2CF
CanonicalForm convertNTLzz_pEX2CF(const zz_pEX &f, const Variable &x, const Variable &alpha)
Definition NTLconvert.cc:1092
convertFacCF2NTLzzpX
zz_pX convertFacCF2NTLzzpX(const CanonicalForm &f)
Definition NTLconvert.cc:105
fac_NTL_char
VAR long fac_NTL_char
Definition NTLconvert.cc:46
reduce
CanonicalForm reduce(const CanonicalForm &f, const CanonicalForm &M)
polynomials in M.mvar() are considered coefficients M univariate monic polynomial the coefficients of...
Definition cf_ops.cc:660
level
int level(const CanonicalForm &f)
Definition canonicalform.h:331
getCharacteristic
int FACTORY_PUBLIC getCharacteristic()
Definition cf_char.cc:70
N
const CanonicalForm CFMap CFMap & N
Definition cfGcdAlgExt.cc:57
tryNewtonInterp
static CanonicalForm tryNewtonInterp(const CanonicalForm &alpha, const CanonicalForm &u, const CanonicalForm &newtonPoly, const CanonicalForm &oldInterPoly, const Variable &x, const CanonicalForm &M, bool &fail)
Definition cfGcdAlgExt.cc:358
topLevel
const CanonicalForm CFMap CFMap bool topLevel
Definition cfGcdAlgExt.cc:58
tryvcontent
static CanonicalForm tryvcontent(const CanonicalForm &f, const Variable &x, const CanonicalForm &M, bool &fail)
Definition cfGcdAlgExt.cc:1048
firstLC
CanonicalForm firstLC(const CanonicalForm &f)
Definition cfGcdAlgExt.cc:956
myCompress
int myCompress(const CanonicalForm &F, const CanonicalForm &G, CFMap &M, CFMap &N, bool topLevel)
compressing two polynomials F and G, M is used for compressing, N to reverse the compression
Definition cfModGcd.cc:92
x
Variable x
Definition cfModGcd.cc:4090
cf
CanonicalForm cf
Definition cfModGcd.cc:4091
cg
CanonicalForm cg
Definition cfModGcd.cc:4091
tryNTLGCD
void tryNTLGCD(zz_pEX &x, const zz_pEX &a, const zz_pEX &b, bool &fail)
compute the GCD x of a and b, fail is set to true if a zero divisor is encountered
Definition cfNTLzzpEXGCD.cc:242
tryFdivides
bool tryFdivides(const CanonicalForm &f, const CanonicalForm &g, const CanonicalForm &M, bool &fail)
CFMap
class CFMap
Definition cf_map.h:85
CanonicalForm::tryDiv
CanonicalForm & tryDiv(const CanonicalForm &, const CanonicalForm &, bool &)
same as divremt but handles zero divisors in case we are in Z_p[x]/(f) where f is not irreducible
Definition canonicalform.cc:905
FFGenerator
generate all elements in F_p starting from 0
Definition cf_generator.h:56
alpha
Variable alpha
Definition facAbsBiFact.cc:52
v
const Variable & v
< [in] a sqrfree bivariate poly
Definition facBivar.h:39
next
ListNode * next
Definition janet.h:31

## trycf_content()

CanonicalForm trycf_content | ( | const CanonicalForm & | f,
|  | const CanonicalForm & | g,
|  | const CanonicalForm & | M,
|  | bool & | fail )

static

Definition at line 1067 of file cfGcdAlgExt.cc.

{ // as cf_content, but takes care of zero divisors
if ( f.inPolyDomain() || ( f.inExtension() && ! getReduce( f.mvar() ) ) )
{
CFIterator i = f;
CanonicalForm tmp = g, result;
while ( i.hasTerms() && ! tmp.isOne() && ! fail )
{
tryBrownGCD( i.coeff(), tmp, M, result, fail );
tmp = result;
i++;
}
return result;
}
return abs( f );
}
tryBrownGCD
void tryBrownGCD(const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &M, CanonicalForm &result, bool &fail, bool topLevel)
modular gcd over F_p[x]/(M) for not necessarily irreducible M. If a zero divisor is encountered fail ...
Definition cfGcdAlgExt.cc:387
getReduce
bool getReduce(const Variable &alpha)
Definition variable.cc:232

## trycontent()

CanonicalForm trycontent | ( | const CanonicalForm & | f,
|  | const Variable & | x,
|  | const CanonicalForm & | M,
|  | bool & | fail )

static

Definition at line 1036 of file cfGcdAlgExt.cc.

{ // as 'content', but takes care of zero divisors
ASSERT( x.level() > 0, "cannot calculate content with respect to algebraic variable" );
Variable y = f.mvar();
if ( y == x )
return trycf_content( f, 0, M, fail );
if ( y < x )
return f;
return swapvar( trycontent( swapvar( f, y, x ), y, M, fail ), y, x );
}
swapvar
CanonicalForm FACTORY_PUBLIC swapvar(const CanonicalForm &, const Variable &, const Variable &)
swapvar() - swap variables x1 and x2 in f.
Definition cf_ops.cc:168
trycontent
static CanonicalForm trycontent(const CanonicalForm &f, const Variable &x, const CanonicalForm &M, bool &fail)
Definition cfGcdAlgExt.cc:1036
trycf_content
static CanonicalForm trycf_content(const CanonicalForm &f, const CanonicalForm &g, const CanonicalForm &M, bool &fail)
Definition cfGcdAlgExt.cc:1067
ASSERT
#define ASSERT(expression, message)
Definition cf_assert.h:99
y
const CanonicalForm int const CFList const Variable & y
Definition facAbsFact.cc:53

## tryInvert()

void tryInvert | ( | const CanonicalForm & | F,
|  | const CanonicalForm & | M,
|  | CanonicalForm & | inv,
|  | bool & | fail )

Definition at line 222 of file cfGcdAlgExt.cc.

{ // F, M are required to be "univariate" polynomials in an algebraic variable
// we try to invert F modulo M
if(F.inBaseDomain())
{
if(F.isZero())
{
fail = true;
return;
}
inv = 1/F;
return;
}
CanonicalForm b;
Variable a = M.mvar();
Variable x = Variable(1);
if(!extgcd( replacevar( F, a, x ), replacevar( M, a, x ), inv, b ).isOne())
fail = true;
else
inv = replacevar( inv, x, a ); // change back to alg var
}
replacevar
CanonicalForm FACTORY_PUBLIC replacevar(const CanonicalForm &, const Variable &, const Variable &)
CanonicalForm replacevar ( const CanonicalForm & f, const Variable & x1, const Variable & x2 ).
Definition cf_ops.cc:271
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

## tryNewtonInterp()

CanonicalForm tryNewtonInterp | ( | const CanonicalForm & | alpha,
|  | const CanonicalForm & | u,
|  | const CanonicalForm & | newtonPoly,
|  | const CanonicalForm & | oldInterPoly,
|  | const Variable & | x,
|  | const CanonicalForm & | M,
|  | bool & | fail )

inline static

Definition at line 358 of file cfGcdAlgExt.cc.

{
CanonicalForm interPoly;
CanonicalForm inv;
tryInvert (newtonPoly (alpha, x), M, inv, fail);
if (fail)
return 0;
interPoly= oldInterPoly+reduce ((u - oldInterPoly (alpha, x))*inv*newtonPoly, M);
return interPoly;
}

## tryvcontent()

CanonicalForm tryvcontent | ( | const CanonicalForm & | f,
|  | const Variable & | x,
|  | const CanonicalForm & | M,
|  | bool & | fail )

static

Definition at line 1048 of file cfGcdAlgExt.cc.

{ // as vcontent, but takes care of zero divisors
ASSERT( x.level() > 0, "cannot calculate vcontent with respect to algebraic variable" );
if ( f.mvar() <= x )
return trycontent( f, x, M, fail );
CFIterator i;
CanonicalForm d = 0, e, ret;
for ( i = f; i.hasTerms() && ! d.isOne() && ! fail; i++ )
{
e = tryvcontent( i.coeff(), x, M, fail );
if(fail)
break;
tryBrownGCD( d, e, M, ret, fail );
d = ret;
}
return d;
}

## Variable Documentation

## both_non_zero

int both_non_zero = 0

Definition at line 69 of file cfGcdAlgExt.cc.

## both_zero

int both_zero = 0

Definition at line 72 of file cfGcdAlgExt.cc.

## degsf

degsf = NEW_ARRAY(int,n + 1)

Definition at line 60 of file cfGcdAlgExt.cc.

## degsg

degsg = NEW_ARRAY(int,n + 1)

Definition at line 61 of file cfGcdAlgExt.cc.

## else

else

Initial value:
{
for (int i= 1; i <= n; i++)
{
if (degsf[i] == 0 && degsg[i] == 0)
{
both_zero++;
continue;
}
else
{
if (both_zero != 0)
{
M.newpair (Variable (i), Variable (i - both_zero));
N.newpair (Variable (i - both_zero), Variable (i));
}
}
}
}
DELETE_ARRAY(degsf)
both_zero
int both_zero
Definition cfGcdAlgExt.cc:72

Definition at line 195 of file cfGcdAlgExt.cc.

{
//arrange Variables such that no gaps occur
for (int i= 1; i <= n; i++)
{
if (degsf[i] == 0 && degsg[i] == 0)
{
both_zero++;
continue;
}
else
{
if (both_zero != 0)
{
M.newpair (Variable (i), Variable (i - both_zero));
N.newpair (Variable (i - both_zero), Variable (i));
}
}
}
}

## f_zero

int f_zero = 0

Definition at line 70 of file cfGcdAlgExt.cc.

## Flevel

int Flevel =F.level()

Definition at line 73 of file cfGcdAlgExt.cc.

## G

const CanonicalForm& G

Definition at line 56 of file cfGcdAlgExt.cc.

## g_zero

int g_zero = 0

Definition at line 71 of file cfGcdAlgExt.cc.

## Glevel

int Glevel =G.level()

Definition at line 74 of file cfGcdAlgExt.cc.

## M

const CanonicalForm CFMap& M

Definition at line 56 of file cfGcdAlgExt.cc.

## N

const CanonicalForm CFMap CFMap& N

Definition at line 57 of file cfGcdAlgExt.cc.

## return

return

Definition at line 219 of file cfGcdAlgExt.cc.

## topLevel

const CanonicalForm CFMap CFMap bool topLevel

Initial value:
{
int n= tmax (F.level(), G.level())

Definition at line 57 of file cfGcdAlgExt.cc.