Changeset 17a710 in git


Ignore:
Timestamp:
Apr 19, 2013, 3:34:15 PM (11 years ago)
Author:
Martin Lee <martinlee84@…>
Branches:
(u'spielwiese', 'fe61d9c35bf7c61f2b6cbf1b56e25e2f08d536cc')
Children:
f7ed7dc2ce876bee21478ae1bf88e2734992bcd8
Parents:
efd410fc3eac114fbca114eb66d8070b61317f92
git-author:
Martin Lee <martinlee84@web.de>2013-04-19 15:34:15+02:00
git-committer:
Martin Lee <martinlee84@web.de>2013-05-02 11:42:37+02:00
Message:
chg: debug output and experimental code
Location:
factory
Files:
3 edited

Legend:

Unmodified
Added
Removed

  • factory/cfNewtonPolygon.cc

    refd410 r17a710  
    2323#include "templates/ftmpl_functions.h"
    2424#include "algext.h"
     25
     26void out_cf(const char *s1,const CanonicalForm &f,const char *s2);
    2527
    2628static
    … …  
    921923  CanonicalForm g= gcd (newtonPolyg[0][0], newtonPolyg[0][1]); //maybe it's better to use plain intgcd
    922924
     925
     926  printf ("newtonPolyg[0]= %d\n", newtonPolyg[0][0]);
     927  printf ("newtonPolyg[0]= %d\n", newtonPolyg[0][1]);
     928  out_cf ("g= ", g, "\n");
     929
    923930  int i= 1;
    924931  while (!g.isOne() && i < sizeOfNewtonPolygon)
    … …  
    926933    g= gcd (g, newtonPolyg[i][0]);
    927934    g= gcd (g, newtonPolyg[i][1]);
     935  printf ("newtonPolyg[i]= %d\n", newtonPolyg[i][0]);
     936  printf ("newtonPolyg[i]= %d\n", newtonPolyg[i][1]);
     937    out_cf ("g= ", g, "\n");
    928938    i++;
    929939  }
    … …  
    10291039  ASSERT (factorize (F).length() <= 2, " expected irreducible polynomial");
    10301040
     1041  out_cf ("F= ", F, "\n");
    10311042  bool isRat= isOn (SW_RATIONAL);
    10321043  if (isRat)
    … …  
    10531064  Fp= Fp (y+E[2], y);
    10541065
     1066  out_cf ("Fp0= ", Fp, "\n");
     1067
    10551068  if (tdeg == totaldegree (Fp))
    10561069  {
    … …  
    10601073      if (factors.length() == 2 && factors.getLast().exp() == 1)
    10611074      {
     1075        printf ("irred1\n");
    10621076        if (isRat)
    10631077          On (SW_RATIONAL);
    … …  
    10731087  Fp= Fp (y+E[2], y);
    10741088
     1089  out_cf ("Fp1= ", Fp, "\n");
     1090
    10751091  if (tdeg == totaldegree (Fp))
    10761092  {
    … …  
    10801096      if (factors.length() == 2 && factors.getLast().exp() == 1)
    10811097      {
     1098        printf ("irred2\n");
    10821099        if (isRat)
    10831100          On (SW_RATIONAL);
    … …  
    11091126          if (factors.length() == 2 && factors.getLast().exp() == 1)
    11101127          {
     1128            out_cf ("Fp2= ", Fp, "\n");
     1129            printf ("tdeg= %d\n", tdeg);
     1130            printf ("totaldegree (Fp)= %d\n", totaldegree (Fp));
     1131            printf ("degree (Fp, 1)= %d\n", degree (Fp,1));
     1132            printf ("degree (Fp,2)= %d\n", degree (Fp,2));
     1133            printf ("degree (F,1)= %d\n", degree (F,1));
     1134            printf ("degree (F,2)= %d\n", degree (F,2));
     1135            out_cf ("factors.getFirst()= ", factors.getFirst().factor(), "\n");
     1136            printf ("char= %d\n", getCharacteristic ());
     1137            printf ("irred3\n");
    11111138            if (isRat)
    11121139              On (SW_RATIONAL);

  • factory/facAbsFact.cc

    refd410 r17a710  
    99/*****************************************************************************/
    1010
     11#include "timing.h"
     12#include "debug.h"
     13
    1114#include "facAbsFact.h"
     15#include "facBivar.h"
     16#include "facFqBivar.h"
    1217#include "cf_reval.h"
    1318#include "cf_primes.h"
    … …  
    2530#ifdef HAVE_NTL
    2631
     32TIMING_DEFINE_PRINT(fac_Qa_factorize)
     33TIMING_DEFINE_PRINT(fac_evalpoint)
     34
     35void out_cf(const char *s1,const CanonicalForm &f,const char *s2);
     36
    2737//TODO optimize choice of p -> choose p as large as possible (better than small p since factorization mod p does not require field extension, also less lifting)
    28 int choosePoint (const CanonicalForm& F, int tdegF, CFArray& eval)
     38int choosePoint (const CanonicalForm& F, int tdegF, CFArray& eval, bool rec)
    2939{
    3040  REvaluation E1 (1, 1, IntRandom (25));
    3141  REvaluation E2 (2, 2, IntRandom (25));
    32   //E1.nextpoint();
    33   //E2.nextpoint();
    34   CanonicalForm f, Fp;
    35   int i;
     42  if (rec)
     43  {
     44    E1.nextpoint();
     45    E2.nextpoint();
     46  }
     47  CanonicalForm f, f1, f2, Fp;
     48  int i, p;
    3649  eval=CFArray (2);
     50  printf ("tdegF= %d\n", tdegF);
     51  out_cf ("F= ", F, "\n");
     52  printf ("getCharacteristic()= %d\n", getCharacteristic());
    3753  while (1)
    3854  {
    39     f= E1(F);
    40     if (!f.isZero() && factorize (f).length() == 2)
     55    f1= E1(F);
     56    if (!f1.isZero() && factorize (f1).length() == 2)
    4157    {
    4258      Off (SW_RATIONAL);
    43       f= E2(f);
    44       if (!f.isZero() && f > cf_getSmallPrime (cf_getNumSmallPrimes()))
     59      f= E2(f1);
     60      f2= E2 (F);
     61      out_cf ("f= ", f, "\n");
     62      printf ("isOn (SW_RATIONAL)= %d\n", isOn (SW_RATIONAL));
     63      printf ("cf_getSmallPrime (cf_getNumSmallPrimes())= %d\n", cf_getSmallPrime (cf_getNumSmallPrimes()-1));
     64      if ((!f.isZero()) && (abs(f) > cf_getSmallPrime (cf_getNumSmallPrimes()-1)))
    4565      {
    46         for (i= cf_getNumPrimes()-1; i > 0; i--)
     66        printf ("hier0\n");
     67        for (i= cf_getNumPrimes()-1; i >= 0; i--)
    4768        {
    4869          if (f % CanonicalForm (cf_getPrime (i)) == 0)
    4970          {
    50             Fp= mod (F,cf_getPrime(i));
    51             if (totaldegree (Fp) == tdegF)
     71            p= cf_getPrime(i);
     72            Fp= mod (F,p);
     73            out_cf ("Fp0= ", Fp, "\n");
     74            if (totaldegree (Fp) == tdegF && degree (mod (f2,p), 1) == degree (F,1) && degree (mod (f1, p),2) == degree (F,2))
    5275            {
    5376              eval[0]= E1[1];
    5477              eval[1]= E2[2];
    55               return cf_getPrime(i);
     78              return p;
    5679            }
    5780          }
    … …  
    6083      else if (!f.isZero())
    6184      {
    62         for (i= cf_getNumSmallPrimes()-1; i > 0; i--)
     85        printf ("hier\n");
     86        for (i= cf_getNumSmallPrimes()-1; i >= 0; i--)
    6387        {
    6488          if (f % CanonicalForm (cf_getSmallPrime (i)) == 0)
    6589          {
    66             Fp= mod (F,cf_getSmallPrime(i));
    67             if (totaldegree (Fp) == tdegF)
     90            p= cf_getSmallPrime (i);
     91            Fp= mod (F,p);
     92            out_cf ("Fp1= ", Fp, "\n");
     93            if (totaldegree (Fp) == tdegF && degree (mod (f2, p),1) == degree (F,1) && degree (mod (f1,p),2) == degree (F,2))
    6894            {
    6995              eval[0]= E1[1];
    7096              eval[1]= E2[2];
    71               return cf_getSmallPrime(i);
     97              return p;
    7298            }
    7399          }
    … …  
    76102      E2.nextpoint();
    77103      On (SW_RATIONAL);
     104      out_cf ("E2= ", E2[2], "\n");
    78105    }
    79106    E1.nextpoint();
     107    out_cf ("E1= ", E1[1], "\n");
    80108  }
    81109  return 0;
    82110}
    83111
     112//TODO sowohl bzgl. x als auch y auswerten und minpoly berechnen
    84113CFAFList absFactorizeMain (const CanonicalForm& G)
    85114{
    86115  //F is assumed to be bivariate, irreducible over Q, primitive wrt x and y, compressed
    87116
     117  out_cf ("F= ", G, "\n");
    88118  CanonicalForm F= bCommonDen (G)*G;
    89119  Off (SW_RATIONAL);
    90120  F /= icontent (F);
     121  out_cf ("F after icontent= ", F, "\n");
    91122  On (SW_RATIONAL);
    92123  CFArray eval;
    … …  
    94125  CanonicalForm Fp, smallestFactor;
    95126  int p;
     127  CFFList factors;
     128  Variable alpha;
     129  bool rec= false;
     130  Variable x= Variable (1);
     131  Variable y= Variable (2);
     132differentevalpoint:
    96133  while (1)
    97134  {
    98     p= choosePoint (F, tdegF, eval);
     135    TIMING_START (fac_evalpoint);
     136    p= choosePoint (F, tdegF, eval, rec);
     137    TIMING_END_AND_PRINT (fac_evalpoint, "time to find eval point: ");
    99138
    100139    setCharacteristic (p);
    101140    Fp=F.mapinto();
    102     CFFList factors= factorize (Fp);
    103 
    104     factors.removeFirst();
     141    factors= factorize (Fp);
     142
     143    for (CFFListIterator iter= factors; iter.hasItem(); iter++)
     144    {
     145      out_cf ("factors= ", iter.getItem().factor(), "\n");
     146      printf ("exp= %d\n", iter.getItem().exp());
     147    }
     148    printf ("p= %d\n", p);
     149    if (factors.getFirst().factor().inCoeffDomain())
     150      factors.removeFirst();
     151    printf ("factors.length()= %d\n", factors.length());
     152    printf ("factors.getFirst().exp()= %d\n", factors.getFirst().exp());
     153    if (factors.length() == 1 && factors.getFirst().exp() == 1)
     154    {
     155      if (absIrredTest (Fp)) //TODO absIrredTest mit shift, modular absIrredTest
     156      {
     157        printf ("irred after test\n");
     158        printf ("absIrred\n");
     159        setCharacteristic(0);
     160        alpha= rootOf (x);
     161    out_cf ("G= ", G, "\n");
     162    out_cf ("getMipo (alpha)= ", getMipo (alpha), "\n");
     163
     164        return CFAFList (CFAFactor (G, getMipo (alpha), 1));
     165      }
     166      else
     167      {
     168        setCharacteristic (0);
     169        if (modularIrredTestWithShift (F))
     170        {
     171          printf ("irred after modular test\n");
     172          alpha= rootOf (x);
     173          return CFAFList (CFAFactor (G, getMipo (alpha), 1));
     174        }
     175        rec= true;
     176        continue;
     177      }
     178    }
    105179    CFFListIterator iter= factors;
    106180    smallestFactor= iter.getItem().factor();
     181    out_cf ("smallestFactor before= ", smallestFactor, "\n");
     182    while (smallestFactor.isUnivariate() && iter.hasItem())
     183    {
     184      iter++;
     185      if (!iter.hasItem())
     186        break;
     187      out_cf ("factors= ", iter.getItem().factor(), "\n");
     188      printf ("exp= %d\n", iter.getItem().exp());
     189      smallestFactor= iter.getItem().factor();
     190    }
     191    //TODO univariate Faktoren rausschmeißen!
    107192    minTdeg= totaldegree (smallestFactor);
    108     iter++;
     193    if (iter.hasItem())
     194      iter++;
    109195    for (; iter.hasItem(); iter++)
    110196    {
    111       if (totaldegree (iter.getItem().factor()) < minTdeg)
     197      out_cf ("factors= ", iter.getItem().factor(), "\n");
     198      printf ("exp= %d\n", iter.getItem().exp());
     199      if (!iter.getItem().factor().isUnivariate() && (totaldegree (iter.getItem().factor()) < minTdeg))
    112200      {
    113201        smallestFactor= iter.getItem().factor();
    … …  
    117205    if (tdegF % minTdeg == 0)
    118206      break;
    119     //TODO else
     207    setCharacteristic(0);
     208    rec=true;
    120209  }
    121210  CanonicalForm Gp= Fp/smallestFactor;
    122   Gp= Gp (eval[0].mapinto(), 1);
    123   Fp= Fp (eval[0].mapinto(), 1);
    124   CanonicalForm smallestFactorEval= smallestFactor (eval[0].mapinto(),1);
     211  out_cf ("Gp before= ", Gp, "\n");
     212  out_cf ("smallestFactor= ", smallestFactor, "\n");
     213  printf ("degree (Gp,1)= %d\n", degree (Gp, 1));
     214  printf ("degree smallestFactor= %d\n", degree (smallestFactor, 1));
     215  printf ("degree (Fp,1)= %d\n", degree (Fp,1));
     216  printf ("degree (F,1)= %d\n", degree (F,1));
     217  out_cf ("eval[1]= ", eval[1], "\n");
     218  out_cf ("eval[0]= ", eval[0], "\n");
     219  //printf ("Gp*smallestFactor==Fp ? %d\n", Gp*smallestFactor == Fp);
     220  Gp= Gp /Lc (Gp);
     221
     222  CanonicalForm Gpy= Gp (eval[0].mapinto(), 1);
     223  CanonicalForm smallestFactorEvaly= smallestFactor (eval[0].mapinto(),1);
     224  CanonicalForm Gpx= Gp (eval[1].mapinto(), 2);
     225  out_cf ("Gp eval= ", Gp, "\n");
     226  CanonicalForm smallestFactorEvalx= smallestFactor (eval[1].mapinto(),2);
     227
     228  out_cf ("smallestFactorEvalx= ", smallestFactorEvalx, "\n");
     229  out_cf ("gcd (Gpx, smallestFactorEvalx)= ", gcd (Gpx, smallestFactorEvalx), "\n");
     230  bool xValid= !(Gpx.inCoeffDomain() || smallestFactorEvalx.inCoeffDomain() || !gcd (Gpx, smallestFactorEvalx).inCoeffDomain());
     231  bool yValid= !(Gpy.inCoeffDomain() || smallestFactorEvaly.inCoeffDomain() || !gcd (Gpy, smallestFactorEvaly).inCoeffDomain());
     232  if (!xValid && !yValid)
     233  {
     234    rec= true;
     235    setCharacteristic (0);
     236    printf ("goto1\n");
     237    goto differentevalpoint;
     238  }
     239
    125240  setCharacteristic (0);
    126   CanonicalForm F1= F(eval[0],1);
    127   int s= tdegF/minTdeg + 1;
    128   CanonicalForm bound= power (maxNorm (F1), 2*(s-1));
    129   bound *= power (CanonicalForm (s),s-1);
    130   bound *= power (CanonicalForm (2), ((s-1)*(s-1))/2); //possible int overflow
    131 
    132   CanonicalForm B = p;
    133   long k = 1L;
    134   while ( B < bound ) {
    135     B *= p;
    136     k++;
    137   }
    138 
    139   //TODO take floor (log2(k))
    140   k= k+1;
    141   fmpz_poly_t FLINTF1;
    142   convertFacCF2Fmpz_poly_t (FLINTF1, F1);
    143   setCharacteristic (p);
    144   nmod_poly_t FLINTFp, FLINTGp;
    145   convertFacCF2nmod_poly_t (FLINTFp, smallestFactorEval/lc (smallestFactorEval));
    146   convertFacCF2nmod_poly_t (FLINTGp, Gp/lc (Gp));
    147   nmod_poly_factor_t nmodFactors;
    148   nmod_poly_factor_init (nmodFactors);
    149   nmod_poly_factor_insert (nmodFactors, FLINTFp, 1L);
    150   nmod_poly_factor_insert (nmodFactors, FLINTGp, 1L);
    151 
    152   long * link= new long [2];
    153   fmpz_poly_t *v= new fmpz_poly_t[2];
    154   fmpz_poly_t *w= new fmpz_poly_t[2];
    155   fmpz_poly_init(v[0]);
    156   fmpz_poly_init(v[1]);
    157   fmpz_poly_init(w[0]);
    158   fmpz_poly_init(w[1]);
    159 
    160   fmpz_poly_factor_t liftedFactors;
    161   fmpz_poly_factor_init (liftedFactors);
    162   _fmpz_poly_hensel_start_lift(liftedFactors, link, v, w, FLINTF1, nmodFactors, k); //lift factors up to p^k
    163 
    164   nmod_poly_factor_clear (nmodFactors);
    165   nmod_poly_clear (FLINTFp);
    166   nmod_poly_clear (FLINTGp);
    167 
    168   setCharacteristic(0);
    169   modpk pk= modpk (p,k);
    170   CanonicalForm liftedSmallestFactor= convertFmpz_poly_t2FacCF ((fmpz_poly_t &)liftedFactors->p[0],Variable (2));
    171 
    172   CanonicalForm otherFactor= convertFmpz_poly_t2FacCF ((fmpz_poly_t &)liftedFactors->p[1],Variable (2));
    173   CanonicalForm test= pk (otherFactor*liftedSmallestFactor);
    174 
    175   Off (SW_SYMMETRIC_FF);
    176   liftedSmallestFactor= pk (liftedSmallestFactor);
    177   liftedSmallestFactor= liftedSmallestFactor (eval[1],2);
    178   On (SW_SYMMETRIC_FF);
    179   CFMatrix M= CFMatrix (s, s);
    180   M(s,s)= power (CanonicalForm (p), k);
    181   for (int i= 1; i < s; i++)
    182   {
    183     M (i,i)= 1;
    184     M (i+1,i)= -liftedSmallestFactor;
    185   }
    186 
    187   mat_ZZ NTLM= *convertFacCFMatrix2NTLmat_ZZ (M);
    188 
    189   ZZ det;
    190 
    191   transpose (NTLM, NTLM);
    192   long r=LLL (det, NTLM, 1L, 1L); //use floating point LLL ?
    193   transpose (NTLM, NTLM);
    194   M= *convertNTLmat_ZZ2FacCFMatrix (NTLM);
    195 
    196   CanonicalForm mipo= 0;
    197   Variable x= Variable (1);
    198   for (int i= 1; i <= s; i++)
    199     mipo += M (i,1)*power (x,s-i);
    200 
    201   CFFList mipoFactors= factorize (mipo);
    202   mipoFactors.removeFirst();
    203   //TODO check if mipoFactors has length 1 and multiplicity 1 - if not choose a new point!
     241
     242  CanonicalForm mipo;
     243
     244  int loop, i;
     245  if (xValid && yValid)
     246  {
     247    loop= 3;
     248    i=1;
     249  }
     250  else if (xValid)
     251  {
     252    loop= 3;
     253    i=2;
     254  }
     255  else
     256  {
     257    loop= 2;
     258    i=1;
     259  }
     260
     261  CFArray mipos= CFArray (loop-i);
     262  printf ("loop= %d\n", loop);
     263  printf ("xValid= %d\n", xValid);
     264  printf ("yValid= %d\n", yValid);
     265
     266  for (; i < loop; i++)
     267  {
     268    CanonicalForm Fi= F(eval[i-1],i);
     269    //CanonicalForm Fx= F(eval[0],1);
     270    //CanonicalForm Fy= F(eval[1],2);
     271
     272    int s= tdegF/minTdeg + 1;
     273    CanonicalForm bound= power (maxNorm (Fi), 2*(s-1));
     274    bound *= power (CanonicalForm (s),s-1);
     275    bound *= power (CanonicalForm (2), ((s-1)*(s-1))/2); //possible int overflow
     276
     277    CanonicalForm B = p;
     278    long k = 1L;
     279    while ( B < bound ) {
     280      B *= p;
     281      k++;
     282    }
     283
     284    //TODO take floor (log2(k))
     285    k= k+1;
     286    fmpz_poly_t FLINTFi;
     287    out_cf ("Fi= ", Fi, "\n");
     288    convertFacCF2Fmpz_poly_t (FLINTFi, Fi);
     289    setCharacteristic (p);
     290    printf ("p= %d\n", p);
     291    nmod_poly_t FLINTFpi, FLINTGpi;
     292    if (i == 2)
     293    {
     294      convertFacCF2nmod_poly_t (FLINTFpi, smallestFactorEvalx/lc (smallestFactorEvalx));
     295      convertFacCF2nmod_poly_t (FLINTGpi, Gpx/lc (Gpx));
     296    }
     297    else
     298    {
     299      convertFacCF2nmod_poly_t (FLINTFpi, smallestFactorEvaly/lc (smallestFactorEvaly));
     300      convertFacCF2nmod_poly_t (FLINTGpi, Gpy/lc (Gpy));
     301    }
     302    nmod_poly_factor_t nmodFactors;
     303    nmod_poly_factor_init (nmodFactors);
     304    nmod_poly_factor_insert (nmodFactors, FLINTFpi, 1L);
     305    nmod_poly_factor_insert (nmodFactors, FLINTGpi, 1L);
     306
     307    //out_cf ("Gpx= ", Gpx, "\n");
     308    //out_cf ("smallestFactorEvalx= ", smallestFactorEvalx, "\n");
     309    long * link= new long [2];
     310    fmpz_poly_t *v= new fmpz_poly_t[2];
     311    fmpz_poly_t *w= new fmpz_poly_t[2];
     312    fmpz_poly_init(v[0]);
     313    fmpz_poly_init(v[1]);
     314    fmpz_poly_init(w[0]);
     315    fmpz_poly_init(w[1]);
     316
     317    printf ("k= %ld\n", k);
     318    fmpz_poly_factor_t liftedFactors;
     319    fmpz_poly_factor_init (liftedFactors);
     320    _fmpz_poly_hensel_start_lift(liftedFactors, link, v, w, FLINTFi, nmodFactors, k); //lift factors up to p^k
     321
     322    nmod_poly_factor_clear (nmodFactors);
     323    nmod_poly_clear (FLINTFpi);
     324    nmod_poly_clear (FLINTGpi);
     325
     326    setCharacteristic(0);
     327    modpk pk= modpk (p,k);
     328    CanonicalForm liftedSmallestFactor= convertFmpz_poly_t2FacCF ((fmpz_poly_t &)liftedFactors->p[0],Variable (1));
     329
     330    CanonicalForm otherFactor= convertFmpz_poly_t2FacCF ((fmpz_poly_t &)liftedFactors->p[1],Variable (1));
     331    CanonicalForm test= pk (otherFactor*liftedSmallestFactor);
     332
     333    Off (SW_SYMMETRIC_FF);
     334    liftedSmallestFactor= pk (liftedSmallestFactor);
     335    if (i == 2)
     336      liftedSmallestFactor= liftedSmallestFactor (eval[0],1);
     337    else
     338      liftedSmallestFactor= liftedSmallestFactor (eval[1],1);
     339
     340    On (SW_SYMMETRIC_FF);
     341    CFMatrix M= CFMatrix (s, s);
     342    M(s,s)= power (CanonicalForm (p), k);
     343    for (int j= 1; j < s; j++)
     344    {
     345      M (j,j)= 1;
     346      M (j+1,j)= -liftedSmallestFactor;
     347    }
     348
     349    mat_ZZ NTLM= *convertFacCFMatrix2NTLmat_ZZ (M);
     350
     351    ZZ det;
     352
     353    transpose (NTLM, NTLM);
     354    (void) LLL (det, NTLM, 1L, 1L); //use floating point LLL ?
     355    transpose (NTLM, NTLM);
     356    M= *convertNTLmat_ZZ2FacCFMatrix (NTLM);
     357
     358    mipo= 0;
     359    for (int j= 1; j <= s; j++)
     360      mipo += M (j,1)*power (x,s-j);
     361
     362    CFFList mipoFactors= factorize (mipo);
     363    mipoFactors.removeFirst();
     364
     365    fmpz_poly_clear (v[0]);
     366    fmpz_poly_clear (v[1]);
     367    fmpz_poly_clear (w[0]);
     368    fmpz_poly_clear (w[1]);
     369    delete [] v;
     370    delete [] w;
     371    delete [] link;
     372    fmpz_poly_factor_clear (liftedFactors);
     373
     374    if (mipoFactors.length() > 1 ||
     375        (mipoFactors.length() == 1 && mipoFactors.getFirst().exp() > 1))
     376    {
     377      if (i+1 >= loop && ((loop-i == 1) || (loop-i==2 && mipos[0].isZero())))
     378      {
     379        rec=true;
     380        printf ("goto2\n");
     381        goto differentevalpoint;
     382        //TODO check if mipoFactors has length 1 and multiplicity 1 - if not choose a new point!
     383      }
     384    }
     385    else
     386      mipos[loop-i-1]= mipo;
     387  }
     388
    204389  On (SW_RATIONAL);
    205   Variable alpha= rootOf (mipo);
    206   CFFList QaFactors= factorize (F1, alpha);
    207 
     390  if (xValid && yValid && !mipos[0].isZero() && !mipos[1].isZero())
     391  {
     392    if (maxNorm (mipos[0]) < maxNorm (mipos[1]))
     393      alpha= rootOf (mipos[0]);
     394    else
     395      alpha= rootOf (mipos[1]);
     396  }
     397  else
     398    alpha= rootOf (mipo);
     399
     400  for (i= 0; i < mipos.size(); i++)
     401  {
     402    out_cf ("mipos= ", mipos [i], "\n");
     403    out_cf ("maxNorm mipo= ", maxNorm (mipos[i]), "\n");
     404  }
     405
     406  CanonicalForm F1;
     407  CFFList QaF1Factors;
     408  int wrongMipo= 0;
     409  if (xValid && yValid)
     410  {
     411    if (degree (F,1) > minTdeg)
     412      F1= F (eval[1], 2);
     413    else
     414      F1= F (eval[0], 1);
     415  }
     416  else if (xValid)
     417    F1= F (eval[1], 2);
     418  else
     419    F1= F (eval[0], 1);
     420
     421  QaF1Factors= factorize (F1, alpha);
     422  if (QaF1Factors.getFirst().factor().inCoeffDomain())
     423    QaF1Factors.removeFirst();
     424  out_cf ("mipo0= ", getMipo (alpha), "\n");
     425  for (CFFListIterator iter= QaF1Factors; iter.hasItem(); iter++)
     426  {
     427    out_cf ("QaF1Factors= ", iter.getItem().factor(), "\n");
     428    if (degree (iter.getItem().factor()) > minTdeg)
     429      wrongMipo++;
     430  }
     431
     432  if (wrongMipo == QaF1Factors.length())
     433  {
     434    if (xValid && yValid)
     435    {
     436      if (mipo==mipos[0])
     437        alpha= rootOf (mipos[1]);
     438      else
     439        alpha= rootOf (mipos[0]);
     440    }
     441
     442    wrongMipo= 0;
     443    out_cf ("mipo1= ", getMipo (alpha), "\n");
     444    QaF1Factors= factorize (F1, alpha);
     445    if (QaF1Factors.getFirst().factor().inCoeffDomain())
     446      QaF1Factors.removeFirst();
     447    for (CFFListIterator iter= QaF1Factors; iter.hasItem(); iter++)
     448    {
     449      out_cf ("QaF1Factors= ", iter.getItem().factor(), "\n");
     450      if (degree (iter.getItem().factor()) > minTdeg)
     451        wrongMipo++;
     452    }
     453    if (wrongMipo == QaF1Factors.length())
     454    {
     455      rec= true;
     456      printf ("goto30\n");
     457      goto differentevalpoint;
     458    }
     459  }
     460
     461  CanonicalForm A= F;
     462  CanonicalForm Aeval= F1;
     463
     464  out_cf ("F1= ", F1, "\n");
     465  A *= bCommonDen (A);
     466  A= A (y + eval[1], y); //TODO find right evaluation point and swap if necessary
     467
     468  out_cf ("A= ", A, "\n");
     469  out_cf ("A[0]= ", A(0,y), "\n");
     470  int liftBound= degree (A,y) + 1;
     471
     472  modpk b= modpk();
     473
     474  //bool mipoHasDen= false;
     475  CanonicalForm den= 1;
     476
     477  mipo= getMipo (alpha);
     478
     479  CFList uniFactors;
     480  for (CFFListIterator iter=QaF1Factors; iter.hasItem(); iter++)
     481  {
     482    uniFactors.append (iter.getItem().factor());
     483    out_cf ("uniFactors.getLast()= ", uniFactors.getLast(), "\n");
     484  }
     485
     486
     487    A /= Lc (Aeval);
     488    //mipoHasDen= !bCommonDen(mipo).isOne();
     489    //mipo *= bCommonDen (mipo);
     490    ZZX NTLmipo= convertFacCF2NTLZZX (mipo);
     491    ZZX NTLLcf= convertFacCF2NTLZZX (Lc (A*bCommonDen (A)));
     492    ZZ NTLf= resultant (NTLmipo, NTLLcf);
     493    ZZ NTLD= discriminant (NTLmipo);
     494    den= abs (convertZZ2CF (NTLD*NTLf));
     495
     496    // make factors elements of Z(a)[x] disable for modularDiophant
     497    CanonicalForm multiplier= 1;
     498    for (CFListIterator i= uniFactors; i.hasItem(); i++)
     499    {
     500      multiplier *= bCommonDen (i.getItem());
     501      i.getItem()= i.getItem()*bCommonDen(i.getItem());
     502    }
     503    A *= multiplier;
     504    A *= bCommonDen (A);
     505
     506    Off (SW_RATIONAL);
     507    int ii= 0;
     508    CanonicalForm discMipo= convertZZ2CF (NTLD);
     509    findGoodPrime (F*discMipo,ii);
     510    findGoodPrime (Aeval*discMipo,ii);
     511    findGoodPrime (A*discMipo,ii);
     512
     513    int pp=cf_getBigPrime(ii);
     514    b = coeffBound( A, pp, mipo );
     515    modpk bb= coeffBound (Aeval, pp, mipo);
     516    if (bb.getk() > b.getk() ) b=bb;
     517      bb= coeffBound (F, pp, mipo);
     518    if (bb.getk() > b.getk() ) b=bb;
     519
     520  ExtensionInfo dummy= ExtensionInfo (alpha, false);
     521  DegreePattern degs= DegreePattern (uniFactors);
     522
     523  bool earlySuccess= false;
     524  CFList earlyFactors;
     525  TIMING_START (fac_bi_hensel_lift);
     526  uniFactors= henselLiftAndEarly
     527              (A, earlySuccess, earlyFactors, degs, liftBound,
     528               uniFactors, dummy, eval[1], b, den);
     529  TIMING_END_AND_PRINT (fac_bi_hensel_lift,
     530                        "time for bivariate hensel lifting over Q: ");
     531  DEBOUTLN (cerr, "lifted factors= " << uniFactors);
     532
     533  CanonicalForm MODl= power (y, liftBound); //TODO
     534
     535  On (SW_RATIONAL);
     536  A *= bCommonDen (A);
     537  Off (SW_RATIONAL);
     538
     539  printf ("earlyFactors.length()= %d\n", earlyFactors.length());
     540  CFList biFactors;
     541
     542  TIMING_START (fac_bi_factor_recombination);
     543  biFactors= factorRecombination (uniFactors, A, MODl, degs, 1,
     544                                uniFactors.length()/2, b, den);
     545  TIMING_END_AND_PRINT (fac_bi_factor_recombination,
     546                        "time for bivariate factor recombination over Q: ");
     547
     548  On (SW_RATIONAL);
     549
     550  if (earlySuccess)
     551    biFactors= Union (earlyFactors, biFactors);
     552  else if (!earlySuccess && degs.getLength() == 1)
     553    biFactors= earlyFactors;
     554
     555  for (CFListIterator i= biFactors; i.hasItem(); i++)
     556    i.getItem()= i.getItem() (y - eval[1], y); //TODO
     557
     558  bool swap= false;
     559  bool swap2= false;
     560  appendSwapDecompress (biFactors, CFList(), CFList(), swap, swap2, CFMap());
     561  if (isOn (SW_RATIONAL))
     562    normalize (biFactors);
     563
     564  CFAFList result;
     565  bool found= false;
     566
     567  out_cf ("mipo= ", mipo, "\n");
     568  printf ("minTdeg= %d\n", minTdeg);
     569  for (CFListIterator iter= biFactors; iter.hasItem(); iter++)
     570  {
     571    out_cf ("biFactors= ", iter.getItem(), "\n");
     572    printf ("totaldegree ()= %d\n", totaldegree (iter.getItem()));
     573    if (totaldegree (iter.getItem()) == minTdeg)
     574    {
     575      result= CFAFList (CFAFactor (iter.getItem(), getMipo (alpha), 1));
     576      found= true;
     577      break;
     578    }
     579  }
     580
     581  if (found)
     582  {
     583    printf ("thisexitexit\n\n");
     584    return result;
     585  }
     586
     587  /*  A *= bCommonDen (A);
     588  A= A (y + evaluation, y);
     589
     590  int liftBound= degree (A, y) + 1;
     591
     592  modpk b= modpk();
     593  bool mipoHasDen= false;
     594  CanonicalForm den= 1;
     595
     596  if (!extension)
     597  {
     598    Off (SW_RATIONAL);
     599    int i= 0;
     600    findGoodPrime(F,i);
     601    findGoodPrime(Aeval,i);
     602    findGoodPrime(A,i);
     603    if (i >= cf_getNumBigPrimes())
     604      printf ("out of primes\n"); //TODO exit
     605
     606    int p=cf_getBigPrime(i);
     607    b = coeffBound( A, p );
     608    modpk bb= coeffBound (Aeval, p);
     609    if (bb.getk() > b.getk() ) b=bb;
     610      bb= coeffBound (F, p);
     611    if (bb.getk() > b.getk() ) b=bb;
     612  }
     613  else
     614  {
     615    A /= Lc (Aeval);
     616    mipoHasDen= !bCommonDen(mipo).isOne();
     617    mipo *= bCommonDen (mipo);
     618    ZZX NTLmipo= convertFacCF2NTLZZX (mipo);
     619    ZZX NTLLcf= convertFacCF2NTLZZX (Lc (A*bCommonDen (A)));
     620    ZZ NTLf= resultant (NTLmipo, NTLLcf);
     621    ZZ NTLD= discriminant (NTLmipo);
     622    den= abs (convertZZ2CF (NTLD*NTLf));
     623
     624    // make factors elements of Z(a)[x] disable for modularDiophant
     625    CanonicalForm multiplier= 1;
     626    for (CFListIterator i= uniFactors; i.hasItem(); i++)
     627    {
     628      multiplier *= bCommonDen (i.getItem());
     629      i.getItem()= i.getItem()*bCommonDen(i.getItem());
     630    }
     631    A *= multiplier;
     632    A *= bCommonDen (A);
     633
     634    Off (SW_RATIONAL);
     635    int i= 0;
     636    CanonicalForm discMipo= convertZZ2CF (NTLD);
     637    findGoodPrime (F*discMipo,i);
     638    findGoodPrime (Aeval*discMipo,i);
     639    findGoodPrime (A*discMipo,i);
     640
     641    int p=cf_getBigPrime(i);
     642    b = coeffBound( A, p, mipo );
     643    modpk bb= coeffBound (Aeval, p, mipo);
     644    if (bb.getk() > b.getk() ) b=bb;
     645      bb= coeffBound (F, p, mipo);
     646    if (bb.getk() > b.getk() ) b=bb;
     647  }
     648
     649  ExtensionInfo dummy= ExtensionInfo (false);
     650  if (extension)
     651    dummy= ExtensionInfo (v, false);
     652  bool earlySuccess= false;
     653  CFList earlyFactors;
     654  TIMING_START (fac_bi_hensel_lift);
     655  uniFactors= henselLiftAndEarly
     656              (A, earlySuccess, earlyFactors, degs, liftBound,
     657               uniFactors, dummy, evaluation, b, den);
     658  TIMING_END_AND_PRINT (fac_bi_hensel_lift,
     659                        "time for bivariate hensel lifting over Q: ");
     660  DEBOUTLN (cerr, "lifted factors= " << uniFactors);
     661
     662  CanonicalForm MODl= power (y, liftBound);
     663
     664  if (mipoHasDen)
     665  {
     666    Variable vv;
     667    for (CFListIterator iter= uniFactors; iter.hasItem(); iter++)
     668      if (hasFirstAlgVar (iter.getItem(), vv))
     669        break;
     670    for (CFListIterator iter= uniFactors; iter.hasItem(); iter++)
     671      iter.getItem()= replacevar (iter.getItem(), vv, v);
     672  }
     673
     674  On (SW_RATIONAL);
     675  A *= bCommonDen (A);
     676  Off (SW_RATIONAL);
     677
     678  TIMING_START (fac_bi_factor_recombination);
     679  factors= factorRecombination (uniFactors, A, MODl, degs, 1,
     680                                uniFactors.length()/2, b, den);
     681  TIMING_END_AND_PRINT (fac_bi_factor_recombination,
     682                        "time for bivariate factor recombination over Q: ");
     683
     684  On (SW_RATIONAL);
     685
     686  if (earlySuccess)
     687    factors= Union (earlyFactors, factors);
     688  else if (!earlySuccess && degs.getLength() == 1)
     689    factors= earlyFactors;
     690
     691  for (CFListIterator i= factors; i.hasItem(); i++)
     692    i.getItem()= i.getItem() (y - evaluation, y);
     693
     694  appendSwapDecompress (factors, conv (contentAxFactors),
     695                        conv (contentAyFactors), swap, swap2, N);
     696  if (isOn (SW_RATIONAL))
     697    normalize (factors);*/
     698
     699  TIMING_START (fac_Qa_factorize);
     700  CFFList QaFactors= factorize (F, alpha); //TODO lift these factors
     701  TIMING_END_AND_PRINT (fac_Qa_factorize, "time to factorize over Qa: ");
     702
     703  /*mipo= getMipo (alpha);
     704  out_cf ("maxNorm (mipo)= ", maxNorm (mipo), "\n");
    208705  QaFactors.append (CFFactor (mipo, 1)); //last factor is the minimal polynomial that defines the extension
    209   fmpz_poly_clear (v[0]);
    210   fmpz_poly_clear (v[1]);
    211   fmpz_poly_clear (w[0]);
    212   fmpz_poly_clear (w[1]);
    213   delete [] v;
    214   delete [] w;
    215   delete [] link;
    216   fmpz_poly_factor_clear (liftedFactors);
    217   return QaFactors;
     706  if (degree (mipo) < 3)
     707    printf ("scheissescheissescheissescheisse\n");*/
     708  printf ("minTdeg= %d\n", minTdeg);
     709  //CFAFList result;
     710  //bool found= false;
     711  out_cf ("mipo= ", getMipo (alpha), "\n");
     712  for (CFFListIterator iter= QaFactors; iter.hasItem(); iter++)
     713  {
     714    out_cf ("QaFactors= ", iter.getItem().factor(), "\n");
     715    printf ("totaldegree ()= %d\n", totaldegree (iter.getItem().factor()));
     716    if (totaldegree (iter.getItem().factor()) == minTdeg)
     717    {
     718      result= CFAFList (CFAFactor (iter.getItem().factor(), getMipo (alpha), 1));
     719      found= true;
     720      break;
     721    }
     722  }
     723  if (!found && xValid && yValid)
     724  {
     725    if (mipo == mipos [0])
     726      mipo= mipos[1];
     727    else
     728      mipo= mipos[0];
     729    alpha= rootOf (mipo);
     730
     731    QaFactors= factorize (F, alpha);
     732
     733    for (CFFListIterator iter= QaFactors; iter.hasItem(); iter++)
     734    {
     735      out_cf ("QaFactors= ", iter.getItem().factor(), "\n");
     736      printf ("totaldegree ()= %d\n", totaldegree (iter.getItem().factor()));
     737      if (totaldegree (iter.getItem().factor()) == minTdeg)
     738      {
     739        result= CFAFList (CFAFactor (iter.getItem().factor(), getMipo (alpha), 1));
     740        found= true;
     741        break;
     742      }
     743    }
     744    if (!found)
     745    {
     746      rec= true;
     747      printf ("goto31\n");
     748      goto differentevalpoint;
     749    }
     750  }
     751  else if (!found)
     752  {
     753    rec= true;
     754    printf ("goto32\n");
     755    goto differentevalpoint;
     756  }
     757
     758  return result;
    218759}
    219760#endif

  • factory/facBivar.h

    refd410 r17a710  
    5151  CFMap N;
    5252  CanonicalForm F= compress (G, N);
    53   CanonicalForm contentX= content (F, 1);
     53  CanonicalForm contentX= content (F, 1); //erwarte hier primitiven input: primitiv ÃŒber Z bzw. Z[a]
    5454  CanonicalForm contentY= content (F, 2);
    5555  F /= (contentX*contentY);
Note: See TracChangeset for help on using the changeset viewer.