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 {
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 
)