Changeset a355723 in git for libpolys

Timestamp:

Oct 22, 2012, 2:45:46 PM (12 years ago)

Author:

Hans Schoenemann <hannes@…>

Branches:

(u'spielwiese', '6e5adcba05493683b94648c659a729c189812c77')

Children:

314f0a2563593c2ed06d4bf520c64eabc6b2c701becbea965e6c5de8e8ab195c7f480cabc295ac0c

Parents:

b79ed5b79606bbec4a3bc2f9ad41e76c5d40137e

Message:

fix: ntInvers/resultant

Location:

libpolys/polys

Files:

• clapsing.cc (modified) (2 diffs)
• ext_fields/transext.cc (modified) (53 diffs)
• monomials/p_polys.cc (modified) (1 diff)

libpolys/polys/clapsing.cc

@@ -177 +177 @@
       CanonicalForm F( convSingTrPFactoryP( f,r ) ), G( convSingTrPFactoryP( g,r ) );
       res= convFactoryPSingTrP( resultant( F, G, X ),r );
-      if ((nf!=NULL)&&(!n_IsOne(nf,r->cf))&&(!n_IsZero(nf,r->cf)))
+      if ((nf!=NULL)&&(!n_IsOne(nf,r->cf)))
       {
         number n=n_Invers(nf,r->cf);
@@ -188 +188 @@
       }
       n_Delete(&nf,r->cf);
-      if ((ng!=NULL)&&(!n_IsOne(ng,r->cf))&&(!n_IsZero(ng,r->cf)))
+      if ((ng!=NULL)&&(!n_IsOne(ng,r->cf)))
       {
         number n=n_Invers(ng,r->cf);

libpolys/polys/ext_fields/transext.cc

@@ -161 +161 @@
 
   const poly den = DEN(t);
-  
+
   if (den != NULL) // !DENIS1(f)
   {
     assume( p_Test(den, ntRing) );
-     
+
     if(p_IsConstant(den, ntRing) && (n_IsOne(pGetCoeff(den), ntRing->cf)))
     {
@@ -171 +171 @@
       return FALSE;
     }
-       
+
     if( !n_GreaterZero(pGetCoeff(den), ntRing->cf) )
     {
@@ -177 +177 @@
       return FALSE;
     }
-    
+
     // test that den is over integers!?
-     
+
   } else
   {  // num != NULL // den == NULL
-   
+
 //    if( COM(t) != 0 )
 //    {
-//      Print("?//NULL with non-zero complexity: %d in %s:%d\n", COM(t), f, l);    
+//      Print("?//NULL with non-zero complexity: %d in %s:%d\n", COM(t), f, l);
 //      return FALSE;
 //    }
@@ -293 +293 @@
   ntTest(a);
   definiteGcdCancellation(a, cf, FALSE);
-  
+
   if (IS0(a)) return NULL;
 
@@ -307 +307 @@
   {
     assume( DEN (f) == NULL );
-    
+
     number g;
     // TODO/NOTE: the following should not be necessary (due to
-    // Hannes!) as NUM (f) should be over Z!!!    
+    // Hannes!) as NUM (f) should be over Z!!!
     CPolyCoeffsEnumerator itr(NUM(f));
 
@@ -324 +324 @@
 
     // g should be a positive integer now!
-    assume( n_GreaterZero(g, ntRing->cf) );  
-    
+    assume( n_GreaterZero(g, ntRing->cf) );
+
     if( !n_IsOne(g, ntRing->cf) )
     {
       DEN (f) = p_NSet(g, ntRing); // update COM(f)???
       COM (f) ++;
-      assume( DEN (f) != NULL );            
+      assume( DEN (f) != NULL );
     }
     else
@@ -339 +339 @@
 
   // Call ntNormalize instead of above?!?
-  
+
   NUM (result) = p_Copy (NUM (f), ntRing); // ???
   DEN (result) = NULL;
   COM (result) = 0;
-  
+
   ntTest((number)result);
   return (number)result;
@@ -361 +361 @@
 
   const BOOLEAN denis1 = DENIS1 (f);
-  
+
   if( IS0(f) || (denis1 && getCoeffType (ntCoeffs) != n_Q) ) // */1 or 0
   {
@@ -388 +388 @@
 //    NUM (result) = p_One(ntRing); // NOTE: just in order to be sure...
   }
-  
+
   // */1 / Q
   assume( getCoeffType (ntCoeffs) == n_Q );
-  assume( DEN (f) == NULL );      
-    
-  number g;    
+  assume( DEN (f) == NULL );
+
+  number g;
 //    poly num= p_Copy (NUM (f), ntRing); // ???
 
@@ -400 +400 @@
   // Hannes!) as NUM (f) should be over Z!!!
   CPolyCoeffsEnumerator itr(NUM(f));
-   
-  n_ClearDenominators(itr, g, ntRing->cf); // may return -1 :((( 
+
+  n_ClearDenominators(itr, g, ntRing->cf); // may return -1 :(((
 //    nlClearDenominators(itr, g, ntRing->cf);
 
-    
+
   if( !n_GreaterZero(g, ntRing->cf) )
   {
@@ -414 +414 @@
 
   // g should be a positive integer now!
-  assume( n_GreaterZero(g, ntRing->cf) );  
+  assume( n_GreaterZero(g, ntRing->cf) );
 
   if( !n_IsOne(g, ntRing->cf) )
@@ -420 +420 @@
     assume( n_GreaterZero(g, ntRing->cf) );
     assume( !n_IsOne(g, ntRing->cf) );
-    
+
     DEN (f) = p_NSet(g, ntRing); // update COM(f)???
     assume( DEN (f) != NULL );
     COM (f) ++;
-      
+
     NUM (result)= p_Copy (DEN (f), ntRing);
   }
@@ -431 +431 @@
     NUM (result)= p_NSet(g, ntRing); // p_Copy (DEN (f), ntRing);
 //  n_Delete(&g, ntRing->cf);
-  }      
-  
+  }
+
 //    if (!p_IsConstant (num, ntRing) && pNext(num) != NULL)
 //    else
@@ -444 +444 @@
 BOOLEAN ntIsOne(number a, const coeffs cf)
 {
-  ntTest(a); // !!!  
+  ntTest(a); // !!!
   definiteGcdCancellation(a, cf, FALSE);
   fraction f = (fraction)a;
@@ -498 +498 @@
     n_Delete(&n, C);
     return NULL;
-  }    
+  }
 
   fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
