My Project
Functions
cfModGcd.h File Reference

modular and sparse modular GCD algorithms over finite fields and Z. More...

#include "cf_assert.h"
#include "cf_factory.h"

Go to the source code of this file.

Functions

CanonicalForm modGCDFq (const CanonicalForm &F, const CanonicalForm &G, Variable &alpha, CFList &l, bool &top_level)
 
static CanonicalForm modGCDFq (const CanonicalForm &A, const CanonicalForm &B, Variable &alpha)
 GCD of A and B over $ F_{p}(\alpha ) $. More...
 
CanonicalForm modGCDFp (const CanonicalForm &F, const CanonicalForm &G, bool &top_level, CFList &l)
 
CanonicalForm modGCDFp (const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &coF, CanonicalForm &coG, bool &topLevel, CFList &l)
 
static CanonicalForm modGCDFp (const CanonicalForm &A, const CanonicalForm &B)
 GCD of A and B over $ F_{p} $. More...
 
static CanonicalForm modGCDFp (const CanonicalForm &A, const CanonicalForm &B, CanonicalForm &coA, CanonicalForm &coB)
 
CanonicalForm modGCDGF (const CanonicalForm &F, const CanonicalForm &G, CFList &l, bool &top_level)
 
static CanonicalForm modGCDGF (const CanonicalForm &A, const CanonicalForm &B)
 GCD of A and B over GF. More...
 
CanonicalForm sparseGCDFp (const CanonicalForm &F, const CanonicalForm &G, bool &topLevel, CFList &l)
 
static CanonicalForm sparseGCDFp (const CanonicalForm &A, const CanonicalForm &B)
 Zippel's sparse GCD over Fp. More...
 
CanonicalForm sparseGCDFq (const CanonicalForm &F, const CanonicalForm &G, const Variable &alpha, CFList &l, bool &topLevel)
 
static CanonicalForm sparseGCDFq (const CanonicalForm &A, const CanonicalForm &B, const Variable &alpha)
 Zippel's sparse GCD over Fq. More...
 
CFArray getMonoms (const CanonicalForm &F)
 extract monomials of F, parts in algebraic variable are considered coefficients More...
 
bool terminationTest (const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &coF, const CanonicalForm &coG, const CanonicalForm &cand)
 
CanonicalForm modGCDZ (const CanonicalForm &FF, const CanonicalForm &GG)
 modular GCD over Z More...
 

Detailed Description

modular and sparse modular GCD algorithms over finite fields and Z.

Author
Martin Lee
Date
22.10.2009
Copyright:
(c) by The SINGULAR Team, see LICENSE file

Definition in file cfModGcd.h.

Function Documentation

◆ getMonoms()

CFArray getMonoms ( const CanonicalForm F)

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

Parameters
[in]Fsome poly

Definition at line 1959 of file cfModGcd.cc.

1960 {
1961  if (F.inCoeffDomain())
1962  {
1963  CFArray result= CFArray (1);
1964  result [0]= 1;
1965  return result;
1966  }
1967  if (F.isUnivariate())
1968  {
1969  CFArray result= CFArray (size(F));
1970  int j= 0;
1971  for (CFIterator i= F; i.hasTerms(); i++, j++)
1972  result[j]= power (F.mvar(), i.exp());
1973  return result;
1974  }
1975  int numMon= size (F);
1976  CFArray result= CFArray (numMon);
1977  int j= 0;
1978  CFArray recResult;
1979  Variable x= F.mvar();
1980  CanonicalForm powX;
1981  for (CFIterator i= F; i.hasTerms(); i++)
1982  {
1983  powX= power (x, i.exp());
1984  recResult= getMonoms (i.coeff());
1985  for (int k= 0; k < recResult.size(); k++)
1986  result[j+k]= powX*recResult[k];
1987  j += recResult.size();
1988  }
1989  return result;
1990 }
CanonicalForm power(const CanonicalForm &f, int n)
exponentiation
int size(const CanonicalForm &f, const Variable &v)
int size ( const CanonicalForm & f, const Variable & v )
Definition: cf_ops.cc:600
Array< CanonicalForm > CFArray
int k
Definition: cfEzgcd.cc:99
CFArray getMonoms(const CanonicalForm &F)
extract monomials of F, parts in algebraic variable are considered coefficients
Definition: cfModGcd.cc:1959
Variable x
Definition: cfModGcd.cc:4084
int i
Definition: cfModGcd.cc:4080
int size() const
Definition: ftmpl_array.cc:92
class to iterate through CanonicalForm's
Definition: cf_iter.h:44
factory's main class
Definition: canonicalform.h:86
Variable mvar() const
mvar() returns the main variable of CO or Variable() if CO is in a base domain.
bool inCoeffDomain() const
bool isUnivariate() const
factory's class for variables
Definition: factory.h:134
return result
Definition: facAbsBiFact.cc:75
int j
Definition: facHensel.cc:110

◆ modGCDFp() [1/4]

static CanonicalForm modGCDFp ( const CanonicalForm A,
const CanonicalForm B 
)
inlinestatic

GCD of A and B over $ F_{p} $.

Parameters
[in]Apoly over F_p
[in]Bpoly over F_p

Definition at line 50 of file cfModGcd.h.

53 {
54  CFList list;
55  bool top_level= true;
56  return modGCDFp (A, B, top_level, list);
57 }
CanonicalForm modGCDFp(const CanonicalForm &F, const CanonicalForm &G, bool &top_level, CFList &l)
Definition: cfModGcd.cc:1214
b *CanonicalForm B
Definition: facBivar.cc:52
#define A
Definition: sirandom.c:24

◆ modGCDFp() [2/4]

static CanonicalForm modGCDFp ( const CanonicalForm A,
const CanonicalForm B,
CanonicalForm coA,
CanonicalForm coB 
)
inlinestatic

Definition at line 60 of file cfModGcd.h.

62 {
63  CFList list;
64  bool top_level= true;
65  return modGCDFp (A, B, coA, coB, top_level, list);
66 }

◆ modGCDFp() [3/4]

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

Definition at line 1214 of file cfModGcd.cc.

1216 {
1217  CanonicalForm dummy1, dummy2;
1218  CanonicalForm result= modGCDFp (F, G, dummy1, dummy2, topLevel, l);
1219  return result;
1220 }
int l
Definition: cfEzgcd.cc:100
const CanonicalForm CFMap CFMap bool topLevel
Definition: cfGcdAlgExt.cc:57
CanonicalForm modGCDFp(const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &coF, CanonicalForm &coG, bool &topLevel, CFList &l)
Definition: cfModGcd.cc:1225
const CanonicalForm & G
Definition: cfModGcd.cc:66

◆ modGCDFp() [4/4]

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

Definition at line 1225 of file cfModGcd.cc.

