Changeset fb4075b in git for polys

Timestamp:

Nov 8, 2010, 5:49:37 PM (13 years ago)

Author:

Hans Schoenemann <hannes@…>

Branches:

(u'spielwiese', 'fe61d9c35bf7c61f2b6cbf1b56e25e2f08d536cc')

Children:

a7ee69c27ecef9708a4ea7677de1715f1ff6fa86

Parents:

f34215cd2606c6fc99bffc25e1a276b416786dca

git-author:

Hans Schoenemann <hannes@mathematik.uni-kl.de>2010-11-08 17:49:37+01:00

git-committer:

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

Message:

moved p_Divide, p_DivideM to p_polys

Location:

polys

Files:

• monomials/p_polys.cc (modified) (1 diff)
• monomials/p_polys.h (modified) (1 diff)
• polys.cc (modified) (2 diffs)
• polys.h (modified) (1 diff)

polys/monomials/p_polys.cc

@@ -1272 +1272 @@
   }
 }
+/*2
+* assumes that the head term of b is a multiple of the head term of a
+* and return the multiplicant *m
+* Frank's observation: If LM(b) = LM(a)*m, then we may actually set
+* negative(!) exponents in the below loop. I suspect that the correct
+* comment should be "assumes that LM(a) = LM(b)*m, for some monomial m..."
+*/
+poly p_Divide(poly a, poly b, const ring r)
+{
+  assume((p_GetComp(a,r)==p_GetComp(b,r)) || (p_GetComp(b,r)==0));
+  int i;
+  poly result = pInit();
+
+  for(i=(int)r->N; i; i--)
+    p_SetExp(result,i, p_GetExp(a,i,r)- p_GetExp(b,i,r),r);
+  p_SetComp(result, p_GetComp(a,r) - p_GetComp(b,r),r);
+  p_Setm(result,r);
+  return result;
+}
+
+/*2
+* divides a by the monomial b, ignores monomials which are not divisible
+* assumes that b is not NULL
+*/
+poly p_DivideM(poly a, poly b, const ring r)
+{
+  if (a==NULL) return NULL;
+  poly result=a;
+  poly prev=NULL;
+  int i;
+#ifdef HAVE_RINGS
+  number inv=pGetCoeff(b);
+#else
+  number inv=n_Invers(pGetCoeff(b),r->cf);
+#endif
+
+  while (a!=NULL)
+  {
+    if (p_DivisibleBy(b,a,r))
+    {
+      for(i=(int)r->N; i; i--)
+         p_SubExp(a,i, p_GetExp(b,i,r),r);
+      p_SubComp(a, p_GetComp(b,r),r);
+      p_Setm(a,r);
+      prev=a;
+      pIter(a);
+    }
+    else
+    {
+      if (prev==NULL)
+      {
+        p_DeleteLm(&result,r);
+        a=result;
+      }
+      else
+      {
+        p_DeleteLm(&pNext(prev),r);
+        a=pNext(prev);
+      }
+    }
+  }
+#ifdef HAVE_RINGS
+  if (n_IsUnit(inv,r->cf))
+  {
+    inv = n_Invers(inv,r->cf);
+    p_Mult_nn(result,inv,r);
+    n_Delete(&inv, r->cf);
+  }
+  else
+  {
+    p_Div_nn(result,inv,r);
+  }
+#else
+  p_Mult_nn(result,inv,r);
+  n_Delete(&inv, r->cf);
+#endif
+  p_Delete(&b, r);
+  return result;
+}
 
 /***************************************************************

polys/monomials/p_polys.h

@@ -1755 +1755 @@
 poly      p_mInit(const char *s, BOOLEAN &ok, const ring r); /* monom s -> poly, interpreter */
 const char *    p_Read(const char *s, poly &p,const ring r); /* monom -> poly */
+poly      p_Divide(poly a, poly b, const ring r);
+poly      p_DivideM(poly a, poly b, const ring r);
 #endif // P_POLYS_H
 

polys/polys.cc

@@ -91 +91 @@
 }
 
-/*2
-* assumes that the head term of b is a multiple of the head term of a
-* and return the multiplicant *m
-* Frank's observation: If LM(b) = LM(a)*m, then we may actually set
-* negative(!) exponents in the below loop. I suspect that the correct
-* comment should be "assumes that LM(a) = LM(b)*m, for some monomial m..."
-*/
-#define pDivide(a,b) p_Divide(a,b,currRing)
-poly p_Divide(poly a, poly b, cont ring r)
-{
-  assume((pGetComp(a)==pGetComp(b)) || (pGetComp(b)==0));
-  int i;
-  poly result = pInit();
-
-  for(i=(int)r->N; i; i--)
-    p_SetExp(result,i, p_GetExp(a,i,r)- p_GetExp(b,i,r),r);
-  p_SetComp(result, p_GetComp(a,r) - p_GetComp(b,r),r);
-  p_Setm(result,r);
-  return result;
-}
-
 #ifdef HAVE_RINGS   //TODO Oliver
 #define pDiv_nn(p, n)              p_Div_nn(p, n, currRing)
@@ -146 +125 @@
 }
 #endif
-
-/*2
-* divides a by the monomial b, ignores monomials which are not divisible
-* assumes that b is not NULL, destroys b
-*/
-poly p_DivideM(poly a, poly b, const ring r)
-{
-  if (a==NULL) { pDelete(&b); return NULL; }
-  poly result=a;
-  poly prev=NULL;
-  int i;
-#ifdef HAVE_RINGS
-  number inv=pGetCoeff(b);
-#else
-  number inv=nInvers(pGetCoeff(b));
-#endif
-
-  while (a!=NULL)
-  {
-    if (p_DivisibleBy(b,a,r))
-    {
-      assume((p_GetComp(a,r)==p_GetComp(b,r)) || (p_GetComp(b,r)==0));
-      for(i=(int)r->N; i; i--)
-         p_SubExp(a,i, p_GetExp(b,i,r),r);
-      p_SubComp(a, p_GetComp(b,r),r);
-      p_Setm(a,r);
-      prev=a;
-      pIter(a);
-    }
-    else
-    {
-      if (prev==NULL)
-      {
-        p_DeleteLm(&result,r);
-        a=result;
-      }
-      else
-      {
-        p_DeleteLm(&pNext(prev),r);
-        a=pNext(prev);
-      }
-    }
-  }
-#ifdef HAVE_RINGS
-  if (n_IsUnit(inv,r->cf))
-  {
-    inv = n_Invers(inv,r->cf);
-    p_Mult_nn(result,inv,r);
-    n_Delete(&inv, r->cf);
-  }
-  else
-  {
-    p_Div_nn(result,inv,r);
-  }
-#else
-  p_Mult_nn(result,inv,r);
-  n_Delete(&inv, r->cf);
-#endif
-  p_Delete(&b, r);
-  return result;
-}
 
 /*2

polys/polys.h

@@ -289 +289 @@
 // ----------------- define to enable new p_procs -----*/
 
-poly      p_Divide(poly a, poly b, const ring r);
-poly      p_DivideM(poly a, poly b, const ring r);
+#define pDivide(a,b) p_Divide(a,b,currRing)
+
 void      pLcm(poly a, poly b, poly m);
 poly      pDiff(poly a, int k);