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 1957 of file cfModGcd.cc.

1958 {
1959  if (F.inCoeffDomain())
1960  {
1961  CFArray result= CFArray (1);
1962  result [0]= 1;
1963  return result;
1964  }
1965  if (F.isUnivariate())
1966  {
1967  CFArray result= CFArray (size(F));
1968  int j= 0;
1969  for (CFIterator i= F; i.hasTerms(); i++, j++)
1970  result[j]= power (F.mvar(), i.exp());
1971  return result;
1972  }
1973  int numMon= size (F);
1974  CFArray result= CFArray (numMon);
1975  int j= 0;
1976  CFArray recResult;
1977  Variable x= F.mvar();
1978  CanonicalForm powX;
1979  for (CFIterator i= F; i.hasTerms(); i++)
1980  {
1981  powX= power (x, i.exp());
1982  recResult= getMonoms (i.coeff());
1983  for (int k= 0; k < recResult.size(); k++)
1984  result[j+k]= powX*recResult[k];
1985  j += recResult.size();
1986  }
1987  return result;
1988 }
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:1957
Variable x
Definition: cfModGcd.cc:4082
int i
Definition: cfModGcd.cc:4078
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:127
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:1212
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 1212 of file cfModGcd.cc.

1214 {
1215  CanonicalForm dummy1, dummy2;
1216  CanonicalForm result= modGCDFp (F, G, dummy1, dummy2, topLevel, l);
1217  return result;
1218 }
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:1223
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 1223 of file cfModGcd.cc.

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

◆ modGCDFq() [2/2]

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

Definition at line 462 of file cfModGcd.cc.

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

◆ 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:858
#define GaloisFieldDomain
Definition: cf_defs.h:18
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 858 of file cfModGcd.cc.

860 {
861  CanonicalForm dummy1, dummy2;
862  CanonicalForm result= modGCDGF (F, G, dummy1, dummy2, l, topLevel);
863  return result;
864 }
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:872

◆ 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:3562
#define FiniteFieldDomain
Definition: cf_defs.h:19

◆ sparseGCDFp() [2/2]

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

Definition at line 3562 of file cfModGcd.cc.

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

◆ 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:3130

◆ sparseGCDFq() [2/2]

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

Definition at line 3130 of file cfModGcd.cc.

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