-   
+
   number den = n_GetDenom(n, C);
-   
+
   assume( n_GreaterZero(den, C) );
-   
+
   if( n_IsOne(den, C) )
   {
      NUM(result) = p_NSet(n, A);
-     DEN(result) = NULL;        
-     n_Delete(&den, C); 
+     DEN(result) = NULL;
+     n_Delete(&den, C);
   } else
   {
-     DEN(result) = p_NSet(den, A);     
-     NUM(result) = p_NSet(n_GetNumerator(n, C), A); 
-     n_Delete(&n, C);     
+     DEN(result) = p_NSet(den, A);
+     NUM(result) = p_NSet(n_GetNumerator(n, C), A);
+     n_Delete(&n, C);
   }
 
@@ -525 +525 @@
 }
 
-   
+
 number ntInit(long i, const coeffs cf)
 {
@@ -545 +545 @@
 {
   if (p == 0) return NULL;
-   
+
     number g;
     // TODO/NOTE: the following should not be necessary (due to
-    // Hannes!) as NUM (f) should be over Z!!!    
+    // Hannes!) as NUM (f) should be over Z!!!
     CPolyCoeffsEnumerator itr(p);
 
@@ -561 +561 @@
 
     // g should be a positive integer now!
-    assume( n_GreaterZero(g, ntRing->cf) );  
+    assume( n_GreaterZero(g, ntRing->cf) );
 
     fraction f = (fraction)omAlloc0Bin(fractionObjectBin);
-     
+
     if( !n_IsOne(g, ntRing->cf) )
     {
-      DEN (f) = p_NSet(g, ntRing); 
+      DEN (f) = p_NSet(g, ntRing);
 //      COM (f) ++; // update COM(f)???
-      assume( DEN (f) != NULL );            
+      assume( DEN (f) != NULL );
     }
     else
     {
-      DEN(f) = NULL;      
+      DEN(f) = NULL;
       n_Delete(&g, ntRing->cf);
-    }    
-   
+    }
+
     NUM(f) = p;
     COM(f) = 0;
@@ -591 +591 @@
   fraction f = (fraction)a;
   if (!DENIS1(f)) return 0;
-  
+
   const poly aAsPoly = NUM(f);
 
   if(aAsPoly == NULL)
     return 0;
-  
+
   if (!p_IsConstant(aAsPoly, ntRing))
     return 0;
 
   assume( aAsPoly != NULL );
-  
+
   return n_Int(p_GetCoeff(aAsPoly, ntRing), ntCoeffs);
 }
