Changeset 2391e4 in git

Timestamp:

May 3, 2013, 8:56:15 PM (10 years ago)

Author:

Oleksandr Motsak <motsak@…>

Branches:

(u'jengelh-datetime', 'ceac47cbc86fe4a15902392bdbb9bd2ae0ea02c6')(u'spielwiese', 'a800fe4b3e9d37a38c5a10cc0ae9dfa0c15a4ee6')

Children:

a8f51f2cbfd5e70ab9d6ada456b3c4f409b753bd

Parents:

552d26572655f537d219086f12ab104e239f7db2ed95e6148ac6926aac6e83629cf4545c5c95a7a7

Message:

Merge pull request #326 from mmklee/factory_absfact_sw

Factory absfact sw

Files:

• Singular/LIB/absfact.lib (modified) (2 diffs)
• Singular/LIB/primdec.lib (modified) (2 diffs)
• Singular/extra.cc (modified) (1 diff)
• Singular/subexpr.cc (modified) (1 diff)
• Tst/Short.lst (modified) (1 diff)
• Tst/Short/absbifact.res.gz.uu (added)
• Tst/Short/absbifact.stat (added)
• Tst/Short/absbifact.tst (added)
• Tst/Short/ok_s.lst (modified) (1 diff)
• factory/Makefile.am (modified) (3 diffs)
• factory/canonicalform.h (modified) (2 diffs)
• factory/cfNewtonPolygon.cc (modified) (5 diffs)
• factory/cfNewtonPolygon.h (modified) (1 diff)
• factory/facAbsFact.cc (added)
• factory/facAbsFact.h (added)
• factory/facBivar.h (modified) (1 diff)
• factory/facFqBivarUtil.cc (modified) (2 diffs)
• factory/factory.template (modified) (2 diffs)
• factory/ftmpl_inst.cc (modified) (3 diffs)
• factory/include/factory/Makefile.am (modified) (1 diff)
• factory/include/factory/templates/ftmpl_afactor.h (added)
• factory/templates/ftmpl_afactor.cc (added)
• libpolys/polys/clapconv.cc (modified) (1 diff)
• libpolys/polys/clapsing.cc (modified) (1 diff)
• libpolys/polys/clapsing.h (modified) (1 diff)

Singular/LIB/absfact.lib

@@ -26 +26 @@
 PROCEDURES:
   absFactorize();        absolute factorization of poly
+  absBiFactorize();      absolute factorization of bivariate poly
 ";
 
@@ -609 +610 @@
 }
 
