Changeset c28ecf in git for libpolys/polys

Timestamp:

May 19, 2011, 2:52:08 PM (13 years ago)

Author:

Frank Seelisch <seelisch@…>

Branches:

(u'fieker-DuVal', '117eb8c30fc9e991c4decca4832b1d19036c4c65')(u'spielwiese', 'b50cc4fc80fdf2a2005655a08d9229aa4b32d4a5')

Children:

4a2260efc22c19ac8eaf58d529af8decbb9420bb

Parents:

f0b01f1bf67efed9ec49d9df0d55f75284f5d968

git-author:

Frank Seelisch <seelisch@mathematik.uni-kl.de>2011-05-19 14:52:08+02:00

git-committer:

Mohamed Barakat <mohamed.barakat@rwth-aachen.de>2011-11-09 12:36:12+01:00

Message:

alg ext implementation now passing all coeffs and polys tests

Location:

libpolys/polys

Files:

• ext_fields/algext.cc (modified) (6 diffs)
• ext_fields/algext.h (modified) (1 diff)
• monomials/p_polys.cc (modified) (5 diffs)
• monomials/p_polys.h (modified) (1 diff)

libpolys/polys/ext_fields/algext.cc

@@ -75 +75 @@
   assume(getCoeffType(cf) == naID);
   if (a == NULL) return TRUE;
-  p_Test((poly)a, naRing); // cannot use omCheckAddrSize(a, sizeof(*a)) here
+  p_Test((poly)a, naRing);
   assume(p_Deg((poly)a, naRing) <= p_Deg(naMinpoly, naRing));
   return TRUE;
@@ -90 +90 @@
 }
 