@@ -720 +720 @@
 
 //  ntTest((number)result);
-   
+
   return (number)result;
 }
@@ -763 +763 @@
   ntTest(a); // !!!?
   ntTest(b); // !!!?
-  
+
   if (IS0(a) || IS0(b)) return NULL;
 
@@ -780 +780 @@
   const poly db = DEN(fb);
 
-  
+
   if (db == NULL)
   {
     // b = ? // NULL
-    
+
     if(da == NULL)
     { // both fa && fb are ?? // NULL!
@@ -791 +791 @@
       COM(result) = 0;
     }
-    else 
+    else
     {
       assume (da != NULL && db == NULL);
@@ -817 +817 @@
 
 //  ntTest((number)result);
-  
+
   return (number)result;
 }
@@ -1046 +1046 @@
     p_Delete(&DEN(f), ntRing); DEN(f) = NULL;
   }
-   
+
   if( DEN(f) != NULL )
     if( !n_GreaterZero(pGetCoeff(DEN(f)), ntCoeffs) )
@@ -1052 +1052 @@
       NUM(f) = p_Neg(NUM(f), ntRing);
       DEN(f) = p_Neg(DEN(f), ntRing);
-    }  
-   
+    }
+
   ntTest((number)f); // TODO!
 }
@@ -1076 +1076 @@
   } else
   {
-    if (COM(f) > BOUND_COMPLEXITY) 
+    if (COM(f) > BOUND_COMPLEXITY)
       definiteGcdCancellation(a, cf, TRUE);
 
@@ -1085 +1085 @@
       NUM(f) = p_Neg(NUM(f), ntRing);
       DEN(f) = p_Neg(DEN(f), ntRing);
-    }  
-  }
-  
-  
+    }
+  }
+
+
   ntTest(a); // !!!!????
 }
@@ -1164 +1164 @@
       NUM(f) = p_Neg(NUM(f), ntRing);
       DEN(f) = p_Neg(DEN(f), ntRing);
-    }   
+    }
 #endif /* HAVE_FACTORY */
-   
+
   ntTest(a); // !!!!
 }
@@ -1192 +1192 @@
       p_String0Long(DEN(f), ntRing, ntRing);
       if (!omitBrackets) StringAppendS(")");
-    }    
+    }
   }
   ntTest(a); // !!!!
@@ -1259 +1259 @@
   }
 
-  return FALSE;  
+  return FALSE;
 }
 
@@ -1287 +1287 @@
 
     ntTest((number)result); // !!!!
-     
+
     return (number)result;
   }
-   
+
 
   /* return pa*pb/gcd */
@@ -1299 +1299 @@
     ntTest((number)result); // !!!!
     return (number)result;
-  
+
 #else
   Print("// factory needed: transext.cc:ntLcm\n");
@@ -1378 +1378 @@
   ntTest(a);
   if (IS0(a))
-  {    
+  {
     WerrorS(nDivBy0);
     return NULL;
-  }  
+  }
   fraction f = (fraction)a;
   assume( f != NULL );
@@ -1389 +1389 @@
   assume( NUM(f) != NULL );
   const poly den = DEN(f);
-  
+
   if (den == NULL)
     NUM(result) = p_One(ntRing);
@@ -1397 +1397 @@
   if( !NUMIS1(f) )
   {
-    DEN(result) = p_Copy(NUM(f), ntRing);
+    poly num_f=NUM(f);
+    BOOLEAN neg= !n_GreaterZero(pGetCoeff(num_f),ntRing->cf);
+    if (neg)
+    {
+      num_f=p_Neg(p_Copy(num_f, ntRing), ntRing);
+      NUM(result)=p_Neg(NUM(result), ntRing);
+    }
+    else
+    {
+      num_f=p_Copy(num_f, ntRing);
+    }
+    DEN(result) = num_f;
     COM(result) = COM(f);
+    if (neg)
+    {
+      if (p_IsOne(num_f, ntRing))
+      {
+        DEN(result)=NULL;
+        COM(result) = 0;
+        p_Delete(&num_f,ntRing);
+      }
+    }
   }
   else
