source: git/factory/facFqFactorize.cc @ 93134c

spielwiese
Last change on this file since 93134c was 93134c, checked in by Martin Lee <martinlee84@…>, 13 years ago
bug fixes for factor recombination in case an extension is needed git-svn-id: file:///usr/local/Singular/svn/trunk@13801 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 51.8 KB
Line 
1/*****************************************************************************\
2 * Computer Algebra System SINGULAR
3\*****************************************************************************/
4/** @file facFqFactorize.cc
5 *
6 * This file implements functions for factoring a multivariate polynomial over
7 * a finite field.
8 *
9 * ABSTRACT: "Efficient Multivariate Factorization over Finite Fields" by
10 * L. Bernardin & M. Monagon.
11 *
12 * @author Martin Lee
13 *
14 * @internal @version \$Id$
15 *
16 **/
17/*****************************************************************************/
18
19#include <config.h>
20
21#include "assert.h"
22#include "debug.h"
23#include "timing.h"
24
25#include "facFqFactorizeUtil.h"
26#include "facFqFactorize.h"
27#include "cf_random.h"
28#include "facHensel.h"
29#include "cf_gcd_smallp.h"
30#include "cf_map_ext.h"
31#include "algext.h"
32
33#ifdef HAVE_NTL
34#include <NTL/ZZ_pEX.h>
35#include "NTLconvert.h"
36
37TIMING_DEFINE_PRINT(fac_bi_factorizer);
38TIMING_DEFINE_PRINT(fac_hensel_lift);
39TIMING_DEFINE_PRINT(fac_factor_recombination);
40
41static inline
42CanonicalForm
43listGCD (const CFList& L);
44
45static inline
46CanonicalForm
47myContent (const CanonicalForm& F)
48{
49  CanonicalForm G= swapvar (F, F.mvar(), Variable (1));
50  CFList L;
51  for (CFIterator i= G; i.hasTerms(); i++)
52    L.append (i.coeff());
53  if (L.length() == 2)
54    return swapvar (gcd (L.getFirst(), L.getLast()), F.mvar(), Variable (1));
55  if (L.length() == 1)
56    return LC (F, Variable (1));
57  return swapvar (listGCD (L), F.mvar(), Variable (1));
58}
59
60static inline
61CanonicalForm
62listGCD (const CFList& L)
63{
64  if (L.length() == 1)
65    return L.getFirst();
66  if (L.length() == 2)
67    return gcd (L.getFirst(), L.getLast());
68  else
69  {
70    CFList lHi, lLo;
71    CanonicalForm resultHi, resultLo;
72    int length= L.length()/2;
73    int j= 0;
74    for (CFListIterator i= L; j < length; i++, j++)
75      lHi.append (i.getItem());
76    lLo= Difference (L, lHi);
77    resultHi= listGCD (lHi);
78    resultLo= listGCD (lLo);
79    if (resultHi.isOne() || resultLo.isOne())
80      return 1;
81    return gcd (resultHi, resultLo);
82  }
83}
84
85static inline
86CanonicalForm
87myContent (const CanonicalForm& F, const Variable& x)
88{
89  if (degree (F, x) <= 0)
90    return 1;
91  CanonicalForm G= F;
92  bool swap= false;
93  if (x != F.mvar())
94  {
95    swap= true;
96    G= swapvar (F, x, F.mvar());
97  }
98  CFList L;
99  Variable alpha;
100  bool algExt= hasFirstAlgVar (G, alpha);
101  for (CFIterator i= G; i.hasTerms(); i++)
102    L.append (i.coeff());
103  if (L.length() == 2)
104  {
105    if (swap)
106      return swapvar (gcd (L.getFirst(), L.getLast()), F.mvar(), x);
107    else
108      return gcd (L.getFirst(), L.getLast());
109  }
110  if (L.length() == 1)
111  {
112    return LC (F, x);
113  }
114  if (swap)
115    return swapvar (listGCD (L), F.mvar(), x);
116  else
117    return listGCD (L);
118}
119
120static inline
121CanonicalForm myCompress (const CanonicalForm& F, CFMap& N)
122{
123  int n= F.level();
124  int * degsf= new int [n + 1];
125  int ** swap;
126  swap= new int* [n + 1];
127  for (int i= 0; i <= n; i++)
128  {
129    degsf[i]= 0;
130    swap [i]= new int [2];
131    swap [i] [0]= 0;
132    swap [i] [1]= 0;
133  }
134  int i= 1;
135  n= 1;
136  degsf= degrees (F, degsf);
137
138  CanonicalForm result= F;
139  while ( i <= F.level() )
140  {
141    while( degsf[i] == 0 ) i++;
142    swap[n][0]= i;
143    swap[n][1]= degsf[i];
144    if (i != n)
145      result= swapvar (result, Variable (n), Variable(i));
146    n++; i++;
147  }
148
149  int buf1, buf2;
150  n--;
151
152  for (i= 1; i < n; i++)
153  {
154    for (int j= 1; j < n - i + 1; j++)
155    {
156      if (swap[j][1] < swap[j + 1][1])
157      {
158        buf1= swap [j + 1] [0];
159        buf2= swap [j + 1] [1];
160        swap[j + 1] [0]= swap[j] [0];
161        swap[j + 1] [1]= swap[j] [1];
162        swap[j][0]= buf1;
163        swap[j][1]= buf2;
164        result= swapvar (result, Variable (j + 1), Variable (j));
165      }
166    }
167  }
168
169  for (i= n; i > 0; i--)
170  {
171    if (i != swap[i] [0])
172      N.newpair (Variable (i), Variable (swap[i] [0]));
173  }
174
175  for (i= 0; i <= n; i++)
176    delete [] swap[i];
177  delete [] swap;
178
179  delete [] degsf;
180
181  return result;
182}
183
184CFList
185monicFactorRecombi (const CFList& factors,const CanonicalForm& F, const
186                    CanonicalForm& M, const DegreePattern& degs)
187{
188  if (degs.getLength() == 1)
189    return CFList(F);
190
191  CFList T, S;
192
193  T= factors;
194
195  int s= 1;
196  CFList result;
197  CanonicalForm buf= F;
198  CanonicalForm LCBuf= LC (buf, Variable (1));
199  CanonicalForm g, gg;
200  int v [T.length()];
201  for (int i= 0; i < T.length(); i++)
202    v[i]= 0;
203  bool noSubset= false;
204  CFArray TT;
205  int subsetDeg;
206  DegreePattern bufDegs1= degs, bufDegs2;
207  TT= copy (factors);
208  while (T.length() >= 2*s)
209  {
210    while (noSubset == false)
211    {
212      if (T.length() == s)
213      {
214        result.append (prodMod (T, M));
215        return result;
216      }
217      S= subset (v, s, TT, noSubset);
218      if (noSubset) break;
219      subsetDeg= subsetDegree (S);
220      if (!degs.find (subsetDeg))
221        continue;
222      else
223      {
224        g= prodMod0 (S, M);
225        gg= mod (g*LCBuf, M); //univariate polys
226        gg /= content (gg);
227        if (fdivides (LC (gg), LCBuf))
228        {
229          g= prodMod (S, M);
230          gg= mulMod2 (LCBuf, g, M);
231          gg /= content (gg, Variable (1));
232          if (fdivides (gg, buf))
233          {
234            result.append (g);
235            buf /= gg;
236            LCBuf= LC (buf, Variable(1));
237            T= Difference (T, S);
238            // compute new possible degree pattern
239            bufDegs2= DegreePattern (T);
240            bufDegs1.intersect (bufDegs2);
241            bufDegs1.refine ();
242            if (T.length() < 2*s || T.length() == s ||
243                bufDegs1.getLength() == 1)
244            {
245              result.append (prodMod (T, M));
246              return result;
247            }
248            TT= copy (T);
249            indexUpdate (v, s, T.length(), noSubset);
250            if (noSubset) break;
251          }
252        }
253      }
254    }
255    s++;
256    if (T.length() < 2*s || T.length() == s)
257    {
258      result.append (prodMod (T, M));
259      return result;
260    }
261    int v [T.length()];
262    for (int i= 0; i < T.length(); i++)
263      v[i]= 0;
264    noSubset= false;
265  }
266  if (T.length() < 2*s)
267    result.append (prodMod (T, M));
268
269  return result;
270}
271
272CFList
273earlyMonicFactorDetect (CanonicalForm& F, CFList& factors,
274                        int& adaptedLiftBound, DegreePattern& degs,
275                        bool& success, int deg, const int bound)
276{
277  DegreePattern bufDegs1= degs;
278  DegreePattern bufDegs2;
279  CFList result;
280  CFList T= factors;
281  CanonicalForm buf= F;
282  CanonicalForm LCBuf= LC (buf, Variable (1));
283  CanonicalForm g, gg;
284  CanonicalForm M= power (F.mvar(), deg);
285  int d= bound;
286  int e= 0;
287  adaptedLiftBound= 0;
288  for (CFListIterator i= factors; i.hasItem(); i++)
289  {
290    if (!bufDegs1.find (degree (i.getItem(), 1)))
291      continue;
292    else
293    {
294      gg= i.getItem() (0, 1);
295      gg *= LCBuf;
296      gg= mod (gg, M);
297      if (fdivides (LC (gg), LCBuf))
298      {
299        g= i.getItem();
300        gg= mulMod2 (g, LCBuf, M);
301        gg /= content (gg, Variable (1));
302        if (fdivides (gg, buf))
303        {
304          result.append (g);
305          CanonicalForm bufbuf= buf;
306          buf /= gg;
307          d -= degree (gg) + degree (LC (gg, 1));
308          e= tmax (e, degree (gg) + degree (LC (gg, 1)));
309          LCBuf= LC (buf, Variable (1));
310          T= Difference (T, CFList (i.getItem()));
311
312          // compute new possible degree pattern
313          bufDegs2= DegreePattern (T);
314          bufDegs1.intersect (bufDegs2);
315          bufDegs1.refine ();
316          if (bufDegs1.getLength() <= 1)
317          {
318            result.append (prodMod (T, M));
319            break;
320          }
321        }
322      }
323    }
324  }
325  adaptedLiftBound= d;
326
327  if (adaptedLiftBound < deg)
328  {
329    factors= T;
330    degs= bufDegs1;
331    if (adaptedLiftBound < degree (F) + 1)
332    {
333      if (d == 1)
334      {
335        if (e + 1 > deg)
336        {
337          success= false;
338          adaptedLiftBound= deg;
339        }
340        else
341        {
342          F= buf;
343          success= true;
344          if (e + 1 < degree (F) + 1)
345            adaptedLiftBound= deg;
346          else
347            adaptedLiftBound= e + 1;
348        }
349      }
350      else
351      {
352        success= true;
353        F= buf;
354        adaptedLiftBound= deg;
355      }
356    }
357    else
358    {
359      F= buf;
360      success= true;
361    }
362  }
363  if (bufDegs1.getLength() <= 1)
364  {
365    degs= bufDegs1;
366    if (!success)
367      adaptedLiftBound= deg;
368  }
369
370  return result;
371}
372
373CFList biFactorizer (const CanonicalForm& F, const Variable& alpha,
374                     CanonicalForm& bivarEval, int& liftBound)
375{
376  CanonicalForm A= F;
377  bool GF= (CFFactory::gettype() == GaloisFieldDomain);
378
379  if (A.isUnivariate())
380    return uniFactorizer (F, alpha, GF);
381
382
383  CFMap N;
384  A= compress (A, N);
385  Variable y= A.mvar();
386  A /= Lc(A);
387
388  if (y.level() > 2) return CFList (F);
389  Variable x= Variable (1);
390
391  //remove content and factorize content
392  CanonicalForm contentAy= content (A, y);
393  CanonicalForm contentAx= content (A, x);
394  A= A/(contentAy*contentAx);
395  CFList contentAyFactors= uniFactorizer (contentAy, alpha, GF);
396
397  //trivial case
398  CFList factors;
399  if (A.inCoeffDomain())
400  {
401    append (factors, contentAyFactors);
402    decompress (factors, N);
403    bivarEval= N (y - bivarEval);
404    return factors;
405  }
406  else if (A.isUnivariate())
407  {
408    if (A.mvar() == x)
409      factors= uniFactorizer (A, alpha, GF);
410    append (factors, contentAyFactors);
411    decompress (factors, N);
412    bivarEval= N (y - bivarEval);
413    return factors;
414  }
415  bool fail;
416  CanonicalForm Aeval, evaluation, bufAeval, bufEvaluation, buf;
417  CFList uniFactors, list, bufUniFactors;
418  DegreePattern degs;
419  DegreePattern bufDegs;
420
421  int factorNums= 5;
422  double logarithm= (double) ilog2 (totaldegree (A));
423  logarithm /= log2exp;
424  logarithm= ceil (logarithm);
425  if (factorNums < (int) logarithm)
426    factorNums= (int) logarithm;
427  for (int i= 0; i < factorNums; i++)
428  {
429    if (i == 0)
430      Aeval= A (bivarEval, y);
431    else if (i > 0 && contentAx.inCoeffDomain())
432    {
433      Aeval= A;
434      evaluation= evalPoint (A, Aeval, alpha, list, GF, fail);
435    }
436
437    if (fail && (i != 0))
438      break;
439
440    // univariate factorization
441    uniFactors= uniFactorizer (Aeval, alpha, GF);
442
443    if (uniFactors.length() == 1)
444    {
445      append (factors, contentAyFactors);
446      if (contentAyFactors.isEmpty())
447        factors.append (F/lc(F));
448      else
449      {
450        A= A (y - bivarEval, y);
451        A= A/lc (A);
452        if (!LC(A, 1).isOne())
453        {
454          CanonicalForm s,t;
455          (void) extgcd (LC (A, 1), power (y, liftBound), s, t);
456          A *= s;
457          A= mod (A, power (y, liftBound));
458        }
459        factors.append (A);
460      }
461      decompress (factors, N);
462      bivarEval= N (bivarEval);
463      return factors;
464    }
465
466    // degree analysis
467    degs= DegreePattern (uniFactors);
468
469    if (i == 0)
470    {
471      bufAeval= Aeval;
472      bufEvaluation= bivarEval;
473      bufUniFactors= uniFactors;
474      bufDegs= degs;
475      if (!contentAx.inCoeffDomain())
476        break;
477    }
478    else
479    {
480      bufDegs.intersect (degs);
481      if (uniFactors.length() < bufUniFactors.length())
482      {
483        bufUniFactors= uniFactors;
484        bufAeval= Aeval;
485        bufEvaluation= evaluation;
486      }
487    }
488
489    if (bufDegs.getLength() == 1)
490    {
491      append (factors, contentAyFactors);
492      if (contentAyFactors.isEmpty())
493        factors.append (F/lc(F));
494      else
495      {
496        A= A (y - bivarEval, y);
497        A= A/lc (A);
498        if (!LC(A, 1).isOne())
499        {
500          CanonicalForm s,t;
501          (void) extgcd (LC (A, 1), power (y, liftBound), s, t);
502          A *= s;
503          A= mod (A, power (y, liftBound));
504        }
505        factors.append (A);
506      }
507      decompress (factors, N);
508      bivarEval= N (bivarEval);
509      return factors;
510    }
511    list.append (evaluation);
512  }
513
514  bivarEval= y - bufEvaluation;
515  A= A (y + bufEvaluation, y);
516  bufUniFactors.insert (LC (A, x));
517
518  // Hensel lifting
519  CFList diophant;
520  CFMatrix M= CFMatrix (liftBound, bufUniFactors.length() - 1);
521  CFArray Pi;
522  bool earlySuccess= false;
523  int newLiftBound= 0;
524  CFList earlyFactors;
525  int smallFactorDeg= 11; //tunable parameter
526  if (smallFactorDeg >= liftBound)
527    henselLift12 (A, bufUniFactors, liftBound, Pi, diophant, M);
528  else if (smallFactorDeg >= degree (A, y) + 1)
529  {
530    henselLift12 (A, bufUniFactors, degree (A, y) + 1, Pi, diophant, M);
531    earlyFactors= earlyMonicFactorDetect (A, bufUniFactors, newLiftBound,
532                                           bufDegs, earlySuccess,
533                                           degree (A, y) + 1, liftBound);
534    if (bufDegs.getLength() > 1 && !earlySuccess)
535    {
536      if (newLiftBound > degree (A, y) + 1)
537      {
538        liftBound= newLiftBound;
539        bufUniFactors.insert (LC(A, x));
540        henselLiftResume12 (A, bufUniFactors, degree (A, y) + 1, liftBound,
541                            Pi, diophant, M);
542      }
543    }
544    else if (earlySuccess)
545      liftBound= newLiftBound;
546  }
547  else if (smallFactorDeg < degree (A, y) + 1)
548  {
549    henselLift12 (A, bufUniFactors, smallFactorDeg, Pi, diophant, M);
550    earlyFactors= earlyMonicFactorDetect (A, bufUniFactors, newLiftBound,
551                                           bufDegs, earlySuccess,
552                                           smallFactorDeg, liftBound);
553    if (bufDegs.getLength() > 1 && !earlySuccess)
554    {
555      bufUniFactors.insert (LC (A, x));
556      henselLiftResume12 (A, bufUniFactors, smallFactorDeg, degree (A, y) +
557                          1, Pi, diophant, M);
558      earlyFactors= earlyMonicFactorDetect (A, bufUniFactors, newLiftBound,
559                                             bufDegs, earlySuccess,
560                                             degree (A, y) + 1, liftBound);
561      if (bufDegs.getLength() > 1 && !earlySuccess)
562      {
563        if (newLiftBound > degree (A, y) + 1)
564        {
565          if (newLiftBound < newLiftBound)
566            liftBound= newLiftBound;
567          bufUniFactors.insert (LC(A, x));
568          henselLiftResume12 (A, bufUniFactors, degree (A, y) + 1, liftBound,
569                              Pi, diophant, M);
570        }
571      }
572      else if (earlySuccess)
573        liftBound= newLiftBound;
574    }
575    else if (earlySuccess)
576      liftBound= newLiftBound;
577  }
578
579  if (newLiftBound > 0)
580    liftBound= newLiftBound;
581
582  CanonicalForm MODl= power (y, liftBound);
583
584  if (bufDegs.getLength() > 1)
585    factors= monicFactorRecombi (bufUniFactors, A, MODl, bufDegs);
586
587  if (earlySuccess)
588    factors= Union (earlyFactors, factors);
589  else if (!earlySuccess && bufDegs.getLength() == 1)
590    factors= earlyFactors;
591
592  decompressAppend (factors, contentAyFactors, N);
593
594  bivarEval= N (bivarEval);
595
596  return factors;
597}
598
599CFList
600extFactorRecombination (const CFList& factors, const CanonicalForm& F,
601                        const CFList& M, const ExtensionInfo& info,
602                        const CFList& evaluation)
603{
604  Variable alpha= info.getAlpha();
605  Variable beta= info.getBeta();
606  CanonicalForm gamma= info.getGamma();
607  CanonicalForm delta= info.getDelta();
608  int k= info.getGFDegree();
609  CFList source, dest;
610  if (factors.length() == 1)
611  {
612    CanonicalForm buf= reverseShift (F, evaluation);
613    return CFList (mapDown (buf, info, source, dest));
614  }
615  if (factors.length() < 1)
616    return CFList();
617
618  int degMipoBeta;
619  if (!k && beta.level() == 1)
620    degMipoBeta= 1;
621  else if (!k && beta.level() != 1)
622    degMipoBeta= degree (getMipo (beta));
623
624  CFList T, S;
625  T= factors;
626
627  int s= 1;
628  CFList result;
629  CanonicalForm buf;
630  if (beta != Variable (1))
631    buf= mapDown (F, gamma, delta, alpha, source, dest);
632  else
633    buf= F;
634
635  CanonicalForm g, LCBuf= LC (buf, Variable (1));
636  CanonicalForm buf2;
637  int v [T.length()];
638  for (int i= 0; i < T.length(); i++)
639    v[i]= 0;
640  bool noSubset= false;
641  CFArray TT;
642  int subsetDeg;
643  TT= copy (factors);
644  bool recombination= false;
645  bool trueFactor= false;
646  while (T.length() >= 2*s)
647  {
648    while (noSubset == false)
649    {
650      if (T.length() == s)
651      {
652        if (recombination)
653        {
654          T.insert (LCBuf);
655          g= prodMod (T, M);
656          T.removeFirst();
657          result.append (g/myContent (g));
658          g= reverseShift (g, evaluation);
659          g /= Lc (g);
660          appendTestMapDown (result, g, info, source, dest);
661          return result;
662        }
663        else
664        {
665          buf= reverseShift (buf, evaluation);
666          return CFList (buf);
667        }
668        return result;
669      }
670
671      S= subset (v, s, TT, noSubset);
672      if (noSubset) break;
673
674      S.insert (LCBuf);
675      g= prodMod (S, M);
676      S.removeFirst();
677      g /= myContent (g);
678      if (fdivides (g, buf))
679      {
680        buf2= reverseShift (g, evaluation);
681        buf2 /= Lc (buf2);
682        if (!k && beta == Variable (1))
683        {
684          if (degree (buf2, alpha) < degMipoBeta)
685          {
686            appendTestMapDown (result, buf2, info, source, dest);
687            buf /= g;
688            LCBuf= LC (buf, Variable (1));
689            recombination= true;
690            trueFactor= true;
691          }
692        }
693        else
694        {
695          if (!isInExtension (buf2, delta, k))
696          {
697            appendTestMapDown (result, buf2, info, source, dest);
698            buf /= g;
699            LCBuf= LC (buf, Variable (1));
700            recombination= true;
701            trueFactor= true;
702          }
703        }
704
705        if (trueFactor)
706        {
707          T= Difference (T, S);
708
709          if (T.length() < 2*s || T.length() == s)
710          {
711            buf= reverseShift (buf, evaluation);
712            buf /= Lc (buf);
713            appendTestMapDown (result, buf, info, source, dest);
714            return result;
715          }
716          trueFactor= false;
717          TT= copy (T);
718          indexUpdate (v, s, T.length(), noSubset);
719          if (noSubset) break;
720        }
721      }
722    }
723    s++;
724    if (T.length() < 2*s || T.length() == s)
725    {
726      buf= reverseShift (buf, evaluation);
727      appendTestMapDown (result, buf, info, source, dest);
728      return result;
729    }
730    int v [T.length()];
731    for (int i= 0; i < T.length(); i++)
732      v[i]= 0;
733    noSubset= false;
734  }
735  if (T.length() < 2*s)
736  {
737    buf= reverseShift (F, evaluation);
738    appendMapDown (result, buf, info, source, dest);
739  }
740
741  return result;
742}
743
744CFList
745factorRecombination (const CanonicalForm& F, const CFList& factors,
746                     const CFList& M)
747{
748  if (factors.length() == 1)
749    return CFList(F);
750  if (factors.length() < 1)
751    return CFList();
752
753  CFList T, S;
754
755  T= factors;
756
757  int s= 1;
758  CFList result;
759  CanonicalForm LCBuf= LC (F, Variable (1));
760  CanonicalForm g, buf= F;
761  int v [T.length()];
762  for (int i= 0; i < T.length(); i++)
763    v[i]= 0;
764  bool noSubset= false;
765  CFArray TT;
766  int subsetDeg;
767  TT= copy (factors);
768  Variable y= F.level() - 1;
769  bool recombination= false;
770  CanonicalForm h;
771  while (T.length() >= 2*s)
772  {
773    while (noSubset == false)
774    {
775      if (T.length() == s)
776      {
777        if (recombination)
778        {
779          T.insert (LC (buf));
780          g= prodMod (T, M);
781          result.append (g/myContent (g));
782          return result;
783        }
784        else
785          return CFList (F);
786      }
787      S= subset (v, s, TT, noSubset);
788      if (noSubset) break;
789      S.insert (LCBuf);
790      g= prodMod (S, M);
791      S.removeFirst();
792      g /= myContent (g);
793      if (fdivides (g, buf))
794      {
795        recombination= true;
796        result.append (g);
797        buf /= g;
798        LCBuf= LC (buf, Variable(1));
799        T= Difference (T, S);
800        if (T.length() < 2*s || T.length() == s)
801        {
802          result.append (buf);
803          return result;
804        }
805        TT= copy (T);
806        indexUpdate (v, s, T.length(), noSubset);
807        if (noSubset) break;
808      }
809    }
810    s++;
811    if (T.length() < 2*s || T.length() == s)
812    {
813      result.append (buf);
814      return result;
815    }
816    int v [T.length()];
817    for (int i= 0; i < T.length(); i++)
818      v[i]= 0;
819    noSubset= false;
820  }
821  if (T.length() < 2*s)
822    result.append (F);
823
824  return result;
825}
826
827int
828liftBoundAdaption (const CanonicalForm& F, const CFList& factors, bool&
829                   success, const int deg, const CFList& MOD, const int bound)
830{
831  int adaptedLiftBound= 0;
832  CanonicalForm buf= F;
833  Variable y= F.mvar();
834  CanonicalForm LCBuf= LC (buf, Variable (1));
835  CanonicalForm g;
836  CFList M= MOD;
837  M.append (power (y, deg));
838  int d= bound;
839  int e= 0;
840  int nBuf;
841  for (CFListIterator i= factors; i.hasItem(); i++)
842  {
843    g= mulMod (i.getItem(), LCBuf, M);
844    g /= myContent (g);
845    if (fdivides (g, buf))
846    {
847      nBuf= degree (g, y) + degree (LC (g, 1), y);
848      d -= nBuf;
849      e= tmax (e, nBuf);
850      buf /= g;
851      LCBuf= LC (buf, Variable (1));
852    }
853  }
854  adaptedLiftBound= d;
855
856  if (adaptedLiftBound < deg)
857  {
858    if (adaptedLiftBound < degree (F) + 1)
859    {
860      if (d == 1)
861      {
862        if (e + 1 > deg)
863        {
864          adaptedLiftBound= deg;
865          success= false;
866        }
867        else
868        {
869          success= true;
870          if (e + 1 < degree (F) + 1)
871            adaptedLiftBound= deg;
872          else
873            adaptedLiftBound= e + 1;
874        }
875      }
876      else
877      {
878        success= true;
879        adaptedLiftBound= deg;
880      }
881    }
882    else
883    {
884      success= true;
885    }
886  }
887  return adaptedLiftBound;
888}
889
890int
891extLiftBoundAdaption (const CanonicalForm& F, const CFList& factors, bool&
892                      success, const ExtensionInfo& info, const CFList& eval,
893                      const int deg, const CFList& MOD, const int bound)
894{
895  Variable alpha= info.getAlpha();
896  Variable beta= info.getBeta();
897  CanonicalForm gamma= info.getGamma();
898  CanonicalForm delta= info.getDelta();
899  int k= info.getGFDegree();
900  int adaptedLiftBound= 0;
901  CanonicalForm buf= F;
902  Variable y= F.mvar();
903  CanonicalForm LCBuf= LC (buf, Variable (1));
904  CanonicalForm g, gg;
905  CFList M= MOD;
906  M.append (power (y, deg));
907  adaptedLiftBound= 0;
908  int d= bound;
909  int e= 0;
910  int nBuf;
911  int degMipoBeta;
912  if (!k && beta.level() == 1)
913    degMipoBeta= 1;
914  else if (!k && beta.level() != 1)
915    degMipoBeta= degree (getMipo (beta));
916
917  for (CFListIterator i= factors; i.hasItem(); i++)
918  {
919    g= mulMod (i.getItem(), LCBuf, M);
920    g /= myContent (g);
921    if (fdivides (g, buf))
922    {
923      gg= reverseShift (g, eval);
924      gg /= Lc (gg);
925      if (!k && beta == Variable (1))
926      {
927        if (degree (gg, alpha) < degMipoBeta)
928        {
929          buf /= g;
930          nBuf= degree (g, y) + degree (LC (g, Variable (1)), y);
931          d -= nBuf;
932          e= tmax (e, nBuf);
933          LCBuf= LC (buf, Variable (1));
934        }
935      }
936      else
937      {
938        if (!isInExtension (gg, delta, k))
939        {
940          buf /= g;
941          nBuf= degree (g, y) + degree (LC (g, Variable (1)), y);
942          d -= nBuf;
943          e= tmax (e, nBuf);
944          LCBuf= LC (buf, Variable (1));
945        }
946      }
947    }
948  }
949  adaptedLiftBound= d;
950
951  if (adaptedLiftBound < deg)
952  {
953    if (adaptedLiftBound < degree (F) + 1)
954    {
955      if (d == 1)
956      {
957        if (e + 1 > deg)
958        {
959          adaptedLiftBound= deg;
960          success= false;
961        }
962        else
963        {
964          success= true;
965          if (e + 1 < degree (F) + 1)
966            adaptedLiftBound= deg;
967          else
968            adaptedLiftBound= e + 1;
969        }
970      }
971      else
972      {
973        success= true;
974        adaptedLiftBound= deg;
975      }
976    }
977    else
978    {
979      success= true;
980    }
981  }
982
983  return adaptedLiftBound;
984}
985
986CFList
987earlyFactorDetect (CanonicalForm& F, CFList& factors, int& adaptedLiftBound,
988                   bool& success, const int deg, const CFList& MOD,
989                   const int bound)
990{
991  CFList result;
992  CFList T= factors;
993  CanonicalForm buf= F;
994  Variable y= F.mvar();
995  CanonicalForm LCBuf= LC (buf, Variable (1));
996  CanonicalForm g;
997  CFList M= MOD;
998  M.append (power (y, deg));
999  adaptedLiftBound= 0;
1000  int d= bound;
1001  int e= 0;
1002  int nBuf;
1003  for (CFListIterator i= factors; i.hasItem(); i++)
1004  {
1005    g= mulMod (i.getItem(), LCBuf, M);
1006    g /= myContent (g);
1007    if (fdivides (g, buf))
1008    {
1009      result.append (g);
1010      nBuf= degree (g, y) + degree (LC (g, Variable (1)), y);
1011      d -= nBuf;
1012      e= tmax (e, nBuf);
1013      buf /= g;
1014      LCBuf= LC (buf, Variable (1));
1015      T= Difference (T, CFList (i.getItem()));
1016    }
1017  }
1018  adaptedLiftBound= d;
1019
1020  if (adaptedLiftBound < deg)
1021  {
1022    if (adaptedLiftBound < degree (F) + 1)
1023    {
1024      if (d == 1)
1025        adaptedLiftBound= tmin (e + 1, deg);
1026      else
1027        adaptedLiftBound= deg;
1028    }
1029    factors= T;
1030    F= buf;
1031    success= true;
1032  }
1033  return result;
1034}
1035
1036CFList
1037extEarlyFactorDetect (CanonicalForm& F, CFList& factors, int& adaptedLiftBound,
1038                      bool& success, const ExtensionInfo& info, const CFList&
1039                      eval, const int deg, const CFList& MOD, const int bound)
1040{
1041  Variable alpha= info.getAlpha();
1042  Variable beta= info.getBeta();
1043  CanonicalForm gamma= info.getGamma();
1044  CanonicalForm delta= info.getDelta();
1045  int k= info.getGFDegree();
1046  CFList result;
1047  CFList T= factors;
1048  CanonicalForm buf= F;
1049  Variable y= F.mvar();
1050  CanonicalForm LCBuf= LC (buf, Variable (1));
1051  CanonicalForm g, gg;
1052  CFList M= MOD;
1053  M.append (power (y, deg));
1054  adaptedLiftBound= 0;
1055  int d= bound;
1056  int e= 0;
1057  int nBuf;
1058  CFList source, dest;
1059
1060  int degMipoBeta;
1061  if (!k && beta.level() == 1)
1062    degMipoBeta= 1;
1063  else if (!k && beta.level() != 1)
1064    degMipoBeta= degree (getMipo (beta));
1065
1066  for (CFListIterator i= factors; i.hasItem(); i++)
1067  {
1068    g= mulMod (i.getItem(), LCBuf, M);
1069    g /= myContent (g);
1070    if (fdivides (g, buf))
1071    {
1072      gg= reverseShift (g, eval);
1073      gg /= Lc (gg);
1074      if (!k && beta == Variable (1))
1075      {
1076        if (degree (gg, alpha) < degMipoBeta)
1077        {
1078          appendTestMapDown (result, gg, info, source, dest);
1079          buf /= g;
1080          nBuf= degree (g, y) + degree (LC (g, Variable (1)), y);
1081          d -= nBuf;
1082          e= tmax (e, nBuf);
1083          LCBuf= LC (buf, Variable (1));
1084          T= Difference (T, CFList (i.getItem()));
1085        }
1086      }
1087      else
1088      {
1089        if (!isInExtension (gg, delta, k))
1090        {
1091          appendTestMapDown (result, gg, info, source, dest);
1092          buf /= g;
1093          nBuf= degree (g, y) + degree (LC (g, Variable (1)), y);
1094          d -= nBuf;
1095          e= tmax (e, nBuf);
1096          LCBuf= LC (buf, Variable (1));
1097          T= Difference (T, CFList (i.getItem()));
1098         }
1099      }
1100    }
1101  }
1102  adaptedLiftBound= d;
1103
1104  if (adaptedLiftBound < deg)
1105  {
1106    if (adaptedLiftBound < degree (F) + 1)
1107    {
1108      if (d == 1)
1109        adaptedLiftBound= tmin (e + 1, deg);
1110      else
1111        adaptedLiftBound= deg;
1112    }
1113    success= true;
1114    factors= T;
1115    F= buf;
1116  }
1117  return result;
1118}
1119
1120CFList
1121evalPoints (const CanonicalForm& F, CFList & eval, const Variable& alpha,
1122            CFList& list, const bool& GF, bool& fail)
1123{
1124  int k= F.level() - 1;
1125  Variable x= Variable (1);
1126  CFList result;
1127  FFRandom genFF;
1128  GFRandom genGF;
1129  int p= getCharacteristic ();
1130  double bound;
1131  if (alpha != Variable (1))
1132  {
1133    bound= pow ((double) p, (double) degree (getMipo(alpha)));
1134    bound= pow ((double) bound, (double) k);
1135  }
1136  else if (GF)
1137  {
1138    bound= pow ((double) p, (double) getGFDegree());
1139    bound= pow ((double) bound, (double) k);
1140  }
1141  else
1142    bound= pow ((double) p, (double) k);
1143
1144  CanonicalForm random;
1145  CanonicalForm deriv_x, gcd_deriv;
1146  do
1147  {
1148    random= 0;
1149    // possible overflow if list.length() does not fit into a int
1150    if (list.length() >= bound)
1151    {
1152      fail= true;
1153      break;
1154    }
1155    for (int i= 0; i < k; i++)
1156    {
1157      if (list.isEmpty())
1158        result.append (0);
1159      else if (GF)
1160      {
1161        result.append (genGF.generate());
1162        random += result.getLast()*power (x, i);
1163      }
1164      else if (alpha.level() != 1)
1165      {
1166        AlgExtRandomF genAlgExt (alpha);
1167        result.append (genAlgExt.generate());
1168        random += result.getLast()*power (x, i);
1169      }
1170      else
1171      {
1172        result.append (genFF.generate());
1173        random += result.getLast()*power (x, i);
1174      }
1175    }
1176    if (find (list, random))
1177    {
1178      result= CFList();
1179      continue;
1180    }
1181    int l= F.level();
1182    eval.insert (F);
1183    for (CFListIterator i= result; i.hasItem(); i++, l--)
1184      eval.insert (eval.getFirst()(i.getItem(), l));
1185
1186    if (degree (eval.getFirst()) != degree (F, 1))
1187    {
1188      if (!find (list, random))
1189        list.append (random);
1190      result= CFList();
1191      eval= CFList();
1192      continue;
1193    }
1194
1195    deriv_x= deriv (eval.getFirst(), x);
1196    gcd_deriv= gcd (eval.getFirst(), deriv_x);
1197    if (degree (gcd_deriv) > 0)
1198    {
1199      if (!find (list, random))
1200        list.append (random);
1201      result= CFList();
1202      eval= CFList();
1203      continue;
1204    }
1205    CFListIterator i= eval;
1206    i++;
1207    CanonicalForm contentx= content (i.getItem(), x);
1208    if (degree (contentx) > 0)
1209    {
1210      if (!find (list, random))
1211        list.append (random);
1212      result= CFList();
1213      eval= CFList();
1214      continue;
1215    }
1216
1217    if (list.length() >= bound)
1218    {
1219      fail= true;
1220      break;
1221    }
1222  } while (find (list, random));
1223
1224  if (!eval.isEmpty())
1225    eval.removeFirst();
1226
1227  return result;
1228}
1229
1230static inline
1231int newMainVariableSearch (CanonicalForm& A, CFList& Aeval, CFList&
1232                           evaluation, const Variable& alpha, const int lev)
1233{
1234  Variable x= Variable (1);
1235  CanonicalForm derivI, buf;
1236  bool GF= (CFFactory::gettype() == GaloisFieldDomain);
1237  int swapLevel= 0;
1238  CFList list;
1239  bool fail= false;
1240  buf= A;
1241  Aeval= CFList();
1242  evaluation= CFList();
1243  for (int i= lev; i <= A.level(); i++)
1244  {
1245    derivI= deriv (buf, Variable (i));
1246    if (!derivI.isZero())
1247    {
1248      buf= swapvar (buf, x, Variable (i));
1249      Aeval= CFList();
1250      evaluation= CFList();
1251      fail= false;
1252      evaluation= evalPoints (buf, Aeval, alpha, list, GF, fail);
1253      if (!fail)
1254      {
1255        A= buf;
1256        swapLevel= i;
1257        break;
1258      }
1259      else
1260        buf= A;
1261    }
1262  }
1263  return swapLevel;
1264}
1265
1266static inline
1267CanonicalForm lcmContent (const CanonicalForm& A, CFList& contentAi)
1268{
1269  int i= A.level();
1270  contentAi.append (myContent (A, i));
1271  contentAi.append (myContent (A, i - 1));
1272  CanonicalForm result= lcm (contentAi.getFirst(), contentAi.getLast());
1273  for (i= i - 2; i > 0; i--)
1274  {
1275    contentAi.append (content (A, i));
1276    result= lcm (result, contentAi.getLast());
1277  }
1278  return result;
1279}
1280
1281CFList
1282henselLiftAndEarly (CanonicalForm& A, CFList& MOD, int*& liftBounds, bool&
1283                    earlySuccess, CFList& earlyFactors, const CFList& Aeval,
1284                    const CFList& biFactors, const CFList& evaluation,
1285                    const ExtensionInfo& info)
1286{
1287  bool extension= info.isInExtension();
1288  CFList bufFactors= biFactors;
1289  bufFactors.insert (LC (Aeval.getFirst(), 1));
1290
1291  sortList (bufFactors, Variable (1));
1292
1293  CFList diophant;
1294  CFArray Pi;
1295  int smallFactorDeg= 11; //tunable parameter
1296  CFList result;
1297  int newLiftBound= 0;
1298  int adaptedLiftBound= 0;
1299  int liftBound= liftBounds[1];
1300
1301  earlySuccess= false;
1302  bool earlyReconst= false;
1303  CFList earlyReconstFactors;
1304  CFListIterator j= Aeval;
1305  j++;
1306  CanonicalForm buf= j.getItem();
1307  CFMatrix Mat= CFMatrix (liftBound, bufFactors.length() - 1);
1308  MOD= CFList (power (Variable (2), liftBounds[0]));
1309  if (smallFactorDeg >= liftBound)
1310  {
1311    result= henselLift23 (Aeval, bufFactors, liftBounds, diophant, Pi, Mat);
1312  }
1313  else if (smallFactorDeg >= degree (buf) + 1)
1314  {
1315    liftBounds[1]= degree (buf) + 1;
1316    result= henselLift23 (Aeval, bufFactors, liftBounds, diophant, Pi, Mat);
1317    if (Aeval.length() == 2)
1318    {
1319      if (!extension)
1320        earlyFactors= earlyFactorDetect
1321                       (buf, result, adaptedLiftBound, earlySuccess,
1322                        degree (buf) + 1, MOD, liftBound);
1323      else
1324        earlyFactors= extEarlyFactorDetect
1325                       (buf, result, adaptedLiftBound, earlySuccess,
1326                        info, evaluation, degree
1327                        (buf) + 1, MOD, liftBound);
1328    }
1329    else
1330    {
1331      if (!extension)
1332        adaptedLiftBound= liftBoundAdaption (buf, result, earlySuccess,
1333                                             degree (buf) + 1, MOD, liftBound);
1334      else
1335        adaptedLiftBound= extLiftBoundAdaption (buf, result, earlySuccess, info,
1336                                                evaluation, degree (buf) + 1,
1337                                                MOD, liftBound);
1338    }
1339    if (!earlySuccess)
1340    {
1341      result.insert (LC (buf, 1));
1342      liftBounds[1]= adaptedLiftBound;
1343      liftBound= adaptedLiftBound;
1344      henselLiftResume (buf, result, degree (buf) + 1, liftBound,
1345                        Pi, diophant, Mat, MOD);
1346    }
1347    else
1348      liftBounds[1]= adaptedLiftBound;
1349  }
1350  else if (smallFactorDeg < degree (buf) + 1)
1351  {
1352    liftBounds[1]= smallFactorDeg;
1353    result= henselLift23 (Aeval, bufFactors, liftBounds, diophant, Pi, Mat);
1354    if (Aeval.length() == 2)
1355    {
1356      if (!extension)
1357        earlyFactors= earlyFactorDetect (buf, result, adaptedLiftBound,
1358                                         earlySuccess, smallFactorDeg, MOD,
1359                                         liftBound);
1360      else
1361        earlyFactors= extEarlyFactorDetect (buf, result, adaptedLiftBound,
1362                                            earlySuccess, info, evaluation,
1363                                            smallFactorDeg, MOD, liftBound);
1364    }
1365    else
1366    {
1367      if (!extension)
1368        adaptedLiftBound= liftBoundAdaption (buf, result, earlySuccess,
1369                                             smallFactorDeg, MOD, liftBound);
1370      else
1371        adaptedLiftBound= extLiftBoundAdaption (buf, result, earlySuccess, info,
1372                                                evaluation, smallFactorDeg, MOD,
1373                                                liftBound);
1374    }
1375
1376    if (!earlySuccess)
1377    {
1378      result.insert (LC (buf, 1));
1379      henselLiftResume (buf, result, smallFactorDeg, degree (buf) + 1,
1380                        Pi, diophant, Mat, MOD);
1381      if (Aeval.length() == 2)
1382      {
1383         if (!extension)
1384           earlyFactors= earlyFactorDetect (buf, result, adaptedLiftBound,
1385                                            earlySuccess, degree (buf) + 1,
1386                                            MOD, liftBound);
1387         else
1388           earlyFactors= extEarlyFactorDetect (buf, result, adaptedLiftBound,
1389                                               earlySuccess, info, evaluation,
1390                                               degree (buf) + 1, MOD,
1391                                               liftBound);
1392      }
1393      else
1394      {
1395        if (!extension)
1396          adaptedLiftBound= liftBoundAdaption (buf, result, earlySuccess,
1397                                               degree (buf) + 1, MOD,liftBound);
1398        else
1399          adaptedLiftBound= extLiftBoundAdaption (buf, result, earlySuccess,
1400                                                  info, evaluation,
1401                                                  degree (buf) + 1, MOD,
1402                                                  liftBound);
1403      }
1404      if (!earlySuccess)
1405      {
1406        result.insert (LC (buf, 1));
1407        liftBounds[1]= adaptedLiftBound;
1408        liftBound= adaptedLiftBound;
1409        henselLiftResume (buf, result, degree (buf) + 1, liftBound,
1410                          Pi, diophant, Mat, MOD);
1411      }
1412      else
1413        liftBounds[1]= adaptedLiftBound;
1414    }
1415    else
1416      liftBounds[1]= adaptedLiftBound;
1417  }
1418
1419  MOD.append (power (Variable (3), liftBounds[1]));
1420
1421  if (Aeval.length() > 2)
1422  {
1423    CFListIterator j= Aeval;
1424    j++;
1425    CFList bufEval;
1426    bufEval.append (j.getItem());
1427    j++;
1428    int liftBoundsLength= Aeval.getLast().level() - 1;
1429    for (int i= 2; i <= liftBoundsLength && j.hasItem(); i++, j++)
1430    {
1431      earlySuccess= false;
1432      result.insert (LC (bufEval.getFirst(), 1));
1433      bufEval.append (j.getItem());
1434      liftBound= liftBounds[i];
1435      Mat= CFMatrix (liftBounds[i], result.length() - 1);
1436
1437      buf= j.getItem();
1438      if (smallFactorDeg >= liftBound)
1439        result= henselLift (bufEval, result, MOD, diophant, Pi, Mat,
1440                            liftBounds[i -  1], liftBounds[i]);
1441      else if (smallFactorDeg >= degree (buf) + 1)
1442      {
1443        result= henselLift (bufEval, result, MOD, diophant, Pi, Mat,
1444                            liftBounds[i -  1], degree (buf) + 1);
1445
1446        if (Aeval.length() == i + 1)
1447        {
1448          if (!extension)
1449            earlyFactors= earlyFactorDetect
1450                           (buf, result, adaptedLiftBound, earlySuccess,
1451                            degree (buf) + 1, MOD, liftBound);
1452          else
1453            earlyFactors= extEarlyFactorDetect
1454                           (buf, result, adaptedLiftBound, earlySuccess,
1455                            info, evaluation, degree (buf) + 1, MOD, liftBound);
1456        }
1457        else
1458        {
1459          if (!extension)
1460            adaptedLiftBound= liftBoundAdaption
1461                                (buf, result, earlySuccess, degree (buf)
1462                                 + 1,  MOD, liftBound);
1463          else
1464            adaptedLiftBound= extLiftBoundAdaption
1465                                (buf, result, earlySuccess, info, evaluation,
1466                                 degree (buf) + 1, MOD, liftBound);
1467        }
1468
1469        if (!earlySuccess)
1470        {
1471          result.insert (LC (buf, 1));
1472          liftBounds[i]= adaptedLiftBound;
1473          liftBound= adaptedLiftBound;
1474          henselLiftResume (buf, result, degree (buf) + 1, liftBound,
1475                            Pi, diophant, Mat, MOD);
1476        }
1477        else
1478        {
1479          liftBounds[i]= adaptedLiftBound;
1480        }
1481      }
1482      else if (smallFactorDeg < degree (buf) + 1)
1483      {
1484        result= henselLift (bufEval, result, MOD, diophant, Pi, Mat,
1485                            liftBounds[i -  1], smallFactorDeg);
1486
1487        if (Aeval.length() == i + 1)
1488        {
1489          if (!extension)
1490            earlyFactors= earlyFactorDetect
1491                           (buf, result, adaptedLiftBound, earlySuccess,
1492                            smallFactorDeg, MOD, liftBound);
1493          else
1494            earlyFactors= extEarlyFactorDetect
1495                           (buf, result, adaptedLiftBound, earlySuccess,
1496                            info, evaluation, smallFactorDeg, MOD, liftBound);
1497        }
1498        else
1499        {
1500          if (!extension)
1501            adaptedLiftBound= liftBoundAdaption
1502                                (buf, result, earlySuccess,
1503                                 smallFactorDeg, MOD, liftBound);
1504          else
1505            adaptedLiftBound= extLiftBoundAdaption
1506                                (buf, result, earlySuccess, info, evaluation,
1507                                 smallFactorDeg, MOD, liftBound);
1508        }
1509
1510        if (!earlySuccess)
1511        {
1512          result.insert (LC (buf, 1));
1513          henselLiftResume (buf, result, smallFactorDeg,
1514                            degree (buf) + 1, Pi, diophant, Mat, MOD);
1515          if (Aeval.length() == i + 1)
1516          {
1517            if (!extension)
1518              earlyFactors= earlyFactorDetect
1519                             (buf, result, adaptedLiftBound, earlySuccess,
1520                              degree (buf) +  1,  MOD, liftBound);
1521            else
1522              earlyFactors= extEarlyFactorDetect
1523                             (buf, result, adaptedLiftBound, earlySuccess,
1524                              info, evaluation, degree (buf) + 1, MOD,
1525                              liftBound);
1526          }
1527          else
1528          {
1529            if (!extension)
1530              adaptedLiftBound= liftBoundAdaption
1531                                  (buf, result, earlySuccess, degree
1532                                   (buf) +  1,  MOD, liftBound);
1533            else
1534              adaptedLiftBound= extLiftBoundAdaption
1535                                  (buf, result, earlySuccess, info, evaluation,
1536                                   degree (buf) + 1,  MOD, liftBound);
1537          }
1538
1539          if (!earlySuccess)
1540          {
1541            result.insert (LC (buf, 1));
1542            liftBounds[i]= adaptedLiftBound;
1543            liftBound= adaptedLiftBound;
1544            henselLiftResume (buf, result, degree (buf) + 1, liftBound,
1545                              Pi, diophant, Mat, MOD);
1546          }
1547          else
1548            liftBounds[i]= adaptedLiftBound;
1549        }
1550        else
1551          liftBounds[i]= adaptedLiftBound;
1552      }
1553      MOD.append (power (Variable (i + 2), liftBounds[i]));
1554      bufEval.removeFirst();
1555    }
1556    bufFactors= result;
1557  }
1558  else
1559    bufFactors= result;
1560
1561  if (earlySuccess)
1562    A= buf;
1563  return result;
1564}
1565
1566CFList
1567extFactorize (const CanonicalForm& F, const ExtensionInfo& info);
1568
1569CFList
1570multiFactorize (const CanonicalForm& F, const ExtensionInfo& info)
1571{
1572
1573  if (F.inCoeffDomain())
1574    return CFList (F);
1575
1576  // compress and find main Variable
1577  CFMap N;
1578  CanonicalForm A= myCompress (F, N);
1579
1580  A /= Lc (A); // make monic
1581
1582  Variable alpha= info.getAlpha();
1583  Variable beta= info.getBeta();
1584  CanonicalForm gamma= info.getGamma();
1585  CanonicalForm delta= info.getDelta();
1586  int k= info.getGFDegree();
1587  bool extension= info.isInExtension();
1588  bool GF= (CFFactory::gettype() == GaloisFieldDomain);
1589  //univariate case
1590  if (F.isUnivariate())
1591  {
1592    if (extension == false)
1593      return uniFactorizer (F, alpha, GF);
1594    else
1595    {
1596      CFList source, dest;
1597      A= mapDown (F, info, source, dest);
1598      return uniFactorizer (A, beta, GF);
1599    }
1600  }
1601
1602  //bivariate case
1603  if (A.level() == 2)
1604  {
1605    CFList buf= biFactorize (F, info);
1606    return buf;
1607  }
1608
1609  Variable x= Variable (1);
1610  Variable y= Variable (2);
1611
1612  // remove content
1613  CFList contentAi;
1614  CanonicalForm lcmCont= lcmContent (A, contentAi);
1615  A /= lcmCont;
1616
1617  // trivial after content removal
1618  CFList contentAFactors;
1619  if (A.inCoeffDomain())
1620  {
1621    for (CFListIterator i= contentAi; i.hasItem(); i++)
1622    {
1623      if (i.getItem().inCoeffDomain())
1624        continue;
1625      else
1626      {
1627        lcmCont /= i.getItem();
1628        contentAFactors=
1629        Union (multiFactorize (lcmCont, info),
1630               multiFactorize (i.getItem(), info));
1631        break;
1632      }
1633    }
1634    decompress (contentAFactors, N);
1635    normalize (contentAFactors);
1636    return contentAFactors;
1637  }
1638
1639  // factorize content
1640  contentAFactors= multiFactorize (lcmCont, info);
1641
1642  // univariate after content removal
1643  CFList factors;
1644  if (A.isUnivariate ())
1645  {
1646    factors= uniFactorizer (A, alpha, GF);
1647    append (factors, contentAFactors);
1648    decompress (factors, N);
1649    return factors;
1650  }
1651
1652  // check main variable
1653  int swapLevel= 0;
1654  CanonicalForm derivZ;
1655  CanonicalForm gcdDerivZ;
1656  CanonicalForm bufA= A;
1657  Variable z;
1658  for (int i= 1; i <= A.level(); i++)
1659  {
1660    z= Variable (i);
1661    derivZ= deriv (bufA, z);
1662    if (derivZ.isZero())
1663    {
1664      if (i == 1)
1665        swapLevel= 1;
1666      else
1667        continue;
1668    }
1669    else
1670    {
1671      if (swapLevel == 1)
1672      {
1673        swapLevel= i;
1674        A= swapvar (A, x, z);
1675      }
1676      gcdDerivZ= gcd (bufA, derivZ);
1677      if (degree (gcdDerivZ) > 0 && !derivZ.isZero())
1678      {
1679        CanonicalForm g= bufA/gcdDerivZ;
1680        CFList factorsG=
1681        Union (multiFactorize (g, info),
1682               multiFactorize (gcdDerivZ, info));
1683        appendSwapDecompress (factorsG, contentAFactors, N, swapLevel, x);
1684        normalize (factorsG);
1685        return factorsG;
1686      }
1687    }
1688  }
1689
1690
1691  CFList Aeval, list, evaluation, bufEvaluation, bufAeval;
1692  bool fail= false;
1693  int swapLevel2= 0;
1694  int level;
1695  int factorNums= 3;
1696  CanonicalForm bivarEval;
1697  CFList biFactors, bufBiFactors;
1698  CanonicalForm evalPoly;
1699  int lift, bufLift;
1700  double logarithm= (double) ilog2 (totaldegree (A));
1701  logarithm /= log2exp;
1702  logarithm= ceil (logarithm);
1703  if (factorNums < (int) logarithm)
1704    factorNums= (int) logarithm;
1705  // several bivariate factorizations
1706  for (int i= 0; i < factorNums; i++)
1707  {
1708    bufA= A;
1709    bufAeval= CFList();
1710    bufEvaluation= evalPoints (bufA, bufAeval, alpha, list, GF, fail);
1711    evalPoly= 0;
1712
1713    if (fail && (i == 0))
1714    {
1715      if (!swapLevel)
1716        level= 2;
1717      else
1718        level= swapLevel + 1;
1719
1720      swapLevel2= newMainVariableSearch (A, Aeval, evaluation, alpha, level);
1721
1722      if (!swapLevel2) // need to pass to an extension
1723      {
1724        factors= extFactorize (A, info);
1725        appendSwapDecompress (factors, contentAFactors, N, swapLevel, x);
1726        normalize (factors);
1727        return factors;
1728      }
1729      else
1730      {
1731        fail= false;
1732        bufAeval= Aeval;
1733        bufA= A;
1734        bufEvaluation= evaluation;
1735      }
1736    }
1737    else if (fail && (i > 0))
1738      break;
1739
1740    bivarEval= bufEvaluation.getLast();
1741    bufEvaluation.removeLast();
1742
1743    bufLift= degree (A, y) + 1 + degree (LC(A, x), y);
1744
1745    TIMING_START (fac_bi_factorizer);
1746    bufBiFactors= biFactorizer (bufAeval.getFirst(), alpha, bivarEval, bufLift);
1747    TIMING_END_AND_PRINT (fac_bi_factorizer,
1748                          "time for bivariate factorization: ");
1749
1750    if (bufBiFactors.length() == 1)
1751    {
1752      if (extension)
1753      {
1754        CFList source, dest;
1755        A= mapDown (A, info, source, dest);
1756      }
1757      factors.append (A);
1758      appendSwapDecompress (factors, contentAFactors, N, swapLevel,
1759                            swapLevel2, x);
1760      normalize (factors);
1761      return factors;
1762    }
1763
1764    bufEvaluation.append (-bivarEval[0]);
1765    if (i == 0)
1766    {
1767      Aeval= bufAeval;
1768      evaluation= bufEvaluation;
1769      biFactors= bufBiFactors;
1770      lift= bufLift;
1771    }
1772    else
1773    {
1774      if (bufBiFactors.length() < biFactors.length() ||
1775          ((bufLift < lift) && (bufBiFactors.length() == biFactors.length())))
1776      {
1777        Aeval= bufAeval;
1778        evaluation= bufEvaluation;
1779        biFactors= bufBiFactors;
1780        lift= bufLift;
1781      }
1782    }
1783    int k= 0;
1784    for (CFListIterator j= bufEvaluation; j.hasItem(); j++, k++)
1785      evalPoly += j.getItem()*power (x, k);
1786    list.append (evalPoly);
1787  }
1788
1789  //shifting to zero
1790  A= shift2Zero (A, Aeval, evaluation);
1791
1792  int* liftBounds;
1793  int liftBoundsLength= F.level() - 1;
1794  liftBounds= liftingBounds (A, lift);
1795
1796  CFList MOD;
1797  bool earlySuccess;
1798  CFList earlyFactors;
1799  TIMING_START (fac_hensel_lift);
1800  CFList liftedFactors= henselLiftAndEarly
1801                        (A, MOD, liftBounds, earlySuccess, earlyFactors,
1802                         Aeval, biFactors, evaluation, info);
1803  TIMING_END_AND_PRINT (fac_hensel_lift, "time for hensel lifting: ");
1804
1805  if (!extension)
1806  {
1807    TIMING_START (fac_factor_recombination);
1808    factors= factorRecombination (A, liftedFactors, MOD);
1809    TIMING_END_AND_PRINT (fac_factor_recombination,
1810                          "time for factor recombination: ");
1811  }
1812  else
1813  {
1814    TIMING_START (fac_factor_recombination);
1815    factors= extFactorRecombination (liftedFactors, A, MOD, info, evaluation);
1816    TIMING_END_AND_PRINT (fac_factor_recombination,
1817                          "time for factor recombination: ");
1818  }
1819
1820  if (earlySuccess)
1821    factors= Union (factors, earlyFactors);
1822
1823  if (!extension)
1824  {
1825    for (CFListIterator i= factors; i.hasItem(); i++)
1826    {
1827      int kk= Aeval.getLast().level();
1828      for (CFListIterator j= evaluation; j.hasItem(); j++, kk--)
1829      {
1830        if (i.getItem().level() < kk)
1831          continue;
1832        i.getItem()= i.getItem() (Variable (kk) - j.getItem(), kk);
1833      }
1834    }
1835  }
1836
1837  swap (factors, swapLevel, swapLevel2, x);
1838  append (factors, contentAFactors);
1839  decompress (factors, N);
1840  normalize (factors);
1841
1842  delete[] liftBounds;
1843
1844  return factors;
1845}
1846
1847/// multivariate factorization over an extension of the initial field
1848CFList
1849extFactorize (const CanonicalForm& F, const ExtensionInfo& info)
1850{
1851  CanonicalForm A= F;
1852
1853  Variable alpha= info.getAlpha();
1854  Variable beta= info.getBeta();
1855  int k= info.getGFDegree();
1856  char cGFName= info.getGFName();
1857  CanonicalForm delta= info.getDelta();
1858  bool GF= (CFFactory::gettype() == GaloisFieldDomain);
1859  Variable w= Variable (1);
1860
1861  CFList factors;
1862  if (!GF && alpha == w)  // we are in F_p
1863  {
1864    CFList factors;
1865    bool extension= true;
1866    int p= getCharacteristic();
1867    if (p*p < (1<<16)) // pass to GF if possible
1868    {
1869      setCharacteristic (getCharacteristic(), 2, 'Z');
1870      ExtensionInfo info= ExtensionInfo (extension);
1871      A= A.mapinto();
1872      factors= multiFactorize (A, info);
1873
1874      Variable vBuf= rootOf (gf_mipo);
1875      setCharacteristic (getCharacteristic());
1876      for (CFListIterator j= factors; j.hasItem(); j++)
1877        j.getItem()= GF2FalphaRep (j.getItem(), vBuf);
1878    }
1879    else  // not able to pass to GF, pass to F_p(\alpha)
1880    {
1881      Variable v= chooseExtension (A, beta);
1882      ExtensionInfo info= ExtensionInfo (v, extension);
1883      factors= multiFactorize (A, info);
1884    }
1885    return factors;
1886  }
1887  else if (!GF && (alpha != w)) // we are in F_p(\alpha)
1888  {
1889    if (k == 1) // need factorization over F_p
1890    {
1891      int extDeg= degree (getMipo (alpha));
1892      extDeg++;
1893      CanonicalForm mipo= randomIrredpoly (extDeg + 1, Variable (1));
1894      Variable v= rootOf (mipo);
1895      ExtensionInfo info= ExtensionInfo (v);
1896      factors= biFactorize (A, info);
1897    }
1898    else
1899    {
1900      if (beta == Variable (1))
1901      {
1902        Variable v= chooseExtension (A, alpha);
1903        CanonicalForm primElem, imPrimElem;
1904        bool primFail= false;
1905        Variable vBuf;
1906        primElem= primitiveElement (alpha, vBuf, primFail);
1907        ASSERT (!primFail, "failure in integer factorizer");
1908        if (primFail)
1909          ; //ERROR
1910        else
1911          imPrimElem= mapPrimElem (primElem, vBuf, v);
1912
1913        CFList source, dest;
1914        CanonicalForm bufA= mapUp (A, alpha, v, primElem, imPrimElem,
1915                                   source, dest);
1916        ExtensionInfo info= ExtensionInfo (v, alpha, imPrimElem, primElem);
1917        factors= biFactorize (bufA, info);
1918      }
1919      else
1920      {
1921        Variable v= chooseExtension (A, alpha);
1922        CanonicalForm primElem, imPrimElem;
1923        bool primFail= false;
1924        Variable vBuf;
1925        ASSERT (!primFail, "failure in integer factorizer");
1926        if (primFail)
1927          ; //ERROR
1928        else
1929          imPrimElem= mapPrimElem (delta, beta, v); //oder mapPrimElem (primElem, vBuf, v);
1930
1931        CFList source, dest;
1932        CanonicalForm bufA= mapDown (A, info, source, dest);
1933        source= CFList();
1934        dest= CFList();
1935        bufA= mapUp (bufA, beta, v, delta, imPrimElem, source, dest);
1936        ExtensionInfo info= ExtensionInfo (v, beta, imPrimElem, delta);
1937        factors= biFactorize (bufA, info);
1938      }
1939    }
1940    return factors;
1941  }
1942  else // we are in GF (p^k)
1943  {
1944    int p= getCharacteristic();
1945    int extensionDeg= getGFDegree();
1946    bool extension= true;
1947    if (k == 1) // need factorization over F_p
1948    {
1949      extensionDeg++;
1950      if (pow ((double) p, (double) extensionDeg) < (1<<16))
1951      // pass to GF(p^k+1)
1952      {
1953        setCharacteristic (p);
1954        Variable vBuf= rootOf (gf_mipo);
1955        A= GF2FalphaRep (A, vBuf);
1956        setCharacteristic (p, extensionDeg, 'Z');
1957        ExtensionInfo info= ExtensionInfo (extension);
1958        factors= multiFactorize (A.mapinto(), info);
1959      }
1960      else // not able to pass to another GF, pass to F_p(\alpha)
1961      {
1962        setCharacteristic (p);
1963        Variable vBuf= rootOf (gf_mipo);
1964        A= GF2FalphaRep (A, vBuf);
1965        Variable v= chooseExtension (A, beta);
1966        ExtensionInfo info= ExtensionInfo (v, extension);
1967        factors= multiFactorize (A, info);
1968      }
1969    }
1970    else // need factorization over GF (p^k)
1971    {
1972      if (pow ((double) p, (double) 2*extensionDeg) < (1<<16))
1973      // pass to GF(p^2k)
1974      {
1975        setCharacteristic (p, 2*extensionDeg, 'Z');
1976        ExtensionInfo info= ExtensionInfo (k, cGFName, extension);
1977        factors= multiFactorize (GFMapUp (A, extensionDeg), info);
1978        setCharacteristic (p, extensionDeg, cGFName);
1979      }
1980      else // not able to pass to GF (p^2k), pass to F_p (\alpha)
1981      {
1982        setCharacteristic (p);
1983        Variable v1= rootOf (gf_mipo);
1984        A= GF2FalphaRep (A, v1);
1985        Variable v2= chooseExtension (A, v1);
1986        CanonicalForm primElem, imPrimElem;
1987        bool primFail= false;
1988        Variable vBuf;
1989        primElem= primitiveElement (v1, vBuf, primFail);
1990        if (primFail)
1991        {
1992          ; //ERROR
1993        }
1994        else
1995        {
1996          imPrimElem= mapPrimElem (primElem, vBuf, v2);
1997        }
1998        CFList source, dest;
1999        CanonicalForm bufA= mapUp (A, alpha, beta, primElem, imPrimElem,
2000                                     source, dest);
2001        ExtensionInfo info= ExtensionInfo (v2, v1, imPrimElem, primElem);
2002        factors= multiFactorize (bufA, info);
2003        setCharacteristic (p, k, cGFName);
2004        for (CFListIterator i= factors; i.hasItem(); i++)
2005          i.getItem()= Falpha2GFRep (i.getItem());
2006      }
2007    }
2008    return factors;
2009  }
2010}
2011
2012#endif
2013/* HAVE_NTL */
2014
Note: See TracBrowser for help on using the repository browser.