+////////////////////////////////////////////////////////////////////////////////
+proc absBiFactorize(poly p, list #)
+"USAGE:  absBiFactorize(p [,s]);   p poly, s string
+ASSUME: coefficient field is the field of rational numbers or a
+        transcendental extension thereof
+RETURN: ring @code{R} which is obtained from the current basering
+        by adding a new parameter (if a string @code{s} is given as a
+        second input, the new parameter gets this string as name). The ring
+        @code{R} comes with a list @code{absolute_factors} with the
+        following entries:
+@format
+    absolute_factors[1]: ideal   (the absolute factors)
+    absolute_factors[2]: intvec  (the multiplicities)
+    absolute_factors[3]: ideal   (the minimal polynomials)
+    absolute_factors[4]: int     (total number of nontriv. absolute factors)
+@end format
+        The entry @code{absolute_factors[1][1]} is a constant, the
+        entry @code{absolute_factors[3][1]} is the parameter added to the
+        current ring.@*
+        Each of the remaining entries @code{absolute_factors[1][j]} stands for
+        a class of conjugated absolute factors. The corresponding entry
+        @code{absolute_factors[3][j]} is the minimal polynomial of the
+        field extension over which the factor is minimally defined (its degree
+        is the number of conjugates in the class). If the entry
+        @code{absolute_factors[3][j]} coincides with @code{absolute_factors[3][1]},
+        no field extension was necessary for the @code{j}th (class of)
+        absolute factor(s).
+NOTE:   The input polynomial p is assumed to be bivariate or to contain one para-
+        meter and one variable. All factors are presented denominator- and 
+        content-free. The constant factor (first entry) is chosen such that the
+        product of all (!) the (denominator- and content-free) absolute factors
+        of @code{p} equals @code{p / absolute_factors[1][1]}.
+SEE ALSO: factorize, absPrimdecGTZ, absFactorize
+EXAMPLE: example absBiFactorize; shows an example
+"
+{
+  int dblevel = printlevel - voice + 2;
+  dbprint(dblevel,"Entering absfact.lib::absBiFactorize with ",p);
+  def MP = basering;
+  int i;
+  if (char(MP) != 0)
+  {
+    ERROR("// absfact.lib::absBiFactorize is only implemented for "+
+          "characteristic 0");
+  }
+  if(minpoly!=0)
+  {
+    ERROR("// absfact.lib::absBiFactorize is not implemented for algebraic "
+          +"extensions");
+  }
+
+  int n = nvars(MP);
+  int pa=npars(MP);
+  list lMP= ringlist(MP);
+  intvec vv,vk;
+  for(i=1;i<=n;i++){vv[i]=1;}
+  vk=vv,1;
+
+  //if the basering has parameters, add the parameters to the variables
+  //takes care about coefficients and possible denominators
+  if(pa>0)
+  {
+    poly qh=cleardenom(p);
+    if (p==0)
+    {
+      number cok=0;
+    }
+    else
+    {
+      number cok=leadcoef(p)/leadcoef(qh);
+    }
+    p=qh;
+    string sp;
+    for(i=1;i<=npars(basering);i++)
+    {
+      sp=string(par(i));
+      sp=sp[2..size(sp)-1];
+      lMP[2][n+i]=sp;
+      vv=vv,1;
+    }
+    lMP[1]=0;
+    n=n+npars(MP);
+  }
+
+  // MPz is obtained by adding the new variable @z to MP
+  // ordering is wp(1...1)
+  // All the above subroutines work in MPz
+  string newvar;
+  if(size(#)>0)
+  {
+    if(typeof(#[1])=="string")
+    {
+      newvar=#[1];
+    }
+    else
+    {
+      newvar = "a";
+    }
+  }
+  else
+  {
+    newvar = "a";
+  }
+  if (newvar=="a")
+  {
+    if(belongTo(newvar, lMP[2])||defined(a)){newvar = "b";}
+    if(belongTo(newvar, lMP[2])||defined(b)){newvar = "c";}
+    if(belongTo(newvar, lMP[2])||defined(c)){newvar = "@c";}
+    while(belongTo(newvar, lMP[2]))
+    {
+       newvar = "@" + newvar;
+    }
+  }
+
+  // create ring with extra parameter `newvar` for output:
+  setring(MP);
+  list Lout=ringlist(MP);
+  if(!pa)
+  {
+    list Lpar=list(char(MP),list(newvar),list(list("lp",intvec(1))),ideal(0));
+  }
+  else
+  {
+    list Lpar=Lout[1];
+    Lpar[2][size(Lpar[2])+1]=newvar;
+    vv=Lpar[3][1][2];
+    vv=vv,1;
+    Lpar[3][1][2]=vv;
+  }
+  Lout[1]=Lpar;
+  def MPo=ring(Lout);
+  setring(MPo);
+
+  poly p=imap(MP,p);
+
+  if (size (variables (p)) > 2)
+  {
+    ERROR("// absfact.lib::absBiFactorize is not implemented for polynomials "
+          +"with more than 2 variables");
+  }
+
+  // special treatment in the homogeneous case, dropping one variable
+  // by substituting the first variable by 1
+  int ho=homog(p);
+  if(ho)
+  {
+    int dh=deg(p);
+    p=subst(p,var(1),1);
+    int di=deg(p);
+  }
+
+  list tmpf=system ("absBiFact", p);
+
+  // the homogeneous case, homogenizing the result
+  // the new variable has to have degree 0
+  // need to change the ring
+  if(ho)
+  {
+    list ll=ringlist(MPo);
+    vv[size(vv)]=0;
+    ll[3][1][2]=vv;
+    def MPhelp=ring(ll);
+    setring(MPhelp);
+    list tmpf=imap(MPo,tmpf);
+    for(i=2;i<=size(tmpf[1]);i++)
+    {
+      tmpf[1][i]=homog(tmpf[1][i],var(1));
+    }
+    if(dh>di)
+    {
+      tmpf[1][size(tmpf[1])+1]=var(1);
+      tmpf[2][size(tmpf[2])+1]=dh-di;
+      tmpf[3][size(tmpf[3])+1]=par(npars(MPo));
+      tmpf[4]= tmpf[4]+dh-di;
+    }
+    setring(MPo);
+    tmpf=imap(MPhelp,tmpf);
+  }
+
+  if (pa)
+  {
+    number cok=imap(MP,cok);
+    tmpf[1][1]=cok*tmpf[1][1];
+  }
+
+  // if we want the output as string
+  if(size(#)>0)
+  {
+    if(typeof(#[1])=="int")
+    {
+      if(#[1]==77)
+      {  // undocumented feature for Gerhard's absPrimdecGTZ
+        if (size(tmpf[1])<2){ list abs_fac = list(var(n+1),poly(1)); }
+        else
+        {
+          list abs_fac= tmpf[3][2];
+          abs_fac= abs_fac, tmpf[1][2];
+          for (i= 3; i <= size(tmpf[1]); i++)
+          {
+            abs_fac=abs_fac,tmpf[3][i];
+            abs_fac=abs_fac,tmpf[1][i];
+          }
+        }
+        abs_fac=abs_fac,newvar;
+        return(string(abs_fac));
+      }
+    }
+  }
+
+  list absolute_factors= tmpf;
+  export absolute_factors;
+  setring(MP);
+
+  dbprint( printlevel-voice+3,"
+// 'absBiFactorize' created a ring, in which a list absolute_factors (the
+// absolute factors) is stored.
+// To access the list of absolute factors, type (if the name S was assigned
+// to the return value):
+        setring(S); absolute_factors; ");
+  return(MPo);
+}
+example
+{
+  "EXAMPLE:"; echo = 2;
+  ring R = (0), (x,y), lp;
+  poly p = (-7*x^2 + 2*x*y^2 + 6*x + y^4 + 14*y^2 + 47)*(5x2+y2)^3*(x-y)^4;
+  def S = absBiFactorize(p) ;
+  setring(S);
+  absolute_factors;
+}
+
+
 /*
   ring r=0,(x,t),dp;

Singular/LIB/primdec.lib

@@ -3769 +3769 @@
     I=interred(I+F);option(set,save);return(I);
    }
+   I=simplify(I,1);
 
    for(i=1;i<=n;i++)             //consider all polynomials over
@@ -4444 +4445 @@
         if(check==0) // H==I => U is a primary decomposition
         {
-          option(set,save);
+          option(set,save);
           return(U);
         }

Singular/extra.cc

@@ -3683 +3683 @@
       }
       else
+  /*================= absBiFact ======================*/
+      if (strcmp(sys_cmd, "absBiFact") == 0)
+      {
+        if (h!=NULL)
+        {
+          res->rtyp=LIST_CMD;
+          if (h->Typ()==POLY_CMD)
+          {
+            intvec *v=NULL;
+            ideal mipos= NULL;
+            int n= 0;
+            ideal f=singclap_absBiFactorize((poly)(h->Data()), mipos, &v, n, currRing);
+            if (f==NULL) return TRUE;
+            ivTest(v);
+            lists l=(lists)omAllocBin(slists_bin);
+            l->Init(4);
+            l->m[0].rtyp=IDEAL_CMD;
+            l->m[0].data=(void *)f;
+            l->m[1].rtyp=INTVEC_CMD;
+            l->m[1].data=(void *)v;
+            l->m[2].rtyp=IDEAL_CMD;
+            l->m[2].data=(void*) mipos;
+            l->m[3].rtyp=INT_CMD;
+            l->m[3].data=(void*) (long) n;
+            res->data=(void *)l;
+            return FALSE;
+          }
+          else return TRUE;
+        }
+        else return TRUE;
+      }
+      else
       #endif
   /*================= probIrredTest ======================*/

Singular/subexpr.cc

@@ -1066 +1066 @@
         return  ((idhdl)h->data.ustring)->data.ustring;
       }
-      case VECHO:      return (void *)si_echo;
-      case VPRINTLEVEL:return (void *)printlevel;
-      case VCOLMAX:    return (void *)colmax;
-      case VTIMER:     return (void *)getTimer();
-      case VRTIMER:    return (void *)getRTimer();
-      case VOICE:      return (void *)(myynest+1);
-      case VMAXDEG:    return (void *)Kstd1_deg;
-      case VMAXMULT:   return (void *)Kstd1_mu;
-      case TRACE:      return (void *)traceit;
-      case VSHORTOUT:  return (void *)(currRing != NULL ? currRing->ShortOut : 0);
+      case VECHO:      return (void *)(long)si_echo;
+      case VPRINTLEVEL:return (void *)(long)printlevel;
+      case VCOLMAX:    return (void *)(long)colmax;
+      case VTIMER:     return (void *)(long)getTimer();
+      case VRTIMER:    return (void *)(long)getRTimer();
+      case VOICE:      return (void *)(long)(myynest+1);
+      case VMAXDEG:    return (void *)(long)Kstd1_deg;
+      case VMAXMULT:   return (void *)(long)Kstd1_mu;
+      case TRACE:      return (void *)(long)traceit;
+      case VSHORTOUT:  return (void *)(long)(currRing != NULL ? currRing->ShortOut : 0);
       case VMINPOLY:
         if ( (currRing != NULL)  && nCoeff_is_algExt(currRing->cf) && !nCoeff_is_GF(currRing->cf))

Tst/Short.lst

@@ +1 @@
+Short/absbifact.tst
 Short/alexpoly.tst
 Short/algebralib.tst

Tst/Short/ok_s.lst

@@ -1 +1 @@
 ;; Test which should pass shifted exponents
+absbifact
 abusalem
 alexpoly

factory/Makefile.am

@@ -62 +62 @@
                 DegreePattern.cc \
                 ExtensionInfo.cc \
+                facAbsFact.cc \
                 facAlgExt.cc \
                 facBivar.cc \
@@ -144 +145 @@
                 DegreePattern.h \
                 ExtensionInfo.h \
+                facAbsFact.h \
                 facAlgExt.h \
                 facBivar.h \
@@ -248 +250 @@
 
 templatesrc =   templates/ftmpl_array.cc \
+                templates/ftmpl_afactor.cc \
                 templates/ftmpl_factor.cc \
                 templates/ftmpl_functions.h \

factory/canonicalform.h

@@ -23 +23 @@
 #include <factory/templates/ftmpl_list.h>
 #include <factory/templates/ftmpl_array.h>
+#include <factory/templates/ftmpl_afactor.h>
 #include <factory/templates/ftmpl_factor.h>
 #include <factory/templates/ftmpl_matrix.h>
@@ -360 +361 @@
 
 //{{{ type definitions
+typedef AFactor<CanonicalForm> CFAFactor;
+typedef List <CFAFactor> CFAFList;
+typedef ListIterator<CFAFactor> CFAFListIterator;
 typedef Factor<CanonicalForm> CFFactor;
 typedef List<CFFactor> CFFList;

factory/cfNewtonPolygon.cc

@@ -13 +13 @@
 
 #include "config.h"
+
+#include "cf_assert.h"
+
 #include <stdlib.h>
 
@@ -18 +21 @@
 #include "cf_iter.h"
 #include "cf_algorithm.h"
+#include "cf_primes.h"
+#include "cf_reval.h"
+#include "cf_factory.h"
+#include "gfops.h"
 #include "cfNewtonPolygon.h"
 #include "templates/ftmpl_functions.h"
@@ -857 +864 @@
 bool irreducibilityTest (const CanonicalForm& F)
 {
+  ASSERT (getNumVars (F) == 2, "expected bivariate polynomial");
+  ASSERT (getCharacteristic() == 0, "expected polynomial over integers or rationals");
+
   int sizeOfNewtonPolygon;
   int ** newtonPolyg= newtonPolygon (F, sizeOfNewtonPolygon);
@@ -872 +882 @@
         if (isRat)
           Off (SW_RATIONAL);
-        CanonicalForm tmp= gcd (newtonPolyg[0][0],newtonPolyg[0][1]);
+        CanonicalForm tmp= gcd (newtonPolyg[0][0],newtonPolyg[0][1]); // maybe it's better to use plain intgcd
         tmp= gcd (tmp, newtonPolyg[1][0]);
         tmp= gcd (tmp, newtonPolyg[1][1]);
@@ -891 +901 @@
   return false;
 }
+
+bool absIrredTest (const CanonicalForm& F)
+{
+  ASSERT (getNumVars (F) == 2, "expected bivariate polynomial");
+  ASSERT (factorize (F).length() <= 2, " expected irreducible polynomial");
+
+  int sizeOfNewtonPolygon;
+  int ** newtonPolyg= newtonPolygon (F, sizeOfNewtonPolygon);
+  bool isRat= isOn (SW_RATIONAL);
+  if (isRat)
+    Off (SW_RATIONAL);
+  int p=getCharacteristic();
+  int d=1;
+  char bufGFName='Z';
+  bool GF= (CFFactory::gettype()==GaloisFieldDomain);
+  if (GF)
+  {
+    d= getGFDegree();
+    bufGFName=gf_name;
+  }
+
+  setCharacteristic(0);
+
+  CanonicalForm g= gcd (newtonPolyg[0][0], newtonPolyg[0][1]); //maybe it's better to use plain intgcd
+
+  int i= 1;
+  while (!g.isOne() && i < sizeOfNewtonPolygon)
+  {
+    g= gcd (g, newtonPolyg[i][0]);
+    g= gcd (g, newtonPolyg[i][1]);
+    i++;
+  }
+
+  bool result= g.isOne();
+
+  if (GF)
+    setCharacteristic (p, d, bufGFName);
+  else
+    setCharacteristic(p);
+
+  if (isRat)
+    On (SW_RATIONAL);
+
+  for (int i= 0; i < sizeOfNewtonPolygon; i++)
+    delete [] newtonPolyg[i];
+
+  delete [] newtonPolyg;
+
+  return result;
+}
+
+bool modularIrredTest (const CanonicalForm& F)
+{
+  ASSERT (getNumVars (F) == 2, "expected bivariate polynomial");
+  ASSERT (factorize (F).length() <= 2, " expected irreducible polynomial");
+
+  bool isRat= isOn (SW_RATIONAL);
+  if (isRat)
+    Off (SW_RATIONAL);
+
+  if (isRat)
+    ASSERT (bCommonDen (F).isOne(), "poly over Z expected");
+
+  CanonicalForm Fp, N= maxNorm (F);
+  int tdeg= totaldegree (F);
+
+  int i= 0;
+  //TODO: maybe it's better to choose the characteristic as large as possible
+  //      as factorization over large finite field will be faster
+  //TODO: handle those cases where our factory primes are not enough
+  //TODO: factorize coefficients corresponding to the vertices of the Newton
+  //      polygon and only try the obtained factors
+  if (N < cf_getSmallPrime (cf_getNumSmallPrimes()-1))
+  {
+    while (i < cf_getNumSmallPrimes() && N > cf_getSmallPrime(i))
+    {
+      setCharacteristic (cf_getSmallPrime (i));
+      Fp= F.mapinto();
+      i++;
+      if (totaldegree (Fp) == tdeg)
+      {
+        if (absIrredTest (Fp))
+        {
+          CFFList factors= factorize (Fp);
+          if (factors.length() == 2 && factors.getLast().exp() == 1)
+          {
+            if (isRat)
+              On (SW_RATIONAL);
+            setCharacteristic (0);
+            return true;
+          }
+        }
+      }
+      setCharacteristic (0);
+    }
+  }
+  else
+  {
+    while (i < cf_getNumPrimes() && N > cf_getPrime (i))
+    {
+      setCharacteristic (cf_getPrime (i));
+      Fp= F.mapinto();
+      i++;
+      if (totaldegree (Fp) == tdeg)
+      {
+        if (absIrredTest (Fp))
+        {
+          CFFList factors= factorize (Fp);
+          if (factors.length() == 2 && factors.getLast().exp() == 1)
+          {
+            if (isRat)
+              On (SW_RATIONAL);
+            setCharacteristic (0);
+            return true;
+          }
+        }
+      }
+      setCharacteristic (0);
+    }
+  }
+
+  if (isRat)
+    On (SW_RATIONAL);
+
+  return false;
+}
+
+bool
+modularIrredTestWithShift (const CanonicalForm& F)
+{
+  ASSERT (getNumVars (F) == 2, "expected bivariate polynomial");
+  ASSERT (factorize (F).length() <= 2, " expected irreducible polynomial");
+
+  bool isRat= isOn (SW_RATIONAL);
+  if (isRat)
+    Off (SW_RATIONAL);
+
+  if (isRat)
+    ASSERT (bCommonDen (F).isOne(), "poly over Z expected");
+
+  Variable x= Variable (1);
+  Variable y= Variable (2);
+  CanonicalForm Fp;
+  int tdeg= totaldegree (F);
+
+  REvaluation E;
+
+  setCharacteristic (2);
+  Fp= F.mapinto();
+
+  E= REvaluation (1,2, FFRandom());
+
+  E.nextpoint();
+
+  Fp= Fp (x+E[1], x);
+  Fp= Fp (y+E[2], y);
+
+  if (tdeg == totaldegree (Fp))
+  {
+    if (absIrredTest (Fp))
+    {
+      CFFList factors= factorize (Fp);
+      if (factors.length() == 2 && factors.getLast().exp() == 1)
+      {
+        if (isRat)
+          On (SW_RATIONAL);
+        setCharacteristic (0);
+        return true;
+      }
+    }
+  }
+
+  E.nextpoint();
+
+  Fp= Fp (x+E[1], x);
+  Fp= Fp (y+E[2], y);
+
+  if (tdeg == totaldegree (Fp))
+  {
+    if (absIrredTest (Fp))
+    {
+      CFFList factors= factorize (Fp);
+      if (factors.length() == 2 && factors.getLast().exp() == 1)
+      {
+        if (isRat)
+          On (SW_RATIONAL);
+        setCharacteristic (0);
+        return true;
+      }
+    }
+  }
+
+  int i= 0;
+  while (cf_getSmallPrime (i) < 102)
+  {
+    setCharacteristic (cf_getSmallPrime (i));
+    i++;
+    E= REvaluation (1, 2, FFRandom());
+
+    for (int j= 0; j < 3; j++)
+    {
+      Fp= F.mapinto();
+      E.nextpoint();
+      Fp= Fp (x+E[1], x);
+      Fp= Fp (y+E[2], y);
+
+      if (tdeg == totaldegree (Fp))
+      {
+        if (absIrredTest (Fp))
+        {
+          CFFList factors= factorize (Fp);
+          if (factors.length() == 2 && factors.getLast().exp() == 1)
+          {
+            if (isRat)
+              On (SW_RATIONAL);
+            setCharacteristic (0);
+            return true;
+          }
+        }
+      }
+    }
+  }
+
+  setCharacteristic (0);
+  if (isRat)
+    On (SW_RATIONAL);
+
+  return false;
+}
+
+

factory/cfNewtonPolygon.h

@@ -74 +74 @@
                    );
 
+/// absolute irreducibility test as described in "Modular Las Vegas Algorithms
+/// for Polynomial Absolute Factorization" by C. Bertone, G. Cheze, A. Galligo
+///
+/// @return true if F satisfies condition (C) from the above paper and thus
+/// is absolutely irreducible, false otherwise
+bool
+absIrredTest (const CanonicalForm& F ///< [in] a bivariate polynomial
+                                     ///< irreducible over ground field
+             );
+
+/// modular absolute irreducibility test as described in "Modular Las Vegas
+/// Algorithms for Polynomial Absolute Factorization" by C. Bertone, G. Cheze,
+/// A. Galligo
+///
+/// @return true if F is absolutely irreducible, false otherwise
+bool
+modularIrredTest (const CanonicalForm& F ///< [in] a bivariate polynomial
+                                         ///< irreducible over Z
+                 );
+
+/// modular absolute irreducibility test with shift as described in "Modular Las
+/// Vegas Algorithms for Polynomial Absolute Factorization" by C. Bertone,
+/// G. Cheze, A. Galligo
+///
+/// @return true if F is absolutely irreducible, false otherwise
+bool
+modularIrredTestWithShift (const CanonicalForm& F ///< [in] a bivariate polynomial
+                                                  ///< irreducible over Z
+                          );
+
+
 #ifdef HAVE_NTL
 /// Algorithm 5 as described in Convex-Dense Bivariate Polynomial Factorization

factory/facBivar.h

@@ -51 +51 @@
   CFMap N;
   CanonicalForm F= compress (G, N);
-  CanonicalForm contentX= content (F, 1);
+  CanonicalForm contentX= content (F, 1); //erwarte hier primitiven input: primitiv ÃŒber Z bzw. Z[a]
   CanonicalForm contentY= content (F, 2);
   F /= (contentX*contentY);

factory/facFqBivarUtil.cc

@@ -680 +680 @@
 {
   n= degree (F, 1);
+  int* result= new int [n];
+  int sizeOfNewtonPolygon;
+  int** newtonPolyg= newtonPolygon (F, sizeOfNewtonPolygon);
+
+  isIrreducible= false;
+  if (sizeOfNewtonPolygon == 3)
+  {
+    bool check1=
+        (newtonPolyg[0][0]==0 || newtonPolyg[1][0]==0 || newtonPolyg[2][0]==0);
+    if (check1)
+    {
+      bool check2=
+        (newtonPolyg[0][1]==0 || newtonPolyg[1][1]==0 || newtonPolyg[2][0]==0);
+      if (check2)
+      {
+        int p=getCharacteristic();
+        int d=1;
+        char bufGFName='Z';
+        bool GF= (CFFactory::gettype()==GaloisFieldDomain);
+        if (GF)
+        {
+          d= getGFDegree();
+          bufGFName=gf_name;
+        }
+        setCharacteristic(0);
+        CanonicalForm tmp= gcd (newtonPolyg[0][0],newtonPolyg[0][1]);  // maybe it's better to use plain intgcd
+        tmp= gcd (tmp, newtonPolyg[1][0]);
+        tmp= gcd (tmp, newtonPolyg[1][1]);
+        tmp= gcd (tmp, newtonPolyg[2][0]);
+        tmp= gcd (tmp, newtonPolyg[2][1]);
+        isIrreducible= (tmp==1);
+        if (GF)
+          setCharacteristic (p, d, bufGFName);
+        else
+          setCharacteristic(p);
+      }
+    }
+  }
+
+  int minX, minY, maxX, maxY;
+  minX= newtonPolyg [0] [0];
+  minY= newtonPolyg [0] [1];
+  maxX= minX;
+  maxY= minY;
+  for (int i= 1; i < sizeOfNewtonPolygon; i++)
+  {
+    if (minX > newtonPolyg [i] [0])
+      minX= newtonPolyg [i] [0];
+    if (maxX < newtonPolyg [i] [0])
+      maxX= newtonPolyg [i] [0];
+    if (minY > newtonPolyg [i] [1])
+      minY= newtonPolyg [i] [1];
+    if (maxY < newtonPolyg [i] [1])
+      maxY= newtonPolyg [i] [1];
+  }
+
+  int k= maxX;
+  for (int i= 0; i < n; i++)
+  {
+    if (i + 1 > maxY || i + 1 < minY)
+    {
+      result [i]= 0;
+      continue;
+    }
+    int* point= new int [2];
+    point [0]= k;
+    point [1]= i + 1;
+    while (!isInPolygon (newtonPolyg, sizeOfNewtonPolygon, point) && k > 0)
+    {
+      k--;
+      point [0]= k;
+    }
+    result [i]= k;
+    k= maxX;
+    delete [] point;
+  }
+
+  return result;
+}
+
+int *
+computeBoundsWrtDiffMainvar (const CanonicalForm& F, int& n,
+                             bool& isIrreducible)
+{
+  n= degree (F, 2);
   int* result= new int [n];
   int sizeOfNewtonPolygon;
@@ -736 +821 @@
   }
 
-  int k= maxX;
-  for (int i= 0; i < n; i++)
-  {
-    if (i + 1 > maxY || i + 1 < minY)
-    {
-      result [i]= 0;
-      continue;
-    }
-    int* point= new int [2];
-    point [0]= k;
-    point [1]= i + 1;
-    while (!isInPolygon (newtonPolyg, sizeOfNewtonPolygon, point) && k > 0)
-    {
-      k--;
-      point [0]= k;
-    }
-    result [i]= k;
-    k= maxX;
-    delete [] point;
-  }
-
-  return result;
-}
-
-int *
-computeBoundsWrtDiffMainvar (const CanonicalForm& F, int& n,
-                             bool& isIrreducible)
-{
-  n= degree (F, 2);
-  int* result= new int [n];
-  int sizeOfNewtonPolygon;
-  int** newtonPolyg= newtonPolygon (F, sizeOfNewtonPolygon);
-
-  isIrreducible= false;
-  if (sizeOfNewtonPolygon == 3)
-  {
-    bool check1=
-        (newtonPolyg[0][0]==0 || newtonPolyg[1][0]==0 || newtonPolyg[2][0]==0);
-    if (check1)
-    {
-      bool check2=
-        (newtonPolyg[0][1]==0 || newtonPolyg[1][1]==0 || newtonPolyg[2][0]==0);
-      if (check2)
-      {
-        int p=getCharacteristic();
-        int d=1;
-        char bufGFName='Z';
-        bool GF= (CFFactory::gettype()==GaloisFieldDomain);
-        if (GF)
-        {
-          d= getGFDegree();
-          bufGFName=gf_name;
-        }
-        setCharacteristic(0);
-        CanonicalForm tmp= gcd (newtonPolyg[0][0],newtonPolyg[0][1]);
-        tmp= gcd (tmp, newtonPolyg[1][0]);
-        tmp= gcd (tmp, newtonPolyg[1][1]);
-        tmp= gcd (tmp, newtonPolyg[2][0]);
-        tmp= gcd (tmp, newtonPolyg[2][1]);
-        isIrreducible= (tmp==1);
-        if (GF)
-          setCharacteristic (p, d, bufGFName);
-        else
-          setCharacteristic(p);
-      }
-    }
-  }
-
-  int minX, minY, maxX, maxY;
-  minX= newtonPolyg [0] [0];
-  minY= newtonPolyg [0] [1];
-  maxX= minX;
-  maxY= minY;
-  for (int i= 1; i < sizeOfNewtonPolygon; i++)
-  {
-    if (minX > newtonPolyg [i] [0])
-      minX= newtonPolyg [i] [0];
-    if (maxX < newtonPolyg [i] [0])
-      maxX= newtonPolyg [i] [0];
-    if (minY > newtonPolyg [i] [1])
-      minY= newtonPolyg [i] [1];
-    if (maxY < newtonPolyg [i] [1])
-      maxY= newtonPolyg [i] [1];
-  }
-
   int k= maxY;
   for (int i= 0; i < n; i++)

factory/factory.template

@@ -38 +38 @@
 
 #include <factory/templates/ftmpl_array.h>
+#include <factory/templates/ftmpl_afactor.h>
 #include <factory/templates/ftmpl_factor.h>
 #include <factory/templates/ftmpl_list.h>
@@ -105 +106 @@
 #include "facIrredTest.h"
 
+/*MAKEHEADER PUBLIC ONLY*/
+#include "facAbsFact.h"
+
 #endif /* ! INCL_FACTORY_H */

factory/ftmpl_inst.cc

@@ -23 +23 @@
 
 #include "templates/ftmpl_array.cc"
+#include "templates/ftmpl_afactor.cc"
 #include "templates/ftmpl_factor.cc"
 #include "templates/ftmpl_list.cc"
@@ -34 +35 @@
 template class ListItem<CFFactor>;
 template class ListIterator<CFFactor>;
+template class AFactor<CanonicalForm>;
+template class List<CFAFactor>;
+template class ListItem<CFAFactor>;
+template class ListIterator<CFAFactor>;
 template class List<CanonicalForm>;
 template class ListItem<CanonicalForm>;
@@ -77 +82 @@
 
 template int operator == ( const Factor<CanonicalForm> &, const Factor<CanonicalForm> & );
+template int operator == ( const AFactor<CanonicalForm> &, const AFactor<CanonicalForm> & );
 
 template List<CFFactor> Union ( const List<CFFactor> &, const List<CFFactor> & );
+template List<CFAFactor> Union ( const List<CFAFactor> &, const List<CFAFactor> & );
 
 #if ! defined(WINNT) || defined(__GNUC__)

factory/include/factory/Makefile.am

@@ -2 +2 @@
 
 templateincl =  templates/ftmpl_array.h \
+                templates/ftmpl_afactor.h \
                 templates/ftmpl_factor.h \
                 templates/ftmpl_list.h \

libpolys/polys/clapconv.cc

@@ -319 +319 @@
   {
     n_Normalize(p_GetCoeff(p, r), r->cf);
-    CanonicalForm term=convSingPFactoryP(NUM(p_GetCoeff(p, r)),r->cf->extRing);
-
-    if ((DEN(p_GetCoeff(p,r))!=NULL)
-    && (!errorreported))
-    {
+
+    // test if denominator is constant
+    if (!p_IsConstantPoly(DEN (p_GetCoeff (p,r)),r) && !errorreported)
       WerrorS("conversion error: denominator!= 1");
+
+    CanonicalForm term=convSingPFactoryP(NUM (p_GetCoeff(p, r)),r->cf->extRing);
+
+    // if denominator is not NULL it should be a constant at this point
+    if (DEN (p_GetCoeff(p,r)) != NULL)
+    {
+      CanonicalForm den= convSingPFactoryP(DEN (p_GetCoeff(p, r)),r->cf->extRing);
+      if (rChar (r) == 0)
+        On (SW_RATIONAL);
+      term /= den;
     }
 

libpolys/polys/clapsing.cc

@@ -1025 +1025 @@
   return res;
 }
+#ifdef HAVE_NTL
+ideal singclap_absBiFactorize ( poly f, ideal & mipos, intvec ** exps, int & numFactors, const ring r)
+{
+  p_Test(f, r);
+
+  ideal res=NULL;
+
+  int offs = rPar(r);
+  if (f==NULL)
+  {
+    res= idInit (1, 1);
+    mipos= idInit (1, 1);
+    mipos->m[0]= convFactoryPSingTrP (Variable (offs), r); //overkill
+    (*exps)=new intvec (1);
+    (**exps)[0]= 1;
+    numFactors= 0;
+    return res;
+  }
+  CanonicalForm F( convSingTrPFactoryP( f, r) );
+
+  if (getNumVars (F) > 2)
+  {
+    WerrorS( feNotImplemented );
+    return res;
+  }
+  CFAFList absFactors= absFactorize (F);
+
+  int n= absFactors.length();
+  *exps=new intvec (n);
+
+  res= idInit (n, 1);
+
+  mipos= idInit (n, 1);
+
+  Variable x= Variable (offs);
+  Variable alpha;
+  int i= 0;
+  numFactors= 0;
+  int count;
+  CFAFListIterator iter= absFactors;
+  CanonicalForm lead= iter.getItem().factor();
+  if (iter.getItem().factor().inCoeffDomain())
+  {
+    i++;
+    iter++;
+  }
+  bool isRat= isOn (SW_RATIONAL);
+  if (!isRat)
+    On (SW_RATIONAL);
+  for (; iter.hasItem(); iter++, i++)
+  {
+    (**exps)[i]= iter.getItem().exp();
+    alpha= iter.getItem().minpoly().mvar();
+    if (iter.getItem().minpoly().isOne())
+      lead /= power (bCommonDen (iter.getItem().factor()), iter.getItem().exp());
+    else
+      lead /= power (power (bCommonDen (iter.getItem().factor()), degree (iter.getItem().minpoly())), iter.getItem().exp());
+    res->m[i]= convFactoryPSingTrP (replacevar (iter.getItem().factor()*bCommonDen (iter.getItem().factor()), alpha, x), r);
+    if (iter.getItem().minpoly().isOne())
+    {
+      count= iter.getItem().exp();
+      mipos->m[i]= convFactoryPSingTrP (x,r);
+    }
+    else
+    {
+      count= iter.getItem().exp()*degree (iter.getItem().minpoly());
+      mipos->m[i]= convFactoryPSingTrP (replacevar (iter.getItem().minpoly(), alpha, x), r);
+    }
+    numFactors += count;
+  }
+  if (!isRat)
+    Off (SW_RATIONAL);
+
+  (**exps)[0]= 1;
+  res->m[0]= convFactoryPSingTrP (lead, r);
+  mipos->m[0]= convFactoryPSingTrP (x, r);
+  return res;
+}
+#endif
 ideal singclap_sqrfree ( poly f, intvec ** v , int with_exps, const ring r)
 {

libpolys/polys/clapsing.h

@@ -51 +51 @@
 matrix singntl_LLL(matrix A, const ring r);
 intvec* singntl_LLL(intvec* A, const ring r);
+ideal singclap_absBiFactorize ( poly f, ideal & mipos, intvec ** exps, int & n, const ring r);
 #endif
 #endif