@@ -1439 +1459 @@
 {
 //  if (n_IsZero(a, cf)) return NULL;
-   
+
   ntTest(a);
 
   if (IS0(a)) return NULL;
-   
+
   const ring rSrc = cf->extRing;
   const ring rDst = dst->extRing;
-  
+
   if( rSrc == rDst )
     return ntCopy(a, dst); // USUALLY WRONG!
-  
+
   fraction f = (fraction)a;
   poly g = prCopyR(NUM(f), rSrc, rDst);
-   
+
   poly h = NULL;
-  
+
   if (!DENIS1(f))
      h = prCopyR(DEN(f), rSrc, rDst);
-   
+
   fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
-   
+
   NUM(result) = g;
   DEN(result) = h;
   COM(result) = COM(f);
   assume(n_Test((number)result, dst));
-  return (number)result;  
+  return (number)result;
 }
 
@@ -1471 +1491 @@
   assume( n_Test(a, cf) );
   if (n_IsZero(a, cf)) return NULL;
-   
-  fraction f = (fraction)omAlloc0Bin(fractionObjectBin); 
+
+  fraction f = (fraction)omAlloc0Bin(fractionObjectBin);
   // DEN(f) = NULL; COM(f) = 0;
   NUM(f) = prCopyR((poly)a, cf->extRing, dst->extRing);
@@ -1492 +1512 @@
     return NULL;
   }
-  
+
   poly g = p_NSet(q, dst->extRing);
 
@@ -1576 +1596 @@
   if (nCoeff_is_Q(bSrc) && nCoeff_is_Q(bDst))
   {
-    if (rVar(src->extRing) > rVar(dst->extRing)) 
+    if (rVar(src->extRing) > rVar(dst->extRing))
        return NULL;
-     
+
     for (int i = 0; i < rVar(src->extRing); i++)
-      if (strcmp(rRingVar(i, src->extRing), rRingVar(i, dst->extRing)) != 0) 
-         return NULL;
-     
+      if (strcmp(rRingVar(i, src->extRing), rRingVar(i, dst->extRing)) != 0)
+         return NULL;
+
     if (src->type==n_transExt)
        return ntCopyMap;          /// Q(T')   --> Q(T)
@@ -1591 +1611 @@
   if (nCoeff_is_Zp(bSrc) && nCoeff_is_Zp(bDst))
   {
-    if (rVar(src->extRing) > rVar(dst->extRing)) 
+    if (rVar(src->extRing) > rVar(dst->extRing))
        return NULL;
-     
+
     for (int i = 0; i < rVar(src->extRing); i++)
-      if (strcmp(rRingVar(i, src->extRing), rRingVar(i, dst->extRing)) != 0) 
-         return NULL;
-     
+      if (strcmp(rRingVar(i, src->extRing), rRingVar(i, dst->extRing)) != 0)
+         return NULL;
+
     if (src->type==n_transExt)
        return ntCopyMap;         /// Z/p(T') --> Z/p(T)
@@ -1716 +1736 @@
   const ring   R = cf->extRing;
   assume(R != NULL);
-  const coeffs Q = R->cf; 
-  assume(Q != NULL); 
-  assume(nCoeff_is_Q(Q));  
-
-  
+  const coeffs Q = R->cf;
+  assume(Q != NULL);
+  assume(nCoeff_is_Q(Q));
+
+
   numberCollectionEnumerator.Reset();
 
@@ -1760 +1780 @@
   while( numberCollectionEnumerator.MoveNext() ) ;
 
-   
+
   // part2: all coeffs = all coeffs * cand
   if( cand != NULL )
@@ -1766 +1786 @@
   if( !p_IsConstant(cand, R) )
   {
-    c = ntInit(cand, cf);   
+    c = ntInit(cand, cf);
     numberCollectionEnumerator.Reset();
     while (numberCollectionEnumerator.MoveNext() )
@@ -1778 +1798 @@
   else { p_Delete(&cand, R); cand = NULL; }
   }
-   
+
   // Quick and dirty fix for constant content clearing: consider numerators???
   CRecursivePolyCoeffsEnumerator<NTNumConverter> itr(numberCollectionEnumerator); // recursively treat the NUM(numbers) as polys!
   number cc;
-   
+
 //  nlClearContentNoPositiveLead(itr, cc, Q); // TODO: get rid of (-LC) normalization!?
   nlClearContent(itr, cc, Q);
   number g = ntInit(p_NSet(cc, R), cf);
-   
+
   if( cand != NULL )
   {
     number gg = ntMult(g, c, cf);
-    ntDelete(&g, cf); 
+    ntDelete(&g, cf);
     ntDelete(&c, cf); c = gg;
-  } else 
+  } else
     c = g;
   ntTest(c);
@@ -1817 +1837 @@
 
   const ring R = cf->extRing;
-  assume(R != NULL); 
-
-  const coeffs Q = R->cf; 
-  assume(Q != NULL); 
-//  assume(nCoeff_is_Q(Q));  
+  assume(R != NULL);
+
+  const coeffs Q = R->cf;
+  assume(Q != NULL);
+//  assume(nCoeff_is_Q(Q));
 
   do
   {
     number &n = numberCollectionEnumerator.Current();
-    
+
     ntNormalize(n, cf);
 
@@ -1847 +1867 @@
       poly gcd = singclap_gcd(p_Copy(cand, R), p_Copy(den, R), R); // gcd(cand, den) is monic no mater leading coeffs! :((((
 //      assume( n_IsOne(pGetCoeff(gcd), Q) ); // TODO: this may be wrong...
-      cand = p_Mult_q(cand, p_Copy(den, R), R); // cand *= den      
+      cand = p_Mult_q(cand, p_Copy(den, R), R); // cand *= den
       const poly t = singclap_pdivide( cand, gcd, R ); // cand' * den / gcd(cand', den)
       p_Delete(&cand, R);
@@ -1860 +1880 @@
     c = ntInit(1, cf);
     return;
-  }  
-
-  c = ntInit(cand, cf); 
+  }
+
+  c = ntInit(cand, cf);
 
   numberCollectionEnumerator.Reset();
