My Project
Loading...
Searching...
No Matches
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 ) $.
 
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} $.
 
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.
 
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.
 
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.
 
CFArray getMonoms (const CanonicalForm &F)
 extract monomials of F, parts in algebraic variable are considered coefficients
 
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
 

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

1966{
1967 if (F.inCoeffDomain())
1968 {
1969 CFArray result= CFArray (1);
1970 result [0]= 1;
1971 return result;
1972 }
1973 if (F.isUnivariate())
1974 {
1975 CFArray result= CFArray (size(F));
1976 int j= 0;
1977 for (CFIterator i= F; i.hasTerms(); i++, j++)
1978 result[j]= power (F.mvar(), i.exp());
1979 return result;
1980 }
1981 int numMon= size (F);
1982 CFArray result= CFArray (numMon);
1983 int j= 0;
1984 CFArray recResult;
1985 Variable x= F.mvar();
1986 CanonicalForm powX;
1987 for (CFIterator i= F; i.hasTerms(); i++)
1988 {
1989 powX= power (x, i.exp());
1990 recResult= getMonoms (i.coeff());
1991 for (int k= 0; k < recResult.size(); k++)
1992 result[j+k]= powX*recResult[k];
1993 j += recResult.size();
1994 }
1995 return result;
1996}
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 i
Definition cfEzgcd.cc:132
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:1965
Variable x
Definition cfModGcd.cc:4090
int size() const
class to iterate through CanonicalForm's
Definition cf_iter.h:44
factory's main class
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
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}
List< CanonicalForm > CFList
CanonicalForm modGCDFp(const CanonicalForm &F, const CanonicalForm &G, bool &top_level, CFList &l)
Definition cfModGcd.cc:1213
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 1213 of file cfModGcd.cc.

1215{
1216 CanonicalForm dummy1, dummy2;
1217 CanonicalForm result= modGCDFp (F, G, dummy1, dummy2, topLevel, l);
1218 return result;
1219}
int l
Definition cfEzgcd.cc:100
const CanonicalForm CFMap CFMap bool topLevel
CanonicalForm modGCDFp(const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &coF, CanonicalForm &coG, bool &topLevel, CFList &l)
Definition cfModGcd.cc:1224
STATIC_VAR TreeM * G
Definition janet.cc:31

◆ modGCDFp() [4/4]

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

Definition at line 1224 of file cfModGcd.cc.

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

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

◆ modGCDFq() [2/2]

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

Definition at line 463 of file cfModGcd.cc.

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

◆ 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:859
#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 859 of file cfModGcd.cc.

861{
862 CanonicalForm dummy1, dummy2;
863 CanonicalForm result= modGCDGF (F, G, dummy1, dummy2, l, topLevel);
864 return result;
865}
CanonicalForm modGCDGF(const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &coF, CanonicalForm &coG, CFList &l, bool &topLevel)
GCD of F and G over GF, based on Alg. 7.2. as described in "Algorithms forComputer Algebra" by Geddes...
Definition cfModGcd.cc:873

◆ modGCDZ()

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:3570
#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 3570 of file cfModGcd.cc.

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

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

◆ sparseGCDFq() [2/2]

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

Definition at line 3138 of file cfModGcd.cc.

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