1228 {
1229  CanonicalForm A= F;
1230  CanonicalForm B= G;
1231  if (F.isZero() && degree(G) > 0)
1232  {
1233  coF= 0;
1234  coG= Lc (G);
1235  return G/Lc(G);
1236  }
1237  else if (G.isZero() && degree (F) > 0)
1238  {
1239  coF= Lc (F);
1240  coG= 0;
1241  return F/Lc(F);
1242  }
1243  else if (F.isZero() && G.isZero())
1244  {
1245  coF= coG= 0;
1246  return F.genOne();
1247  }
1248  if (F.inBaseDomain() || G.inBaseDomain())
1249  {
1250  coF= F;
1251  coG= G;
1252  return F.genOne();
1253  }
1254  if (F.isUnivariate() && fdivides(F, G, coG))
1255  {
1256  coF= Lc (F);
1257  return F/Lc(F);
1258  }
1259  if (G.isUnivariate() && fdivides(G, F, coF))
1260  {
1261  coG= Lc (G);
1262  return G/Lc(G);
1263  }
1264  if (F == G)
1265  {
1266  coF= coG= Lc (F);
1267  return F/Lc(F);
1268  }
1269  CFMap M,N;
1270  int best_level= myCompress (A, B, M, N, topLevel);
1271 
1272  if (best_level == 0)
1273  {
1274  coF= F;
1275  coG= G;
1276  return B.genOne();
1277  }
1278 
1279  A= M(A);
1280  B= M(B);
1281 
1282  Variable x= Variable (1);
1283 
1284  //univariate case
1285  if (A.isUnivariate() && B.isUnivariate())
1286  {
1287  CanonicalForm result= gcd (A, B);
1288  coF= N (A/result);
1289  coG= N (B/result);
1290  return N (result);
1291  }
1292 
1293  CanonicalForm cA, cB; // content of A and B
1294  CanonicalForm ppA, ppB; // primitive part of A and B
1295  CanonicalForm gcdcAcB;
1296 
1297  cA = uni_content (A);
1298  cB = uni_content (B);
1299  gcdcAcB= gcd (cA, cB);
1300  ppA= A/cA;
1301  ppB= B/cB;
1302 
1303  CanonicalForm lcA, lcB; // leading coefficients of A and B
1304  CanonicalForm gcdlcAlcB;
1305  lcA= uni_lcoeff (ppA);
1306  lcB= uni_lcoeff (ppB);
1307 
1308  gcdlcAlcB= gcd (lcA, lcB);
1309 
1310  int d= totaldegree (ppA, Variable (2), Variable (ppA.level()));
1311  int d0;
1312 
1313  if (d == 0)
1314  {
1315  coF= N (ppA*(cA/gcdcAcB));
1316  coG= N (ppB*(cB/gcdcAcB));
1317  return N(gcdcAcB);
1318  }
1319 
1320  d0= totaldegree (ppB, Variable (2), Variable (ppB.level()));
1321 
1322  if (d0 < d)
1323  d= d0;
1324 
1325  if (d == 0)
1326  {
1327  coF= N (ppA*(cA/gcdcAcB));
1328  coG= N (ppB*(cB/gcdcAcB));
1329  return N(gcdcAcB);
1330  }
1331 
1332  CanonicalForm m, random_element, G_m, G_random_element, H, cH, ppH;
1333  CanonicalForm newtonPoly, coF_random_element, coG_random_element,
1334  coF_m, coG_m, ppCoF, ppCoG;
1335 
1336  newtonPoly= 1;
1337  m= gcdlcAlcB;
1338  H= 0;
1339  coF= 0;
1340  coG= 0;
1341  G_m= 0;
1342  Variable alpha, V_buf, cleanUp;
1343  bool fail= false;
1344  bool inextension= false;
1345  topLevel= false;
1346  CFList source, dest;
1347  int bound1= degree (ppA, 1);
1348  int bound2= degree (ppB, 1);
1349  do
1350  {
1351  if (inextension)
1352  random_element= randomElement (m*lcA*lcB, V_buf, l, fail);
1353  else
1354  random_element= FpRandomElement (m*lcA*lcB, l, fail);
1355 
1356  if (!fail && !inextension)
1357  {
1358  CFList list;
1359  TIMING_START (gcd_recursion);
1360  G_random_element=
1361  modGCDFp (ppA (random_element,x), ppB (random_element,x),
1362  coF_random_element, coG_random_element, topLevel,
1363  list);
1364  TIMING_END_AND_PRINT (gcd_recursion,
1365  "time for recursive call: ");
1366  DEBOUTLN (cerr, "G_random_element= " << G_random_element);
1367  }
1368  else if (!fail && inextension)
1369  {
1370  CFList list;
1371  TIMING_START (gcd_recursion);
1372  G_random_element=
1373  modGCDFq (ppA (random_element, x), ppB (random_element, x),
1374  coF_random_element, coG_random_element, V_buf,
1375  list, topLevel);
1376  TIMING_END_AND_PRINT (gcd_recursion,
1377  "time for recursive call: ");
1378  DEBOUTLN (cerr, "G_random_element= " << G_random_element);
1379  }
1380  else if (fail && !inextension)
1381  {
1382  source= CFList();
1383  dest= CFList();
1384  CFList list;
1386  int deg= 2;
1387  bool initialized= false;
1388  do
1389  {
1390  mipo= randomIrredpoly (deg, x);
1391  if (initialized)
1392  setMipo (alpha, mipo);
1393  else
1394  alpha= rootOf (mipo);
1395  inextension= true;
1396  initialized= true;
1397  fail= false;
1398  random_element= randomElement (m*lcA*lcB, alpha, l, fail);
1399  deg++;
1400  } while (fail);
1401  list= CFList();
1402  V_buf= alpha;
1403  cleanUp= alpha;
1404  TIMING_START (gcd_recursion);
1405  G_random_element=
1406  modGCDFq (ppA (random_element, x), ppB (random_element, x),
1407  coF_random_element, coG_random_element, alpha,
1408  list, topLevel);
1409  TIMING_END_AND_PRINT (gcd_recursion,
1410  "time for recursive call: ");
1411  DEBOUTLN (cerr, "G_random_element= " << G_random_element);
1412  }
1413  else if (fail && inextension)
1414  {
1415  source= CFList();
1416  dest= CFList();
1417 
1418  Variable V_buf3= V_buf;
1419  V_buf= chooseExtension (V_buf);
1420  bool prim_fail= false;
1421  Variable V_buf2;
1422  CanonicalForm prim_elem, im_prim_elem;
1423  prim_elem= primitiveElement (alpha, V_buf2, prim_fail);
1424 
1425  if (V_buf3 != alpha)
1426  {
1427  m= mapDown (m, prim_elem, im_prim_elem, alpha, source, dest);
1428  G_m= mapDown (G_m, prim_elem, im_prim_elem, alpha, source, dest);
1429  coF_m= mapDown (coF_m, prim_elem, im_prim_elem, alpha, source, dest);
1430  coG_m= mapDown (coG_m, prim_elem, im_prim_elem, alpha, source, dest);
1431  newtonPoly= mapDown (newtonPoly, prim_elem, im_prim_elem, alpha,
1432  source, dest);
1433  ppA= mapDown (ppA, prim_elem, im_prim_elem, alpha, source, dest);
1434  ppB= mapDown (ppB, prim_elem, im_prim_elem, alpha, source, dest);
1435  gcdlcAlcB= mapDown (gcdlcAlcB, prim_elem, im_prim_elem, alpha, source,
1436  dest);
1437  lcA= mapDown (lcA, prim_elem, im_prim_elem, alpha, source, dest);
1438  lcB= mapDown (lcB, prim_elem, im_prim_elem, alpha, source, dest);
1439  for (CFListIterator i= l; i.hasItem(); i++)
1440  i.getItem()= mapDown (i.getItem(), prim_elem, im_prim_elem, alpha,
1441  source, dest);
1442  }
1443 
1444  ASSERT (!prim_fail, "failure in integer factorizer");
1445  if (prim_fail)
1446  ; //ERROR
1447  else
1448  im_prim_elem= mapPrimElem (prim_elem, alpha, V_buf);
1449 
1450  DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (alpha));
1451  DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (V_buf2));
1452 
1453  for (CFListIterator i= l; i.hasItem(); i++)
1454  i.getItem()= mapUp (i.getItem(), alpha, V_buf, prim_elem,
1455  im_prim_elem, source, dest);
1456  m= mapUp (m, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
1457  G_m= mapUp (G_m, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
1458  coF_m= mapUp (coF_m, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
1459  coG_m= mapUp (coG_m, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
1460  newtonPoly= mapUp (newtonPoly, alpha, V_buf, prim_elem, im_prim_elem,
1461  source, dest);
1462  ppA= mapUp (ppA, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
1463  ppB= mapUp (ppB, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
1464  gcdlcAlcB= mapUp (gcdlcAlcB, alpha, V_buf, prim_elem, im_prim_elem,
1465  source, dest);
1466  lcA= mapUp (lcA, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
1467  lcB= mapUp (lcB, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
1468  fail= false;
1469  random_element= randomElement (m*lcA*lcB, V_buf, l, fail );
1470  DEBOUTLN (cerr, "fail= " << fail);
1471  CFList list;
1472  TIMING_START (gcd_recursion);
1473  G_random_element=
1474  modGCDFq (ppA (random_element, x), ppB (random_element, x),
1475  coF_random_element, coG_random_element, V_buf,
1476  list, topLevel);
1477  TIMING_END_AND_PRINT (gcd_recursion,
1478  "time for recursive call: ");
1479  DEBOUTLN (cerr, "G_random_element= " << G_random_element);
1480  alpha= V_buf;
1481  }
1482 
1483  if (!G_random_element.inCoeffDomain())
1484  d0= totaldegree (G_random_element, Variable(2),
1485  Variable (G_random_element.level()));
1486  else
1487  d0= 0;
1488 
1489  if (d0 == 0)
1490  {
1491  if (inextension)
1492  prune (cleanUp);
1493  coF= N (ppA*(cA/gcdcAcB));
1494  coG= N (ppB*(cB/gcdcAcB));
1495  return N(gcdcAcB);
1496  }
1497 
1498  if (d0 > d)
1499  {
1500  if (!find (l, random_element))
1501  l.append (random_element);
1502  continue;
1503  }
1504 
1505  G_random_element= (gcdlcAlcB(random_element,x)/uni_lcoeff(G_random_element))
1506  *G_random_element;
1507 
1508  coF_random_element= (lcA(random_element,x)/uni_lcoeff(coF_random_element))
1509  *coF_random_element;
1510  coG_random_element= (lcB(random_element,x)/uni_lcoeff(coG_random_element))
1511  *coG_random_element;
1512 
1513  if (!G_random_element.inCoeffDomain())
1514  d0= totaldegree (G_random_element, Variable(2),
1515  Variable (G_random_element.level()));
1516  else
1517  d0= 0;
1518 
1519  if (d0 < d)
1520  {
1521  m= gcdlcAlcB;
1522  newtonPoly= 1;
1523  G_m= 0;
1524  d= d0;
1525  coF_m= 0;
1526  coG_m= 0;
1527  }
1528 
1529  TIMING_START (newton_interpolation);
1530  H= newtonInterp (random_element, G_random_element, newtonPoly, G_m, x);
1531  coF= newtonInterp (random_element, coF_random_element, newtonPoly, coF_m,x);
1532  coG= newtonInterp (random_element, coG_random_element, newtonPoly, coG_m,x);
1533  TIMING_END_AND_PRINT (newton_interpolation,
1534  "time for newton_interpolation: ");
1535 
1536  //termination test
1537  if ((uni_lcoeff (H) == gcdlcAlcB) || (G_m == H))
1538  {
1539  TIMING_START (termination_test);
1540  if (gcdlcAlcB.isOne())
1541  cH= 1;
1542  else
1543  cH= uni_content (H);
1544  ppH= H/cH;
1545  ppH /= Lc (ppH);
1546  CanonicalForm lcppH= gcdlcAlcB/cH;
1547  CanonicalForm ccoF= lcppH/Lc (lcppH);
1548  CanonicalForm ccoG= lcppH/Lc (lcppH);
1549  ppCoF= coF/ccoF;
1550  ppCoG= coG/ccoG;
1551  DEBOUTLN (cerr, "ppH= " << ppH);
1552  if (((degree (ppCoF,1)+degree (ppH,1) == bound1) &&
1553  (degree (ppCoG,1)+degree (ppH,1) == bound2) &&
1554  terminationTest (ppA, ppB, ppCoF, ppCoG, ppH)) ||
1555  (fdivides (ppH, ppA, ppCoF) && fdivides (ppH, ppB, ppCoG)))
1556  {
1557  if (inextension)
1558  prune (cleanUp);
1559  coF= N ((cA/gcdcAcB)*ppCoF);
1560  coG= N ((cB/gcdcAcB)*ppCoG);
1561  TIMING_END_AND_PRINT (termination_test,
1562  "time for successful termination Fp: ");
1563  return N(gcdcAcB*ppH);
1564  }
1565  TIMING_END_AND_PRINT (termination_test,
1566  "time for unsuccessful termination Fp: ");
1567  }
1568 
1569  G_m= H;
1570  coF_m= coF;
1571  coG_m= coG;
1572  newtonPoly= newtonPoly*(x - random_element);
1573  m= m*(x - random_element);
1574  if (!find (l, random_element))
1575  l.append (random_element);
1576  } while (1);
1577 }
int degree(const CanonicalForm &f)
int totaldegree(const CanonicalForm &f)
int totaldegree ( const CanonicalForm & f )
Definition: cf_ops.cc:523
CanonicalForm Lc(const CanonicalForm &f)
List< CanonicalForm > CFList
const CanonicalForm CFMap CFMap & N
Definition: cfEzgcd.cc:56
int m
Definition: cfEzgcd.cc:128
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:480
const CanonicalForm const CanonicalForm const CanonicalForm & coG
Definition: cfModGcd.cc:67
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:93
static Variable chooseExtension(const Variable &alpha)
Definition: cfModGcd.cc:422
static CanonicalForm uni_content(const CanonicalForm &F)
compute the content of F, where F is considered as an element of
Definition: cfModGcd.cc:283
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:341
const CanonicalForm const CanonicalForm & coF
Definition: cfModGcd.cc:67
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:381
static CanonicalForm FpRandomElement(const CanonicalForm &F, CFList &list, bool &fail)
Definition: cfModGcd.cc:1173
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:366
bool terminationTest(const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &coF, const CanonicalForm &coG, const CanonicalForm &cand)
bool fdivides(const CanonicalForm &f, const CanonicalForm &g)
bool fdivides ( const CanonicalForm & f, const CanonicalForm & g )
#define ASSERT(expression, message)
Definition: cf_assert.h:99
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
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:450
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:342
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:123
static CanonicalForm mapUp(const Variable &alpha, const Variable &beta)
and is a primitive element, returns the image of
Definition: cf_map_ext.cc:70
class CFMap
Definition: cf_map.h:85
CanonicalForm genOne() const
CF_NO_INLINE bool isZero() const
CF_NO_INLINE bool isOne() const
int level() const
level() returns the level of CO.
bool inBaseDomain() const
#define DEBOUTLN(stream, objects)
Definition: debug.h:49
Variable alpha
Definition: facAbsBiFact.cc:51
CanonicalForm H
Definition: facAbsFact.cc:60
CanonicalForm mipo
Definition: facAlgExt.cc:57
TIMING_END_AND_PRINT(fac_alg_resultant, "time to compute resultant0: ")
TIMING_START(fac_alg_resultant)
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
CanonicalForm getMipo(const Variable &alpha, const Variable &x)
Definition: variable.cc:207
void FACTORY_PUBLIC prune(Variable &alpha)
Definition: variable.cc:261
void setMipo(const Variable &alpha, const CanonicalForm &mipo)
Definition: variable.cc:219
template bool find(const List< CanonicalForm > &, const CanonicalForm &)
#define M
Definition: sirandom.c:25
int gcd(int a, int b)
Definition: walkSupport.cc:836

◆ modGCDFq() [1/2]

static CanonicalForm modGCDFq ( const CanonicalForm A,
const CanonicalForm B,
Variable alpha 
)
inlinestatic

GCD of A and B over $ F_{p}(\alpha ) $.

Parameters
[in]Apoly over F_q
[in]Bpoly over F_q
[in]alphaalgebraic variable

Definition at line 28 of file cfModGcd.h.

32 {
33  CFList list;
34  bool top_level= true;
35  return modGCDFq (A, B, alpha, list, top_level);
36 }
CanonicalForm modGCDFq(const CanonicalForm &F, const CanonicalForm &G, Variable &alpha, CFList &l, bool &top_level)
Definition: cfModGcd.cc:464

◆ modGCDFq() [2/2]

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

Definition at line 464 of file cfModGcd.cc.

466 {
467  CanonicalForm dummy1, dummy2;
468  CanonicalForm result= modGCDFq (F, G, dummy1, dummy2, alpha, l,
469  topLevel);
470  return result;
471 }

◆ modGCDGF() [1/2]

static CanonicalForm modGCDGF ( const CanonicalForm A,
const CanonicalForm B 
)
inlinestatic

GCD of A and B over GF.

Parameters
[in]Apoly over GF
[in]Bpoly over GF

Definition at line 74 of file cfModGcd.h.

77 {
79  "GF as base field expected");
80  CFList list;
81  bool top_level= true;
82  return modGCDGF (A, B, list, top_level);
83 }
CanonicalForm modGCDGF(const CanonicalForm &F, const CanonicalForm &G, CFList &l, bool &top_level)
Definition: cfModGcd.cc:860
#define GaloisFieldDomain
Definition: cf_defs.h:24
static int gettype()
Definition: cf_factory.h:28

◆ modGCDGF() [2/2]

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

Definition at line 860 of file cfModGcd.cc.

862 {
863  CanonicalForm dummy1, dummy2;
864  CanonicalForm result= modGCDGF (F, G, dummy1, dummy2, l, topLevel);
865  return result;
866 }
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 Gedde...
Definition: cfModGcd.cc:874

◆ modGCDZ()

CanonicalForm modGCDZ ( const CanonicalForm FF,
const CanonicalForm GG 
)

modular GCD over Z

Parameters
[in]FFpoly over Z
[in]GGpoly over Z

◆ sparseGCDFp() [1/2]

static CanonicalForm sparseGCDFp ( const CanonicalForm A,
const CanonicalForm B 
)
inlinestatic

Zippel's sparse GCD over Fp.

Parameters
[in]Apoly over F_p
[in]Bpoly over F_p

Definition at line 90 of file cfModGcd.h.

93 {
95  "Fp as base field expected");
96  CFList list;
97  bool topLevel= true;
98  return sparseGCDFp (A, B, topLevel, list);
99 }
CanonicalForm sparseGCDFp(const CanonicalForm &F, const CanonicalForm &G, bool &topLevel, CFList &l)
Definition: cfModGcd.cc:3564
#define FiniteFieldDomain
Definition: cf_defs.h:25

◆ sparseGCDFp() [2/2]

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

Definition at line 3564 of file cfModGcd.cc.

3566 {
3567  CanonicalForm A= F;
3568  CanonicalForm B= G;
3569  if (F.isZero() && degree(G) > 0) return G/Lc(G);
3570  else if (G.isZero() && degree (F) > 0) return F/Lc(F);
3571  else if (F.isZero() && G.isZero()) return F.genOne();
3572  if (F.inBaseDomain() || G.inBaseDomain()) return F.genOne();
3573  if (F.isUnivariate() && fdivides(F, G)) return F/Lc(F);
3574  if (G.isUnivariate() && fdivides(G, F)) return G/Lc(G);
3575  if (F == G) return F/Lc(F);
3576 
3577  CFMap M,N;
3578  int best_level= myCompress (A, B, M, N, topLevel);
3579 
3580  if (best_level == 0) return B.genOne();
3581 
3582  A= M(A);
3583  B= M(B);
3584 
3585  Variable x= Variable (1);
3586 
3587  //univariate case
3588  if (A.isUnivariate() && B.isUnivariate())
3589  return N (gcd (A, B));
3590 
3591  CanonicalForm cA, cB; // content of A and B
3592  CanonicalForm ppA, ppB; // primitive part of A and B
3593  CanonicalForm gcdcAcB;
3594 
3595  cA = uni_content (A);
3596  cB = uni_content (B);
3597  gcdcAcB= gcd (cA, cB);
3598  ppA= A/cA;
3599  ppB= B/cB;
3600 
3601  CanonicalForm lcA, lcB; // leading coefficients of A and B
3602  CanonicalForm gcdlcAlcB;
3603  lcA= uni_lcoeff (ppA);
3604  lcB= uni_lcoeff (ppB);
3605 
3606  if (fdivides (lcA, lcB))
3607  {
3608  if (fdivides (A, B))
3609  return F/Lc(F);
3610  }
3611  if (fdivides (lcB, lcA))
3612  {
3613  if (fdivides (B, A))
3614  return G/Lc(G);
3615  }
3616 
3617  gcdlcAlcB= gcd (lcA, lcB);
3618 
3619  int d= totaldegree (ppA, Variable (2), Variable (ppA.level()));
3620  int d0;
3621 
3622  if (d == 0)
3623  return N(gcdcAcB);
3624 
3625  d0= totaldegree (ppB, Variable (2), Variable (ppB.level()));
3626 
3627  if (d0 < d)
3628  d= d0;
3629 
3630  if (d == 0)
3631  return N(gcdcAcB);
3632 
3633  CanonicalForm m, random_element, G_m, G_random_element, H, cH, ppH, skeleton;
3634  CanonicalForm newtonPoly= 1;
3635  m= gcdlcAlcB;
3636  G_m= 0;
3637  H= 0;
3638  bool fail= false;
3639  topLevel= false;
3640  bool inextension= false;
3641  Variable V_buf, alpha, cleanUp;
3642  CanonicalForm prim_elem, im_prim_elem;
3643  CFList source, dest;
3644  do //first do
3645  {
3646  if (inextension)
3647  random_element= randomElement (m, V_buf, l, fail);
3648  else
3649  random_element= FpRandomElement (m, l, fail);
3650  if (random_element == 0 && !fail)
3651  {
3652  if (!find (l, random_element))
3653  l.append (random_element);
3654  continue;
3655  }
3656 
3657  if (!fail && !inextension)
3658  {
3659  CFList list;
3660  TIMING_START (gcd_recursion);
3661  G_random_element=
3662  sparseGCDFp (ppA (random_element,x), ppB (random_element,x), topLevel,
3663  list);
3664  TIMING_END_AND_PRINT (gcd_recursion,
3665  "time for recursive call: ");
3666  DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3667  }
3668  else if (!fail && inextension)
3669  {
3670  CFList list;
3671  TIMING_START (gcd_recursion);
3672  G_random_element=
3673  sparseGCDFq (ppA (random_element, x), ppB (random_element, x), alpha,
3674  list, topLevel);
3675  TIMING_END_AND_PRINT (gcd_recursion,
3676  "time for recursive call: ");
3677  DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3678  }
3679  else if (fail && !inextension)
3680  {
3681  source= CFList();
3682  dest= CFList();
3683  CFList list;
3685  int deg= 2;
3686  bool initialized= false;
3687  do
3688  {
3689  mipo= randomIrredpoly (deg, x);
3690  if (initialized)
3691  setMipo (alpha, mipo);
3692  else
3693  alpha= rootOf (mipo);
3694  inextension= true;
3695  fail= false;
3696  initialized= true;
3697  random_element= randomElement (m, alpha, l, fail);
3698  deg++;
3699  } while (fail);
3700  cleanUp= alpha;
3701  V_buf= alpha;
3702  list= CFList();
3703  TIMING_START (gcd_recursion);
3704  G_random_element=
3705  sparseGCDFq (ppA (random_element, x), ppB (random_element, x), alpha,
3706  list, topLevel);
3707  TIMING_END_AND_PRINT (gcd_recursion,
3708  "time for recursive call: ");
3709  DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3710  }
3711  else if (fail && inextension)
3712  {
3713  source= CFList();
3714  dest= CFList();
3715 
3716  Variable V_buf3= V_buf;
3717  V_buf= chooseExtension (V_buf);
3718  bool prim_fail= false;
3719  Variable V_buf2;
3720  CanonicalForm prim_elem, im_prim_elem;
3721  prim_elem= primitiveElement (alpha, V_buf2, prim_fail);
3722 
3723  if (V_buf3 != alpha)
3724  {
3725  m= mapDown (m, prim_elem, im_prim_elem, alpha, source, dest);
3726  G_m= mapDown (m, prim_elem, im_prim_elem, alpha, source, dest);
3727  newtonPoly= mapDown (newtonPoly, prim_elem, im_prim_elem, alpha, source,
3728  dest);
3729  ppA= mapDown (ppA, prim_elem, im_prim_elem, alpha, source, dest);
3730  ppB= mapDown (ppB, prim_elem, im_prim_elem, alpha, source, dest);
3731  gcdlcAlcB= mapDown (gcdlcAlcB, prim_elem, im_prim_elem, alpha, source,
3732  dest);
3733  for (CFListIterator i= l; i.hasItem(); i++)
3734  i.getItem()= mapDown (i.getItem(), prim_elem, im_prim_elem, alpha,
3735  source, dest);
3736  }
3737 
3738  ASSERT (!prim_fail, "failure in integer factorizer");
3739  if (prim_fail)
3740  ; //ERROR
3741  else
3742  im_prim_elem= mapPrimElem (prim_elem, alpha, V_buf);
3743 
3744  DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (alpha));
3745  DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (V_buf2));
3746 
3747  for (CFListIterator i= l; i.hasItem(); i++)
3748  i.getItem()= mapUp (i.getItem(), alpha, V_buf, prim_elem,
3749  im_prim_elem, source, dest);
3750  m= mapUp (m, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
3751  G_m= mapUp (G_m, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
3752  newtonPoly= mapUp (newtonPoly, alpha, V_buf, prim_elem, im_prim_elem,
3753  source, dest);
3754  ppA= mapUp (ppA, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
3755  ppB= mapUp (ppB, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
3756  gcdlcAlcB= mapUp (gcdlcAlcB, alpha, V_buf, prim_elem, im_prim_elem,
3757  source, dest);
3758  fail= false;
3759  random_element= randomElement (m, V_buf, l, fail );
3760  DEBOUTLN (cerr, "fail= " << fail);
3761  CFList list;
3762  TIMING_START (gcd_recursion);
3763  G_random_element=
3764  sparseGCDFq (ppA (random_element, x), ppB (random_element, x), V_buf,
3765  list, topLevel);
3766  TIMING_END_AND_PRINT (gcd_recursion,
3767  "time for recursive call: ");
3768  DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3769  alpha= V_buf;
3770  }
3771 
3772  if (!G_random_element.inCoeffDomain())
3773  d0= totaldegree (G_random_element, Variable(2),
3774  Variable (G_random_element.level()));
3775  else
3776  d0= 0;
3777 
3778  if (d0 == 0)
3779  {
3780  if (inextension)
3781  prune (cleanUp);
3782  return N(gcdcAcB);
3783  }
3784  if (d0 > d)
3785  {
3786  if (!find (l, random_element))
3787  l.append (random_element);
3788  continue;
3789  }
3790 
3791  G_random_element=
3792  (gcdlcAlcB(random_element, x)/uni_lcoeff (G_random_element))
3793  * G_random_element;
3794 
3795  skeleton= G_random_element;
3796 
3797  if (!G_random_element.inCoeffDomain())
3798  d0= totaldegree (G_random_element, Variable(2),
3799  Variable (G_random_element.level()));
3800  else
3801  d0= 0;
3802 
3803  if (d0 < d)
3804  {
3805  m= gcdlcAlcB;
3806  newtonPoly= 1;
3807  G_m= 0;
3808  d= d0;
3809  }
3810 
3811  H= newtonInterp (random_element, G_random_element, newtonPoly, G_m, x);
3812 
3813  if (uni_lcoeff (H) == gcdlcAlcB)
3814  {
3815  cH= uni_content (H);
3816  ppH= H/cH;
3817  ppH /= Lc (ppH);
3818  DEBOUTLN (cerr, "ppH= " << ppH);
3819 
3820  if (fdivides (ppH, ppA) && fdivides (ppH, ppB))
3821  {
3822  if (inextension)
3823  prune (cleanUp);
3824  return N(gcdcAcB*ppH);
3825  }
3826  }
3827  G_m= H;
3828  newtonPoly= newtonPoly*(x - random_element);
3829  m= m*(x - random_element);
3830  if (!find (l, random_element))
3831  l.append (random_element);
3832 
3833  if ((getCharacteristic() > 3 && size (skeleton) < 200))
3834  {
3835  CFArray Monoms;
3836  CFArray* coeffMonoms;
3837 
3838  do //second do
3839  {
3840  if (inextension)
3841  random_element= randomElement (m, alpha, l, fail);
3842  else
3843  random_element= FpRandomElement (m, l, fail);
3844  if (random_element == 0 && !fail)
3845  {
3846  if (!find (l, random_element))
3847  l.append (random_element);
3848  continue;
3849  }
3850 
3851  bool sparseFail= false;
3852  if (!fail && !inextension)
3853  {
3854  CFList list;
3855  TIMING_START (gcd_recursion);
3856  if (LC (skeleton).inCoeffDomain())
3857  G_random_element=
3858  monicSparseInterpol(ppA(random_element, x), ppB (random_element, x),
3859  skeleton, x, sparseFail, coeffMonoms,
3860  Monoms);
3861  else
3862  G_random_element=
3863  nonMonicSparseInterpol(ppA(random_element,x), ppB(random_element,x),
3864  skeleton, x, sparseFail,
3865  coeffMonoms, Monoms);
3866  TIMING_END_AND_PRINT (gcd_recursion,
3867  "time for recursive call: ");
3868  DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3869  }
3870  else if (!fail && inextension)
3871  {
3872  CFList list;
3873  TIMING_START (gcd_recursion);
3874  if (LC (skeleton).inCoeffDomain())
3875  G_random_element=
3876  monicSparseInterpol(ppA (random_element,x), ppB (random_element, x),
3877  skeleton, alpha, sparseFail, coeffMonoms,
3878  Monoms);
3879  else
3880  G_random_element=
3881  nonMonicSparseInterpol(ppA(random_element,x), ppB(random_element,x),
3882  skeleton, alpha, sparseFail, coeffMonoms,
3883  Monoms);
3884  TIMING_END_AND_PRINT (gcd_recursion,
3885  "time for recursive call: ");
3886  DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3887  }
3888  else if (fail && !inextension)
3889  {
3890  source= CFList();
3891  dest= CFList();
3892  CFList list;
3894  int deg= 2;
3895  bool initialized= false;
3896  do
3897  {
3898  mipo= randomIrredpoly (deg, x);
3899  if (initialized)
3900  setMipo (alpha, mipo);
3901  else
3902  alpha= rootOf (mipo);
3903  inextension= true;
3904  fail= false;
3905  initialized= true;
3906  random_element= randomElement (m, alpha, l, fail);
3907  deg++;
3908  } while (fail);
3909  cleanUp= alpha;
3910  V_buf= alpha;
3911  list= CFList();
3912  TIMING_START (gcd_recursion);
3913  if (LC (skeleton).inCoeffDomain())
3914  G_random_element=
3915  monicSparseInterpol (ppA (random_element,x), ppB (random_element,x),
3916  skeleton, alpha, sparseFail, coeffMonoms,
3917  Monoms);
3918  else
3919  G_random_element=
3920  nonMonicSparseInterpol(ppA(random_element,x), ppB(random_element,x),
3921  skeleton, alpha, sparseFail, coeffMonoms,
3922  Monoms);
3923  TIMING_END_AND_PRINT (gcd_recursion,
3924  "time for recursive call: ");
3925  DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3926  }
3927  else if (fail && inextension)
3928  {
3929  source= CFList();
3930  dest= CFList();
3931 
3932  Variable V_buf3= V_buf;
3933  V_buf= chooseExtension (V_buf);
3934  bool prim_fail= false;
3935  Variable V_buf2;
3936  CanonicalForm prim_elem, im_prim_elem;
3937  prim_elem= primitiveElement (alpha, V_buf2, prim_fail);
3938 
3939  if (V_buf3 != alpha)
3940  {
3941  m= mapDown (m, prim_elem, im_prim_elem, alpha, source, dest);
3942  G_m= mapDown (m, prim_elem, im_prim_elem, alpha, source, dest);
3943  newtonPoly= mapDown (newtonPoly, prim_elem, im_prim_elem, alpha,
3944  source, dest);
3945  ppA= mapDown (ppA, prim_elem, im_prim_elem, alpha, source, dest);
3946  ppB= mapDown (ppB, prim_elem, im_prim_elem, alpha, source, dest);
3947  gcdlcAlcB= mapDown (gcdlcAlcB, prim_elem, im_prim_elem, alpha,
3948  source, dest);
3949  for (CFListIterator i= l; i.hasItem(); i++)
3950  i.getItem()= mapDown (i.getItem(), prim_elem, im_prim_elem, alpha,
3951  source, dest);
3952  }
3953 
3954  ASSERT (!prim_fail, "failure in integer factorizer");
3955  if (prim_fail)
3956  ; //ERROR
3957  else
3958  im_prim_elem= mapPrimElem (prim_elem, alpha, V_buf);
3959 
3960  DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (alpha));
3961  DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (V_buf2));
3962 
3963  for (CFListIterator i= l; i.hasItem(); i++)
3964  i.getItem()= mapUp (i.getItem(), alpha, V_buf, prim_elem,
3965  im_prim_elem, source, dest);
3966  m= mapUp (m, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
3967  G_m= mapUp (G_m, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
3968  newtonPoly= mapUp (newtonPoly, alpha, V_buf, prim_elem, im_prim_elem,
3969  source, dest);
3970  ppA= mapUp (ppA, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
3971  ppB= mapUp (ppB, alpha, V_buf, prim_elem, im_prim_elem, source, dest);
3972  gcdlcAlcB= mapUp (gcdlcAlcB, alpha, V_buf, prim_elem, im_prim_elem,
3973  source, dest);
3974  fail= false;
3975  random_element= randomElement (m, V_buf, l, fail );
3976  DEBOUTLN (cerr, "fail= " << fail);
3977  CFList list;
3978  TIMING_START (gcd_recursion);
3979  if (LC (skeleton).inCoeffDomain())
3980  G_random_element=
3981  monicSparseInterpol (ppA (random_element, x), ppB (random_element, x),
3982  skeleton, V_buf, sparseFail, coeffMonoms,
3983  Monoms);
3984  else
3985  G_random_element=
3986  nonMonicSparseInterpol (ppA(random_element,x), ppB(random_element,x),
3987  skeleton, V_buf, sparseFail,
3988  coeffMonoms, Monoms);
3989  TIMING_END_AND_PRINT (gcd_recursion,
3990  "time for recursive call: ");
3991  DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3992  alpha= V_buf;
3993  }
3994 
3995  if (sparseFail)
3996  break;
3997 
3998  if (!G_random_element.inCoeffDomain())
3999  d0= totaldegree (G_random_element, Variable(2),
4000  Variable (G_random_element.level()));
4001  else
4002  d0= 0;
4003 
4004  if (d0 == 0)
4005  {
4006  if (inextension)
4007  prune (cleanUp);
4008  return N(gcdcAcB);
4009  }
4010  if (d0 > d)
4011  {
4012  if (!find (l, random_element))
4013  l.append (random_element);
4014  continue;
4015  }
4016 
4017  G_random_element=
4018  (gcdlcAlcB(random_element, x)/uni_lcoeff (G_random_element))
4019  * G_random_element;
4020 
4021  if (!G_random_element.inCoeffDomain())
4022  d0= totaldegree (G_random_element, Variable(2),
4023  Variable (G_random_element.level()));
4024  else
4025  d0= 0;
4026 
4027  if (d0 < d)
4028  {
4029  m= gcdlcAlcB;
4030  newtonPoly= 1;
4031  G_m= 0;
4032  d= d0;
4033  }
4034 
4035  TIMING_START (newton_interpolation);
4036  H= newtonInterp (random_element, G_random_element, newtonPoly, G_m, x);
4037  TIMING_END_AND_PRINT (newton_interpolation,
4038  "time for newton interpolation: ");
4039 
4040  //termination test
4041  if (uni_lcoeff (H) == gcdlcAlcB)
4042  {
4043  cH= uni_content (H);
4044  ppH= H/cH;
4045  ppH /= Lc (ppH);
4046  DEBOUTLN (cerr, "ppH= " << ppH);
4047  if (fdivides (ppH, ppA) && fdivides (ppH, ppB))
4048  {
4049  if (inextension)
4050  prune (cleanUp);
4051  return N(gcdcAcB*ppH);
4052  }
4053  }
4054 
4055  G_m= H;
4056  newtonPoly= newtonPoly*(x - random_element);
4057  m= m*(x - random_element);
4058  if (!find (l, random_element))
4059  l.append (random_element);
4060 
4061  } while (1); //end of second do
4062  }
4063  else
4064  {
4065  if (inextension)
4066  prune (cleanUp);
4067  return N(gcdcAcB*modGCDFp (ppA, ppB));
4068  }
4069  } while (1); //end of first do
4070 }
CanonicalForm LC(const CanonicalForm &f)
int FACTORY_PUBLIC getCharacteristic()
Definition: cf_char.cc:70
CanonicalForm sparseGCDFq(const CanonicalForm &F, const CanonicalForm &G, const Variable &alpha, CFList &l, bool &topLevel)
Definition: cfModGcd.cc:3132
CanonicalForm monicSparseInterpol(const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &skeleton, const Variable &alpha, bool &fail, CFArray *&coeffMonoms, CFArray &Monoms)
Definition: cfModGcd.cc:2201
CanonicalForm nonMonicSparseInterpol(const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &skeleton, const Variable &alpha, bool &fail, CFArray *&coeffMonoms, CFArray &Monoms)
Definition: cfModGcd.cc:2485
CanonicalForm sparseGCDFp(const CanonicalForm &F, const CanonicalForm &G, bool &topLevel, CFList &l)
Definition: cfModGcd.cc:3564

◆ sparseGCDFq() [1/2]

static CanonicalForm sparseGCDFq ( const CanonicalForm A,
const CanonicalForm B,
const Variable alpha 
)
inlinestatic

Zippel's sparse GCD over Fq.

Parameters
[in]Apoly over F_q
[in]Bpoly over F_q
[in]alphaalgebraic variable

Definition at line 108 of file cfModGcd.h.

112 {
113  CFList list;
114  bool topLevel= true;
115  return sparseGCDFq (A, B, alpha, list, topLevel);
116 }
CanonicalForm sparseGCDFq(const CanonicalForm &F, const CanonicalForm &G, const Variable &alpha, CFList &l, bool &topLevel)
Definition: cfModGcd.cc:3132

◆ sparseGCDFq() [2/2]

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

Definition at line 3132 of file cfModGcd.cc.

3134 {
3135  CanonicalForm A= F;
3136  CanonicalForm B= G;
3137  if (F.isZero() && degree(G) > 0) return G/Lc(G);
3138  else if (G.isZero() && degree (F) > 0) return F/Lc(F);
3139  else if (F.isZero() && G.isZero()) return F.genOne();
3140  if (F.inBaseDomain() || G.inBaseDomain()) return F.genOne();
3141  if (F.isUnivariate() && fdivides(F, G)) return F/Lc(F);
3142  if (G.isUnivariate() && fdivides(G, F)) return G/Lc(G);
3143  if (F == G) return F/Lc(F);
3144 
3145  CFMap M,N;
3146  int best_level= myCompress (A, B, M, N, topLevel);
3147 
3148  if (best_level == 0) return B.genOne();
3149 
3150  A= M(A);
3151  B= M(B);
3152 
3153  Variable x= Variable (1);
3154 
3155  //univariate case
3156  if (A.isUnivariate() && B.isUnivariate())
3157  return N (gcd (A, B));
3158 
3159  CanonicalForm cA, cB; // content of A and B
3160  CanonicalForm ppA, ppB; // primitive part of A and B
3161  CanonicalForm gcdcAcB;
3162 
3163  cA = uni_content (A);
3164  cB = uni_content (B);
3165  gcdcAcB= gcd (cA, cB);
3166  ppA= A/cA;
3167  ppB= B/cB;
3168 
3169  CanonicalForm lcA, lcB; // leading coefficients of A and B
3170  CanonicalForm gcdlcAlcB;
3171  lcA= uni_lcoeff (ppA);
3172  lcB= uni_lcoeff (ppB);
3173 
3174  if (fdivides (lcA, lcB))
3175  {
3176  if (fdivides (A, B))
3177  return F/Lc(F);
3178  }
3179  if (fdivides (lcB, lcA))
3180  {
3181  if (fdivides (B, A))
3182  return G/Lc(G);
3183  }
3184 
3185  gcdlcAlcB= gcd (lcA, lcB);
3186 
3187  int d= totaldegree (ppA, Variable (2), Variable (ppA.level()));
3188  int d0;
3189 
3190  if (d == 0)
3191  return N(gcdcAcB);
3192  d0= totaldegree (ppB, Variable (2), Variable (ppB.level()));
3193 
3194  if (d0 < d)
3195  d= d0;
3196 
3197  if (d == 0)
3198  return N(gcdcAcB);
3199 
3200  CanonicalForm m, random_element, G_m, G_random_element, H, cH, ppH, skeleton;
3201  CanonicalForm newtonPoly= 1;
3202  m= gcdlcAlcB;
3203  G_m= 0;
3204  H= 0;
3205  bool fail= false;
3206  topLevel= false;
3207  bool inextension= false;
3208  Variable V_buf= alpha, V_buf4= alpha;
3209  CanonicalForm prim_elem, im_prim_elem;
3210  CanonicalForm prim_elem_alpha, im_prim_elem_alpha;
3211  CFList source, dest;
3212  do // first do
3213  {
3214  random_element= randomElement (m, V_buf, l, fail);
3215  if (random_element == 0 && !fail)
3216  {
3217  if (!find (l, random_element))
3218  l.append (random_element);
3219  continue;
3220  }
3221  if (fail)
3222  {
3223  source= CFList();
3224  dest= CFList();
3225 
3226  Variable V_buf3= V_buf;
3227  V_buf= chooseExtension (V_buf);
3228  bool prim_fail= false;
3229  Variable V_buf2;
3230  prim_elem= primitiveElement (V_buf4, V_buf2, prim_fail);
3231  if (V_buf4 == alpha)
3232  prim_elem_alpha= prim_elem;
3233 
3234  if (V_buf3 != V_buf4)
3235  {
3236  m= mapDown (m, prim_elem, im_prim_elem, V_buf4, source, dest);
3237  G_m= mapDown (m, prim_elem, im_prim_elem, V_buf4, source, dest);
3238  newtonPoly= mapDown (newtonPoly, prim_elem, im_prim_elem, V_buf4,
3239  source, dest);
3240  ppA= mapDown (ppA, prim_elem, im_prim_elem, V_buf4, source, dest);
3241  ppB= mapDown (ppB, prim_elem, im_prim_elem, V_buf4, source, dest);
3242  gcdlcAlcB= mapDown (gcdlcAlcB, prim_elem, im_prim_elem, V_buf4, source,
3243  dest);
3244  for (CFListIterator i= l; i.hasItem(); i++)
3245  i.getItem()= mapDown (i.getItem(), prim_elem, im_prim_elem, V_buf4,
3246  source, dest);
3247  }
3248 
3249  ASSERT (!prim_fail, "failure in integer factorizer");
3250  if (prim_fail)
3251  ; //ERROR
3252  else
3253  im_prim_elem= mapPrimElem (prim_elem, V_buf4, V_buf);
3254 
3255  if (V_buf4 == alpha)
3256  im_prim_elem_alpha= im_prim_elem;
3257  else
3258  im_prim_elem_alpha= mapUp (im_prim_elem_alpha, V_buf4, V_buf, prim_elem,
3259  im_prim_elem, source, dest);
3260 
3261  DEBOUTLN (cerr, "getMipo (V_buf4)= " << getMipo (V_buf4));
3262  DEBOUTLN (cerr, "getMipo (V_buf2)= " << getMipo (V_buf2));
3263  inextension= true;
3264  for (CFListIterator i= l; i.hasItem(); i++)
3265  i.getItem()= mapUp (i.getItem(), V_buf4, V_buf, prim_elem,
3266  im_prim_elem, source, dest);
3267  m= mapUp (m, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
3268  G_m= mapUp (G_m, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
3269  newtonPoly= mapUp (newtonPoly, V_buf4, V_buf, prim_elem, im_prim_elem,
3270  source, dest);
3271  ppA= mapUp (ppA, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
3272  ppB= mapUp (ppB, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
3273  gcdlcAlcB= mapUp (gcdlcAlcB, V_buf4, V_buf, prim_elem, im_prim_elem,
3274  source, dest);
3275 
3276  fail= false;
3277  random_element= randomElement (m, V_buf, l, fail );
3278  DEBOUTLN (cerr, "fail= " << fail);
3279  CFList list;
3280  TIMING_START (gcd_recursion);
3281  G_random_element=
3282  sparseGCDFq (ppA (random_element, x), ppB (random_element, x), V_buf,
3283  list, topLevel);
3284  TIMING_END_AND_PRINT (gcd_recursion,
3285  "time for recursive call: ");
3286  DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3287  V_buf4= V_buf;
3288  }
3289  else
3290  {
3291  CFList list;
3292  TIMING_START (gcd_recursion);
3293  G_random_element=
3294  sparseGCDFq (ppA(random_element, x), ppB(random_element, x), V_buf,
3295  list, topLevel);
3296  TIMING_END_AND_PRINT (gcd_recursion,
3297  "time for recursive call: ");
3298  DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3299  }
3300 
3301  if (!G_random_element.inCoeffDomain())
3302  d0= totaldegree (G_random_element, Variable(2),
3303  Variable (G_random_element.level()));
3304  else
3305  d0= 0;
3306 
3307  if (d0 == 0)
3308  {
3309  if (inextension)
3310  prune1 (alpha);
3311  return N(gcdcAcB);
3312  }
3313  if (d0 > d)
3314  {
3315  if (!find (l, random_element))
3316  l.append (random_element);
3317  continue;
3318  }
3319 
3320  G_random_element=
3321  (gcdlcAlcB(random_element, x)/uni_lcoeff (G_random_element))
3322  * G_random_element;
3323 
3324  skeleton= G_random_element;
3325  if (!G_random_element.inCoeffDomain())
3326  d0= totaldegree (G_random_element, Variable(2),
3327  Variable (G_random_element.level()));
3328  else
3329  d0= 0;
3330 
3331  if (d0 < d)
3332  {
3333  m= gcdlcAlcB;
3334  newtonPoly= 1;
3335  G_m= 0;
3336  d= d0;
3337  }
3338 
3339  H= newtonInterp (random_element, G_random_element, newtonPoly, G_m, x);
3340  if (uni_lcoeff (H) == gcdlcAlcB)
3341  {
3342  cH= uni_content (H);
3343  ppH= H/cH;
3344  if (inextension)
3345  {
3346  CFList u, v;
3347  //maybe it's better to test if ppH is an element of F(\alpha) before
3348  //mapping down
3349  if (fdivides (ppH, ppA) && fdivides (ppH, ppB))
3350  {
3351  DEBOUTLN (cerr, "ppH before mapDown= " << ppH);
3352  ppH= mapDown (ppH, prim_elem_alpha, im_prim_elem_alpha, alpha, u, v);
3353  ppH /= Lc(ppH);
3354  DEBOUTLN (cerr, "ppH after mapDown= " << ppH);
3355  prune1 (alpha);
3356  return N(gcdcAcB*ppH);
3357  }
3358  }
3359  else if (fdivides (ppH, ppA) && fdivides (ppH, ppB))
3360  return N(gcdcAcB*ppH);
3361  }
3362  G_m= H;
3363  newtonPoly= newtonPoly*(x - random_element);
3364  m= m*(x - random_element);
3365  if (!find (l, random_element))
3366  l.append (random_element);
3367 
3368  if (getCharacteristic () > 3 && size (skeleton) < 100)
3369  {
3370  CFArray Monoms;
3371  CFArray *coeffMonoms;
3372  do //second do
3373  {
3374  random_element= randomElement (m, V_buf, l, fail);
3375  if (random_element == 0 && !fail)
3376  {
3377  if (!find (l, random_element))
3378  l.append (random_element);
3379  continue;
3380  }
3381  if (fail)
3382  {
3383  source= CFList();
3384  dest= CFList();
3385 
3386  Variable V_buf3= V_buf;
3387  V_buf= chooseExtension (V_buf);
3388  bool prim_fail= false;
3389  Variable V_buf2;
3390  prim_elem= primitiveElement (V_buf4, V_buf2, prim_fail);
3391  if (V_buf4 == alpha)
3392  prim_elem_alpha= prim_elem;
3393 
3394  if (V_buf3 != V_buf4)
3395  {
3396  m= mapDown (m, prim_elem, im_prim_elem, V_buf4, source, dest);
3397  G_m= mapDown (m, prim_elem, im_prim_elem, V_buf4, source, dest);
3398  newtonPoly= mapDown (newtonPoly, prim_elem, im_prim_elem, V_buf4,
3399  source, dest);
3400  ppA= mapDown (ppA, prim_elem, im_prim_elem, V_buf4, source, dest);
3401  ppB= mapDown (ppB, prim_elem, im_prim_elem, V_buf4, source, dest);
3402  gcdlcAlcB= mapDown (gcdlcAlcB, prim_elem, im_prim_elem, V_buf4,
3403  source, dest);
3404  for (CFListIterator i= l; i.hasItem(); i++)
3405  i.getItem()= mapDown (i.getItem(), prim_elem, im_prim_elem, V_buf4,
3406  source, dest);
3407  }
3408 
3409  ASSERT (!prim_fail, "failure in integer factorizer");
3410  if (prim_fail)
3411  ; //ERROR
3412  else
3413  im_prim_elem= mapPrimElem (prim_elem, V_buf4, V_buf);
3414 
3415  if (V_buf4 == alpha)
3416  im_prim_elem_alpha= im_prim_elem;
3417  else
3418  im_prim_elem_alpha= mapUp (im_prim_elem_alpha, V_buf4, V_buf,
3419  prim_elem, im_prim_elem, source, dest);
3420 
3421  DEBOUTLN (cerr, "getMipo (V_buf4)= " << getMipo (V_buf4));
3422  DEBOUTLN (cerr, "getMipo (V_buf2)= " << getMipo (V_buf2));
3423  inextension= true;
3424  for (CFListIterator i= l; i.hasItem(); i++)
3425  i.getItem()= mapUp (i.getItem(), V_buf4, V_buf, prim_elem,
3426  im_prim_elem, source, dest);
3427  m= mapUp (m, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
3428  G_m= mapUp (G_m, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
3429  newtonPoly= mapUp (newtonPoly, V_buf4, V_buf, prim_elem, im_prim_elem,
3430  source, dest);
3431  ppA= mapUp (ppA, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
3432  ppB= mapUp (ppB, V_buf4, V_buf, prim_elem, im_prim_elem, source, dest);
3433 
3434  gcdlcAlcB= mapUp (gcdlcAlcB, V_buf4, V_buf, prim_elem, im_prim_elem,
3435  source, dest);
3436 
3437  fail= false;
3438  random_element= randomElement (m, V_buf, l, fail);
3439  DEBOUTLN (cerr, "fail= " << fail);
3440  CFList list;
3441  TIMING_START (gcd_recursion);
3442 
3443  V_buf4= V_buf;
3444 
3445  //sparseInterpolation
3446  bool sparseFail= false;
3447  if (LC (skeleton).inCoeffDomain())
3448  G_random_element=
3449  monicSparseInterpol (ppA (random_element, x), ppB(random_element,x),
3450  skeleton,V_buf, sparseFail, coeffMonoms,Monoms);
3451  else
3452  G_random_element=
3453  nonMonicSparseInterpol (ppA(random_element,x),ppB(random_element,x),
3454  skeleton, V_buf, sparseFail, coeffMonoms,
3455  Monoms);
3456  if (sparseFail)
3457  break;
3458 
3459  TIMING_END_AND_PRINT (gcd_recursion,
3460  "time for recursive call: ");
3461  DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3462  }
3463  else
3464  {
3465  CFList list;
3466  TIMING_START (gcd_recursion);
3467  bool sparseFail= false;
3468  if (LC (skeleton).inCoeffDomain())
3469  G_random_element=
3470  monicSparseInterpol (ppA (random_element, x),ppB(random_element, x),
3471  skeleton,V_buf, sparseFail,coeffMonoms, Monoms);
3472  else
3473  G_random_element=
3474  nonMonicSparseInterpol (ppA(random_element,x),ppB(random_element,x),
3475  skeleton, V_buf, sparseFail, coeffMonoms,
3476  Monoms);
3477  if (sparseFail)
3478  break;
3479 
3480  TIMING_END_AND_PRINT (gcd_recursion,
3481  "time for recursive call: ");
3482  DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3483  }
3484 
3485  if (!G_random_element.inCoeffDomain())
3486  d0= totaldegree (G_random_element, Variable(2),
3487  Variable (G_random_element.level()));
3488  else
3489  d0= 0;
3490 
3491  if (d0 == 0)
3492  {
3493  if (inextension)
3494  prune1 (alpha);
3495  return N(gcdcAcB);
3496  }
3497  if (d0 > d)
3498  {
3499  if (!find (l, random_element))
3500  l.append (random_element);
3501  continue;
3502  }
3503 
3504  G_random_element=
3505  (gcdlcAlcB(random_element, x)/uni_lcoeff (G_random_element))
3506  * G_random_element;
3507 
3508  if (!G_random_element.inCoeffDomain())
3509  d0= totaldegree (G_random_element, Variable(2),
3510  Variable (G_random_element.level()));
3511  else
3512  d0= 0;
3513 
3514  if (d0 < d)
3515  {
3516  m= gcdlcAlcB;
3517  newtonPoly= 1;
3518  G_m= 0;
3519  d= d0;
3520  }
3521 
3522  TIMING_START (newton_interpolation);
3523  H= newtonInterp (random_element, G_random_element, newtonPoly, G_m, x);
3524  TIMING_END_AND_PRINT (newton_interpolation,
3525  "time for newton interpolation: ");
3526 
3527  //termination test
3528  if (uni_lcoeff (H) == gcdlcAlcB)
3529  {
3530  cH= uni_content (H);
3531  ppH= H/cH;
3532  if (inextension)
3533  {
3534  CFList u, v;
3535  //maybe it's better to test if ppH is an element of F(\alpha) before
3536  //mapping down
3537  if (fdivides (ppH, ppA) && fdivides (ppH, ppB))
3538  {
3539  DEBOUTLN (cerr, "ppH before mapDown= " << ppH);
3540  ppH= mapDown (ppH, prim_elem_alpha, im_prim_elem_alpha, alpha, u, v);
3541  ppH /= Lc(ppH);
3542  DEBOUTLN (cerr, "ppH after mapDown= " << ppH);
3543  prune1 (alpha);
3544  return N(gcdcAcB*ppH);
3545  }
3546  }
3547  else if (fdivides (ppH, ppA) && fdivides (ppH, ppB))
3548  {
3549  return N(gcdcAcB*ppH);
3550  }
3551  }
3552 
3553  G_m= H;
3554  newtonPoly= newtonPoly*(x - random_element);
3555  m= m*(x - random_element);
3556  if (!find (l, random_element))
3557  l.append (random_element);
3558 
3559  } while (1);
3560  }
3561  } while (1);
3562 }
const Variable & v
< [in] a sqrfree bivariate poly
Definition: facBivar.h:39
void prune1(const Variable &alpha)
Definition: variable.cc:291

◆ terminationTest()

bool terminationTest ( const CanonicalForm F,
const CanonicalForm G,
const CanonicalForm coF,
const CanonicalForm coG,
const CanonicalForm cand 
)