-   
+
   number d = NULL;
 
@@ -1876 +1896 @@
     ntNormalize(t, cf); // TODO: needed?
     n = t;
-    
+
     fraction f = (fraction)t;
     assume( f != NULL );
@@ -1886 +1906 @@
       assume( p_IsConstant(den, R) );
       assume( pNext(den) == NULL );
-          
+
       if( d == NULL )
-        d = n_Copy(pGetCoeff(den), Q);
+        d = n_Copy(pGetCoeff(den), Q);
       else
       {
-        number g = n_Lcm(d, pGetCoeff(den), Q);
+        number g = n_Lcm(d, pGetCoeff(den), Q);
         n_Delete(&d, Q); d = g;
       }
     }
   }
-   
+
   if( d != NULL )
   {
@@ -1915 +1935 @@
         assume( p_IsConstant(den, R) );
         assume( pNext(den) == NULL );
-         
-        number ddd = n_Div(d, pGetCoeff(den), Q); // but be an integer now!!!
+
+        number ddd = n_Div(d, pGetCoeff(den), Q); // but be an integer now!!!
         NUM(f) = p_Mult_nn(NUM(f), ddd, R);
-        n_Delete(&ddd, Q);
-        
-        p_Delete(&DEN(f), R); 
-        DEN(f) = NULL; // TODO: check if this is needed!?
-      }      
-       
+        n_Delete(&ddd, Q);
+
+        p_Delete(&DEN(f), R);
+        DEN(f) = NULL; // TODO: check if this is needed!?
+      }
+
       assume( DEN(f) == NULL );
     }
-     
+
     NUM(c) = p_Mult_nn(NUM(c), d, R);
     n_Delete(&d, Q);
   }
-   
-   
+
+
   ntTest(c);
 }
@@ -1951 +1971 @@
   ring R = e->r;
   assume(R != NULL);
-  
+
   R->ref ++; // increase the ref.counter for the ground poly. ring!
 
@@ -1959 +1979 @@
   cf->ch = R->cf->ch;
   cf->factoryVarOffset = R->cf->factoryVarOffset + rVar(R);
-  
+
   cf->cfGreaterZero  = ntGreaterZero;
   cf->cfGreater      = ntGreater;

libpolys/polys/monomials/p_polys.cc

@@ -2712 +2712 @@
   }
 
-  assume( n_GreaterZero(pGetCoeff(ph),C) );
   if(!n_GreaterZero(pGetCoeff(ph),C))
   {