-void naDelete(number *a, const coeffs cf)
-{
-  if (a == NULL) return;
-  p_Delete((poly *)a, naRing);
+void naDelete(number * a, const coeffs cf)
+{
+  if (*a == NULL) return;
+  poly aAsPoly = (poly)(*a);
+  p_Delete(&aAsPoly, naRing);
   *a = NULL;
 }
@@ -339 +340 @@
 {
   #ifdef LDEBUG
-  omCheckAddr(p); omCheckAddr(reducer);
+  p_Test((poly)p, naRing);
+  p_Test((poly)reducer, naRing);
   #endif
   if (p_Deg(p, naRing) > 10 * p_Deg(reducer, naRing))
@@ -427 +429 @@
 {
   #ifdef LDEBUG
-  omCheckAddr(p); omCheckAddr(reducer);
+  p_Test((poly)p, naRing);
+  p_Test((poly)reducer, naRing);
   #endif
   p_PolyDiv(p, reducer, FALSE, naRing);
@@ -443 +446 @@
   naTest(a);
   if (a == NULL) WerrorS(nDivBy0);
-  poly aFactor; poly mFactor;
+  poly aFactor = NULL; poly mFactor = NULL;
   poly theGcd = p_ExtGcd((poly)a, aFactor, naMinpoly, mFactor, naRing);
-  /* the gcd must be one since naMinpoly is irreducible and a != NULL: */
+  /* the gcd must be 1 since naMinpoly is irreducible and a != NULL: */
   assume(naIsOne((number)theGcd, cf));      
   p_Delete(&theGcd, naRing);
   p_Delete(&mFactor, naRing);
+  /* printf("naInvers\n");
+     p_Write((poly)a, naRing);
+     p_Write(aFactor, naRing); */
   return (number)(aFactor);
 }
@@ -566 +572 @@
   
   #ifdef LDEBUG
-  omCheckAddr(naMinpoly);
+  p_Test((poly)naMinpoly, naRing);
   #endif
   

libpolys/polys/ext_fields/algext.h

@@ -69 +69 @@
 number   naLcm(number a, number b, const coeffs cf);
 number   naSize(number a, const coeffs cf);
-void     naDelete(number *a, const coeffs cf);
+void     naDelete(number * a, const coeffs cf);
 void     naCoeffWrite(const coeffs cf);
 number   naIntDiv(number a, number b, const coeffs cf);

libpolys/polys/monomials/p_polys.cc

@@ -1453 +1453 @@
   int divisorLE = p_GetExp(divisor, 1, r);
 //p_Write(p, r); p_Write(divisor, r);
-  while ((p != NULL) && (p_Deg(p, r) > p_Deg(divisor, r)))
+  while ((p != NULL) && (p_Deg(p, r) >= p_Deg(divisor, r)))
   {
     /* determine t = LT(p) / LT(divisor) */
@@ -1462 +1462 @@
     p_SetExp(t, 1, e, r);
     p_Setm(t, r);
-//printf("t\n"); p_Write(t, r);
     if (needResult) result = p_Add_q(result, p_Copy(t, r), r);
-//t = p_Mult_q(t, p_Copy(divisor, r), r); p_Write(t, r);
-//t = p_Neg(t, r); p_Write(t, r);
-//printf("p\n"); p_Write(p, r);
-//t = p_Add_q(p, t, r); p_Write(t, r);
-//p = t;
     p = p_Add_q(p, p_Neg(p_Mult_q(t, p_Copy(divisor, r), r), r), r);
-//printf("EXECUTION STOPPED\n"); while (1) { };
   }
   n_Delete(&divisorLC, r->cf);
+//printf("p_PolyDiv result:\n");
+//p_Write(result, r);
   return result;
+}
+
+/* returns NULL if p == NULL, otherwise makes p monic by dividing
+   by its leading coefficient (only done if this is not already 1);
+   this assumes that we are over a ground field so that division
+   is well-defined;
+   modifies p */
+void p_Monic(poly &p, ring r)
+{
+  if (p == NULL) return;
+  poly pp = p;
+  number lc = p_GetCoeff(p, r);
+  if (n_IsOne(lc, r->cf)) return;
+  number n = n_Init(1, r->cf);
+  p_SetCoeff(p, n, r);
+  p = pIter(p);
+  while (p != NULL)
+  {
+    number c = p_GetCoeff(p, r);
+    number n = n_Div(c, lc, r->cf);
+    n_Delete(&c, r->cf);
+    p_SetCoeff(p, n, r);
+    p = pIter(p);
+  }
+  n_Delete(&lc, r->cf);
+  p = pp;
 }
 
@@ -1481 +1502 @@
 poly p_GcdHelper(poly &p, poly &q, ring r)
 {
-  if (q == NULL) return p;
+  if (q == NULL)
+  {
+    /* We have to make p monic before we return it, so that if the
+       gcd is a unit in the ground field, we will actually return 1. */
+    p_Monic(p, r);
+    return p;
+  }
   else
   {
@@ -1511 +1538 @@
                     ring r)
 {
-//printf("p_ExtGcdHelper:\n");
   if (q == NULL)
   {
-    qFactor = NULL; pFactor = p_ISet(1, r); return p;
+    qFactor = NULL;
+    pFactor = p_ISet(1, r);
+    number n = p_GetCoeff(pFactor, r);
+    p_SetCoeff(pFactor, n_Invers(p_GetCoeff(p, r), r->cf), r);
+    n_Delete(&n, r->cf);
+    p_Monic(p, r);
+    return p;
   }
   else
   {
+//printf("p_ExtGcdHelper1:\n");
 //p_Write(p, r); p_Write(q, r);
     poly pDivQ = p_PolyDiv(p, q, TRUE, r);
-    poly ppFactor; poly qqFactor;
-    poly theGcd = p_ExtGcdHelper(q, ppFactor, p, qqFactor, r);
-    pFactor = qqFactor;
-    qFactor = p_Add_q(ppFactor,
-                      p_Neg(p_Mult_q(pDivQ, p_Copy(qqFactor, r), r), r),
+    poly ppFactor = NULL; poly qqFactor = NULL;
+//printf("p_ExtGcdHelper2:\n");
+//p_Write(p, r); p_Write(ppFactor, r);
+//p_Write(q, r); p_Write(qqFactor, r);
+    poly theGcd = p_ExtGcdHelper(q, qqFactor, p, ppFactor, r);
+//printf("p_ExtGcdHelper3:\n");
+//p_Write(q, r); p_Write(qqFactor, r);
+//p_Write(p, r); p_Write(ppFactor, r);
+//p_Write(theGcd, r);
+    pFactor = ppFactor;
+    qFactor = p_Add_q(qqFactor,
+                      p_Neg(p_Mult_q(pDivQ, p_Copy(ppFactor, r), r), r),
                       r);
+//printf("p_ExtGcdHelper4:\n");
+//p_Write(pFactor, r); p_Write(qFactor, r);
     return theGcd;
   }
@@ -1540 +1582 @@
 poly p_ExtGcd(poly p, poly &pFactor, poly q, poly &qFactor, ring r)
 {
-//printf("p_ExtGcd:\n");
+//printf("p_ExtGcd1:\n");
 //p_Write(p, r); p_Write(q, r);
   assume((p != NULL) || (q != NULL));  
-  poly a = p; poly b = q;
-  poly aFactor = pFactor; poly bFactor = qFactor; 
+  poly a = p; poly b = q; BOOLEAN aCorrespondsToP = TRUE;
   if (p_Deg(a, r) < p_Deg(b, r))
-  { a = q; b = p; aFactor = qFactor, bFactor = pFactor; }
+  { a = q; b = p; aCorrespondsToP = FALSE; }
   a = p_Copy(a, r); b = p_Copy(b, r);
-  return p_ExtGcdHelper(a, aFactor, b, bFactor, r);
+  poly aFactor = NULL; poly bFactor = NULL;
+  poly theGcd = p_ExtGcdHelper(a, aFactor, b, bFactor, r);
+  if (aCorrespondsToP) { pFactor = aFactor; qFactor = bFactor; }
+  else                 { pFactor = bFactor; qFactor = aFactor; }
+  return theGcd;
 }
 

libpolys/polys/monomials/p_polys.h

@@ -1839 +1839 @@
 poly      p_PolyDiv(poly &p, poly divisor, BOOLEAN needResult, ring r);
 
+/* returns NULL if p == NULL, otherwise makes p monic by dividing
+   by its leading coefficient (only done if this is not already 1);
+   this assumes that we are over a ground field so that division
+   is well-defined;
+   modifies p */
+void      p_Monic(poly &p, ring r);
+
 /* assumes that p and q are univariate polynomials in r,
    mentioning the same variable;