Loading...
Searching...
No Matches
cfModGcd.cc File Reference

This file implements the GCD of two polynomials over , , GF or Z based on Alg. More...

`#include "config.h"`
`#include <math.h>`
`#include "cf_assert.h"`
`#include "debug.h"`
`#include "timing.h"`
`#include "canonicalform.h"`
`#include "cfGcdUtil.h"`
`#include "cf_map.h"`
`#include "cf_util.h"`
`#include "cf_irred.h"`
`#include "templates/ftmpl_functions.h"`
`#include "cf_random.h"`
`#include "cf_reval.h"`
`#include "facHensel.h"`
`#include "cf_iter.h"`
`#include "cfNewtonPolygon.h"`
`#include "cf_algorithm.h"`
`#include "cf_primes.h"`
`#include "cf_map_ext.h"`
`#include "NTLconvert.h"`
`#include "FLINTconvert.h"`
`#include "cfModGcd.h"`

Go to the source code of this file.

## Functions

| TIMING_DEFINE_PRINT (gcd_recursion) TIMING_DEFINE_PRINT(newton_interpolation) TIMING_DEFINE_PRINT(termination_test) TIMING_DEFINE_PRINT(ez_p_compress) TIMING_DEFINE_PRINT(ez_p_hensel_lift) TIMING_DEFINE_PRINT(ez_p_eval) TIMING_DEFINE_PRINT(ez_p_content) bool terminationTest(const CanonicalForm &F
| if (LCCand *abs(LC(coF))==abs(LC(F)))
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
static CanonicalForm | uni_content (const CanonicalForm &F)
| compute the content of F, where F is considered as an element of
CanonicalForm | uni_content (const CanonicalForm &F, const Variable &x)
CanonicalForm | extractContents (const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &contentF, CanonicalForm &contentG, CanonicalForm &ppF, CanonicalForm &ppG, const int d)
static CanonicalForm | uni_lcoeff (const CanonicalForm &F)
| compute the leading coefficient of F, where F is considered as an element of , order on is dp.
static CanonicalForm | newtonInterp (const CanonicalForm &alpha, const CanonicalForm &u, const CanonicalForm &newtonPoly, const CanonicalForm &oldInterPoly, const Variable &x)

Newton interpolation - Incremental algorithm. Given a list of values alpha_i and a list of polynomials u_i, 1 <= i <= n, computes the interpolation polynomial assuming that the polynomials in u are the results of evaluating the variabe x of the unknown polynomial at the alpha values. This incremental version receives only the values of alpha_n and u_n and the previous interpolation polynomial for points 1 <= i <= n-1, and computes the polynomial interpolating in all the points. newtonPoly must be equal to (x - alpha_1) * ... * (x - alpha_{n-1}).

static CanonicalForm | randomElement (const CanonicalForm &F, const Variable &alpha, CFList &list, bool &fail)
| compute a random element a of , s.t. F(a) , F is a univariate polynomial, returns fail if there are no field elements left which have not been used before
static Variable | chooseExtension (const Variable &alpha)
CanonicalForm | modGCDFq (const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &coF, CanonicalForm &coG, Variable &alpha, CFList &l, bool &topLevel)
| GCD of F and G over , l and topLevel are only used internally, output is monic based on Alg. 7.2. as described in "Algorithms for Computer Algebra" by Geddes, Czapor, Labahn.
CanonicalForm | modGCDFq (const CanonicalForm &F, const CanonicalForm &G, Variable &alpha, CFList &l, bool &topLevel)
static CanonicalForm | GFRandomElement (const CanonicalForm &F, CFList &list, bool &fail)
| compute a random element a of GF, s.t. F(a) , F is a univariate polynomial, returns fail if there are no field elements left which have not been used before
CanonicalForm | modGCDGF (const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &coF, CanonicalForm &coG, CFList &l, bool &topLevel)
| GCD of F and G over GF, based on Alg. 7.2. as described in "Algorithms for Computer Algebra" by Geddes, Czapor, Labahn Usually this algorithm will be faster than modGCDFq since GF has faster field arithmetics, however it might fail if the input is large since the size of the base field is bounded by 2^16, output is monic.
CanonicalForm | modGCDGF (const CanonicalForm &F, const CanonicalForm &G, CFList &l, bool &topLevel)
static CanonicalForm | FpRandomElement (const CanonicalForm &F, CFList &list, bool &fail)
CanonicalForm | modGCDFp (const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &coF, CanonicalForm &coG, bool &topLevel, CFList &l)
CanonicalForm | modGCDFp (const CanonicalForm &F, const CanonicalForm &G, bool &topLevel, CFList &l)
CFArray | solveVandermonde (const CFArray &M, const CFArray &A)
CFArray | solveGeneralVandermonde (const CFArray &M, const CFArray &A)
CFArray | readOffSolution (const CFMatrix &M, const long rk)
| M in row echolon form, rk rank of M.
CFArray | readOffSolution (const CFMatrix &M, const CFArray &L, const CFArray &partialSol)
long | gaussianElimFp (CFMatrix &M, CFArray &L)
long | gaussianElimFq (CFMatrix &M, CFArray &L, const Variable &alpha)
CFArray | solveSystemFp (const CFMatrix &M, const CFArray &L)
CFArray | solveSystemFq (const CFMatrix &M, const CFArray &L, const Variable &alpha)
CFArray | getMonoms (const CanonicalForm &F)
| extract monomials of F, parts in algebraic variable are considered coefficients
CFArray | evaluateMonom (const CanonicalForm &F, const CFList &evalPoints)
CFArray | evaluate (const CFArray &A, const CFList &evalPoints)
CFList | evaluationPoints (const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &Feval, CanonicalForm &Geval, const CanonicalForm &LCF, const bool &GF, const Variable &alpha, bool &fail, CFList &list)
void | mult (CFList &L1, const CFList &L2)
| multiply two lists componentwise
void | eval (const CanonicalForm &A, const CanonicalForm &B, CanonicalForm &Aeval, CanonicalForm &Beval, const CFList &L)
CanonicalForm | monicSparseInterpol (const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &skeleton, const Variable &alpha, bool &fail, CFArray *&coeffMonoms, CFArray &Monoms)
CanonicalForm | nonMonicSparseInterpol (const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &skeleton, const Variable &alpha, bool &fail, CFArray *&coeffMonoms, CFArray &Monoms)
CanonicalForm | sparseGCDFq (const CanonicalForm &F, const CanonicalForm &G, const Variable &alpha, CFList &l, bool &topLevel)
CanonicalForm | sparseGCDFp (const CanonicalForm &F, const CanonicalForm &G, bool &topLevel, CFList &l)
| TIMING_DEFINE_PRINT (modZ_termination) TIMING_DEFINE_PRINT(modZ_recursion) CanonicalForm modGCDZ(const CanonicalForm &FF
| modular gcd algorithm, see Keith, Czapor, Geddes "Algorithms for Computer Algebra", Algorithm 7.1
| for (i=tmax(f.level(), g.level());i > 0;i--)
| if (i==0) return gcdcfcg
| for (;i > 0;i--)
| while (true)

## Variables

const CanonicalForm & | G
const CanonicalForm const CanonicalForm & | coF
const CanonicalForm const CanonicalForm const CanonicalForm & | coG
const CanonicalForm const CanonicalForm const CanonicalForm const CanonicalForm & | cand
return | false
const CanonicalForm & | GG
int | p
int | i = cf_getNumBigPrimes() - 1
int | dp_deg
int | d_deg =-1
CanonicalForm | cd (bCommonDen(FF)) = bCommonDen( GG )
| f =cd*FF
Variable | x = Variable (1)
CanonicalForm | cf = icontent (f)
CanonicalForm | cg = icontent (g)
| g =cd*GG
CanonicalForm | Dn
CanonicalForm | test = 0
CanonicalForm | lcf = f.lc()
CanonicalForm | lcg = g.lc()
| cl = gcd (f.lc(),g.lc())
CanonicalForm | gcdcfcg = gcd (cf, cg)
CanonicalForm | fp
CanonicalForm | gp
CanonicalForm | b = 1
int | minCommonDeg = 0
bool | equal = false
CanonicalForm | cof
CanonicalForm | cog
CanonicalForm | cofp
CanonicalForm | cogp
CanonicalForm | newCof
CanonicalForm | newCog
CanonicalForm | cofn
CanonicalForm | cogn
CanonicalForm | cDn
int | maxNumVars = tmax (getNumVars (f), getNumVars (g))

## Detailed Description

This file implements the GCD of two polynomials over , , GF or Z based on Alg.

7.1. and 7.2. as described in "Algorithms for Computer Algebra" by Geddes, Czapor, Labahn via modular computations. And sparse modular variants as described in Zippel "Effective Polynomial Computation", deKleine, Monagan, Wittkopf "Algorithms for the non-monic case of the sparse modular GCD algorithm" and Javadi "A new solution to the normalization problem"

Author
  Martin Lee

Date
  22.10.2009

Copyright:
  (c) by The SINGULAR Team, see LICENSE file

Definition in file cfModGcd.cc.

## Function Documentation

## chooseExtension()

Variable chooseExtension | ( | const Variable & | alpha | ) |

inline static

Definition at line 421 of file cfModGcd.cc.

{
int i, m;
// extension of F_p needed
if (alpha.level() == 1)
{
i= 1;
m= 2;
} //extension of F_p(alpha)
if (alpha.level() != 1)
{
i= 4;
m= degree (getMipo (alpha));
}
#ifdef HAVE_FLINT
nmod_poly_t Irredpoly;
nmod_poly_init(Irredpoly,getCharacteristic());
nmod_poly_randtest_monic_irreducible(Irredpoly, FLINTrandom, i*m+1);
CanonicalForm newMipo=convertnmod_poly_t2FacCF(Irredpoly,Variable(1));
nmod_poly_clear(Irredpoly);
#else
if (fac_NTL_char != getCharacteristic())
{
fac_NTL_char= getCharacteristic();
zz_p::init (getCharacteristic());
}
zz_pX NTLIrredpoly;
BuildIrred (NTLIrredpoly, i*m);
CanonicalForm newMipo= convertNTLzzpX2CF (NTLIrredpoly, Variable (1));
#endif
return rootOf (newMipo);
}
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
convertNTLzzpX2CF
CanonicalForm convertNTLzzpX2CF(const zz_pX &poly, const Variable &x)
Definition NTLconvert.cc:255
fac_NTL_char
VAR long fac_NTL_char
Definition NTLconvert.cc:46
degree
int degree(const CanonicalForm &f)
Definition canonicalform.h:316
getCharacteristic
int FACTORY_PUBLIC getCharacteristic()
Definition cf_char.cc:70
m
int m
Definition cfEzgcd.cc:128
i
int i
Definition cfEzgcd.cc:132
FLINTrandom
GLOBAL_VAR flint_rand_t FLINTrandom
Definition cf_random.cc:23
CanonicalForm
factory's main class
Definition canonicalform.h:86
Variable
factory's class for variables
Definition factory.h:127
alpha
Variable alpha
Definition facAbsBiFact.cc:52
nmod_poly_init
nmod_poly_init(FLINTmipo, getCharacteristic())
nmod_poly_clear
nmod_poly_clear(FLINTmipo)
rootOf
Variable FACTORY_PUBLIC rootOf(const CanonicalForm &, char name='@')
returns a symbolic root of polynomial with name name Use it to define algebraic variables
Definition variable.cc:162
getMipo
CanonicalForm getMipo(const Variable &alpha, const Variable &x)
Definition variable.cc:207

## eval()

void eval | ( | const CanonicalForm & | A,
|  | const CanonicalForm & | B,
|  | CanonicalForm & | Aeval,
|  | CanonicalForm & | Beval,
|  | const CFList & | L )

Definition at line 2193 of file cfModGcd.cc.

{
Aeval= A;
Beval= B;
int j= 1;
for (CFListIterator i= L; i.hasItem(); i++, j++)
{
Aeval= Aeval (i.getItem(), j);
Beval= Beval (i.getItem(), j);
}
}
CFListIterator
ListIterator< CanonicalForm > CFListIterator
Definition canonicalform.h:396
B
b *CanonicalForm B
Definition facBivar.cc:52
Aeval
CFList *& Aeval
<[in] poly
Definition facFactorize.h:31
j
int j
Definition facHensel.cc:110
A
#define A
Definition sirandom.c:24

## evaluate()

CFArray evaluate | ( | const CFArray & | A,
|  | const CFList & | evalPoints )

Definition at line 2039 of file cfModGcd.cc.

{
CFArray result= A.size();
CanonicalForm tmp;
int k;
for (int i= 0; i < A.size(); i++)
{
tmp= A[i];
k= 1;
for (CFListIterator j= evalPoints; j.hasItem(); j++, k++)
tmp= tmp (j.getItem(), k);
result[i]= tmp;
}
return result;
}
CFArray
Array< CanonicalForm > CFArray
Definition canonicalform.h:397
k
int k
Definition cfEzgcd.cc:99
result
return result
Definition facAbsBiFact.cc:76
evalPoints
CFList evalPoints(const CanonicalForm &F, CFList &eval, const Variable &alpha, CFList &list, const bool &GF, bool &fail)
evaluation point search for multivariate factorization, looks for a (F.level() - 1)-tuple such that t...
Definition facFqFactorize.cc:751

## evaluateMonom()

CFArray evaluateMonom | ( | const CanonicalForm & | F,
|  | const CFList & | evalPoints )

Definition at line 2000 of file cfModGcd.cc.

{
if (F.inCoeffDomain())
{
CFArray result= CFArray (1);
result [0]= F;
return result;
}
if (F.isUnivariate())
{
ASSERT (evalPoints.length() == 1,
"expected an eval point with only one component");
CFArray result= CFArray (size(F));
int j= 0;
CanonicalForm evalPoint= evalPoints.getLast();
for (CFIterator i= F; i.hasTerms(); i++, j++)
result[j]= power (evalPoint, i.exp());
return result;
}
int numMon= size (F);
CFArray result= CFArray (numMon);
int j= 0;
CanonicalForm evalPoint= evalPoints.getLast();
CFList buf= evalPoints;
buf.removeLast();
CFArray recResult;
CanonicalForm powEvalPoint;
for (CFIterator i= F; i.hasTerms(); i++)
{
powEvalPoint= power (evalPoint, i.exp());
recResult= evaluateMonom (i.coeff(), buf);
for (int k= 0; k < recResult.size(); k++)
result[j+k]= powEvalPoint*recResult[k];
j += recResult.size();
}
return result;
}
power
CanonicalForm power(const CanonicalForm &f, int n)
exponentiation
Definition canonicalform.cc:1896
size
int size(const CanonicalForm &f, const Variable &v)
int size ( const CanonicalForm & f, const Variable & v )
Definition cf_ops.cc:600
CFList
List< CanonicalForm > CFList
Definition canonicalform.h:395
evaluateMonom
CFArray evaluateMonom(const CanonicalForm &F, const CFList &evalPoints)
Definition cfModGcd.cc:2000
evalPoint
static void evalPoint(const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &FEval, CanonicalForm &GEval, CFGenerator &evalPoint)
Definition cfModResultant.cc:311
ASSERT
#define ASSERT(expression, message)
Definition cf_assert.h:99
Array::size
int size() const
Definition ftmpl_array.cc:92
CFIterator
class to iterate through CanonicalForm's
Definition cf_iter.h:44
CanonicalForm::inCoeffDomain
bool inCoeffDomain() const
Definition canonicalform.cc:122
CanonicalForm::isUnivariate
bool isUnivariate() const
Definition canonicalform.cc:155
buf
int status int void * buf
Definition si_signals.h:69

## evaluationPoints()

CFList evaluationPoints | ( | const CanonicalForm & | F,
|  | const CanonicalForm & | G,
|  | CanonicalForm & | Feval,
|  | CanonicalForm & | Geval,
|  | const CanonicalForm & | LCF,
|  | const bool & | GF,
|  | const Variable & | alpha,
|  | bool & | fail,
|  | CFList & | list )

Definition at line 2056 of file cfModGcd.cc.

{
int k= tmax (F.level(), G.level()) - 1;
Variable x= Variable (1);
CFList result;
FFRandom genFF;
GFRandom genGF;
int p= getCharacteristic ();
double bound;
if (alpha != Variable (1))
{
bound= pow ((double) p, (double) degree (getMipo(alpha)));
bound= pow (bound, (double) k);
}
else if (GF)
{
bound= pow ((double) p, (double) getGFDegree());
bound= pow ((double) bound, (double) k);
}
else
bound= pow ((double) p, (double) k);
CanonicalForm random;
int j;
bool zeroOneOccured= false;
bool allEqual= false;
CanonicalForm buf;
do
{
random= 0;
// possible overflow if list.length() does not fit into a int
if (list.length() >= bound)
{
fail= true;
break;
}
for (int i= 0; i < k; i++)
{
if (GF)
{
result.append (genGF.generate());
random += result.getLast()*power (x, i);
}
else if (alpha.level() != 1)
{
AlgExtRandomF genAlgExt (alpha);
result.append (genAlgExt.generate());
random += result.getLast()*power (x, i);
}
else
{
result.append (genFF.generate());
random += result.getLast()*power (x, i);
}
if (result.getLast().isOne() || result.getLast().isZero())
zeroOneOccured= true;
}
if (find (list, random))
{
zeroOneOccured= false;
allEqual= false;
result= CFList();
continue;
}
if (zeroOneOccured)
{
list.append (random);
zeroOneOccured= false;
allEqual= false;
result= CFList();
continue;
}
// no zero at this point
if (k > 1)
{
allEqual= true;
CFIterator iter= random;
buf= iter.coeff();
iter++;
for (; iter.hasTerms(); iter++)
if (buf != iter.coeff())
allEqual= false;
}
if (allEqual)
{
list.append (random);
allEqual= false;
zeroOneOccured= false;
result= CFList();
continue;
}
Feval= F;
Geval= G;
CanonicalForm LCeval= LCF;
j= 1;
for (CFListIterator i= result; i.hasItem(); i++, j++)
{
Feval= Feval (i.getItem(), j);
Geval= Geval (i.getItem(), j);
LCeval= LCeval (i.getItem(), j);
}
if (LCeval.isZero())
{
if (!find (list, random))
list.append (random);
zeroOneOccured= false;
allEqual= false;
result= CFList();
continue;
}
if (list.length() >= bound)
{
fail= true;
break;
}
} while (find (list, random));
return result;
}
pow
Rational pow(const Rational &a, int e)
Definition GMPrat.cc:411
getGFDegree
int getGFDegree()
Definition cf_char.cc:75
x
Variable x
Definition cfModGcd.cc:4090
p
int p
Definition cfModGcd.cc:4086
bound
static CanonicalForm bound(const CFMatrix &M)
Definition cf_linsys.cc:460
AlgExtRandomF
generate random elements in F_p(alpha)
Definition cf_random.h:70
CanonicalForm::isZero
CF_NO_INLINE bool isZero() const
CanonicalForm::level
int level() const
level() returns the level of CO.
Definition canonicalform.cc:576
FFRandom
generate random elements in F_p
Definition cf_random.h:44
FFRandom::generate
CanonicalForm generate() const
Definition cf_random.cc:65
GFRandom
generate random elements in GF
Definition cf_random.h:32
GFRandom::generate
CanonicalForm generate() const
Definition cf_random.cc:75
List::length
int length() const
Definition ftmpl_list.cc:273
List::append
void append(const T &)
Definition ftmpl_list.cc:256
iter
CFFListIterator iter
Definition facAbsBiFact.cc:54
Feval
CanonicalForm Feval
Definition facAbsFact.cc:60
LCF
CanonicalForm LCF
Definition facFactorize.cc:52
tmax
template CanonicalForm tmax(const CanonicalForm &, const CanonicalForm &)
find
template bool find(const List< CanonicalForm > &, const CanonicalForm &)
G
STATIC_VAR TreeM * G
Definition janet.cc:31

## extractContents()

CanonicalForm extractContents | ( | const CanonicalForm & | F,
|  | const CanonicalForm & | G,
|  | CanonicalForm & | contentF,
|  | CanonicalForm & | contentG,
|  | CanonicalForm & | ppF,
|  | CanonicalForm & | ppG,
|  | const int | d )

Definition at line 312 of file cfModGcd.cc.

{
CanonicalForm uniContentF, uniContentG, gcdcFcG;
contentF= 1;
contentG= 1;
ppF= F;
ppG= G;
CanonicalForm result= 1;
for (int i= 1; i <= d; i++)
{
uniContentF= uni_content (F, Variable (i));
uniContentG= uni_content (G, Variable (i));
gcdcFcG= gcd (uniContentF, uniContentG);
contentF *= uniContentF;
contentG *= uniContentG;
ppF /= uniContentF;
ppG /= uniContentG;
result *= gcdcFcG;
}
return result;
}
uni_content
static CanonicalForm uni_content(const CanonicalForm &F)
compute the content of F, where F is considered as an element of
Definition cfModGcd.cc:282
gcd
int gcd(int a, int b)
Definition walkSupport.cc:836

## for() [1/2]

for | ( | ; | i,
|  | 0;i-- | )

Definition at line 4125 of file cfModGcd.cc.

{
if (degree (f, i) <= 0 || degree (g, i) <= 0)
continue;
else
minCommonDeg= tmin (minCommonDeg, tmin (degree (g, i), degree (f, i)));
}
g
g
Definition cfModGcd.cc:4098
minCommonDeg
int minCommonDeg
Definition cfModGcd.cc:4112
f
FILE * f
Definition checklibs.c:9
tmin
template CanonicalForm tmin(const CanonicalForm &, const CanonicalForm &)

## for() [2/2]

for | ( | i | = tmax (f.level(), g.level()); i,
|  | 0;i-- | )

Definition at line 4113 of file cfModGcd.cc.

{
if (degree (f, i) <= 0 || degree (g, i) <= 0)
continue;
else
{
minCommonDeg= tmin (degree (g, i), degree (f, i));
break;
}
}

## FpRandomElement()

CanonicalForm FpRandomElement | ( | const CanonicalForm & | F,
|  | CFList & | list,
|  | bool & | fail )

inline static

Definition at line 1172 of file cfModGcd.cc.

{
fail= false;
Variable x= F.mvar();
FFRandom genFF;
CanonicalForm random;
int p= getCharacteristic();
int bound= p;
do
{
if (list.length() == bound)
{
fail= true;
break;
}
if (list.length() < 1)
random= 0;
else
{
random= genFF.generate();
while (find (list, random))
random= genFF.generate();
}
if (F (random, x) == 0)
{
list.append (random);
continue;
}
} while (find (list, random));
return random;
}
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

## gaussianElimFp()

long gaussianElimFp | ( | CFMatrix & | M,
|  | CFArray & | L )

Definition at line 1740 of file cfModGcd.cc.

{
ASSERT (L.size() <= M.rows(), "dimension exceeded");
CFMatrix *N;
N= new CFMatrix (M.rows(), M.columns() + 1);
for (int i= 1; i <= M.rows(); i++)
for (int j= 1; j <= M.columns(); j++)
(*N) (i, j)= M (i, j);
int j= 1;
for (int i= 0; i < L.size(); i++, j++)
(*N) (j, M.columns() + 1)= L[i];
#ifdef HAVE_FLINT
nmod_mat_t FLINTN;
convertFacCFMatrix2nmod_mat_t (FLINTN, *N);
long rk= nmod_mat_rref (FLINTN);
delete N;
N= convertNmod_mat_t2FacCFMatrix (FLINTN);
nmod_mat_clear (FLINTN);
#else
int p= getCharacteristic ();
if (fac_NTL_char != p)
{
fac_NTL_char= p;
zz_p::init (p);
}
mat_zz_p *NTLN= convertFacCFMatrix2NTLmat_zz_p(*N);
delete N;
long rk= gauss (*NTLN);
N= convertNTLmat_zz_p2FacCFMatrix (*NTLN);
delete NTLN;
#endif
L= CFArray (M.rows());
for (int i= 0; i < M.rows(); i++)
L[i]= (*N) (i + 1, M.columns() + 1);
M= (*N) (1, M.rows(), 1, M.columns());
delete N;
return rk;
}
convertNmod_mat_t2FacCFMatrix
CFMatrix * convertNmod_mat_t2FacCFMatrix(const nmod_mat_t m)
conversion of a FLINT matrix over Z/p to a factory matrix
Definition FLINTconvert.cc:630
convertFacCFMatrix2nmod_mat_t
void convertFacCFMatrix2nmod_mat_t(nmod_mat_t M, const CFMatrix &m)
conversion of a factory matrix over Z/p to a nmod_mat_t
Definition FLINTconvert.cc:612
convertNTLmat_zz_p2FacCFMatrix
CFMatrix * convertNTLmat_zz_p2FacCFMatrix(const mat_zz_p &m)
Definition NTLconvert.cc:1183
convertFacCFMatrix2NTLmat_zz_p
mat_zz_p * convertFacCFMatrix2NTLmat_zz_p(const CFMatrix &m)
Definition NTLconvert.cc:1167
CFMatrix
Matrix< CanonicalForm > CFMatrix
Definition canonicalform.h:398
N
const CanonicalForm CFMap CFMap & N
Definition cfEzgcd.cc:56
M
#define M
Definition sirandom.c:25

## gaussianElimFq()

long gaussianElimFq | ( | CFMatrix & | M,
|  | CFArray & | L,
|  | const Variable & | alpha )

Definition at line 1785 of file cfModGcd.cc.

{
ASSERT (L.size() <= M.rows(), "dimension exceeded");
CFMatrix *N;
N= new CFMatrix (M.rows(), M.columns() + 1);
for (int i= 1; i <= M.rows(); i++)
for (int j= 1; j <= M.columns(); j++)
(*N) (i, j)= M (i, j);
int j= 1;
for (int i= 0; i < L.size(); i++, j++)
(*N) (j, M.columns() + 1)= L[i];
#ifdef HAVE_FLINT
// convert mipo
nmod_poly_t mipo1;
convertFacCF2nmod_poly_t(mipo1,getMipo(alpha));
fq_nmod_ctx_t ctx;
fq_nmod_ctx_init_modulus(ctx,mipo1,"t");
nmod_poly_clear(mipo1);
// convert matrix
fq_nmod_mat_t FLINTN;
convertFacCFMatrix2Fq_nmod_mat_t (FLINTN, ctx, *N);
// rank
#if __FLINT_RELEASE >= 30100
long rk= fq_nmod_mat_rref (FLINTN,FLINTN,ctx);
#else
long rk= fq_nmod_mat_rref (FLINTN,ctx);
#endif
// clean up
fq_nmod_mat_clear (FLINTN,ctx);
fq_nmod_ctx_clear(ctx);
#elif defined(HAVE_NTL)
int p= getCharacteristic ();
if (fac_NTL_char != p)
{
fac_NTL_char= p;
zz_p::init (p);
}
zz_pX NTLMipo= convertFacCF2NTLzzpX (getMipo (alpha));
zz_pE::init (NTLMipo);
mat_zz_pE *NTLN= convertFacCFMatrix2NTLmat_zz_pE(*N);
long rk= gauss (*NTLN);
N= convertNTLmat_zz_pE2FacCFMatrix (*NTLN, alpha);
delete NTLN;
#else
factoryError("NTL/FLINT missing: gaussianElimFq");
#endif
M= (*N) (1, M.rows(), 1, M.columns());
L= CFArray (M.rows());
for (int i= 0; i < M.rows(); i++)
L[i]= (*N) (i + 1, M.columns() + 1);
delete N;
return rk;
}
convertFacCFMatrix2Fq_nmod_mat_t
void convertFacCFMatrix2Fq_nmod_mat_t(fq_nmod_mat_t M, const fq_nmod_ctx_t fq_con, const CFMatrix &m)
conversion of a factory matrix over F_q to a fq_nmod_mat_t
Definition FLINTconvert.cc:646
convertNTLmat_zz_pE2FacCFMatrix
CFMatrix * convertNTLmat_zz_pE2FacCFMatrix(const mat_zz_pE &m, const Variable &alpha)
Definition NTLconvert.cc:1212
convertFacCFMatrix2NTLmat_zz_pE
mat_zz_pE * convertFacCFMatrix2NTLmat_zz_pE(const CFMatrix &m)
Definition NTLconvert.cc:1196
convertFacCF2NTLzzpX
zz_pX convertFacCF2NTLzzpX(const CanonicalForm &f)
Definition NTLconvert.cc:105
factoryError
VAR void(* factoryError)(const char *s)
Definition cf_util.cc:80
fq_nmod_ctx_clear
fq_nmod_ctx_clear(fq_con)
fq_nmod_ctx_init_modulus
fq_nmod_ctx_init_modulus(fq_con, FLINTmipo, "Z")
convertFacCF2nmod_poly_t
convertFacCF2nmod_poly_t(FLINTmipo, M)

## getMonoms()

CFArray getMonoms | ( | const CanonicalForm & | F | ) |

extract monomials of F, parts in algebraic variable are considered coefficients

Parameters

  [in] | F | some poly

Definition at line 1965 of file cfModGcd.cc.

{
if (F.inCoeffDomain())
{
CFArray result= CFArray (1);
result [0]= 1;
return result;
}
if (F.isUnivariate())
{
CFArray result= CFArray (size(F));
int j= 0;
for (CFIterator i= F; i.hasTerms(); i++, j++)
result[j]= power (F.mvar(), i.exp());
return result;
}
int numMon= size (F);
CFArray result= CFArray (numMon);
int j= 0;
CFArray recResult;
Variable x= F.mvar();
CanonicalForm powX;
for (CFIterator i= F; i.hasTerms(); i++)
{
powX= power (x, i.exp());
recResult= getMonoms (i.coeff());
for (int k= 0; k < recResult.size(); k++)
result[j+k]= powX*recResult[k];
j += recResult.size();
}
return result;
}
getMonoms
CFArray getMonoms(const CanonicalForm &F)
extract monomials of F, parts in algebraic variable are considered coefficients
Definition cfModGcd.cc:1965

## GFRandomElement()

CanonicalForm GFRandomElement | ( | const CanonicalForm & | F,
|  | CFList & | list,
|  | bool & | fail )

inline static

compute a random element a of GF, s.t. F(a) , F is a univariate polynomial, returns fail if there are no field elements left which have not been used before

Definition at line 820 of file cfModGcd.cc.

{
fail= false;
Variable x= F.mvar();
GFRandom genGF;
CanonicalForm random;
int p= getCharacteristic();
int d= getGFDegree();
int bound= ipower (p, d);
do
{
if (list.length() == bound)
{
fail= true;
break;
}
if (list.length() < 1)
random= 0;
else
{
random= genGF.generate();
while (find (list, random))
random= genGF.generate();
}
if (F (random, x) == 0)
{
list.append (random);
continue;
}
} while (find (list, random));
return random;
}
ipower
int ipower(int b, int m)
int ipower ( int b, int m )
Definition cf_util.cc:27

## if() [1/2]

if | ( | i | = =0 | ) |

## if() [2/2]

if | ( | LCCand * | absLC(coF) = = abs (LC (F)) | ) |

Definition at line 72 of file cfModGcd.cc.

{
if (LCCand*abs (LC (coG)) == abs (LC (G)))
{
if (abs (cand)*abs (coF) == abs (F))
{
if (abs (cand)*abs (coG) == abs (G))
return true;
}
return false;
}
return false;
}
abs
Rational abs(const Rational &a)
Definition GMPrat.cc:436
LC
CanonicalForm LC(const CanonicalForm &f)
Definition canonicalform.h:310
coG
const CanonicalForm const CanonicalForm const CanonicalForm & coG
Definition cfModGcd.cc:68
cand
const CanonicalForm const CanonicalForm const CanonicalForm const CanonicalForm & cand
Definition cfModGcd.cc:70
coF
const CanonicalForm const CanonicalForm & coF
Definition cfModGcd.cc:68

## modGCDFp() [1/2]

CanonicalForm modGCDFp | ( | const CanonicalForm & | F,
|  | const CanonicalForm & | G,
|  | bool & | topLevel,
|  | CFList & | l )

Definition at line 1213 of file cfModGcd.cc.

{
CanonicalForm dummy1, dummy2;
CanonicalForm result= modGCDFp (F, G, dummy1, dummy2, topLevel, l);
return result;
}
l
int l
Definition cfEzgcd.cc:100
topLevel
const CanonicalForm CFMap CFMap bool topLevel
Definition cfGcdAlgExt.cc:58
modGCDFp
CanonicalForm modGCDFp(const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &coF, CanonicalForm &coG, bool &topLevel, CFList &l)
Definition cfModGcd.cc:1224

## modGCDFp() [2/2]

CanonicalForm modGCDFp | ( | const CanonicalForm & | F,
|  | const CanonicalForm & | G,
|  | CanonicalForm & | coF,
|  | CanonicalForm & | coG,
|  | bool & | topLevel,
|  | CFList & | l )

Definition at line 1224 of file cfModGcd.cc.

{
CanonicalForm A= F;
CanonicalForm B= G;
if (F.isZero() && degree(G) > 0)
{
coF= 0;
coG= Lc (G);
return G/Lc(G);
}
else if (G.isZero() && degree (F) > 0)
{
coF= Lc (F);
coG= 0;
return F/Lc(F);
}
else if (F.isZero() && G.isZero())
{
coF= coG= 0;
return F.genOne();
}
if (F.inBaseDomain() || G.inBaseDomain())
{
coF= F;
coG= G;
return F.genOne();
}
if (F.isUnivariate() && fdivides(F, G, coG))
{
coF= Lc (F);
return F/Lc(F);
}
if (G.isUnivariate() && fdivides(G, F, coF))
{
coG= Lc (G);
return G/Lc(G);
}
if (F == G)
{
coF= coG= Lc (F);
return F/Lc(F);
}
CFMap M,N;
int best_level= myCompress (A, B, M, N, topLevel);
if (best_level == 0)
{
coF= F;
coG= G;
return B.genOne();
}
A= M(A);
B= M(B);
Variable x= Variable (1);
//univariate case
if (A.isUnivariate() && B.isUnivariate())
{
CanonicalForm result= gcd (A, B);
coF= N (A/result);
coG= N (B/result);
return N (result);
}
CanonicalForm cA, cB; // content of A and B
CanonicalForm ppA, ppB; // primitive part of A and B
CanonicalForm gcdcAcB;
cA = uni_content (A);
cB = uni_content (B);
gcdcAcB= gcd (cA, cB);
ppA= A/cA;
ppB= B/cB;
CanonicalForm lcA, lcB; // leading coefficients of A and B
CanonicalForm gcdlcAlcB;
lcA= uni_lcoeff (ppA);
lcB= uni_lcoeff (ppB);
gcdlcAlcB= gcd (lcA, lcB);
int d= totaldegree (ppA, Variable (2), Variable (ppA.level()));
int d0;
if (d == 0)
{
coF= N (ppA*(cA/gcdcAcB));
coG= N (ppB*(cB/gcdcAcB));
return N(gcdcAcB);
}
d0= totaldegree (ppB, Variable (2), Variable (ppB.level()));
if (d0 < d)
d= d0;
if (d == 0)
{
coF= N (ppA*(cA/gcdcAcB));
coG= N (ppB*(cB/gcdcAcB));
return N(gcdcAcB);
}
CanonicalForm m, random_element, G_m, G_random_element, H, cH, ppH;
CanonicalForm newtonPoly, coF_random_element, coG_random_element,
coF_m, coG_m, ppCoF, ppCoG;
newtonPoly= 1;
m= gcdlcAlcB;
H= 0;
coF= 0;
coG= 0;
G_m= 0;
Variable alpha, V_buf, cleanUp;
bool fail= false;
bool inextension= false;
topLevel= false;
CFList source, dest;
int bound1= degree (ppA, 1);
int bound2= degree (ppB, 1);
do
{
if (inextension)
random_element= randomElement (m*lcA*lcB, V_buf, l, fail);
else
random_element= FpRandomElement (m*lcA*lcB, l, fail);
if (!fail && !inextension)
{
CFList list;
TIMING_START (gcd_recursion);
G_random_element=
modGCDFp (ppA (random_element,x), ppB (random_element,x),
coF_random_element, coG_random_element, topLevel,
list);
TIMING_END_AND_PRINT (gcd_recursion,
"time for recursive call: ");
DEBOUTLN (cerr, "G_random_element= " << G_random_element);
}
else if (!fail && inextension)
{
CFList list;
TIMING_START (gcd_recursion);
G_random_element=
modGCDFq (ppA (random_element, x), ppB (random_element, x),
coF_random_element, coG_random_element, V_buf,
list, topLevel);
TIMING_END_AND_PRINT (gcd_recursion,
"time for recursive call: ");
DEBOUTLN (cerr, "G_random_element= " << G_random_element);
}
else if (fail && !inextension)
{
source= CFList();
dest= CFList();
CFList list;
CanonicalForm mipo;
int deg= 2;
bool initialized= false;
do
{
mipo= randomIrredpoly (deg, x);
if (initialized)
setMipo (alpha, mipo);
else
alpha= rootOf (mipo);
inextension= true;
initialized= true;
fail= false;
random_element= randomElement (m*lcA*lcB, alpha, l, fail);
deg++;
} while (fail);
list= CFList();
V_buf= alpha;
cleanUp= alpha;
TIMING_START (gcd_recursion);
G_random_element=
modGCDFq (ppA (random_element, x), ppB (random_element, x),
coF_random_element, coG_random_element, alpha,
list, topLevel);
TIMING_END_AND_PRINT (gcd_recursion,
"time for recursive call: ");
DEBOUTLN (cerr, "G_random_element= " << G_random_element);
}
else if (fail && inextension)
{
source= CFList();
dest= CFList();
Variable V_buf3= V_buf;
V_buf= chooseExtension (V_buf);
bool prim_fail= false;
Variable V_buf2;
CanonicalForm prim_elem, im_prim_elem;
prim_elem= primitiveElement (alpha, V_buf2, prim_fail);
if (V_buf3 != alpha)
{
m= mapDown (m, prim_elem, im_prim_elem, alpha, source, dest);
G_m= mapDown (G_m, prim_elem, im_prim_elem, alpha, source, dest);
coF_m= mapDown (coF_m, prim_elem, im_prim_elem, alpha, source, dest);
coG_m= mapDown (coG_m, prim_elem, im_prim_elem, alpha, source, dest);
newtonPoly= mapDown (newtonPoly, prim_elem, im_prim_elem, alpha,
source, dest);
ppA= mapDown (ppA, prim_elem, im_prim_elem, alpha, source, dest);
ppB= mapDown (ppB, prim_elem, im_prim_elem, alpha, source, dest);
gcdlcAlcB= mapDown (gcdlcAlcB, prim_elem, im_prim_elem, alpha, source,
dest);
lcA= mapDown (lcA, prim_elem, im_prim_elem, alpha, source, dest);
lcB= mapDown (lcB, prim_elem, im_prim_elem, alpha, source, dest);
for (CFListIterator i= l; i.hasItem(); i++)
i.getItem()= mapDown (i.getItem(), prim_elem, im_prim_elem, alpha,
source, dest);
}
ASSERT (!prim_fail, "failure in integer factorizer");
if (prim_fail)
; //ERROR
else
im_prim_elem= mapPrimElem (prim_elem, alpha, V_buf);
DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (alpha));
DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (V_buf2));
for (CFListIterator i= l; i.hasItem(); i++)
i.getItem()= mapUp (i.getItem(), alpha, V_buf, prim_elem,
im_prim_elem, source, dest);
m= mapUp (m, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
G_m= mapUp (G_m, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
coF_m= mapUp (coF_m, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
coG_m= mapUp (coG_m, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
newtonPoly= mapUp (newtonPoly, alpha, V_buf, prim_elem, im_prim_elem,
source, dest);
ppA= mapUp (ppA, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
ppB= mapUp (ppB, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
gcdlcAlcB= mapUp (gcdlcAlcB, alpha, V_buf, prim_elem, im_prim_elem,
source, dest);
lcA= mapUp (lcA, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
lcB= mapUp (lcB, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
fail= false;
random_element= randomElement (m*lcA*lcB, V_buf, l, fail );
DEBOUTLN (cerr, "fail= " << fail);
CFList list;
TIMING_START (gcd_recursion);
G_random_element=
modGCDFq (ppA (random_element, x), ppB (random_element, x),
coF_random_element, coG_random_element, V_buf,
list, topLevel);
TIMING_END_AND_PRINT (gcd_recursion,
"time for recursive call: ");
DEBOUTLN (cerr, "G_random_element= " << G_random_element);
alpha= V_buf;
}
if (!G_random_element.inCoeffDomain())
d0= totaldegree (G_random_element, Variable(2),
Variable (G_random_element.level()));
else
d0= 0;
if (d0 == 0)
{
if (inextension)
prune (cleanUp);
coF= N (ppA*(cA/gcdcAcB));
coG= N (ppB*(cB/gcdcAcB));
return N(gcdcAcB);
}
if (d0 > d)
{
if (!find (l, random_element))
l.append (random_element);
continue;
}
G_random_element= (gcdlcAlcB(random_element,x)/uni_lcoeff(G_random_element))
*G_random_element;
coF_random_element= (lcA(random_element,x)/uni_lcoeff(coF_random_element))
*coF_random_element;
coG_random_element= (lcB(random_element,x)/uni_lcoeff(coG_random_element))
*coG_random_element;
if (!G_random_element.inCoeffDomain())
d0= totaldegree (G_random_element, Variable(2),
Variable (G_random_element.level()));
else
d0= 0;
if (d0 < d)
{
m= gcdlcAlcB;
newtonPoly= 1;
G_m= 0;
d= d0;
coF_m= 0;
coG_m= 0;
}
TIMING_START (newton_interpolation);
H= newtonInterp (random_element, G_random_element, newtonPoly, G_m, x);
coF= newtonInterp (random_element, coF_random_element, newtonPoly, coF_m,x);
coG= newtonInterp (random_element, coG_random_element, newtonPoly, coG_m,x);
TIMING_END_AND_PRINT (newton_interpolation,
"time for newton_interpolation: ");
//termination test
if ((uni_lcoeff (H) == gcdlcAlcB) || (G_m == H))
{
TIMING_START (termination_test);
if (gcdlcAlcB.isOne())
cH= 1;
else
cH= uni_content (H);
ppH= H/cH;
ppH /= Lc (ppH);
CanonicalForm lcppH= gcdlcAlcB/cH;
CanonicalForm ccoF= lcppH/Lc (lcppH);
CanonicalForm ccoG= lcppH/Lc (lcppH);
ppCoF= coF/ccoF;
ppCoG= coG/ccoG;
DEBOUTLN (cerr, "ppH= " << ppH);
if (((degree (ppCoF,1)+degree (ppH,1) == bound1) &&
(degree (ppCoG,1)+degree (ppH,1) == bound2) &&
terminationTest (ppA, ppB, ppCoF, ppCoG, ppH)) ||
(fdivides (ppH, ppA, ppCoF) && fdivides (ppH, ppB, ppCoG)))
{
if (inextension)
prune (cleanUp);
coF= N ((cA/gcdcAcB)*ppCoF);
coG= N ((cB/gcdcAcB)*ppCoG);
TIMING_END_AND_PRINT (termination_test,
"time for successful termination Fp: ");
return N(gcdcAcB*ppH);
}
TIMING_END_AND_PRINT (termination_test,
"time for unsuccessful termination Fp: ");
}
G_m= H;
coF_m= coF;
coG_m= coG;
newtonPoly= newtonPoly*(x - random_element);
m= m*(x - random_element);
if (!find (l, random_element))
l.append (random_element);
} while (1);
}
totaldegree
int totaldegree(const CanonicalForm &f)
int totaldegree ( const CanonicalForm & f )
Definition cf_ops.cc:523
Lc
CanonicalForm Lc(const CanonicalForm &f)
Definition canonicalform.h:307
modGCDFq
CanonicalForm modGCDFq(const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &coF, CanonicalForm &coG, Variable &alpha, CFList &l, bool &topLevel)
GCD of F and G over , l and topLevel are only used internally, output is monic based on Alg....
Definition cfModGcd.cc:479
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
chooseExtension
static Variable chooseExtension(const Variable &alpha)
Definition cfModGcd.cc:421
uni_lcoeff
static CanonicalForm uni_lcoeff(const CanonicalForm &F)
compute the leading coefficient of F, where F is considered as an element of , order on is dp.
Definition cfModGcd.cc:340
randomElement
static CanonicalForm randomElement(const CanonicalForm &F, const Variable &alpha, CFList &list, bool &fail)
compute a random element a of , s.t. F(a) , F is a univariate polynomial, returns fail if there are...
Definition cfModGcd.cc:380
FpRandomElement
static CanonicalForm FpRandomElement(const CanonicalForm &F, CFList &list, bool &fail)
Definition cfModGcd.cc:1172
newtonInterp
static CanonicalForm newtonInterp(const CanonicalForm &alpha, const CanonicalForm &u, const CanonicalForm &newtonPoly, const CanonicalForm &oldInterPoly, const Variable &x)
Newton interpolation - Incremental algorithm. Given a list of values alpha_i and a list of polynomial...
Definition cfModGcd.cc:365
terminationTest
bool terminationTest(const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &coF, const CanonicalForm &coG, const CanonicalForm &cand)
fdivides
bool fdivides(const CanonicalForm &f, const CanonicalForm &g)
randomIrredpoly
CanonicalForm randomIrredpoly(int i, const Variable &x)
computes a random monic irreducible univariate polynomial in x over Fp of degree i via NTL/FLINT
Definition cf_irred.cc:26
mapPrimElem
CanonicalForm mapPrimElem(const CanonicalForm &primElem, const Variable &alpha, const Variable &beta)
compute the image of a primitive element of in . We assume .
Definition cf_map_ext.cc:451
primitiveElement
CanonicalForm primitiveElement(const Variable &alpha, Variable &beta, bool &fail)
determine a primitive element of , is a primitive element of a field which is isomorphic to
Definition cf_map_ext.cc:343
mapDown
static CanonicalForm mapDown(const CanonicalForm &F, const Variable &alpha, const CanonicalForm &G, CFList &source, CFList &dest)
the CanonicalForm G is the output of map_up, returns F considered as an element over ,...
Definition cf_map_ext.cc:124
mapUp
static CanonicalForm mapUp(const Variable &alpha, const Variable &beta)
and is a primitive element, returns the image of
Definition cf_map_ext.cc:71
CFMap
class CFMap
Definition cf_map.h:85
CanonicalForm::genOne
CanonicalForm genOne() const
Definition canonicalform.cc:1881
CanonicalForm::isOne
CF_NO_INLINE bool isOne() const
CanonicalForm::inBaseDomain
bool inBaseDomain() const
Definition canonicalform.cc:104
DEBOUTLN
#define DEBOUTLN(stream, objects)
Definition debug.h:49
H
CanonicalForm H
Definition facAbsFact.cc:60
mipo
CanonicalForm mipo
Definition facAlgExt.cc:57
prune
void FACTORY_PUBLIC prune(Variable &alpha)
Definition variable.cc:261
setMipo
void setMipo(const Variable &alpha, const CanonicalForm &mipo)
Definition variable.cc:219
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

## modGCDFq() [1/2]

CanonicalForm modGCDFq | ( | const CanonicalForm & | F,
|  | const CanonicalForm & | G,
|  | CanonicalForm & | coF,
|  | CanonicalForm & | coG,
|  | Variable & | alpha,
|  | CFList & | l,
|  | bool & | topLevel )

GCD of F and G over , l and topLevel are only used internally, output is monic based on Alg. 7.2. as described in "Algorithms for Computer Algebra" by Geddes, Czapor, Labahn.

Definition at line 479 of file cfModGcd.cc.

{
CanonicalForm A= F;
CanonicalForm B= G;
if (F.isZero() && degree(G) > 0)
{
coF= 0;
coG= Lc (G);
return G/Lc(G);
}
else if (G.isZero() && degree (F) > 0)
{
coF= Lc (F);
coG= 0;
return F/Lc(F);
}
else if (F.isZero() && G.isZero())
{
coF= coG= 0;
return F.genOne();
}
if (F.inBaseDomain() || G.inBaseDomain())
{
coF= F;
coG= G;
return F.genOne();
}
if (F.isUnivariate() && fdivides(F, G, coG))
{
coF= Lc (F);
return F/Lc(F);
}
if (G.isUnivariate() && fdivides(G, F, coF))
{
coG= Lc (G);
return G/Lc(G);
}
if (F == G)
{
coF= coG= Lc (F);
return F/Lc(F);
}
CFMap M,N;
int best_level= myCompress (A, B, M, N, topLevel);
if (best_level == 0)
{
coF= F;
coG= G;
return B.genOne();
}
A= M(A);
B= M(B);
Variable x= Variable(1);
//univariate case
if (A.isUnivariate() && B.isUnivariate())
{
CanonicalForm result= gcd (A, B);
coF= N (A/result);
coG= N (B/result);
return N (result);
}
CanonicalForm cA, cB; // content of A and B
CanonicalForm ppA, ppB; // primitive part of A and B
CanonicalForm gcdcAcB;
cA = uni_content (A);
cB = uni_content (B);
gcdcAcB= gcd (cA, cB);
ppA= A/cA;
ppB= B/cB;
CanonicalForm lcA, lcB; // leading coefficients of A and B
CanonicalForm gcdlcAlcB;
lcA= uni_lcoeff (ppA);
lcB= uni_lcoeff (ppB);
gcdlcAlcB= gcd (lcA, lcB);
int d= totaldegree (ppA, Variable(2), Variable (ppA.level()));
if (d == 0)
{
coF= N (ppA*(cA/gcdcAcB));
coG= N (ppB*(cB/gcdcAcB));
return N(gcdcAcB);
}
int d0= totaldegree (ppB, Variable(2), Variable (ppB.level()));
if (d0 < d)
d= d0;
if (d == 0)
{
coF= N (ppA*(cA/gcdcAcB));
coG= N (ppB*(cB/gcdcAcB));
return N(gcdcAcB);
}
CanonicalForm m, random_element, G_m, G_random_element, H, cH, ppH;
CanonicalForm newtonPoly, coF_random_element, coG_random_element, coF_m,
coG_m, ppCoF, ppCoG;
newtonPoly= 1;
m= gcdlcAlcB;
G_m= 0;
coF= 0;
coG= 0;
H= 0;
bool fail= false;
topLevel= false;
bool inextension= false;
Variable V_buf= alpha, V_buf4= alpha;
CanonicalForm prim_elem, im_prim_elem;
CanonicalForm prim_elem_alpha, im_prim_elem_alpha;
CFList source, dest;
int bound1= degree (ppA, 1);
int bound2= degree (ppB, 1);
do
{
random_element= randomElement (m*lcA*lcB, V_buf, l, fail);
if (fail)
{
source= CFList();
dest= CFList();
Variable V_buf3= V_buf;
V_buf= chooseExtension (V_buf);
bool prim_fail= false;
Variable V_buf2;
prim_elem= primitiveElement (V_buf4, V_buf2, prim_fail);
if (V_buf4 == alpha)
prim_elem_alpha= prim_elem;
if (V_buf3 != V_buf4)
{
m= mapDown (m, prim_elem, im_prim_elem, V_buf4, source, dest);
G_m= mapDown (G_m, prim_elem, im_prim_elem, V_buf4, source, dest);
coF_m= mapDown (coF_m, prim_elem, im_prim_elem, V_buf4, source, dest);
coG_m= mapDown (coG_m, prim_elem, im_prim_elem, V_buf4, source, dest);
newtonPoly= mapDown (newtonPoly, prim_elem, im_prim_elem, V_buf4,
source, dest);
ppA= mapDown (ppA, prim_elem, im_prim_elem, V_buf4, source, dest);
ppB= mapDown (ppB, prim_elem, im_prim_elem, V_buf4, source, dest);
gcdlcAlcB= mapDown (gcdlcAlcB, prim_elem, im_prim_elem, V_buf4,
source, dest);
lcA= mapDown (lcA, prim_elem, im_prim_elem, V_buf4, source, dest);
lcB= mapDown (lcB, prim_elem, im_prim_elem, V_buf4, source, dest);
for (CFListIterator i= l; i.hasItem(); i++)
i.getItem()= mapDown (i.getItem(), prim_elem, im_prim_elem, V_buf4,
source, dest);
}
ASSERT (!prim_fail, "failure in integer factorizer");
if (prim_fail)
; //ERROR
else
im_prim_elem= mapPrimElem (prim_elem, V_buf4, V_buf);
if (V_buf4 == alpha)
im_prim_elem_alpha= im_prim_elem;
else
im_prim_elem_alpha= mapUp (im_prim_elem_alpha, V_buf4, V_buf, prim_elem,
im_prim_elem, source, dest);
DEBOUTLN (cerr, "getMipo (V_buf4)= " << getMipo (V_buf4));
DEBOUTLN (cerr, "getMipo (V_buf2)= " << getMipo (V_buf2));
inextension= true;
for (CFListIterator i= l; i.hasItem(); i++)
i.getItem()= mapUp (i.getItem(), V_buf4, V_buf, prim_elem,
im_prim_elem, source, dest);
m= mapUp (m, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
G_m= mapUp (G_m, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
coF_m= mapUp (coF_m, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
coG_m= mapUp (coG_m, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
newtonPoly= mapUp (newtonPoly, V_buf4, V_buf, prim_elem, im_prim_elem,
source, dest);
ppA= mapUp (ppA, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
ppB= mapUp (ppB, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
gcdlcAlcB= mapUp (gcdlcAlcB, V_buf4, V_buf, prim_elem, im_prim_elem,
source, dest);
lcA= mapUp (lcA, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
lcB= mapUp (lcB, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
fail= false;
random_element= randomElement (m*lcA*lcB, V_buf, l, fail );
DEBOUTLN (cerr, "fail= " << fail);
CFList list;
TIMING_START (gcd_recursion);
G_random_element=
modGCDFq (ppA (random_element, x), ppB (random_element, x),
coF_random_element, coG_random_element, V_buf,
list, topLevel);
TIMING_END_AND_PRINT (gcd_recursion,
"time for recursive call: ");
DEBOUTLN (cerr, "G_random_element= " << G_random_element);
V_buf4= V_buf;
}
else
{
CFList list;
TIMING_START (gcd_recursion);
G_random_element=
modGCDFq (ppA(random_element, x), ppB(random_element, x),
coF_random_element, coG_random_element, V_buf,
list, topLevel);
TIMING_END_AND_PRINT (gcd_recursion,
"time for recursive call: ");
DEBOUTLN (cerr, "G_random_element= " << G_random_element);
}
if (!G_random_element.inCoeffDomain())
d0= totaldegree (G_random_element, Variable(2),
Variable (G_random_element.level()));
else
d0= 0;
if (d0 == 0)
{
if (inextension)
{
CFList u, v;
ppA= mapDown (ppA, prim_elem_alpha, im_prim_elem_alpha, alpha, u, v);
ppB= mapDown (ppB, prim_elem_alpha, im_prim_elem_alpha, alpha, u, v);
prune1 (alpha);
}
coF= N (ppA*(cA/gcdcAcB));
coG= N (ppB*(cB/gcdcAcB));
return N(gcdcAcB);
}
if (d0 > d)
{
if (!find (l, random_element))
l.append (random_element);
continue;
}
G_random_element=
(gcdlcAlcB(random_element, x)/uni_lcoeff (G_random_element))
* G_random_element;
coF_random_element= (lcA(random_element,x)/uni_lcoeff(coF_random_element))
*coF_random_element;
coG_random_element= (lcB(random_element,x)/uni_lcoeff(coG_random_element))
*coG_random_element;
if (!G_random_element.inCoeffDomain())
d0= totaldegree (G_random_element, Variable(2),
Variable (G_random_element.level()));
else
d0= 0;
if (d0 < d)
{
m= gcdlcAlcB;
newtonPoly= 1;
G_m= 0;
d= d0;
coF_m= 0;
coG_m= 0;
}
TIMING_START (newton_interpolation);
H= newtonInterp (random_element, G_random_element, newtonPoly, G_m, x);
coF= newtonInterp (random_element, coF_random_element, newtonPoly, coF_m,x);
coG= newtonInterp (random_element, coG_random_element, newtonPoly, coG_m,x);
TIMING_END_AND_PRINT (newton_interpolation,
"time for newton interpolation: ");
//termination test
if ((uni_lcoeff (H) == gcdlcAlcB) || (G_m == H))
{
TIMING_START (termination_test);
if (gcdlcAlcB.isOne())
cH= 1;
else
cH= uni_content (H);
ppH= H/cH;
ppH /= Lc (ppH);
CanonicalForm lcppH= gcdlcAlcB/cH;
CanonicalForm ccoF= lcppH/Lc (lcppH);
CanonicalForm ccoG= lcppH/Lc (lcppH);
ppCoF= coF/ccoF;
ppCoG= coG/ccoG;
if (inextension)
{
if (((degree (ppCoF,1)+degree (ppH,1) == bound1) &&
(degree (ppCoG,1)+degree (ppH,1) == bound2) &&
terminationTest (ppA, ppB, ppCoF, ppCoG, ppH)) ||
(fdivides (ppH, ppA, ppCoF) && fdivides (ppH, ppB, ppCoG)))
{
CFList u, v;
DEBOUTLN (cerr, "ppH before mapDown= " << ppH);
ppH= mapDown (ppH, prim_elem_alpha, im_prim_elem_alpha, alpha, u, v);
ppCoF= mapDown (ppCoF, prim_elem_alpha, im_prim_elem_alpha, alpha, u, v);
ppCoG= mapDown (ppCoG, prim_elem_alpha, im_prim_elem_alpha, alpha, u, v);
DEBOUTLN (cerr, "ppH after mapDown= " << ppH);
coF= N ((cA/gcdcAcB)*ppCoF);
coG= N ((cB/gcdcAcB)*ppCoG);
TIMING_END_AND_PRINT (termination_test,
"time for successful termination test Fq: ");
prune1 (alpha);
return N(gcdcAcB*ppH);
}
}
else if (((degree (ppCoF,1)+degree (ppH,1) == bound1) &&
(degree (ppCoG,1)+degree (ppH,1) == bound2) &&
terminationTest (ppA, ppB, ppCoF, ppCoG, ppH)) ||
(fdivides (ppH, ppA, ppCoF) && fdivides (ppH, ppB, ppCoG)))
{
coF= N ((cA/gcdcAcB)*ppCoF);
coG= N ((cB/gcdcAcB)*ppCoG);
TIMING_END_AND_PRINT (termination_test,
"time for successful termination test Fq: ");
return N(gcdcAcB*ppH);
}
TIMING_END_AND_PRINT (termination_test,
"time for unsuccessful termination test Fq: ");
}
G_m= H;
coF_m= coF;
coG_m= coG;
newtonPoly= newtonPoly*(x - random_element);
m= m*(x - random_element);
if (!find (l, random_element))
l.append (random_element);
} while (1);
}
v
const Variable & v
< [in] a sqrfree bivariate poly
Definition facBivar.h:39
prune1
void prune1(const Variable &alpha)
Definition variable.cc:291

## modGCDFq() [2/2]

CanonicalForm modGCDFq | ( | const CanonicalForm & | F,
|  | const CanonicalForm & | G,
|  | Variable & | alpha,
|  | CFList & | l,
|  | bool & | topLevel )

Definition at line 463 of file cfModGcd.cc.

{
CanonicalForm dummy1, dummy2;
CanonicalForm result= modGCDFq (F, G, dummy1, dummy2, alpha, l,
topLevel);
return result;
}

## modGCDGF() [1/2]

CanonicalForm modGCDGF | ( | const CanonicalForm & | F,
|  | const CanonicalForm & | G,
|  | CanonicalForm & | coF,
|  | CanonicalForm & | coG,
|  | CFList & | l,
|  | bool & | topLevel )

GCD of F and G over GF, based on Alg. 7.2. as described in "Algorithms for Computer Algebra" by Geddes, Czapor, Labahn Usually this algorithm will be faster than modGCDFq since GF has faster field arithmetics, however it might fail if the input is large since the size of the base field is bounded by 2^16, output is monic.

Definition at line 873 of file cfModGcd.cc.

{
CanonicalForm A= F;
CanonicalForm B= G;
if (F.isZero() && degree(G) > 0)
{
coF= 0;
coG= Lc (G);
return G/Lc(G);
}
else if (G.isZero() && degree (F) > 0)
{
coF= Lc (F);
coG= 0;
return F/Lc(F);
}
else if (F.isZero() && G.isZero())
{
coF= coG= 0;
return F.genOne();
}
if (F.inBaseDomain() || G.inBaseDomain())
{
coF= F;
coG= G;
return F.genOne();
}
if (F.isUnivariate() && fdivides(F, G, coG))
{
coF= Lc (F);
return F/Lc(F);
}
if (G.isUnivariate() && fdivides(G, F, coF))
{
coG= Lc (G);
return G/Lc(G);
}
if (F == G)
{
coF= coG= Lc (F);
return F/Lc(F);
}
CFMap M,N;
int best_level= myCompress (A, B, M, N, topLevel);
if (best_level == 0)
{
coF= F;
coG= G;
return B.genOne();
}
A= M(A);
B= M(B);
Variable x= Variable(1);
//univariate case
if (A.isUnivariate() && B.isUnivariate())
{
CanonicalForm result= gcd (A, B);
coF= N (A/result);
coG= N (B/result);
return N (result);
}
CanonicalForm cA, cB; // content of A and B
CanonicalForm ppA, ppB; // primitive part of A and B
CanonicalForm gcdcAcB;
cA = uni_content (A);
cB = uni_content (B);
gcdcAcB= gcd (cA, cB);
ppA= A/cA;
ppB= B/cB;
CanonicalForm lcA, lcB; // leading coefficients of A and B
CanonicalForm gcdlcAlcB;
lcA= uni_lcoeff (ppA);
lcB= uni_lcoeff (ppB);
gcdlcAlcB= gcd (lcA, lcB);
int d= totaldegree (ppA, Variable(2), Variable (ppA.level()));
if (d == 0)
{
coF= N (ppA*(cA/gcdcAcB));
coG= N (ppB*(cB/gcdcAcB));
return N(gcdcAcB);
}
int d0= totaldegree (ppB, Variable(2), Variable (ppB.level()));
if (d0 < d)
d= d0;
if (d == 0)
{
coF= N (ppA*(cA/gcdcAcB));
coG= N (ppB*(cB/gcdcAcB));
return N(gcdcAcB);
}
CanonicalForm m, random_element, G_m, G_random_element, H, cH, ppH;
CanonicalForm newtonPoly, coF_random_element, coG_random_element, coF_m,
coG_m, ppCoF, ppCoG;
newtonPoly= 1;
m= gcdlcAlcB;
G_m= 0;
coF= 0;
coG= 0;
H= 0;
bool fail= false;
topLevel= false;
bool inextension= false;
int p=-1;
int k= getGFDegree();
int kk;
int expon;
char gf_name_buf= gf_name;
int bound1= degree (ppA, 1);
int bound2= degree (ppB, 1);
do
{
random_element= GFRandomElement (m*lcA*lcB, l, fail);
if (fail)
{
p= getCharacteristic();
expon= 2;
kk= getGFDegree();
if (ipower (p, kk*expon) < (1 << 16))
setCharacteristic (p, kk*(int)expon, 'b');
else
{
expon= (int) floor((log ((double)((1<<16) - 1)))/(log((double)p)*kk));
ASSERT (expon >= 2, "not enough points in modGCDGF");
setCharacteristic (p, (int)(kk*expon), 'b');
}
inextension= true;
fail= false;
for (CFListIterator i= l; i.hasItem(); i++)
i.getItem()= GFMapUp (i.getItem(), kk);
m= GFMapUp (m, kk);
G_m= GFMapUp (G_m, kk);
newtonPoly= GFMapUp (newtonPoly, kk);
coF_m= GFMapUp (coF_m, kk);
coG_m= GFMapUp (coG_m, kk);
ppA= GFMapUp (ppA, kk);
ppB= GFMapUp (ppB, kk);
gcdlcAlcB= GFMapUp (gcdlcAlcB, kk);
lcA= GFMapUp (lcA, kk);
lcB= GFMapUp (lcB, kk);
random_element= GFRandomElement (m*lcA*lcB, l, fail);
DEBOUTLN (cerr, "fail= " << fail);
CFList list;
TIMING_START (gcd_recursion);
G_random_element= modGCDGF (ppA(random_element, x), ppB(random_element, x),
coF_random_element, coG_random_element,
list, topLevel);
TIMING_END_AND_PRINT (gcd_recursion,
"time for recursive call: ");
DEBOUTLN (cerr, "G_random_element= " << G_random_element);
}
else
{
CFList list;
TIMING_START (gcd_recursion);
G_random_element= modGCDGF (ppA(random_element, x), ppB(random_element, x),
coF_random_element, coG_random_element,
list, topLevel);
TIMING_END_AND_PRINT (gcd_recursion,
"time for recursive call: ");
DEBOUTLN (cerr, "G_random_element= " << G_random_element);
}
if (!G_random_element.inCoeffDomain())
d0= totaldegree (G_random_element, Variable(2),
Variable (G_random_element.level()));
else
d0= 0;
if (d0 == 0)
{
if (inextension)
{
ppA= GFMapDown (ppA, k);
ppB= GFMapDown (ppB, k);
setCharacteristic (p, k, gf_name_buf);
}
coF= N (ppA*(cA/gcdcAcB));
coG= N (ppB*(cB/gcdcAcB));
return N(gcdcAcB);
}
if (d0 > d)
{
if (!find (l, random_element))
l.append (random_element);
continue;
}
G_random_element=
(gcdlcAlcB(random_element, x)/uni_lcoeff(G_random_element)) *
G_random_element;
coF_random_element= (lcA(random_element,x)/uni_lcoeff(coF_random_element))
*coF_random_element;
coG_random_element= (lcB(random_element,x)/uni_lcoeff(coG_random_element))
*coG_random_element;
if (!G_random_element.inCoeffDomain())
d0= totaldegree (G_random_element, Variable(2),
Variable (G_random_element.level()));
else
d0= 0;
if (d0 < d)
{
m= gcdlcAlcB;
newtonPoly= 1;
G_m= 0;
d= d0;
coF_m= 0;
coG_m= 0;
}
TIMING_START (newton_interpolation);
H= newtonInterp (random_element, G_random_element, newtonPoly, G_m, x);
coF= newtonInterp (random_element, coF_random_element, newtonPoly, coF_m,x);
coG= newtonInterp (random_element, coG_random_element, newtonPoly, coG_m,x);
TIMING_END_AND_PRINT (newton_interpolation,
"time for newton interpolation: ");
//termination test
if ((uni_lcoeff (H) == gcdlcAlcB) || (G_m == H))
{
TIMING_START (termination_test);
if (gcdlcAlcB.isOne())
cH= 1;
else
cH= uni_content (H);
ppH= H/cH;
ppH /= Lc (ppH);
CanonicalForm lcppH= gcdlcAlcB/cH;
CanonicalForm ccoF= lcppH/Lc (lcppH);
CanonicalForm ccoG= lcppH/Lc (lcppH);
ppCoF= coF/ccoF;
ppCoG= coG/ccoG;
if (inextension)
{
if (((degree (ppCoF,1)+degree (ppH,1) == bound1) &&
(degree (ppCoG,1)+degree (ppH,1) == bound2) &&
terminationTest (ppA, ppB, ppCoF, ppCoG, ppH)) ||
(fdivides (ppH, ppA, ppCoF) && fdivides (ppH, ppB, ppCoG)))
{
DEBOUTLN (cerr, "ppH before mapDown= " << ppH);
ppH= GFMapDown (ppH, k);
ppCoF= GFMapDown (ppCoF, k);
ppCoG= GFMapDown (ppCoG, k);
DEBOUTLN (cerr, "ppH after mapDown= " << ppH);
coF= N ((cA/gcdcAcB)*ppCoF);
coG= N ((cB/gcdcAcB)*ppCoG);
setCharacteristic (p, k, gf_name_buf);
TIMING_END_AND_PRINT (termination_test,
"time for successful termination GF: ");
return N(gcdcAcB*ppH);
}
}
else
{
if (((degree (ppCoF,1)+degree (ppH,1) == bound1) &&
(degree (ppCoG,1)+degree (ppH,1) == bound2) &&
terminationTest (ppA, ppB, ppCoF, ppCoG, ppH)) ||
(fdivides (ppH, ppA, ppCoF) && fdivides (ppH, ppB, ppCoG)))
{
coF= N ((cA/gcdcAcB)*ppCoF);
coG= N ((cB/gcdcAcB)*ppCoG);
TIMING_END_AND_PRINT (termination_test,
"time for successful termination GF: ");
return N(gcdcAcB*ppH);
}
}
TIMING_END_AND_PRINT (termination_test,
"time for unsuccessful termination GF: ");
}
G_m= H;
coF_m= coF;
coG_m= coG;
newtonPoly= newtonPoly*(x - random_element);
m= m*(x - random_element);
if (!find (l, random_element))
l.append (random_element);
} while (1);
}
setCharacteristic
void FACTORY_PUBLIC setCharacteristic(int c)
Definition cf_char.cc:28
GFRandomElement
static CanonicalForm GFRandomElement(const CanonicalForm &F, CFList &list, bool &fail)
compute a random element a of GF, s.t. F(a) , F is a univariate polynomial, returns fail if there ar...
Definition cfModGcd.cc:820
modGCDGF
CanonicalForm modGCDGF(const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &coF, CanonicalForm &coG, CFList &l, bool &topLevel)
GCD of F and G over GF, based on Alg. 7.2. as described in "Algorithms forComputer Algebra" by Geddes...
Definition cfModGcd.cc:873
GFMapDown
CanonicalForm GFMapDown(const CanonicalForm &F, int k)
maps a polynomial over to a polynomial over , d needs to be a multiple of k
Definition cf_map_ext.cc:277
GFMapUp
CanonicalForm GFMapUp(const CanonicalForm &F, int k)
maps a polynomial over to a polynomial over , d needs to be a multiple of k
Definition cf_map_ext.cc:241
gf_name
VAR char gf_name
Definition gfops.cc:52
log
gmp_float log(const gmp_float &a)
Definition mpr_complex.cc:335
amp::floor
const signed long floor(const ampf< Precision > &x)
Definition amp.h:773

## modGCDGF() [2/2]

CanonicalForm modGCDGF | ( | const CanonicalForm & | F,
|  | const CanonicalForm & | G,
|  | CFList & | l,
|  | bool & | topLevel )

Definition at line 859 of file cfModGcd.cc.

{
CanonicalForm dummy1, dummy2;
CanonicalForm result= modGCDGF (F, G, dummy1, dummy2, l, topLevel);
return result;
}

## monicSparseInterpol()

CanonicalForm monicSparseInterpol | ( | const CanonicalForm & | F,
|  | const CanonicalForm & | G,
|  | const CanonicalForm & | skeleton,
|  | const Variable & | alpha,
|  | bool & | fail,
|  | CFArray *& | coeffMonoms,
|  | CFArray & | Monoms )

Definition at line 2207 of file cfModGcd.cc.

{
CanonicalForm A= F;
CanonicalForm B= G;
if (F.isZero() && degree(G) > 0) return G/Lc(G);
else if (G.isZero() && degree (F) > 0) return F/Lc(F);
else if (F.isZero() && G.isZero()) return F.genOne();
if (F.inBaseDomain() || G.inBaseDomain()) return F.genOne();
if (F.isUnivariate() && fdivides(F, G)) return F/Lc(F);
if (G.isUnivariate() && fdivides(G, F)) return G/Lc(G);
if (F == G) return F/Lc(F);
ASSERT (degree (A, 1) == 0, "expected degree (F, 1) == 0");
ASSERT (degree (B, 1) == 0, "expected degree (G, 1) == 0");
CFMap M,N;
int best_level= myCompress (A, B, M, N, false);
if (best_level == 0)
return B.genOne();
A= M(A);
B= M(B);
Variable x= Variable (1);
//univariate case
if (A.isUnivariate() && B.isUnivariate())
return N (gcd (A, B));
CanonicalForm skel= M(skeleton);
CanonicalForm cA, cB; // content of A and B
CanonicalForm ppA, ppB; // primitive part of A and B
CanonicalForm gcdcAcB;
cA = uni_content (A);
cB = uni_content (B);
gcdcAcB= gcd (cA, cB);
ppA= A/cA;
ppB= B/cB;
CanonicalForm lcA, lcB; // leading coefficients of A and B
CanonicalForm gcdlcAlcB;
lcA= uni_lcoeff (ppA);
lcB= uni_lcoeff (ppB);
if (fdivides (lcA, lcB))
{
if (fdivides (A, B))
return F/Lc(F);
}
if (fdivides (lcB, lcA))
{
if (fdivides (B, A))
return G/Lc(G);
}
gcdlcAlcB= gcd (lcA, lcB);
int skelSize= size (skel, skel.mvar());
int j= 0;
int biggestSize= 0;
for (CFIterator i= skel; i.hasTerms(); i++, j++)
biggestSize= tmax (biggestSize, size (i.coeff()));
CanonicalForm g, Aeval, Beval;
CFList evalPoints;
bool evalFail= false;
CFList list;
bool GF= false;
CanonicalForm LCA= LC (A);
CanonicalForm tmp;
CFArray gcds= CFArray (biggestSize);
CFList * pEvalPoints= new CFList [biggestSize];
Variable V_buf= alpha, V_buf4= alpha;
CFList source, dest;
CanonicalForm prim_elem, im_prim_elem;
CanonicalForm prim_elem_alpha, im_prim_elem_alpha;
for (int i= 0; i < biggestSize; i++)
{
if (i == 0)
evalPoints= evaluationPoints (A, B, Aeval, Beval, LCA, GF, V_buf, evalFail,
list);
else
{
mult (evalPoints, pEvalPoints [0]);
eval (A, B, Aeval, Beval, evalPoints);
}
if (evalFail)
{
if (V_buf.level() != 1)
{
do
{
Variable V_buf3= V_buf;
V_buf= chooseExtension (V_buf);
source= CFList();
dest= CFList();
bool prim_fail= false;
Variable V_buf2;
prim_elem= primitiveElement (V_buf4, V_buf2, prim_fail);
if (V_buf4 == alpha && alpha.level() != 1)
prim_elem_alpha= prim_elem;
ASSERT (!prim_fail, "failure in integer factorizer");
if (prim_fail)
; //ERROR
else
im_prim_elem= mapPrimElem (prim_elem, V_buf4, V_buf);
DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (V_buf));
DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (V_buf2));
if (V_buf4 == alpha && alpha.level() != 1)
im_prim_elem_alpha= im_prim_elem;
else if (alpha.level() != 1)
im_prim_elem_alpha= mapUp (im_prim_elem_alpha, V_buf4, V_buf,
prim_elem, im_prim_elem, source, dest);
for (CFListIterator j= list; j.hasItem(); j++)
j.getItem()= mapUp (j.getItem(), V_buf4, V_buf, prim_elem,
im_prim_elem, source, dest);
for (int k= 0; k < i; k++)
{
for (CFListIterator j= pEvalPoints[k]; j.hasItem(); j++)
j.getItem()= mapUp (j.getItem(), V_buf4, V_buf, prim_elem,
im_prim_elem, source, dest);
gcds[k]= mapUp (gcds[k], V_buf4, V_buf, prim_elem, im_prim_elem,
source, dest);
}
if (alpha.level() != 1)
{
A= mapUp (A, V_buf4, V_buf, prim_elem, im_prim_elem, source,dest);
B= mapUp (B, V_buf4, V_buf, prim_elem, im_prim_elem, source,dest);
}
V_buf4= V_buf;
evalFail= false;
evalPoints= evaluationPoints (A, B, Aeval, Beval, LCA, GF, V_buf,
evalFail, list);
} while (evalFail);
}
else
{
CanonicalForm mipo;
int deg= 2;
bool initialized= false;
do
{
mipo= randomIrredpoly (deg, x);
if (initialized)
setMipo (V_buf, mipo);
else
V_buf= rootOf (mipo);
evalFail= false;
initialized= true;
evalPoints= evaluationPoints (A, B, Aeval, Beval, LCA, GF, V_buf,
evalFail, list);
deg++;
} while (evalFail);
V_buf4= V_buf;
}
}
g= gcd (Aeval, Beval);
g /= Lc (g);
if (degree (g) != degree (skel) || (skelSize != size (g)))
{
delete[] pEvalPoints;
fail= true;
if (alpha.level() != 1 && V_buf != alpha)
prune1 (alpha);
return 0;
}
CFIterator l= skel;
for (CFIterator k= g; k.hasTerms(); k++, l++)
{
if (k.exp() != l.exp())
{
delete[] pEvalPoints;
fail= true;
if (alpha.level() != 1 && V_buf != alpha)
prune1 (alpha);
return 0;
}
}
pEvalPoints[i]= evalPoints;
gcds[i]= g;
tmp= 0;
int j= 0;
for (CFListIterator k= evalPoints; k.hasItem(); k++, j++)
tmp += k.getItem()*power (x, j);
list.append (tmp);
}
if (Monoms.size() == 0)
Monoms= getMonoms (skel);
coeffMonoms= new CFArray [skelSize];
j= 0;
for (CFIterator i= skel; i.hasTerms(); i++, j++)
coeffMonoms[j]= getMonoms (i.coeff());
CFArray* pL= new CFArray [skelSize];
CFArray* pM= new CFArray [skelSize];
for (int i= 0; i < biggestSize; i++)
{
CFIterator l= gcds [i];
evalPoints= pEvalPoints [i];
for (int k= 0; k < skelSize; k++, l++)
{
if (i == 0)
pL[k]= CFArray (biggestSize);
pL[k] [i]= l.coeff();
if (i == 0)
pM[k]= evaluate (coeffMonoms [k], evalPoints);
}
}
CFArray solution;
CanonicalForm result= 0;
int ind= 0;
CFArray bufArray;
CFMatrix Mat;
for (int k= 0; k < skelSize; k++)
{
if (biggestSize != coeffMonoms[k].size())
{
bufArray= CFArray (coeffMonoms[k].size());
for (int i= 0; i < coeffMonoms[k].size(); i++)
bufArray [i]= pL[k] [i];
solution= solveGeneralVandermonde (pM [k], bufArray);
}
else
solution= solveGeneralVandermonde (pM [k], pL[k]);
if (solution.size() == 0)
{
delete[] pEvalPoints;
delete[] pM;
delete[] pL;
delete[] coeffMonoms;
fail= true;
if (alpha.level() != 1 && V_buf != alpha)
prune1 (alpha);
return 0;
}
for (int l= 0; l < solution.size(); l++)
result += solution[l]*Monoms [ind + l];
ind += solution.size();
}
delete[] pEvalPoints;
delete[] pM;
delete[] pL;
delete[] coeffMonoms;
if (alpha.level() != 1 && V_buf != alpha)
{
CFList u, v;
result= mapDown (result, prim_elem_alpha, im_prim_elem_alpha, alpha, u, v);
prune1 (alpha);
}
result= N(result);
if (fdivides (result, F) && fdivides (result, G))
return result;
else
{
fail= true;
return 0;
}
}
mult
void mult(CFList &L1, const CFList &L2)
multiply two lists componentwise
Definition cfModGcd.cc:2184
evaluate
CFArray evaluate(const CFArray &A, const CFList &evalPoints)
Definition cfModGcd.cc:2039
solveGeneralVandermonde
CFArray solveGeneralVandermonde(const CFArray &M, const CFArray &A)
Definition cfModGcd.cc:1633
evaluationPoints
CFList evaluationPoints(const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &Feval, CanonicalForm &Geval, const CanonicalForm &LCF, const bool &GF, const Variable &alpha, bool &fail, CFList &list)
Definition cfModGcd.cc:2056
Variable::level
int level() const
Definition factory.h:143
eval
CFList & eval
Definition facFactorize.cc:47

## mult()

void mult | ( | CFList & | L1,
|  | const CFList & | L2 )

multiply two lists componentwise

Definition at line 2184 of file cfModGcd.cc.

{
ASSERT (L1.length() == L2.length(), "lists of the same size expected");
CFListIterator j= L2;
for (CFListIterator i= L1; i.hasItem(); i++, j++)
i.getItem() *= j.getItem();
}

## 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 at line 92 of file cfModGcd.cc.

{
int n= tmax (F.level(), G.level());
int * degsf= NEW_ARRAY(int,n + 1);
int * degsg= NEW_ARRAY(int,n + 1);
for (int i = n; i >= 0; i--)
degsf[i]= degsg[i]= 0;
degsf= degrees (F, degsf);
degsg= degrees (G, degsg);
int both_non_zero= 0;
int f_zero= 0;
int g_zero= 0;
int both_zero= 0;
if (topLevel)
{
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 <= G.level())
{
f_zero++;
continue;
}
if (degsg[i] == 0 && degsf[i] && i <= F.level())
{
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 <= F.level())
{
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 <= G.level())
{
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 (F.level(), G.level());
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++;
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 (degsf[j] != 0 && degsg [j] != 0 &&
tmax (degsf[j],degsg[j]) <= min_max_deg)
{
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--;
}
}
else
{
//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));
}
}
}
}
DELETE_ARRAY(degsf);
DELETE_ARRAY(degsg);
return 1;
}
degrees
int * degrees(const CanonicalForm &f, int *degs=0)
int * degrees ( const CanonicalForm & f, int * degs )
Definition cf_ops.cc:493
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
g_zero
int g_zero
Definition cfEzgcd.cc:70
degsg
int * degsg
Definition cfEzgcd.cc:60
both_zero
int both_zero
Definition cfGcdAlgExt.cc:72
DELETE_ARRAY
#define DELETE_ARRAY(P)
Definition cf_defs.h:65
NEW_ARRAY
#define NEW_ARRAY(T, N)
Definition cf_defs.h:64

## newtonInterp()

CanonicalForm newtonInterp | ( | const CanonicalForm & | alpha,
|  | const CanonicalForm & | u,
|  | const CanonicalForm & | newtonPoly,
|  | const CanonicalForm & | oldInterPoly,
|  | const Variable & | x )

inline static

Newton interpolation - Incremental algorithm. Given a list of values alpha_i and a list of polynomials u_i, 1 <= i <= n, computes the interpolation polynomial assuming that the polynomials in u are the results of evaluating the variabe x of the unknown polynomial at the alpha values. This incremental version receives only the values of alpha_n and u_n and the previous interpolation polynomial for points 1 <= i <= n-1, and computes the polynomial interpolating in all the points. newtonPoly must be equal to (x - alpha_1) * ... * (x - alpha_{n-1}).

Definition at line 365 of file cfModGcd.cc.

{
CanonicalForm interPoly;
interPoly= oldInterPoly + ((u - oldInterPoly(alpha, x))/newtonPoly(alpha, x))
*newtonPoly;
return interPoly;
}

## nonMonicSparseInterpol()

CanonicalForm nonMonicSparseInterpol | ( | const CanonicalForm & | F,
|  | const CanonicalForm & | G,
|  | const CanonicalForm & | skeleton,
|  | const Variable & | alpha,
|  | bool & | fail,
|  | CFArray *& | coeffMonoms,
|  | CFArray & | Monoms )

Definition at line 2491 of file cfModGcd.cc.

{
CanonicalForm A= F;
CanonicalForm B= G;
if (F.isZero() && degree(G) > 0) return G/Lc(G);
else if (G.isZero() && degree (F) > 0) return F/Lc(F);
else if (F.isZero() && G.isZero()) return F.genOne();
if (F.inBaseDomain() || G.inBaseDomain()) return F.genOne();
if (F.isUnivariate() && fdivides(F, G)) return F/Lc(F);
if (G.isUnivariate() && fdivides(G, F)) return G/Lc(G);
if (F == G) return F/Lc(F);
ASSERT (degree (A, 1) == 0, "expected degree (F, 1) == 0");
ASSERT (degree (B, 1) == 0, "expected degree (G, 1) == 0");
CFMap M,N;
int best_level= myCompress (A, B, M, N, false);
if (best_level == 0)
return B.genOne();
A= M(A);
B= M(B);
Variable x= Variable (1);
//univariate case
if (A.isUnivariate() && B.isUnivariate())
return N (gcd (A, B));
CanonicalForm skel= M(skeleton);
CanonicalForm cA, cB; // content of A and B
CanonicalForm ppA, ppB; // primitive part of A and B
CanonicalForm gcdcAcB;
cA = uni_content (A);
cB = uni_content (B);
gcdcAcB= gcd (cA, cB);
ppA= A/cA;
ppB= B/cB;
CanonicalForm lcA, lcB; // leading coefficients of A and B
CanonicalForm gcdlcAlcB;
lcA= uni_lcoeff (ppA);
lcB= uni_lcoeff (ppB);
if (fdivides (lcA, lcB))
{
if (fdivides (A, B))
return F/Lc(F);
}
if (fdivides (lcB, lcA))
{
if (fdivides (B, A))
return G/Lc(G);
}
gcdlcAlcB= gcd (lcA, lcB);
int skelSize= size (skel, skel.mvar());
int j= 0;
int biggestSize= 0;
int bufSize;
int numberUni= 0;
for (CFIterator i= skel; i.hasTerms(); i++, j++)
{
bufSize= size (i.coeff());
biggestSize= tmax (biggestSize, bufSize);
numberUni += bufSize;
}
numberUni--;
numberUni= (int) ceil ( (double) numberUni / (skelSize - 1));
biggestSize= tmax (biggestSize , numberUni);
numberUni= biggestSize + size (LC(skel)) - 1;
int biggestSize2= tmax (numberUni, biggestSize);
CanonicalForm g, Aeval, Beval;
CFList evalPoints;
CFArray coeffEval;
bool evalFail= false;
CFList list;
bool GF= false;
CanonicalForm LCA= LC (A);
CanonicalForm tmp;
CFArray gcds= CFArray (biggestSize);
CFList * pEvalPoints= new CFList [biggestSize];
Variable V_buf= alpha, V_buf4= alpha;
CFList source, dest;
CanonicalForm prim_elem, im_prim_elem;
CanonicalForm prim_elem_alpha, im_prim_elem_alpha;
for (int i= 0; i < biggestSize; i++)
{
if (i == 0)
{
if (getCharacteristic() > 3)
evalPoints= evaluationPoints (A, B, Aeval, Beval, LCA, GF, V_buf,
evalFail, list);
else
evalFail= true;
if (evalFail)
{
if (V_buf.level() != 1)
{
do
{
Variable V_buf3= V_buf;
V_buf= chooseExtension (V_buf);
source= CFList();
dest= CFList();
bool prim_fail= false;
Variable V_buf2;
prim_elem= primitiveElement (V_buf4, V_buf2, prim_fail);
if (V_buf4 == alpha && alpha.level() != 1)
prim_elem_alpha= prim_elem;
ASSERT (!prim_fail, "failure in integer factorizer");
if (prim_fail)
; //ERROR
else
im_prim_elem= mapPrimElem (prim_elem, V_buf4, V_buf);
DEBOUTLN (cerr, "getMipo (V_buf)= " << getMipo (V_buf));
DEBOUTLN (cerr, "getMipo (V_buf2)= " << getMipo (V_buf2));
if (V_buf4 == alpha && alpha.level() != 1)
im_prim_elem_alpha= im_prim_elem;
else if (alpha.level() != 1)
im_prim_elem_alpha= mapUp (im_prim_elem_alpha, V_buf4, V_buf,
prim_elem, im_prim_elem, source, dest);
for (CFListIterator i= list; i.hasItem(); i++)
i.getItem()= mapUp (i.getItem(), V_buf4, V_buf, prim_elem,
im_prim_elem, source, dest);
if (alpha.level() != 1)
{
A= mapUp (A, V_buf4, V_buf, prim_elem, im_prim_elem, source,dest);
B= mapUp (B, V_buf4, V_buf, prim_elem, im_prim_elem, source,dest);
}
evalFail= false;
V_buf4= V_buf;
evalPoints= evaluationPoints (A, B, Aeval, Beval, LCA, GF, V_buf,
evalFail, list);
} while (evalFail);
}
else
{
CanonicalForm mipo;
int deg= 2;
bool initialized= false;
do
{
mipo= randomIrredpoly (deg, x);
if (initialized)
setMipo (V_buf, mipo);
else
V_buf= rootOf (mipo);
evalFail= false;
initialized= true;
evalPoints= evaluationPoints (A, B, Aeval, Beval, LCA, GF, V_buf,
evalFail, list);
deg++;
} while (evalFail);
V_buf4= V_buf;
}
}
}
else
{
mult (evalPoints, pEvalPoints[0]);
eval (A, B, Aeval, Beval, evalPoints);
}
g= gcd (Aeval, Beval);
g /= Lc (g);
if (degree (g) != degree (skel) || (skelSize != size (g)))
{
delete[] pEvalPoints;
fail= true;
if (alpha.level() != 1 && V_buf != alpha)
prune1 (alpha);
return 0;
}
CFIterator l= skel;
for (CFIterator k= g; k.hasTerms(); k++, l++)
{
if (k.exp() != l.exp())
{
delete[] pEvalPoints;
fail= true;
if (alpha.level() != 1 && V_buf != alpha)
prune1 (alpha);
return 0;
}
}
pEvalPoints[i]= evalPoints;
gcds[i]= g;
tmp= 0;
int j= 0;
for (CFListIterator k= evalPoints; k.hasItem(); k++, j++)
tmp += k.getItem()*power (x, j);
list.append (tmp);
}
if (Monoms.size() == 0)
Monoms= getMonoms (skel);
coeffMonoms= new CFArray [skelSize];
j= 0;
for (CFIterator i= skel; i.hasTerms(); i++, j++)
coeffMonoms[j]= getMonoms (i.coeff());
int minimalColumnsIndex;
if (skelSize > 1)
minimalColumnsIndex= 1;
else
minimalColumnsIndex= 0;
int minimalColumns=-1;
CFArray* pM= new CFArray [skelSize];
CFMatrix Mat;
// find the Matrix with minimal number of columns
for (int i= 0; i < skelSize; i++)
{
pM[i]= CFArray (coeffMonoms[i].size());
if (i == 1)
minimalColumns= coeffMonoms[i].size();
if (i > 1)
{
if (minimalColumns > coeffMonoms[i].size())
{
minimalColumns= coeffMonoms[i].size();
minimalColumnsIndex= i;
}
}
}
CFMatrix* pMat= new CFMatrix [2];
pMat[0]= CFMatrix (biggestSize, coeffMonoms[0].size());
pMat[1]= CFMatrix (biggestSize, coeffMonoms[minimalColumnsIndex].size());
CFArray* pL= new CFArray [skelSize];
for (int i= 0; i < biggestSize; i++)
{
CFIterator l= gcds [i];
evalPoints= pEvalPoints [i];
for (int k= 0; k < skelSize; k++, l++)
{
if (i == 0)
pL[k]= CFArray (biggestSize);
pL[k] [i]= l.coeff();
if (i == 0 && k != 0 && k != minimalColumnsIndex)
pM[k]= evaluate (coeffMonoms[k], evalPoints);
else if (i == 0 && (k == 0 || k == minimalColumnsIndex))
coeffEval= evaluate (coeffMonoms[k], evalPoints);
if (i > 0 && (k == 0 || k == minimalColumnsIndex))
coeffEval= evaluate (coeffMonoms[k], evalPoints);
if (k == 0)
{
if (pMat[k].rows() >= i + 1)
{
for (int ind= 1; ind < coeffEval.size() + 1; ind++)
pMat[k] (i + 1, ind)= coeffEval[ind - 1];
}
}
if (k == minimalColumnsIndex)
{
if (pMat[1].rows() >= i + 1)
{
for (int ind= 1; ind < coeffEval.size() + 1; ind++)
pMat[1] (i + 1, ind)= coeffEval[ind - 1];
}
}
}
}
CFArray solution;
CanonicalForm result= 0;
int ind= 1;
int matRows, matColumns;
matRows= pMat[1].rows();
matColumns= pMat[0].columns() - 1;
matColumns += pMat[1].columns();
Mat= CFMatrix (matRows, matColumns);
for (int i= 1; i <= matRows; i++)
for (int j= 1; j <= pMat[1].columns(); j++)
Mat (i, j)= pMat[1] (i, j);
for (int j= pMat[1].columns() + 1; j <= matColumns;
j++, ind++)
{
for (int i= 1; i <= matRows; i++)
Mat (i, j)= (-pMat [0] (i, ind + 1))*pL[minimalColumnsIndex] [i - 1];
}
CFArray firstColumn= CFArray (pMat[0].rows());
for (int i= 0; i < pMat[0].rows(); i++)
firstColumn [i]= pMat[0] (i + 1, 1);
CFArray bufArray= pL[minimalColumnsIndex];
for (int i= 0; i < biggestSize; i++)
pL[minimalColumnsIndex] [i] *= firstColumn[i];
CFMatrix bufMat= pMat[1];
pMat[1]= Mat;
if (V_buf.level() != 1)
solution= solveSystemFq (pMat[1],
pL[minimalColumnsIndex], V_buf);
else
solution= solveSystemFp (pMat[1],
pL[minimalColumnsIndex]);
if (solution.size() == 0)
{ //Javadi's method failed, try deKleine, Monagan, Wittkopf instead
CFMatrix bufMat0= pMat[0];
delete [] pMat;
pMat= new CFMatrix [skelSize];
pL[minimalColumnsIndex]= bufArray;
CFList* bufpEvalPoints= NULL;
CFArray bufGcds;
if (biggestSize != biggestSize2)
{
bufpEvalPoints= pEvalPoints;
pEvalPoints= new CFList [biggestSize2];
bufGcds= gcds;
gcds= CFArray (biggestSize2);
for (int i= 0; i < biggestSize; i++)
{
pEvalPoints[i]= bufpEvalPoints [i];
gcds[i]= bufGcds[i];
}
for (int i= 0; i < biggestSize2 - biggestSize; i++)
{
mult (evalPoints, pEvalPoints[0]);
eval (A, B, Aeval, Beval, evalPoints);
g= gcd (Aeval, Beval);
g /= Lc (g);
if (degree (g) != degree (skel) || (skelSize != size (g)))
{
delete[] pEvalPoints;
delete[] pMat;
delete[] pL;
delete[] coeffMonoms;
delete[] pM;
if (bufpEvalPoints != NULL)
delete [] bufpEvalPoints;
fail= true;
if (alpha.level() != 1 && V_buf != alpha)
prune1 (alpha);
return 0;
}
CFIterator l= skel;
for (CFIterator k= g; k.hasTerms(); k++, l++)
{
if (k.exp() != l.exp())
{
delete[] pEvalPoints;
delete[] pMat;
delete[] pL;
delete[] coeffMonoms;
delete[] pM;
if (bufpEvalPoints != NULL)
delete [] bufpEvalPoints;
fail= true;
if (alpha.level() != 1 && V_buf != alpha)
prune1 (alpha);
return 0;
}
}
pEvalPoints[i + biggestSize]= evalPoints;
gcds[i + biggestSize]= g;
}
}
for (int i= 0; i < biggestSize; i++)
{
CFIterator l= gcds [i];
evalPoints= pEvalPoints [i];
for (int k= 1; k < skelSize; k++, l++)
{
if (i == 0)
pMat[k]= CFMatrix (biggestSize2,coeffMonoms[k].size()+biggestSize2-1);
if (k == minimalColumnsIndex)
continue;
coeffEval= evaluate (coeffMonoms[k], evalPoints);
if (pMat[k].rows() >= i + 1)
{
for (int ind= 1; ind < coeffEval.size() + 1; ind++)
pMat[k] (i + 1, ind)= coeffEval[ind - 1];
}
}
}
Mat= bufMat0;
pMat[0]= CFMatrix (biggestSize2, coeffMonoms[0].size() + biggestSize2 - 1);
for (int j= 1; j <= Mat.rows(); j++)
for (int k= 1; k <= Mat.columns(); k++)
pMat [0] (j,k)= Mat (j,k);
Mat= bufMat;
for (int j= 1; j <= Mat.rows(); j++)
for (int k= 1; k <= Mat.columns(); k++)
pMat [minimalColumnsIndex] (j,k)= Mat (j,k);
// write old matrix entries into new matrices
for (int i= 0; i < skelSize; i++)
{
bufArray= pL[i];
pL[i]= CFArray (biggestSize2);
for (int j= 0; j < bufArray.size(); j++)
pL[i] [j]= bufArray [j];
}
//write old vector entries into new and add new entries to old matrices
for (int i= 0; i < biggestSize2 - biggestSize; i++)
{
CFIterator l= gcds [i + biggestSize];
evalPoints= pEvalPoints [i + biggestSize];
for (int k= 0; k < skelSize; k++, l++)
{
pL[k] [i + biggestSize]= l.coeff();
coeffEval= evaluate (coeffMonoms[k], evalPoints);
if (pMat[k].rows() >= i + biggestSize + 1)
{
for (int ind= 1; ind < coeffEval.size() + 1; ind++)
pMat[k] (i + biggestSize + 1, ind)= coeffEval[ind - 1];
}
}
}
// begin new
for (int i= 0; i < skelSize; i++)
{
if (pL[i].size() > 1)
{
for (int j= 2; j <= pMat[i].rows(); j++)
pMat[i] (j, coeffMonoms[i].size() + j - 1)=
-pL[i] [j - 1];
}
}
matColumns= biggestSize2 - 1;
matRows= 0;
for (int i= 0; i < skelSize; i++)
{
if (V_buf.level() == 1)
(void) gaussianElimFp (pMat[i], pL[i]);
else
(void) gaussianElimFq (pMat[i], pL[i], V_buf);
if (pMat[i] (coeffMonoms[i].size(), coeffMonoms[i].size()) == 0)
{
delete[] pEvalPoints;
delete[] pMat;
delete[] pL;
delete[] coeffMonoms;
delete[] pM;
if (bufpEvalPoints != NULL)
delete [] bufpEvalPoints;
fail= true;
if (alpha.level() != 1 && V_buf != alpha)
prune1 (alpha);
return 0;
}
matRows += pMat[i].rows() - coeffMonoms[i].size();
}
CFMatrix bufMat;
Mat= CFMatrix (matRows, matColumns);
ind= 0;
bufArray= CFArray (matRows);
CFArray bufArray2;
for (int i= 0; i < skelSize; i++)
{
if (coeffMonoms[i].size() + 1 >= pMat[i].rows() || coeffMonoms[i].size() + 1 >= pMat[i].columns())
{
delete[] pEvalPoints;
delete[] pMat;
delete[] pL;
delete[] coeffMonoms;
delete[] pM;
if (bufpEvalPoints != NULL)
delete [] bufpEvalPoints;
fail= true;
if (alpha.level() != 1 && V_buf != alpha)
prune1 (alpha);
return 0;
}
bufMat= pMat[i] (coeffMonoms[i].size() + 1, pMat[i].rows(),
coeffMonoms[i].size() + 1, pMat[i].columns());
for (int j= 1; j <= bufMat.rows(); j++)
for (int k= 1; k <= bufMat.columns(); k++)
Mat (j + ind, k)= bufMat(j, k);
bufArray2= coeffMonoms[i].size();
for (int j= 1; j <= pMat[i].rows(); j++)
{
if (j > coeffMonoms[i].size())
bufArray [j-coeffMonoms[i].size() + ind - 1]= pL[i] [j - 1];
else
bufArray2 [j - 1]= pL[i] [j - 1];
}
pL[i]= bufArray2;
ind += bufMat.rows();
pMat[i]= pMat[i] (1, coeffMonoms[i].size(), 1, pMat[i].columns());
}
if (V_buf.level() != 1)
solution= solveSystemFq (Mat, bufArray, V_buf);
else
solution= solveSystemFp (Mat, bufArray);
if (solution.size() == 0)
{
delete[] pEvalPoints;
delete[] pMat;
delete[] pL;
delete[] coeffMonoms;
delete[] pM;
if (bufpEvalPoints != NULL)
delete [] bufpEvalPoints;
fail= true;
if (alpha.level() != 1 && V_buf != alpha)
prune1 (alpha);
return 0;
}
ind= 0;
result= 0;
CFArray bufSolution;
for (int i= 0; i < skelSize; i++)
{
bufSolution= readOffSolution (pMat[i], pL[i], solution);
for (int i= 0; i < bufSolution.size(); i++, ind++)
result += Monoms [ind]*bufSolution[i];
}
if (alpha.level() != 1 && V_buf != alpha)
{
CFList u, v;
result= mapDown (result,prim_elem_alpha, im_prim_elem_alpha, alpha, u, v);
prune1 (alpha);
}
result= N(result);
delete[] pEvalPoints;
delete[] pMat;
delete[] pL;
delete[] coeffMonoms;
delete[] pM;
if (bufpEvalPoints != NULL)
delete [] bufpEvalPoints;
if (fdivides (result, F) && fdivides (result, G))
return result;
else
{
fail= true;
return 0;
}
} // end of deKleine, Monagan & Wittkopf
result += Monoms[0];
int ind2= 0, ind3= 2;
ind= 0;
for (int l= 0; l < minimalColumnsIndex; l++)
ind += coeffMonoms[l].size();
for (int l= solution.size() - pMat[0].columns() + 1; l < solution.size();
l++, ind2++, ind3++)
{
result += solution[l]*Monoms [1 + ind2];
for (int i= 0; i < pMat[0].rows(); i++)
firstColumn[i] += solution[l]*pMat[0] (i + 1, ind3);
}
for (int l= 0; l < solution.size() + 1 - pMat[0].columns(); l++)
result += solution[l]*Monoms [ind + l];
ind= coeffMonoms[0].size();
for (int k= 1; k < skelSize; k++)
{
if (k == minimalColumnsIndex)
{
ind += coeffMonoms[k].size();
continue;
}
if (k != minimalColumnsIndex)
{
for (int i= 0; i < biggestSize; i++)
pL[k] [i] *= firstColumn [i];
}
if (biggestSize != coeffMonoms[k].size() && k != minimalColumnsIndex)
{
bufArray= CFArray (coeffMonoms[k].size());
for (int i= 0; i < bufArray.size(); i++)
bufArray [i]= pL[k] [i];
solution= solveGeneralVandermonde (pM[k], bufArray);
}
else
solution= solveGeneralVandermonde (pM[k], pL[k]);
if (solution.size() == 0)
{
delete[] pEvalPoints;
delete[] pMat;
delete[] pL;
delete[] coeffMonoms;
delete[] pM;
fail= true;
if (alpha.level() != 1 && V_buf != alpha)
prune1 (alpha);
return 0;
}
if (k != minimalColumnsIndex)
{
for (int l= 0; l < solution.size(); l++)
result += solution[l]*Monoms [ind + l];
ind += solution.size();
}
}
delete[] pEvalPoints;
delete[] pMat;
delete[] pL;
delete[] pM;
delete[] coeffMonoms;
if (alpha.level() != 1 && V_buf != alpha)
{
CFList u, v;
result= mapDown (result, prim_elem, im_prim_elem, alpha, u, v);
prune1 (alpha);
}
result= N(result);
if (fdivides (result, F) && fdivides (result, G))
return result;
else
{
fail= true;
return 0;
}
}
readOffSolution
CFArray readOffSolution(const CFMatrix &M, const long rk)
M in row echolon form, rk rank of M.
Definition cfModGcd.cc:1689
gaussianElimFq
long gaussianElimFq(CFMatrix &M, CFArray &L, const Variable &alpha)
Definition cfModGcd.cc:1785
gaussianElimFp
long gaussianElimFp(CFMatrix &M, CFArray &L)
Definition cfModGcd.cc:1740
solveSystemFp
CFArray solveSystemFp(const CFMatrix &M, const CFArray &L)
Definition cfModGcd.cc:1844
solveSystemFq
CFArray solveSystemFq(const CFMatrix &M, const CFArray &L, const Variable &alpha)
Definition cfModGcd.cc:1896
Matrix::rows
int rows() const
Definition ftmpl_matrix.h:43
Matrix::columns
int columns() const
Definition ftmpl_matrix.h:44
NULL
#define NULL
Definition omList.c:12

## randomElement()

CanonicalForm randomElement | ( | const CanonicalForm & | F,
|  | const Variable & | alpha,
|  | CFList & | list,
|  | bool & | fail )

inline static

compute a random element a of , s.t. F(a) , F is a univariate polynomial, returns fail if there are no field elements left which have not been used before

Definition at line 380 of file cfModGcd.cc.

{
fail= false;
Variable x= F.mvar();
AlgExtRandomF genAlgExt (alpha);
FFRandom genFF;
CanonicalForm random, mipo;
mipo= getMipo (alpha);
int p= getCharacteristic ();
int d= degree (mipo);
double bound= pow ((double) p, (double) d);
do
{
if (list.length() == bound)
{
fail= true;
break;
}
if (list.length() < p)
{
random= genFF.generate();
while (find (list, random))
random= genFF.generate();
}
else
{
random= genAlgExt.generate();
while (find (list, random))
random= genAlgExt.generate();
}
if (F (random, x) == 0)
{
list.append (random);
continue;
}
} while (find (list, random));
return random;
}

## readOffSolution() [1/2]

CFArray readOffSolution | ( | const CFMatrix & | M,
|  | const CFArray & | L,
|  | const CFArray & | partialSol )

Definition at line 1711 of file cfModGcd.cc.

{
CFArray result= CFArray (M.rows());
CanonicalForm tmp1, tmp2, tmp3;
int k;
for (int i= M.rows(); i >= 1; i--)
{
tmp3= 0;
tmp1= L[i - 1];
k= 0;
for (int j= M.columns(); j >= 1; j--, k++)
{
tmp2= M (i, j);
if (j == i)
break;
else
{
if (k > partialSol.size() - 1)
tmp3 += tmp2*result[j - 1];
else
tmp3 += tmp2*partialSol[partialSol.size() - k - 1];
}
}
result [i - 1]= (tmp1 - tmp3)/tmp2;
}
return result;
}
tmp1
CFList tmp1
Definition facFqBivar.cc:75
tmp2
CFList tmp2
Definition facFqBivar.cc:75

## readOffSolution() [2/2]

CFArray readOffSolution | ( | const CFMatrix & | M,
|  | const long | rk )

M in row echolon form, rk rank of M.

Definition at line 1689 of file cfModGcd.cc.

{
CFArray result= CFArray (rk);
CanonicalForm tmp1, tmp2, tmp3;
for (int i= rk; i >= 1; i--)
{
tmp3= 0;
tmp1= M (i, M.columns());
for (int j= M.columns() - 1; j >= 1; j--)
{
tmp2= M (i, j);
if (j == i)
break;
else
tmp3 += tmp2*result[j - 1];
}
result[i - 1]= (tmp1 - tmp3)/tmp2;
}
return result;
}

## solveGeneralVandermonde()

CFArray solveGeneralVandermonde | ( | const CFArray & | M,
|  | const CFArray & | A )

Definition at line 1633 of file cfModGcd.cc.

{
int r= M.size();
ASSERT (A.size() == r, "vector does not have right size");
if (r == 1)
{
CFArray result= CFArray (1);
result [0]= A[0] / M [0];
return result;
}
// check solvability
bool notDistinct= false;
for (int i= 0; i < r - 1; i++)
{
for (int j= i + 1; j < r; j++)
{
if (M [i] == M [j])
{
notDistinct= true;
break;
}
}
}
if (notDistinct)
return CFArray();
CanonicalForm master= 1;
Variable x= Variable (1);
for (int i= 0; i < r; i++)
master *= x - M [i];
master *= x;
CFList Pj;
CanonicalForm tmp;
for (int i= 0; i < r; i++)
{
tmp= master/(x - M [i]);
tmp /= tmp (M [i], 1);
Pj.append (tmp);
}
CFArray result= CFArray (r);
CFListIterator j= Pj;
for (int i= 1; i <= r; i++, j++)
{
tmp= 0;
for (int l= 1; l <= A.size(); l++)
tmp += A[l - 1]*j.getItem()[l];
result[i - 1]= tmp;
}
return result;
}

## solveSystemFp()

CFArray solveSystemFp | ( | const CFMatrix & | M,
|  | const CFArray & | L )

Definition at line 1844 of file cfModGcd.cc.

{
ASSERT (L.size() <= M.rows(), "dimension exceeded");
CFMatrix *N;
N= new CFMatrix (M.rows(), M.columns() + 1);
for (int i= 1; i <= M.rows(); i++)
for (int j= 1; j <= M.columns(); j++)
(*N) (i, j)= M (i, j);
int j= 1;
for (int i= 0; i < L.size(); i++, j++)
(*N) (j, M.columns() + 1)= L[i];
#ifdef HAVE_FLINT
nmod_mat_t FLINTN;
convertFacCFMatrix2nmod_mat_t (FLINTN, *N);
long rk= nmod_mat_rref (FLINTN);
#else
int p= getCharacteristic ();
if (fac_NTL_char != p)
{
fac_NTL_char= p;
zz_p::init (p);
}
mat_zz_p *NTLN= convertFacCFMatrix2NTLmat_zz_p(*N);
long rk= gauss (*NTLN);
#endif
delete N;
if (rk != M.columns())
{
#ifdef HAVE_FLINT
nmod_mat_clear (FLINTN);
#else
delete NTLN;
#endif
return CFArray();
}
#ifdef HAVE_FLINT
N= convertNmod_mat_t2FacCFMatrix (FLINTN);
nmod_mat_clear (FLINTN);
#else
N= convertNTLmat_zz_p2FacCFMatrix (*NTLN);
delete NTLN;
#endif
CFArray A= readOffSolution (*N, rk);
delete N;
return A;
}

## solveSystemFq()

CFArray solveSystemFq | ( | const CFMatrix & | M,
|  | const CFArray & | L,
|  | const Variable & | alpha )

Definition at line 1896 of file cfModGcd.cc.

{
ASSERT (L.size() <= M.rows(), "dimension exceeded");
CFMatrix *N;
N= new CFMatrix (M.rows(), M.columns() + 1);
for (int i= 1; i <= M.rows(); i++)
for (int j= 1; j <= M.columns(); j++)
(*N) (i, j)= M (i, j);
int j= 1;
for (int i= 0; i < L.size(); i++, j++)
(*N) (j, M.columns() + 1)= L[i];
#ifdef HAVE_FLINT
// convert mipo
nmod_poly_t mipo1;
convertFacCF2nmod_poly_t(mipo1,getMipo(alpha));
fq_nmod_ctx_t ctx;
fq_nmod_ctx_init_modulus(ctx,mipo1,"t");
nmod_poly_clear(mipo1);
// convert matrix
fq_nmod_mat_t FLINTN;
convertFacCFMatrix2Fq_nmod_mat_t (FLINTN, ctx, *N);
// rank
#if __FLINT_RELEASE >= 30100
long rk= fq_nmod_mat_rref (FLINTN,FLINTN,ctx);
#else
long rk= fq_nmod_mat_rref (FLINTN,ctx);
#endif
#elif defined(HAVE_NTL)
int p= getCharacteristic ();
if (fac_NTL_char != p)
{
fac_NTL_char= p;
zz_p::init (p);
}
zz_pX NTLMipo= convertFacCF2NTLzzpX (getMipo (alpha));
zz_pE::init (NTLMipo);
mat_zz_pE *NTLN= convertFacCFMatrix2NTLmat_zz_pE(*N);
long rk= gauss (*NTLN);
#else
factoryError("NTL/FLINT missing: solveSystemFq");
#endif
delete N;
if (rk != M.columns())
{
#if defined(HAVE_NTL) && !defined(HAVE_FLINT)
delete NTLN;
#endif
return CFArray();
}
#ifdef HAVE_FLINT
// convert and clean up
N=convertFq_nmod_mat_t2FacCFMatrix(FLINTN,ctx,alpha);
fq_nmod_mat_clear (FLINTN,ctx);
fq_nmod_ctx_clear(ctx);
#elif defined(HAVE_NTL)
N= convertNTLmat_zz_pE2FacCFMatrix (*NTLN, alpha);
delete NTLN;
#endif
CFArray A= readOffSolution (*N, rk);
delete N;
return A;
}
convertFq_nmod_mat_t2FacCFMatrix
CFMatrix * convertFq_nmod_mat_t2FacCFMatrix(const fq_nmod_mat_t m, const fq_nmod_ctx_t &fq_con, const Variable &alpha)
conversion of a FLINT matrix over F_q to a factory matrix
Definition FLINTconvert.cc:661

## solveVandermonde()

CFArray solveVandermonde | ( | const CFArray & | M,
|  | const CFArray & | A )

Definition at line 1580 of file cfModGcd.cc.

{
int r= M.size();
ASSERT (A.size() == r, "vector does not have right size");
if (r == 1)
{
CFArray result= CFArray (1);
result [0]= A [0] / M [0];
return result;
}
// check solvability
bool notDistinct= false;
for (int i= 0; i < r - 1; i++)
{
for (int j= i + 1; j < r; j++)
{
if (M [i] == M [j])
{
notDistinct= true;
break;
}
}
}
if (notDistinct)
return CFArray();
CanonicalForm master= 1;
Variable x= Variable (1);
for (int i= 0; i < r; i++)
master *= x - M [i];
CFList Pj;
CanonicalForm tmp;
for (int i= 0; i < r; i++)
{
tmp= master/(x - M [i]);
tmp /= tmp (M [i], 1);
Pj.append (tmp);
}
CFArray result= CFArray (r);
CFListIterator j= Pj;
for (int i= 1; i <= r; i++, j++)
{
tmp= 0;
for (int l= 0; l < A.size(); l++)
tmp += A[l]*j.getItem()[l];
result[i - 1]= tmp;
}
return result;
}

## sparseGCDFp()

CanonicalForm sparseGCDFp | ( | const CanonicalForm & | F,
|  | const CanonicalForm & | G,
|  | bool & | topLevel,
|  | CFList & | l )

Definition at line 3570 of file cfModGcd.cc.

{
CanonicalForm A= F;
CanonicalForm B= G;
if (F.isZero() && degree(G) > 0) return G/Lc(G);
else if (G.isZero() && degree (F) > 0) return F/Lc(F);
else if (F.isZero() && G.isZero()) return F.genOne();
if (F.inBaseDomain() || G.inBaseDomain()) return F.genOne();
if (F.isUnivariate() && fdivides(F, G)) return F/Lc(F);
if (G.isUnivariate() && fdivides(G, F)) return G/Lc(G);
if (F == G) return F/Lc(F);
CFMap M,N;
int best_level= myCompress (A, B, M, N, topLevel);
if (best_level == 0) return B.genOne();
A= M(A);
B= M(B);
Variable x= Variable (1);
//univariate case
if (A.isUnivariate() && B.isUnivariate())
return N (gcd (A, B));
CanonicalForm cA, cB; // content of A and B
CanonicalForm ppA, ppB; // primitive part of A and B
CanonicalForm gcdcAcB;
cA = uni_content (A);
cB = uni_content (B);
gcdcAcB= gcd (cA, cB);
ppA= A/cA;
ppB= B/cB;
CanonicalForm lcA, lcB; // leading coefficients of A and B
CanonicalForm gcdlcAlcB;
lcA= uni_lcoeff (ppA);
lcB= uni_lcoeff (ppB);
if (fdivides (lcA, lcB))
{
if (fdivides (A, B))
return F/Lc(F);
}
if (fdivides (lcB, lcA))
{
if (fdivides (B, A))
return G/Lc(G);
}
gcdlcAlcB= gcd (lcA, lcB);
int d= totaldegree (ppA, Variable (2), Variable (ppA.level()));
int d0;
if (d == 0)
return N(gcdcAcB);
d0= totaldegree (ppB, Variable (2), Variable (ppB.level()));
if (d0 < d)
d= d0;
if (d == 0)
return N(gcdcAcB);
CanonicalForm m, random_element, G_m, G_random_element, H, cH, ppH, skeleton;
CanonicalForm newtonPoly= 1;
m= gcdlcAlcB;
G_m= 0;
H= 0;
bool fail= false;
topLevel= false;
bool inextension= false;
Variable V_buf, alpha, cleanUp;
CanonicalForm prim_elem, im_prim_elem;
CFList source, dest;
do //first do
{
if (inextension)
random_element= randomElement (m, V_buf, l, fail);
else
random_element= FpRandomElement (m, l, fail);
if (random_element == 0 && !fail)
{
if (!find (l, random_element))
l.append (random_element);
continue;
}
if (!fail && !inextension)
{
CFList list;
TIMING_START (gcd_recursion);
G_random_element=
sparseGCDFp (ppA (random_element,x), ppB (random_element,x), topLevel,
list);
TIMING_END_AND_PRINT (gcd_recursion,
"time for recursive call: ");
DEBOUTLN (cerr, "G_random_element= " << G_random_element);
}
else if (!fail && inextension)
{
CFList list;
TIMING_START (gcd_recursion);
G_random_element=
sparseGCDFq (ppA (random_element, x), ppB (random_element, x), alpha,
list, topLevel);
TIMING_END_AND_PRINT (gcd_recursion,
"time for recursive call: ");
DEBOUTLN (cerr, "G_random_element= " << G_random_element);
}
else if (fail && !inextension)
{
source= CFList();
dest= CFList();
CFList list;
CanonicalForm mipo;
int deg= 2;
bool initialized= false;
do
{
mipo= randomIrredpoly (deg, x);
if (initialized)
setMipo (alpha, mipo);
else
alpha= rootOf (mipo);
inextension= true;
fail= false;
initialized= true;
random_element= randomElement (m, alpha, l, fail);
deg++;
} while (fail);
cleanUp= alpha;
V_buf= alpha;
list= CFList();
TIMING_START (gcd_recursion);
G_random_element=
sparseGCDFq (ppA (random_element, x), ppB (random_element, x), alpha,
list, topLevel);
TIMING_END_AND_PRINT (gcd_recursion,
"time for recursive call: ");
DEBOUTLN (cerr, "G_random_element= " << G_random_element);
}
else if (fail && inextension)
{
source= CFList();
dest= CFList();
Variable V_buf3= V_buf;
V_buf= chooseExtension (V_buf);
bool prim_fail= false;
Variable V_buf2;
CanonicalForm prim_elem, im_prim_elem;
prim_elem= primitiveElement (alpha, V_buf2, prim_fail);
if (V_buf3 != alpha)
{
m= mapDown (m, prim_elem, im_prim_elem, alpha, source, dest);
G_m= mapDown (m, prim_elem, im_prim_elem, alpha, source, dest);
newtonPoly= mapDown (newtonPoly, prim_elem, im_prim_elem, alpha, source,
dest);
ppA= mapDown (ppA, prim_elem, im_prim_elem, alpha, source, dest);
ppB= mapDown (ppB, prim_elem, im_prim_elem, alpha, source, dest);
gcdlcAlcB= mapDown (gcdlcAlcB, prim_elem, im_prim_elem, alpha, source,
dest);
for (CFListIterator i= l; i.hasItem(); i++)
i.getItem()= mapDown (i.getItem(), prim_elem, im_prim_elem, alpha,
source, dest);
}
ASSERT (!prim_fail, "failure in integer factorizer");
if (prim_fail)
; //ERROR
else
im_prim_elem= mapPrimElem (prim_elem, alpha, V_buf);
DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (alpha));
DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (V_buf2));
for (CFListIterator i= l; i.hasItem(); i++)
i.getItem()= mapUp (i.getItem(), alpha, V_buf, prim_elem,
im_prim_elem, source, dest);
m= mapUp (m, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
G_m= mapUp (G_m, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
newtonPoly= mapUp (newtonPoly, alpha, V_buf, prim_elem, im_prim_elem,
source, dest);
ppA= mapUp (ppA, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
ppB= mapUp (ppB, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
gcdlcAlcB= mapUp (gcdlcAlcB, alpha, V_buf, prim_elem, im_prim_elem,
source, dest);
fail= false;
random_element= randomElement (m, V_buf, l, fail );
DEBOUTLN (cerr, "fail= " << fail);
CFList list;
TIMING_START (gcd_recursion);
G_random_element=
sparseGCDFq (ppA (random_element, x), ppB (random_element, x), V_buf,
list, topLevel);
TIMING_END_AND_PRINT (gcd_recursion,
"time for recursive call: ");
DEBOUTLN (cerr, "G_random_element= " << G_random_element);
alpha= V_buf;
}
if (!G_random_element.inCoeffDomain())
d0= totaldegree (G_random_element, Variable(2),
Variable (G_random_element.level()));
else
d0= 0;
if (d0 == 0)
{
if (inextension)
prune (cleanUp);
return N(gcdcAcB);
}
if (d0 > d)
{
if (!find (l, random_element))
l.append (random_element);
continue;
}
G_random_element=
(gcdlcAlcB(random_element, x)/uni_lcoeff (G_random_element))
* G_random_element;
skeleton= G_random_element;
if (!G_random_element.inCoeffDomain())
d0= totaldegree (G_random_element, Variable(2),
Variable (G_random_element.level()));
else
d0= 0;
if (d0 < d)
{
m= gcdlcAlcB;
newtonPoly= 1;
G_m= 0;
d= d0;
}
H= newtonInterp (random_element, G_random_element, newtonPoly, G_m, x);
if (uni_lcoeff (H) == gcdlcAlcB)
{
cH= uni_content (H);
ppH= H/cH;
ppH /= Lc (ppH);
DEBOUTLN (cerr, "ppH= " << ppH);
if (fdivides (ppH, ppA) && fdivides (ppH, ppB))
{
if (inextension)
prune (cleanUp);
return N(gcdcAcB*ppH);
}
}
G_m= H;
newtonPoly= newtonPoly*(x - random_element);
m= m*(x - random_element);
if (!find (l, random_element))
l.append (random_element);
if ((getCharacteristic() > 3 && size (skeleton) < 200))
{
CFArray Monoms;
CFArray* coeffMonoms;
do //second do
{
if (inextension)
random_element= randomElement (m, alpha, l, fail);
else
random_element= FpRandomElement (m, l, fail);
if (random_element == 0 && !fail)
{
if (!find (l, random_element))
l.append (random_element);
continue;
}
bool sparseFail= false;
if (!fail && !inextension)
{
CFList list;
TIMING_START (gcd_recursion);
if (LC (skeleton).inCoeffDomain())
G_random_element=
monicSparseInterpol(ppA(random_element, x), ppB (random_element, x),
skeleton, x, sparseFail, coeffMonoms,
Monoms);
else
G_random_element=
nonMonicSparseInterpol(ppA(random_element,x), ppB(random_element,x),
skeleton, x, sparseFail,
coeffMonoms, Monoms);
TIMING_END_AND_PRINT (gcd_recursion,
"time for recursive call: ");
DEBOUTLN (cerr, "G_random_element= " << G_random_element);
}
else if (!fail && inextension)
{
CFList list;
TIMING_START (gcd_recursion);
if (LC (skeleton).inCoeffDomain())
G_random_element=
monicSparseInterpol(ppA (random_element,x), ppB (random_element, x),
skeleton, alpha, sparseFail, coeffMonoms,
Monoms);
else
G_random_element=
nonMonicSparseInterpol(ppA(random_element,x), ppB(random_element,x),
skeleton, alpha, sparseFail, coeffMonoms,
Monoms);
TIMING_END_AND_PRINT (gcd_recursion,
"time for recursive call: ");
DEBOUTLN (cerr, "G_random_element= " << G_random_element);
}
else if (fail && !inextension)
{
source= CFList();
dest= CFList();
CFList list;
CanonicalForm mipo;
int deg= 2;
bool initialized= false;
do
{
mipo= randomIrredpoly (deg, x);
if (initialized)
setMipo (alpha, mipo);
else
alpha= rootOf (mipo);
inextension= true;
fail= false;
initialized= true;
random_element= randomElement (m, alpha, l, fail);
deg++;
} while (fail);
cleanUp= alpha;
V_buf= alpha;
list= CFList();
TIMING_START (gcd_recursion);
if (LC (skeleton).inCoeffDomain())
G_random_element=
monicSparseInterpol (ppA (random_element,x), ppB (random_element,x),
skeleton, alpha, sparseFail, coeffMonoms,
Monoms);
else
G_random_element=
nonMonicSparseInterpol(ppA(random_element,x), ppB(random_element,x),
skeleton, alpha, sparseFail, coeffMonoms,
Monoms);
TIMING_END_AND_PRINT (gcd_recursion,
"time for recursive call: ");
DEBOUTLN (cerr, "G_random_element= " << G_random_element);
}
else if (fail && inextension)
{
source= CFList();
dest= CFList();
Variable V_buf3= V_buf;
V_buf= chooseExtension (V_buf);
bool prim_fail= false;
Variable V_buf2;
CanonicalForm prim_elem, im_prim_elem;
prim_elem= primitiveElement (alpha, V_buf2, prim_fail);
if (V_buf3 != alpha)
{
m= mapDown (m, prim_elem, im_prim_elem, alpha, source, dest);
G_m= mapDown (m, prim_elem, im_prim_elem, alpha, source, dest);
newtonPoly= mapDown (newtonPoly, prim_elem, im_prim_elem, alpha,
source, dest);
ppA= mapDown (ppA, prim_elem, im_prim_elem, alpha, source, dest);
ppB= mapDown (ppB, prim_elem, im_prim_elem, alpha, source, dest);
gcdlcAlcB= mapDown (gcdlcAlcB, prim_elem, im_prim_elem, alpha,
source, dest);
for (CFListIterator i= l; i.hasItem(); i++)
i.getItem()= mapDown (i.getItem(), prim_elem, im_prim_elem, alpha,
source, dest);
}
ASSERT (!prim_fail, "failure in integer factorizer");
if (prim_fail)
; //ERROR
else
im_prim_elem= mapPrimElem (prim_elem, alpha, V_buf);
DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (alpha));
DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (V_buf2));
for (CFListIterator i= l; i.hasItem(); i++)
i.getItem()= mapUp (i.getItem(), alpha, V_buf, prim_elem,
im_prim_elem, source, dest);
m= mapUp (m, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
G_m= mapUp (G_m, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
newtonPoly= mapUp (newtonPoly, alpha, V_buf, prim_elem, im_prim_elem,
source, dest);
ppA= mapUp (ppA, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
ppB= mapUp (ppB, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
gcdlcAlcB= mapUp (gcdlcAlcB, alpha, V_buf, prim_elem, im_prim_elem,
source, dest);
fail= false;
random_element= randomElement (m, V_buf, l, fail );
DEBOUTLN (cerr, "fail= " << fail);
CFList list;
TIMING_START (gcd_recursion);
if (LC (skeleton).inCoeffDomain())
G_random_element=
monicSparseInterpol (ppA (random_element, x), ppB (random_element, x),
skeleton, V_buf, sparseFail, coeffMonoms,
Monoms);
else
G_random_element=
nonMonicSparseInterpol (ppA(random_element,x), ppB(random_element,x),
skeleton, V_buf, sparseFail,
coeffMonoms, Monoms);
TIMING_END_AND_PRINT (gcd_recursion,
"time for recursive call: ");
DEBOUTLN (cerr, "G_random_element= " << G_random_element);
alpha= V_buf;
}
if (sparseFail)
break;
if (!G_random_element.inCoeffDomain())
d0= totaldegree (G_random_element, Variable(2),
Variable (G_random_element.level()));
else
d0= 0;
if (d0 == 0)
{
if (inextension)
prune (cleanUp);
return N(gcdcAcB);
}
if (d0 > d)
{
if (!find (l, random_element))
l.append (random_element);
continue;
}
G_random_element=
(gcdlcAlcB(random_element, x)/uni_lcoeff (G_random_element))
* G_random_element;
if (!G_random_element.inCoeffDomain())
d0= totaldegree (G_random_element, Variable(2),
Variable (G_random_element.level()));
else
d0= 0;
if (d0 < d)
{
m= gcdlcAlcB;
newtonPoly= 1;
G_m= 0;
d= d0;
}
TIMING_START (newton_interpolation);
H= newtonInterp (random_element, G_random_element, newtonPoly, G_m, x);
TIMING_END_AND_PRINT (newton_interpolation,
"time for newton interpolation: ");
//termination test
if (uni_lcoeff (H) == gcdlcAlcB)
{
cH= uni_content (H);
ppH= H/cH;
ppH /= Lc (ppH);
DEBOUTLN (cerr, "ppH= " << ppH);
if (fdivides (ppH, ppA) && fdivides (ppH, ppB))
{
if (inextension)
prune (cleanUp);
return N(gcdcAcB*ppH);
}
}
G_m= H;
newtonPoly= newtonPoly*(x - random_element);
m= m*(x - random_element);
if (!find (l, random_element))
l.append (random_element);
} while (1); //end of second do
}
else
{
if (inextension)
prune (cleanUp);
return N(gcdcAcB*modGCDFp (ppA, ppB));
}
} while (1); //end of first do
}
sparseGCDFq
CanonicalForm sparseGCDFq(const CanonicalForm &F, const CanonicalForm &G, const Variable &alpha, CFList &l, bool &topLevel)
Definition cfModGcd.cc:3138
monicSparseInterpol
CanonicalForm monicSparseInterpol(const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &skeleton, const Variable &alpha, bool &fail, CFArray *&coeffMonoms, CFArray &Monoms)
Definition cfModGcd.cc:2207
nonMonicSparseInterpol
CanonicalForm nonMonicSparseInterpol(const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &skeleton, const Variable &alpha, bool &fail, CFArray *&coeffMonoms, CFArray &Monoms)
Definition cfModGcd.cc:2491
sparseGCDFp
CanonicalForm sparseGCDFp(const CanonicalForm &F, const CanonicalForm &G, bool &topLevel, CFList &l)
Definition cfModGcd.cc:3570

## sparseGCDFq()

CanonicalForm sparseGCDFq | ( | const CanonicalForm & | F,
|  | const CanonicalForm & | G,
|  | const Variable & | alpha,
|  | CFList & | l,
|  | bool & | topLevel )

Definition at line 3138 of file cfModGcd.cc.

{
CanonicalForm A= F;
CanonicalForm B= G;
if (F.isZero() && degree(G) > 0) return G/Lc(G);
else if (G.isZero() && degree (F) > 0) return F/Lc(F);
else if (F.isZero() && G.isZero()) return F.genOne();
if (F.inBaseDomain() || G.inBaseDomain()) return F.genOne();
if (F.isUnivariate() && fdivides(F, G)) return F/Lc(F);
if (G.isUnivariate() && fdivides(G, F)) return G/Lc(G);
if (F == G) return F/Lc(F);
CFMap M,N;
int best_level= myCompress (A, B, M, N, topLevel);
if (best_level == 0) return B.genOne();
A= M(A);
B= M(B);
Variable x= Variable (1);
//univariate case
if (A.isUnivariate() && B.isUnivariate())
return N (gcd (A, B));
CanonicalForm cA, cB; // content of A and B
CanonicalForm ppA, ppB; // primitive part of A and B
CanonicalForm gcdcAcB;
cA = uni_content (A);
cB = uni_content (B);
gcdcAcB= gcd (cA, cB);
ppA= A/cA;
ppB= B/cB;
CanonicalForm lcA, lcB; // leading coefficients of A and B
CanonicalForm gcdlcAlcB;
lcA= uni_lcoeff (ppA);
lcB= uni_lcoeff (ppB);
if (fdivides (lcA, lcB))
{
if (fdivides (A, B))
return F/Lc(F);
}
if (fdivides (lcB, lcA))
{
if (fdivides (B, A))
return G/Lc(G);
}
gcdlcAlcB= gcd (lcA, lcB);
int d= totaldegree (ppA, Variable (2), Variable (ppA.level()));
int d0;
if (d == 0)
return N(gcdcAcB);
d0= totaldegree (ppB, Variable (2), Variable (ppB.level()));
if (d0 < d)
d= d0;
if (d == 0)
return N(gcdcAcB);
CanonicalForm m, random_element, G_m, G_random_element, H, cH, ppH, skeleton;
CanonicalForm newtonPoly= 1;
m= gcdlcAlcB;
G_m= 0;
H= 0;
bool fail= false;
topLevel= false;
bool inextension= false;
Variable V_buf= alpha, V_buf4= alpha;
CanonicalForm prim_elem, im_prim_elem;
CanonicalForm prim_elem_alpha, im_prim_elem_alpha;
CFList source, dest;
do // first do
{
random_element= randomElement (m, V_buf, l, fail);
if (random_element == 0 && !fail)
{
if (!find (l, random_element))
l.append (random_element);
continue;
}
if (fail)
{
source= CFList();
dest= CFList();
Variable V_buf3= V_buf;
V_buf= chooseExtension (V_buf);
bool prim_fail= false;
Variable V_buf2;
prim_elem= primitiveElement (V_buf4, V_buf2, prim_fail);
if (V_buf4 == alpha)
prim_elem_alpha= prim_elem;
if (V_buf3 != V_buf4)
{
m= mapDown (m, prim_elem, im_prim_elem, V_buf4, source, dest);
G_m= mapDown (m, prim_elem, im_prim_elem, V_buf4, source, dest);
newtonPoly= mapDown (newtonPoly, prim_elem, im_prim_elem, V_buf4,
source, dest);
ppA= mapDown (ppA, prim_elem, im_prim_elem, V_buf4, source, dest);
ppB= mapDown (ppB, prim_elem, im_prim_elem, V_buf4, source, dest);
gcdlcAlcB= mapDown (gcdlcAlcB, prim_elem, im_prim_elem, V_buf4, source,
dest);
for (CFListIterator i= l; i.hasItem(); i++)
i.getItem()= mapDown (i.getItem(), prim_elem, im_prim_elem, V_buf4,
source, dest);
}
ASSERT (!prim_fail, "failure in integer factorizer");
if (prim_fail)
; //ERROR
else
im_prim_elem= mapPrimElem (prim_elem, V_buf4, V_buf);
if (V_buf4 == alpha)
im_prim_elem_alpha= im_prim_elem;
else
im_prim_elem_alpha= mapUp (im_prim_elem_alpha, V_buf4, V_buf, prim_elem,
im_prim_elem, source, dest);
DEBOUTLN (cerr, "getMipo (V_buf4)= " << getMipo (V_buf4));
DEBOUTLN (cerr, "getMipo (V_buf2)= " << getMipo (V_buf2));
inextension= true;
for (CFListIterator i= l; i.hasItem(); i++)
i.getItem()= mapUp (i.getItem(), V_buf4, V_buf, prim_elem,
im_prim_elem, source, dest);
m= mapUp (m, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
G_m= mapUp (G_m, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
newtonPoly= mapUp (newtonPoly, V_buf4, V_buf, prim_elem, im_prim_elem,
source, dest);
ppA= mapUp (ppA, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
ppB= mapUp (ppB, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
gcdlcAlcB= mapUp (gcdlcAlcB, V_buf4, V_buf, prim_elem, im_prim_elem,
source, dest);
fail= false;
random_element= randomElement (m, V_buf, l, fail );
DEBOUTLN (cerr, "fail= " << fail);
CFList list;
TIMING_START (gcd_recursion);
G_random_element=
sparseGCDFq (ppA (random_element, x), ppB (random_element, x), V_buf,
list, topLevel);
TIMING_END_AND_PRINT (gcd_recursion,
"time for recursive call: ");
DEBOUTLN (cerr, "G_random_element= " << G_random_element);
V_buf4= V_buf;
}
else
{
CFList list;
TIMING_START (gcd_recursion);
G_random_element=
sparseGCDFq (ppA(random_element, x), ppB(random_element, x), V_buf,
list, topLevel);
TIMING_END_AND_PRINT (gcd_recursion,
"time for recursive call: ");
DEBOUTLN (cerr, "G_random_element= " << G_random_element);
}
if (!G_random_element.inCoeffDomain())
d0= totaldegree (G_random_element, Variable(2),
Variable (G_random_element.level()));
else
d0= 0;
if (d0 == 0)
{
if (inextension)
prune1 (alpha);
return N(gcdcAcB);
}
if (d0 > d)
{
if (!find (l, random_element))
l.append (random_element);
continue;
}
G_random_element=
(gcdlcAlcB(random_element, x)/uni_lcoeff (G_random_element))
* G_random_element;
skeleton= G_random_element;
if (!G_random_element.inCoeffDomain())
d0= totaldegree (G_random_element, Variable(2),
Variable (G_random_element.level()));
else
d0= 0;
if (d0 < d)
{
m= gcdlcAlcB;
newtonPoly= 1;
G_m= 0;
d= d0;
}
H= newtonInterp (random_element, G_random_element, newtonPoly, G_m, x);
if (uni_lcoeff (H) == gcdlcAlcB)
{
cH= uni_content (H);
ppH= H/cH;
if (inextension)
{
CFList u, v;
//maybe it's better to test if ppH is an element of F(\alpha) before
//mapping down
if (fdivides (ppH, ppA) && fdivides (ppH, ppB))
{
DEBOUTLN (cerr, "ppH before mapDown= " << ppH);
ppH= mapDown (ppH, prim_elem_alpha, im_prim_elem_alpha, alpha, u, v);
ppH /= Lc(ppH);
DEBOUTLN (cerr, "ppH after mapDown= " << ppH);
prune1 (alpha);
return N(gcdcAcB*ppH);
}
}
else if (fdivides (ppH, ppA) && fdivides (ppH, ppB))
return N(gcdcAcB*ppH);
}
G_m= H;
newtonPoly= newtonPoly*(x - random_element);
m= m*(x - random_element);
if (!find (l, random_element))
l.append (random_element);
if (getCharacteristic () > 3 && size (skeleton) < 100)
{
CFArray Monoms;
CFArray *coeffMonoms;
do //second do
{
random_element= randomElement (m, V_buf, l, fail);
if (random_element == 0 && !fail)
{
if (!find (l, random_element))
l.append (random_element);
continue;
}
if (fail)
{
source= CFList();
dest= CFList();
Variable V_buf3= V_buf;
V_buf= chooseExtension (V_buf);
bool prim_fail= false;
Variable V_buf2;
prim_elem= primitiveElement (V_buf4, V_buf2, prim_fail);
if (V_buf4 == alpha)
prim_elem_alpha= prim_elem;
if (V_buf3 != V_buf4)
{
m= mapDown (m, prim_elem, im_prim_elem, V_buf4, source, dest);
G_m= mapDown (m, prim_elem, im_prim_elem, V_buf4, source, dest);
newtonPoly= mapDown (newtonPoly, prim_elem, im_prim_elem, V_buf4,
source, dest);
ppA= mapDown (ppA, prim_elem, im_prim_elem, V_buf4, source, dest);
ppB= mapDown (ppB, prim_elem, im_prim_elem, V_buf4, source, dest);
gcdlcAlcB= mapDown (gcdlcAlcB, prim_elem, im_prim_elem, V_buf4,
source, dest);
for (CFListIterator i= l; i.hasItem(); i++)
i.getItem()= mapDown (i.getItem(), prim_elem, im_prim_elem, V_buf4,
source, dest);
}
ASSERT (!prim_fail, "failure in integer factorizer");
if (prim_fail)
; //ERROR
else
im_prim_elem= mapPrimElem (prim_elem, V_buf4, V_buf);
if (V_buf4 == alpha)
im_prim_elem_alpha= im_prim_elem;
else
im_prim_elem_alpha= mapUp (im_prim_elem_alpha, V_buf4, V_buf,
prim_elem, im_prim_elem, source, dest);
DEBOUTLN (cerr, "getMipo (V_buf4)= " << getMipo (V_buf4));
DEBOUTLN (cerr, "getMipo (V_buf2)= " << getMipo (V_buf2));
inextension= true;
for (CFListIterator i= l; i.hasItem(); i++)
i.getItem()= mapUp (i.getItem(), V_buf4, V_buf, prim_elem,
im_prim_elem, source, dest);
m= mapUp (m, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
G_m= mapUp (G_m, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
newtonPoly= mapUp (newtonPoly, V_buf4, V_buf, prim_elem, im_prim_elem,
source, dest);
ppA= mapUp (ppA, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
ppB= mapUp (ppB, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
gcdlcAlcB= mapUp (gcdlcAlcB, V_buf4, V_buf, prim_elem, im_prim_elem,
source, dest);
fail= false;
random_element= randomElement (m, V_buf, l, fail);
DEBOUTLN (cerr, "fail= " << fail);
CFList list;
TIMING_START (gcd_recursion);
V_buf4= V_buf;
//sparseInterpolation
bool sparseFail= false;
if (LC (skeleton).inCoeffDomain())
G_random_element=
monicSparseInterpol (ppA (random_element, x), ppB(random_element,x),
skeleton,V_buf, sparseFail, coeffMonoms,Monoms);
else
G_random_element=
nonMonicSparseInterpol (ppA(random_element,x),ppB(random_element,x),
skeleton, V_buf, sparseFail, coeffMonoms,
Monoms);
if (sparseFail)
break;
TIMING_END_AND_PRINT (gcd_recursion,
"time for recursive call: ");
DEBOUTLN (cerr, "G_random_element= " << G_random_element);
}
else
{
CFList list;
TIMING_START (gcd_recursion);
bool sparseFail= false;
if (LC (skeleton).inCoeffDomain())
G_random_element=
monicSparseInterpol (ppA (random_element, x),ppB(random_element, x),
skeleton,V_buf, sparseFail,coeffMonoms, Monoms);
else
G_random_element=
nonMonicSparseInterpol (ppA(random_element,x),ppB(random_element,x),
skeleton, V_buf, sparseFail, coeffMonoms,
Monoms);
if (sparseFail)
break;
TIMING_END_AND_PRINT (gcd_recursion,
"time for recursive call: ");
DEBOUTLN (cerr, "G_random_element= " << G_random_element);
}
if (!G_random_element.inCoeffDomain())
d0= totaldegree (G_random_element, Variable(2),
Variable (G_random_element.level()));
else
d0= 0;
if (d0 == 0)
{
if (inextension)
prune1 (alpha);
return N(gcdcAcB);
}
if (d0 > d)
{
if (!find (l, random_element))
l.append (random_element);
continue;
}
G_random_element=
(gcdlcAlcB(random_element, x)/uni_lcoeff (G_random_element))
* G_random_element;
if (!G_random_element.inCoeffDomain())
d0= totaldegree (G_random_element, Variable(2),
Variable (G_random_element.level()));
else
d0= 0;
if (d0 < d)
{
m= gcdlcAlcB;
newtonPoly= 1;
G_m= 0;
d= d0;
}
TIMING_START (newton_interpolation);
H= newtonInterp (random_element, G_random_element, newtonPoly, G_m, x);
TIMING_END_AND_PRINT (newton_interpolation,
"time for newton interpolation: ");
//termination test
if (uni_lcoeff (H) == gcdlcAlcB)
{
cH= uni_content (H);
ppH= H/cH;
if (inextension)
{
CFList u, v;
//maybe it's better to test if ppH is an element of F(\alpha) before
//mapping down
if (fdivides (ppH, ppA) && fdivides (ppH, ppB))
{
DEBOUTLN (cerr, "ppH before mapDown= " << ppH);
ppH= mapDown (ppH, prim_elem_alpha, im_prim_elem_alpha, alpha, u, v);
ppH /= Lc(ppH);
DEBOUTLN (cerr, "ppH after mapDown= " << ppH);
prune1 (alpha);
return N(gcdcAcB*ppH);
}
}
else if (fdivides (ppH, ppA) && fdivides (ppH, ppB))
{
return N(gcdcAcB*ppH);
}
}
G_m= H;
newtonPoly= newtonPoly*(x - random_element);
m= m*(x - random_element);
if (!find (l, random_element))
l.append (random_element);
} while (1);
}
} while (1);
}

## TIMING_DEFINE_PRINT() [1/2]

TIMING_DEFINE_PRINT | ( | gcd_recursion |  | ) | const &

## TIMING_DEFINE_PRINT() [2/2]

TIMING_DEFINE_PRINT | ( | modZ_termination |  | ) | const &

modular gcd algorithm, see Keith, Czapor, Geddes "Algorithms for Computer Algebra", Algorithm 7.1

## uni_content() [1/2]

CanonicalForm uni_content | ( | const CanonicalForm & | F | ) |

inline static

compute the content of F, where F is considered as an element of

Definition at line 282 of file cfModGcd.cc.

{
if (F.inBaseDomain())
return F.genOne();
if (F.level() == 1 && F.isUnivariate())
return F;
if (F.level() != 1 && F.isUnivariate())
return F.genOne();
if (degree (F,1) == 0) return F.genOne();
int l= F.level();
if (l == 2)
return content(F);
else
{
CanonicalForm pol, c = 0;
CFIterator i = F;
for (; i.hasTerms(); i++)
{
pol= i.coeff();
pol= uni_content (pol);
c= gcd (c, pol);
if (c.isOne())
return c;
}
return c;
}
}
content
CanonicalForm FACTORY_PUBLIC content(const CanonicalForm &)
CanonicalForm content ( const CanonicalForm & f ).
Definition cf_gcd.cc:603

## uni_content() [2/2]

CanonicalForm uni_content | ( | const CanonicalForm & | F,
|  | const Variable & | x )

Definition at line 260 of file cfModGcd.cc.

{
if (F.inCoeffDomain())
return F.genOne();
if (F.level() == x.level() && F.isUnivariate())
return F;
if (F.level() != x.level() && F.isUnivariate())
return F.genOne();
if (x.level() != 1)
{
CanonicalForm f= swapvar (F, x, Variable (1));
CanonicalForm result= uni_content (f);
return swapvar (result, x, Variable (1));
}
else
return uni_content (F);
}
swapvar
CanonicalForm FACTORY_PUBLIC swapvar(const CanonicalForm &, const Variable &, const Variable &)
swapvar() - swap variables x1 and x2 in f.
Definition cf_ops.cc:168

## uni_lcoeff()

CanonicalForm uni_lcoeff | ( | const CanonicalForm & | F | ) |

inline static

compute the leading coefficient of F, where F is considered as an element of , order on is dp.

Definition at line 340 of file cfModGcd.cc.

{
if (F.level() > 1)
{
Variable x= Variable (2);
int deg= totaldegree (F, x, F.mvar());
for (CFIterator i= F; i.hasTerms(); i++)
{
if (i.exp() + totaldegree (i.coeff(), x, i.coeff().mvar()) == deg)
return uni_lcoeff (i.coeff());
}
}
return F;
}

## while()

while | ( | true |  | ) |

Definition at line 4140 of file cfModGcd.cc.

{
p = cf_getBigPrime( i );
i--;
while ( i >= 0 && mod( cl*(lc(f)/cl)*(lc(g)/cl), p ) == 0 )
{
p = cf_getBigPrime( i );
i--;
}
//printf("try p=%d\n",p);
setCharacteristic( p );
fp= mapinto (f);
gp= mapinto (g);
TIMING_START (modZ_recursion)
#if defined(HAVE_NTL) || defined(HAVE_FLINT)
if (size (fp)/maxNumVars > 500 && size (gp)/maxNumVars > 500)
Dp = modGCDFp (fp, gp, cofp, cogp);
else
{
Dp= gcd_poly (fp, gp);
cofp= fp/Dp;
cogp= gp/Dp;
}
#else
Dp= gcd_poly (fp, gp);
cofp= fp/Dp;
cogp= gp/Dp;
#endif
TIMING_END_AND_PRINT (modZ_recursion,
"time for gcd mod p in modular gcd: ");
Dp /=Dp.lc();
Dp *= mapinto (cl);
cofp /= lc (cofp);
cofp *= mapinto (lcf);
cogp /= lc (cogp);
cogp *= mapinto (lcg);
setCharacteristic( 0 );
dp_deg=totaldegree(Dp);
if ( dp_deg == 0 )
{
//printf(" -> 1\n");
return CanonicalForm(gcdcfcg);
}
if ( q.isZero() )
{
D = mapinto( Dp );
cof= mapinto (cofp);
cog= mapinto (cogp);
d_deg=dp_deg;
q = p;
Dn= balance_p (D, p);
cofn= balance_p (cof, p);
cogn= balance_p (cog, p);
}
else
{
if ( dp_deg == d_deg )
{
chineseRemainder( D, q, mapinto( Dp ), p, newD, newq );
chineseRemainder( cof, q, mapinto (cofp), p, newCof, newq);
chineseRemainder( cog, q, mapinto (cogp), p, newCog, newq);
cof= newCof;
cog= newCog;
newqh= newq/2;
Dn= balance_p (newD, newq, newqh);
cofn= balance_p (newCof, newq, newqh);
cogn= balance_p (newCog, newq, newqh);
if (test != Dn) //balance_p (newD, newq))
test= balance_p (newD, newq);
else
equal= true;
q = newq;
D = newD;
}
else if ( dp_deg < d_deg )
{
//printf(" were all bad, try more\n");
// all previous p's are bad primes
q = p;
D = mapinto( Dp );
cof= mapinto (cof);
cog= mapinto (cog);
d_deg=dp_deg;
test= 0;
equal= false;
Dn= balance_p (D, p);
cofn= balance_p (cof, p);
cogn= balance_p (cog, p);
}
else
{
//printf(" was bad, try more\n");
}
//else dp_deg > d_deg: bad prime
}
if ( i >= 0 )
{
cDn= icontent (Dn);
Dn /= cDn;
cofn /= cl/cDn;
//cofn /= icontent (cofn);
cogn /= cl/cDn;
//cogn /= icontent (cogn);
TIMING_START (modZ_termination);
if ((terminationTest (f,g, cofn, cogn, Dn)) ||
((equal || q > b) && fdivides (Dn, f) && fdivides (Dn, g)))
{
TIMING_END_AND_PRINT (modZ_termination,
"time for successful termination in modular gcd: ");
//printf(" -> success\n");
return Dn*gcdcfcg;
}
TIMING_END_AND_PRINT (modZ_termination,
"time for unsuccessful termination in modular gcd: ");
equal= false;
//else: try more primes
}
else
{ // try other method
//printf("try other gcd\n");
Off(SW_USE_CHINREM_GCD);
D=gcd_poly( f, g );
On(SW_USE_CHINREM_GCD);
return D*gcdcfcg;
}
}
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
icontent
CanonicalForm FACTORY_PUBLIC icontent(const CanonicalForm &f)
CanonicalForm icontent ( const CanonicalForm & f ).
Definition cf_gcd.cc:74
mod
CF_NO_INLINE FACTORY_PUBLIC CanonicalForm mod(const CanonicalForm &, const CanonicalForm &)
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
cofp
CanonicalForm cofp
Definition cfModGcd.cc:4137
lcg
CanonicalForm lcg
Definition cfModGcd.cc:4105
dp_deg
int dp_deg
Definition cfModGcd.cc:4086
newCog
CanonicalForm newCog
Definition cfModGcd.cc:4137
cogn
CanonicalForm cogn
Definition cfModGcd.cc:4137
cl
cl
Definition cfModGcd.cc:4108
lcf
CanonicalForm lcf
Definition cfModGcd.cc:4105
maxNumVars
int maxNumVars
Definition cfModGcd.cc:4138
fp
CanonicalForm fp
Definition cfModGcd.cc:4110
cogp
CanonicalForm cogp
Definition cfModGcd.cc:4137
cofn
CanonicalForm cofn
Definition cfModGcd.cc:4137
cof
CanonicalForm cof
Definition cfModGcd.cc:4137
cog
CanonicalForm cog
Definition cfModGcd.cc:4137
gcdcfcg
CanonicalForm gcdcfcg
Definition cfModGcd.cc:4109
Dn
CanonicalForm Dn
Definition cfModGcd.cc:4104
newCof
CanonicalForm newCof
Definition cfModGcd.cc:4137
gp
CanonicalForm gp
Definition cfModGcd.cc:4110
equal
bool equal
Definition cfModGcd.cc:4134
test
CanonicalForm test
Definition cfModGcd.cc:4104
b
CanonicalForm b
Definition cfModGcd.cc:4111
cDn
CanonicalForm cDn
Definition cfModGcd.cc:4137
d_deg
int d_deg
Definition cfModGcd.cc:4086
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_CHINREM_GCD
static const int SW_USE_CHINREM_GCD
set to 1 to use modular gcd over Z
Definition cf_defs.h:41
balance_p
static CanonicalForm balance_p(const CanonicalForm &f, const CanonicalForm &q, const CanonicalForm &qh)
Definition cf_gcd.cc:149
cf_getBigPrime
int cf_getBigPrime(int i)
Definition cf_primes.cc:39
D
#define D(A)
Definition gentable.cc:128

## Variable Documentation

## b

else L b = 1

Definition at line 4111 of file cfModGcd.cc.

## cand

const CanonicalForm const CanonicalForm const CanonicalForm const CanonicalForm& cand

Initial value:
{
CanonicalForm LCCand= abs (LC (cand))

Definition at line 69 of file cfModGcd.cc.

## cd

cd | ( | bCommonDen(FF) |  | ) | = bCommonDen( GG )

Definition at line 4097 of file cfModGcd.cc.

## cDn

CanonicalForm cDn

Definition at line 4137 of file cfModGcd.cc.

## cf

cf = icontent (f)

Definition at line 4091 of file cfModGcd.cc.

## cg

cg = icontent (g)

Definition at line 4091 of file cfModGcd.cc.

## cl

cl = gcd (f.lc(),g.lc())

Definition at line 4108 of file cfModGcd.cc.

## coF

const CanonicalForm const CanonicalForm& coF

Definition at line 68 of file cfModGcd.cc.

## cof

CanonicalForm cof

Definition at line 4137 of file cfModGcd.cc.

## cofn

CanonicalForm cofn

Definition at line 4137 of file cfModGcd.cc.

## cofp

CanonicalForm cofp

Definition at line 4137 of file cfModGcd.cc.

## coG

const CanonicalForm const CanonicalForm const CanonicalForm& coG

Definition at line 68 of file cfModGcd.cc.

## cog

CanonicalForm cog

Definition at line 4137 of file cfModGcd.cc.

## cogn

CanonicalForm cogn

Definition at line 4137 of file cfModGcd.cc.

## cogp

CanonicalForm cogp

Definition at line 4137 of file cfModGcd.cc.

## d_deg

int d_deg =-1

Definition at line 4086 of file cfModGcd.cc.

## Dn

CanonicalForm Dn

Definition at line 4104 of file cfModGcd.cc.

## dp_deg

int dp_deg

Definition at line 4086 of file cfModGcd.cc.

## equal

bool equal = false

Definition at line 4134 of file cfModGcd.cc.

## f

f =cd*FF

Definition at line 4089 of file cfModGcd.cc.

## false

return false

Definition at line 85 of file cfModGcd.cc.

## fp

CanonicalForm fp

Definition at line 4110 of file cfModGcd.cc.

## G

const CanonicalForm& G

Definition at line 67 of file cfModGcd.cc.

## g

g =cd*GG

Definition at line 4098 of file cfModGcd.cc.

## gcdcfcg

CanonicalForm gcdcfcg = gcd (cf, cg)

Definition at line 4109 of file cfModGcd.cc.

## GG

const CanonicalForm& GG

Initial value:
{
CanonicalForm f, g, cl, q(0), Dp, newD, D, newq, newqh

Definition at line 4083 of file cfModGcd.cc.

## gp

CanonicalForm gp

Definition at line 4110 of file cfModGcd.cc.

## i

i = cf_getNumBigPrimes() - 1

Definition at line 4086 of file cfModGcd.cc.

## lcf

lcf = f.lc()

Definition at line 4105 of file cfModGcd.cc.

## lcg

lcg = g.lc()

Definition at line 4105 of file cfModGcd.cc.

## maxNumVars

int maxNumVars = tmax (getNumVars (f), getNumVars (g))

Definition at line 4138 of file cfModGcd.cc.

## minCommonDeg

int minCommonDeg = 0

Definition at line 4112 of file cfModGcd.cc.

## newCof

CanonicalForm newCof

Definition at line 4137 of file cfModGcd.cc.

## newCog

CanonicalForm newCog

Definition at line 4137 of file cfModGcd.cc.

## p

int p

Definition at line 4086 of file cfModGcd.cc.

## test

CanonicalForm test = 0

Definition at line 4104 of file cfModGcd.cc.

## x

Variable x = Variable (1)

Definition at line 4090 of file cfModGcd.cc.