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 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 i
Definition cfModGcd.cc:4086
int size() const
class to iterate through CanonicalForm's
Definition cf_iter.h:44
CF_NO_INLINE int hasTerms() const
check if iterator has reached the end of CanonicalForm
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}
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
const CanonicalForm & G
Definition cfModGcd.cc:67

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