source: git/factory/facFqFactorize.cc @ 7a1151

spielwiese
Last change on this file since 7a1151 was 7a1151, checked in by Martin Lee <martinlee84@…>, 12 years ago
fix: compiler warnings
  • Property mode set to 100644
File size: 75.6 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. Precomputation of leading coefficients is
11 * described in "Sparse Hensel lifting" by E. Kaltofen
12 *
13 * @author Martin Lee
14 *
15 * @internal @version \$Id$
16 *
17 **/
18/*****************************************************************************/
19
20#include "config.h"
21
22#include "cf_assert.h"
23#include "debug.h"
24#include "timing.h"
25
26#include "facFqFactorizeUtil.h"
27#include "facFqFactorize.h"
28#include "cf_random.h"
29#include "facHensel.h"
30#include "cf_gcd_smallp.h"
31#include "cf_map_ext.h"
32#include "algext.h"
33#include "facSparseHensel.h"
34#include "facMul.h"
35
36#ifdef HAVE_NTL
37#include <NTL/ZZ_pEX.h>
38#include "NTLconvert.h"
39
40TIMING_DEFINE_PRINT(fac_bi_factorizer)
41TIMING_DEFINE_PRINT(fac_hensel_lift)
42TIMING_DEFINE_PRINT(fac_factor_recombination)
43
44static inline
45CanonicalForm
46listGCD (const CFList& L);
47
48static inline
49CanonicalForm
50myContent (const CanonicalForm& F)
51{
52  Variable x= Variable (1);
53  CanonicalForm G= swapvar (F, F.mvar(), x);
54  CFList L;
55  for (CFIterator i= G; i.hasTerms(); i++)
56    L.append (i.coeff());
57  if (L.length() == 2)
58    return swapvar (gcd (L.getFirst(), L.getLast()), F.mvar(), x);
59  if (L.length() == 1)
60    return LC (F, x);
61  return swapvar (listGCD (L), F.mvar(), x);
62}
63
64static inline
65CanonicalForm
66listGCD (const CFList& L)
67{
68  if (L.length() == 1)
69    return L.getFirst();
70  if (L.length() == 2)
71    return gcd (L.getFirst(), L.getLast());
72  else
73  {
74    CFList lHi, lLo;
75    CanonicalForm resultHi, resultLo;
76    int length= L.length()/2;
77    int j= 0;
78    for (CFListIterator i= L; j < length; i++, j++)
79      lHi.append (i.getItem());
80    lLo= Difference (L, lHi);
81    resultHi= listGCD (lHi);
82    resultLo= listGCD (lLo);
83    if (resultHi.isOne() || resultLo.isOne())
84      return 1;
85    return gcd (resultHi, resultLo);
86  }
87}
88
89static inline
90CanonicalForm
91myContent (const CanonicalForm& F, const Variable& x)
92{
93  if (degree (F, x) <= 0)
94    return 1;
95  CanonicalForm G= F;
96  bool swap= false;
97  if (x != F.mvar())
98  {
99    swap= true;
100    G= swapvar (F, x, F.mvar());
101  }
102  CFList L;
103  Variable alpha;
104  for (CFIterator i= G; i.hasTerms(); i++)
105    L.append (i.coeff());
106  if (L.length() == 2)
107  {
108    if (swap)
109      return swapvar (gcd (L.getFirst(), L.getLast()), F.mvar(), x);
110    else
111      return gcd (L.getFirst(), L.getLast());
112  }
113  if (L.length() == 1)
114  {
115    return LC (F, x);
116  }
117  if (swap)
118    return swapvar (listGCD (L), F.mvar(), x);
119  else
120    return listGCD (L);
121}
122
123CanonicalForm myCompress (const CanonicalForm& F, CFMap& N)
124{
125  int n= F.level();
126  int * degsf= new int [n + 1];
127  int ** swap;
128  swap= new int* [n + 1];
129  for (int i= 0; i <= n; i++)
130  {
131    degsf[i]= 0;
132    swap [i]= new int [2];
133    swap [i] [0]= 0;
134    swap [i] [1]= 0;
135  }
136  int i= 1;
137  n= 1;
138  degsf= degrees (F, degsf);
139
140  CanonicalForm result= F;
141  while ( i <= F.level() )
142  {
143    while( degsf[i] == 0 ) i++;
144    swap[n][0]= i;
145    swap[n][1]= degsf[i];
146    if (i != n)
147      result= swapvar (result, Variable (n), Variable(i));
148    n++; i++;
149  }
150
151  int buf1, buf2;
152  n--;
153
154  for (i= 1; i < n; i++)
155  {
156    for (int j= 1; j < n - i + 1; j++)
157    {
158      if (swap[j][1] < swap[j + 1][1])
159      {
160        buf1= swap [j + 1] [0];
161        buf2= swap [j + 1] [1];
162        swap[j + 1] [0]= swap[j] [0];
163        swap[j + 1] [1]= swap[j] [1];
164        swap[j][0]= buf1;
165        swap[j][1]= buf2;
166        result= swapvar (result, Variable (j + 1), Variable (j));
167      }
168    }
169  }
170
171  for (i= n; i > 0; i--)
172  {
173    if (i != swap[i] [0])
174      N.newpair (Variable (i), Variable (swap[i] [0]));
175  }
176
177  for (i= 0; i <= n; i++)
178    delete [] swap[i];
179  delete [] swap;
180
181  delete [] degsf;
182
183  return result;
184}
185
186CFList
187extFactorRecombination (const CFList& factors, const CanonicalForm& F,
188                        const CFList& M, const ExtensionInfo& info,
189                        const CFList& evaluation)
190{
191  Variable alpha= info.getAlpha();
192  Variable beta= info.getBeta();
193  CanonicalForm gamma= info.getGamma();
194  CanonicalForm delta= info.getDelta();
195  int k= info.getGFDegree();
196  CFList source, dest;
197  if (factors.length() == 1)
198  {
199    CanonicalForm buf= reverseShift (F, evaluation);
200    return CFList (mapDown (buf, info, source, dest));
201  }
202  if (factors.length() < 1)
203    return CFList();
204
205  int degMipoBeta= 1;
206  if (!k && beta.level() != 1)
207    degMipoBeta= degree (getMipo (beta));
208
209  CFList T, S;
210  T= factors;
211
212  int s= 1;
213  CFList result;
214  CanonicalForm buf;
215
216  buf= F;
217
218  Variable x= Variable (1);
219  CanonicalForm g, LCBuf= LC (buf, x);
220  CanonicalForm buf2, quot;
221  int * v= new int [T.length()];
222  for (int i= 0; i < T.length(); i++)
223    v[i]= 0;
224  bool noSubset= false;
225  CFArray TT;
226  TT= copy (factors);
227  bool recombination= false;
228  bool trueFactor= false;
229  while (T.length() >= 2*s)
230  {
231    while (noSubset == false)
232    {
233      if (T.length() == s)
234      {
235        delete [] v;
236        if (recombination)
237        {
238          T.insert (LCBuf);
239          g= prodMod (T, M);
240          T.removeFirst();
241          result.append (g/myContent (g));
242          g= reverseShift (g, evaluation);
243          g /= Lc (g);
244          appendTestMapDown (result, g, info, source, dest);
245          return result;
246        }
247        else
248        {
249          buf= reverseShift (buf, evaluation);
250          return CFList (buf);
251        }
252      }
253
254      S= subset (v, s, TT, noSubset);
255      if (noSubset) break;
256
257      S.insert (LCBuf);
258      g= prodMod (S, M);
259      S.removeFirst();
260      g /= myContent (g);
261      if (fdivides (g, buf, quot))
262      {
263        buf2= reverseShift (g, evaluation);
264        buf2 /= Lc (buf2);
265        if (!k && beta == x)
266        {
267          if (degree (buf2, alpha) < degMipoBeta)
268          {
269            appendTestMapDown (result, buf2, info, source, dest);
270            buf= quot;
271            LCBuf= LC (buf, x);
272            recombination= true;
273            trueFactor= true;
274          }
275        }
276        else
277        {
278          if (!isInExtension (buf2, gamma, k, delta, source, dest))
279          {
280            appendTestMapDown (result, buf2, info, source, dest);
281            buf /= g;
282            LCBuf= LC (buf, x);
283            recombination= true;
284            trueFactor= true;
285          }
286        }
287
288        if (trueFactor)
289        {
290          T= Difference (T, S);
291
292          if (T.length() < 2*s || T.length() == s)
293          {
294            buf= reverseShift (buf, evaluation);
295            buf /= Lc (buf);
296            appendTestMapDown (result, buf, info, source, dest);
297            delete [] v;
298            return result;
299          }
300          trueFactor= false;
301          TT= copy (T);
302          indexUpdate (v, s, T.length(), noSubset);
303          if (noSubset) break;
304        }
305      }
306    }
307    s++;
308    if (T.length() < 2*s || T.length() == s)
309    {
310      buf= reverseShift (buf, evaluation);
311      appendTestMapDown (result, buf, info, source, dest);
312      delete [] v;
313      return result;
314    }
315    for (int i= 0; i < T.length(); i++)
316      v[i]= 0;
317    noSubset= false;
318  }
319  if (T.length() < 2*s)
320  {
321    buf= reverseShift (F, evaluation);
322    appendMapDown (result, buf, info, source, dest);
323  }
324
325  delete [] v;
326  return result;
327}
328
329CFList
330factorRecombination (const CanonicalForm& F, const CFList& factors,
331                     const CFList& M)
332{
333  if (factors.length() == 1)
334    return CFList(F);
335  if (factors.length() < 1)
336    return CFList();
337
338  CFList T, S;
339
340  T= factors;
341
342  int s= 1;
343  CFList result;
344  CanonicalForm LCBuf= LC (F, Variable (1));
345  CanonicalForm g, buf= F;
346  int * v= new int [T.length()];
347  for (int i= 0; i < T.length(); i++)
348    v[i]= 0;
349  bool noSubset= false;
350  CFArray TT;
351  TT= copy (factors);
352  Variable y= F.level() - 1;
353  bool recombination= false;
354  CanonicalForm h, quot;
355  while (T.length() >= 2*s)
356  {
357    while (noSubset == false)
358    {
359      if (T.length() == s)
360      {
361        delete [] v;
362        if (recombination)
363        {
364          T.insert (LC (buf));
365          g= prodMod (T, M);
366          result.append (g/myContent (g));
367          return result;
368        }
369        else
370          return CFList (F);
371      }
372      S= subset (v, s, TT, noSubset);
373      if (noSubset) break;
374      S.insert (LCBuf);
375      g= prodMod (S, M);
376      S.removeFirst();
377      g /= myContent (g);
378      if (fdivides (g, buf, quot))
379      {
380        recombination= true;
381        result.append (g);
382        buf= quot;
383        LCBuf= LC (buf, Variable(1));
384        T= Difference (T, S);
385        if (T.length() < 2*s || T.length() == s)
386        {
387          result.append (buf);
388          delete [] v;
389          return result;
390        }
391        TT= copy (T);
392        indexUpdate (v, s, T.length(), noSubset);
393        if (noSubset) break;
394      }
395    }
396    s++;
397    if (T.length() < 2*s || T.length() == s)
398    {
399      result.append (buf);
400      delete [] v;
401      return result;
402    }
403    for (int i= 0; i < T.length(); i++)
404      v[i]= 0;
405    noSubset= false;
406  }
407  if (T.length() < 2*s)
408    result.append (F);
409
410  delete [] v;
411  return result;
412}
413
414int
415liftBoundAdaption (const CanonicalForm& F, const CFList& factors, bool&
416                   success, const int deg, const CFList& MOD, const int bound)
417{
418  int adaptedLiftBound= 0;
419  CanonicalForm buf= F;
420  Variable y= F.mvar();
421  Variable x= Variable (1);
422  CanonicalForm LCBuf= LC (buf, x);
423  CanonicalForm g, quot;
424  CFList M= MOD;
425  M.append (power (y, deg));
426  int d= bound;
427  int e= 0;
428  int nBuf;
429  for (CFListIterator i= factors; i.hasItem(); i++)
430  {
431    g= mulMod (i.getItem(), LCBuf, M);
432    g /= myContent (g);
433    if (fdivides (g, buf, quot))
434    {
435      nBuf= degree (g, y) + degree (LC (g, 1), y);
436      d -= nBuf;
437      e= tmax (e, nBuf);
438      buf= quot;
439      LCBuf= LC (buf, x);
440    }
441  }
442  adaptedLiftBound= d;
443
444  if (adaptedLiftBound < deg)
445  {
446    if (adaptedLiftBound < degree (F) + 1)
447    {
448      if (d == 1)
449      {
450        if (e + 1 > deg)
451        {
452          adaptedLiftBound= deg;
453          success= false;
454        }
455        else
456        {
457          success= true;
458          if (e + 1 < degree (F) + 1)
459            adaptedLiftBound= deg;
460          else
461            adaptedLiftBound= e + 1;
462        }
463      }
464      else
465      {
466        success= true;
467        adaptedLiftBound= deg;
468      }
469    }
470    else
471    {
472      success= true;
473    }
474  }
475  return adaptedLiftBound;
476}
477
478int
479extLiftBoundAdaption (const CanonicalForm& F, const CFList& factors, bool&
480                      success, const ExtensionInfo& info, const CFList& eval,
481                      const int deg, const CFList& MOD, const int bound)
482{
483  Variable alpha= info.getAlpha();
484  Variable beta= info.getBeta();
485  CanonicalForm gamma= info.getGamma();
486  CanonicalForm delta= info.getDelta();
487  int k= info.getGFDegree();
488  int adaptedLiftBound= 0;
489  CanonicalForm buf= F;
490  Variable y= F.mvar();
491  Variable x= Variable (1);
492  CanonicalForm LCBuf= LC (buf, x);
493  CanonicalForm g, gg, quot;
494  CFList M= MOD;
495  M.append (power (y, deg));
496  adaptedLiftBound= 0;
497  int d= bound;
498  int e= 0;
499  int nBuf;
500  int degMipoBeta= 1;
501  if (!k && beta.level() != 1)
502    degMipoBeta= degree (getMipo (beta));
503
504  CFList source, dest;
505  for (CFListIterator i= factors; i.hasItem(); i++)
506  {
507    g= mulMod (i.getItem(), LCBuf, M);
508    g /= myContent (g);
509    if (fdivides (g, buf, quot))
510    {
511      gg= reverseShift (g, eval);
512      gg /= Lc (gg);
513      if (!k && beta == x)
514      {
515        if (degree (gg, alpha) < degMipoBeta)
516        {
517          buf= quot;
518          nBuf= degree (g, y) + degree (LC (g, x), y);
519          d -= nBuf;
520          e= tmax (e, nBuf);
521          LCBuf= LC (buf, x);
522        }
523      }
524      else
525      {
526        if (!isInExtension (gg, gamma, k, delta, source, dest))
527        {
528          buf= quot;
529          nBuf= degree (g, y) + degree (LC (g, x), y);
530          d -= nBuf;
531          e= tmax (e, nBuf);
532          LCBuf= LC (buf, x);
533        }
534      }
535    }
536  }
537  adaptedLiftBound= d;
538
539  if (adaptedLiftBound < deg)
540  {
541    if (adaptedLiftBound < degree (F) + 1)
542    {
543      if (d == 1)
544      {
545        if (e + 1 > deg)
546        {
547          adaptedLiftBound= deg;
548          success= false;
549        }
550        else
551        {
552          success= true;
553          if (e + 1 < degree (F) + 1)
554            adaptedLiftBound= deg;
555          else
556            adaptedLiftBound= e + 1;
557        }
558      }
559      else
560      {
561        success= true;
562        adaptedLiftBound= deg;
563      }
564    }
565    else
566    {
567      success= true;
568    }
569  }
570
571  return adaptedLiftBound;
572}
573
574CFList
575earlyFactorDetect (CanonicalForm& F, CFList& factors, int& adaptedLiftBound,
576                   bool& success, const int deg, const CFList& MOD,
577                   const int bound)
578{
579  CFList result;
580  CFList T= factors;
581  CanonicalForm buf= F;
582  Variable y= F.mvar();
583  Variable x= Variable (1);
584  CanonicalForm LCBuf= LC (buf, x);
585  CanonicalForm g, quot;
586  CFList M= MOD;
587  M.append (power (y, deg));
588  adaptedLiftBound= 0;
589  int d= bound;
590  int e= 0;
591  int nBuf;
592  for (CFListIterator i= factors; i.hasItem(); i++)
593  {
594    g= mulMod (i.getItem(), LCBuf, M);
595    g /= myContent (g);
596    if (fdivides (g, buf, quot))
597    {
598      result.append (g);
599      nBuf= degree (g, y) + degree (LC (g, x), y);
600      d -= nBuf;
601      e= tmax (e, nBuf);
602      buf= quot;
603      LCBuf= LC (buf, x);
604      T= Difference (T, CFList (i.getItem()));
605    }
606  }
607  adaptedLiftBound= d;
608
609  if (adaptedLiftBound < deg)
610  {
611    if (adaptedLiftBound < degree (F) + 1)
612    {
613      if (d == 1)
614        adaptedLiftBound= tmin (e + 1, deg);
615      else
616        adaptedLiftBound= deg;
617    }
618    factors= T;
619    F= buf;
620    success= true;
621  }
622  return result;
623}
624
625CFList
626extEarlyFactorDetect (CanonicalForm& F, CFList& factors, int& adaptedLiftBound,
627                      bool& success, const ExtensionInfo& info, const CFList&
628                      eval, const int deg, const CFList& MOD, const int bound)
629{
630  Variable alpha= info.getAlpha();
631  Variable beta= info.getBeta();
632  CanonicalForm gamma= info.getGamma();
633  CanonicalForm delta= info.getDelta();
634  int k= info.getGFDegree();
635  CFList result;
636  CFList T= factors;
637  CanonicalForm buf= F;
638  Variable y= F.mvar();
639  Variable x= Variable (1);
640  CanonicalForm LCBuf= LC (buf, x);
641  CanonicalForm g, gg, quot;
642  CFList M= MOD;
643  M.append (power (y, deg));
644  adaptedLiftBound= 0;
645  int d= bound;
646  int e= 0;
647  int nBuf;
648  CFList source, dest;
649
650  int degMipoBeta= 1;
651  if (!k && beta.level() != 1)
652    degMipoBeta= degree (getMipo (beta));
653
654  for (CFListIterator i= factors; i.hasItem(); i++)
655  {
656    g= mulMod (i.getItem(), LCBuf, M);
657    g /= myContent (g);
658    if (fdivides (g, buf, quot))
659    {
660      gg= reverseShift (g, eval);
661      gg /= Lc (gg);
662      if (!k && beta == x)
663      {
664        if (degree (gg, alpha) < degMipoBeta)
665        {
666          appendTestMapDown (result, gg, info, source, dest);
667          buf= quot;
668          nBuf= degree (g, y) + degree (LC (g, x), y);
669          d -= nBuf;
670          e= tmax (e, nBuf);
671          LCBuf= LC (buf, x);
672          T= Difference (T, CFList (i.getItem()));
673        }
674      }
675      else
676      {
677        if (!isInExtension (gg, gamma, k, delta, source, dest))
678        {
679          appendTestMapDown (result, gg, info, source, dest);
680          buf= quot;
681          nBuf= degree (g, y) + degree (LC (g, x), y);
682          d -= nBuf;
683          e= tmax (e, nBuf);
684          LCBuf= LC (buf, x);
685          T= Difference (T, CFList (i.getItem()));
686         }
687      }
688    }
689  }
690  adaptedLiftBound= d;
691
692  if (adaptedLiftBound < deg)
693  {
694    if (adaptedLiftBound < degree (F) + 1)
695    {
696      if (d == 1)
697        adaptedLiftBound= tmin (e + 1, deg);
698      else
699        adaptedLiftBound= deg;
700    }
701    success= true;
702    factors= T;
703    F= buf;
704  }
705  return result;
706}
707
708CFList
709evalPoints (const CanonicalForm& F, CFList & eval, const Variable& alpha,
710            CFList& list, const bool& GF, bool& fail)
711{
712  int k= F.level() - 1;
713  Variable x= Variable (1);
714  CFList result;
715  FFRandom genFF;
716  GFRandom genGF;
717  int p= getCharacteristic ();
718  double bound;
719  if (alpha != x)
720  {
721    bound= pow ((double) p, (double) degree (getMipo(alpha)));
722    bound= pow ((double) bound, (double) k);
723  }
724  else if (GF)
725  {
726    bound= pow ((double) p, (double) getGFDegree());
727    bound= pow ((double) bound, (double) k);
728  }
729  else
730    bound= pow ((double) p, (double) k);
731
732  CanonicalForm random;
733  CanonicalForm deriv_x, gcd_deriv;
734  do
735  {
736    random= 0;
737    // possible overflow if list.length() does not fit into a int
738    if (list.length() >= bound)
739    {
740      fail= true;
741      break;
742    }
743    for (int i= 0; i < k; i++)
744    {
745      if (list.isEmpty())
746        result.append (0);
747      else if (GF)
748      {
749        result.append (genGF.generate());
750        random += result.getLast()*power (x, i);
751      }
752      else if (alpha.level() != 1)
753      {
754        AlgExtRandomF genAlgExt (alpha);
755        result.append (genAlgExt.generate());
756        random += result.getLast()*power (x, i);
757      }
758      else
759      {
760        result.append (genFF.generate());
761        random += result.getLast()*power (x, i);
762      }
763    }
764    if (find (list, random))
765    {
766      result= CFList();
767      continue;
768    }
769    int l= F.level();
770    eval.insert (F);
771    bool bad= false;
772    for (CFListIterator i= result; i.hasItem(); i++, l--)
773    {
774      eval.insert (eval.getFirst()(i.getItem(), l));
775      if (degree (eval.getFirst(), l - 1) != degree (F, l - 1))
776      {
777        if (!find (list, random))
778          list.append (random);
779        result= CFList();
780        eval= CFList();
781        bad= true;
782        break;
783      }
784    }
785
786    if (bad)
787      continue;
788
789    if (degree (eval.getFirst()) != degree (F, 1))
790    {
791      if (!find (list, random))
792        list.append (random);
793      result= CFList();
794      eval= CFList();
795      continue;
796    }
797
798    deriv_x= deriv (eval.getFirst(), x);
799    gcd_deriv= gcd (eval.getFirst(), deriv_x);
800    if (degree (gcd_deriv) > 0)
801    {
802      if (!find (list, random))
803        list.append (random);
804      result= CFList();
805      eval= CFList();
806      continue;
807    }
808    CFListIterator i= eval;
809    i++;
810    CanonicalForm contentx= content (i.getItem(), x);
811    if (degree (contentx) > 0)
812    {
813      if (!find (list, random))
814        list.append (random);
815      result= CFList();
816      eval= CFList();
817      continue;
818    }
819
820    if (list.length() >= bound)
821    {
822      fail= true;
823      break;
824    }
825  } while (find (list, random));
826
827  if (!eval.isEmpty())
828    eval.removeFirst();
829
830  return result;
831}
832
833static inline
834int newMainVariableSearch (CanonicalForm& A, CFList& Aeval, CFList&
835                           evaluation, const Variable& alpha, const int lev,
836                           CanonicalForm& g
837                          )
838{
839  Variable x= Variable (1);
840  CanonicalForm derivI, buf;
841  bool GF= (CFFactory::gettype() == GaloisFieldDomain);
842  int swapLevel= 0;
843  CFList list;
844  bool fail= false;
845  buf= A;
846  Aeval= CFList();
847  evaluation= CFList();
848  for (int i= lev; i <= A.level(); i++)
849  {
850    derivI= deriv (buf, Variable (i));
851    if (!derivI.isZero())
852    {
853      g= gcd (buf, derivI);
854      if (degree (g) > 0)
855        return -1;
856
857      buf= swapvar (buf, x, Variable (i));
858      Aeval= CFList();
859      evaluation= CFList();
860      fail= false;
861      evaluation= evalPoints (buf, Aeval, alpha, list, GF, fail);
862      if (!fail)
863      {
864        A= buf;
865        swapLevel= i;
866        break;
867      }
868      else
869        buf= A;
870    }
871  }
872  return swapLevel;
873}
874
875CanonicalForm lcmContent (const CanonicalForm& A, CFList& contentAi)
876{
877  int i= A.level();
878  contentAi.append (myContent (A, i));
879  contentAi.append (myContent (A, i - 1));
880  CanonicalForm result= lcm (contentAi.getFirst(), contentAi.getLast());
881  for (i= i - 2; i > 0; i--)
882  {
883    contentAi.append (content (A, i));
884    result= lcm (result, contentAi.getLast());
885  }
886  return result;
887}
888
889CFList
890henselLiftAndEarly (CanonicalForm& A, CFList& MOD, int*& liftBounds, bool&
891                    earlySuccess, CFList& earlyFactors, const CFList& Aeval,
892                    const CFList& biFactors, const CFList& evaluation,
893                    const ExtensionInfo& info)
894{
895  bool extension= info.isInExtension();
896  CFList bufFactors= biFactors;
897  bufFactors.insert (LC (Aeval.getFirst(), 1));
898
899  sortList (bufFactors, Variable (1));
900
901  CFList diophant;
902  CFArray Pi;
903  int smallFactorDeg= 11; //tunable parameter
904  CFList result;
905  int adaptedLiftBound= 0;
906  int liftBound= liftBounds[1];
907
908  earlySuccess= false;
909  CFList earlyReconstFactors;
910  CFListIterator j= Aeval;
911  j++;
912  CanonicalForm buf= j.getItem();
913  CFMatrix Mat= CFMatrix (liftBound, bufFactors.length() - 1);
914  MOD= CFList (power (Variable (2), liftBounds[0]));
915  if (smallFactorDeg >= liftBound)
916  {
917    result= henselLift23 (Aeval, bufFactors, liftBounds, diophant, Pi, Mat);
918  }
919  else if (smallFactorDeg >= degree (buf) + 1)
920  {
921    liftBounds[1]= degree (buf) + 1;
922    result= henselLift23 (Aeval, bufFactors, liftBounds, diophant, Pi, Mat);
923    if (Aeval.length() == 2)
924    {
925      if (!extension)
926        earlyFactors= earlyFactorDetect
927                       (buf, result, adaptedLiftBound, earlySuccess,
928                        degree (buf) + 1, MOD, liftBound);
929      else
930        earlyFactors= extEarlyFactorDetect
931                       (buf, result, adaptedLiftBound, earlySuccess,
932                        info, evaluation, degree
933                        (buf) + 1, MOD, liftBound);
934    }
935    else
936    {
937      if (!extension)
938        adaptedLiftBound= liftBoundAdaption (buf, result, earlySuccess,
939                                             degree (buf) + 1, MOD, liftBound);
940      else
941        adaptedLiftBound= extLiftBoundAdaption (buf, result, earlySuccess, info,
942                                                evaluation, degree (buf) + 1,
943                                                MOD, liftBound);
944    }
945    if (!earlySuccess)
946    {
947      result.insert (LC (buf, 1));
948      liftBounds[1]= adaptedLiftBound;
949      liftBound= adaptedLiftBound;
950      henselLiftResume (buf, result, degree (buf) + 1, liftBound,
951                        Pi, diophant, Mat, MOD);
952    }
953    else
954      liftBounds[1]= adaptedLiftBound;
955  }
956  else if (smallFactorDeg < degree (buf) + 1)
957  {
958    liftBounds[1]= smallFactorDeg;
959    result= henselLift23 (Aeval, bufFactors, liftBounds, diophant, Pi, Mat);
960    if (Aeval.length() == 2)
961    {
962      if (!extension)
963        earlyFactors= earlyFactorDetect (buf, result, adaptedLiftBound,
964                                         earlySuccess, smallFactorDeg, MOD,
965                                         liftBound);
966      else
967        earlyFactors= extEarlyFactorDetect (buf, result, adaptedLiftBound,
968                                            earlySuccess, info, evaluation,
969                                            smallFactorDeg, MOD, liftBound);
970    }
971    else
972    {
973      if (!extension)
974        adaptedLiftBound= liftBoundAdaption (buf, result, earlySuccess,
975                                             smallFactorDeg, MOD, liftBound);
976      else
977        adaptedLiftBound= extLiftBoundAdaption (buf, result, earlySuccess, info,
978                                                evaluation, smallFactorDeg, MOD,
979                                                liftBound);
980    }
981
982    if (!earlySuccess)
983    {
984      result.insert (LC (buf, 1));
985      henselLiftResume (buf, result, smallFactorDeg, degree (buf) + 1,
986                        Pi, diophant, Mat, MOD);
987      if (Aeval.length() == 2)
988      {
989         if (!extension)
990           earlyFactors= earlyFactorDetect (buf, result, adaptedLiftBound,
991                                            earlySuccess, degree (buf) + 1,
992                                            MOD, liftBound);
993         else
994           earlyFactors= extEarlyFactorDetect (buf, result, adaptedLiftBound,
995                                               earlySuccess, info, evaluation,
996                                               degree (buf) + 1, MOD,
997                                               liftBound);
998      }
999      else
1000      {
1001        if (!extension)
1002          adaptedLiftBound= liftBoundAdaption (buf, result, earlySuccess,
1003                                               degree (buf) + 1, MOD,liftBound);
1004        else
1005          adaptedLiftBound= extLiftBoundAdaption (buf, result, earlySuccess,
1006                                                  info, evaluation,
1007                                                  degree (buf) + 1, MOD,
1008                                                  liftBound);
1009      }
1010      if (!earlySuccess)
1011      {
1012        result.insert (LC (buf, 1));
1013        liftBounds[1]= adaptedLiftBound;
1014        liftBound= adaptedLiftBound;
1015        henselLiftResume (buf, result, degree (buf) + 1, liftBound,
1016                          Pi, diophant, Mat, MOD);
1017      }
1018      else
1019        liftBounds[1]= adaptedLiftBound;
1020    }
1021    else
1022      liftBounds[1]= adaptedLiftBound;
1023  }
1024
1025  MOD.append (power (Variable (3), liftBounds[1]));
1026
1027  if (Aeval.length() > 2)
1028  {
1029    CFListIterator j= Aeval;
1030    j++;
1031    CFList bufEval;
1032    bufEval.append (j.getItem());
1033    j++;
1034    int liftBoundsLength= Aeval.getLast().level() - 1;
1035    for (int i= 2; i <= liftBoundsLength && j.hasItem(); i++, j++)
1036    {
1037      earlySuccess= false;
1038      result.insert (LC (bufEval.getFirst(), 1));
1039      bufEval.append (j.getItem());
1040      liftBound= liftBounds[i];
1041      Mat= CFMatrix (liftBounds[i], result.length() - 1);
1042
1043      buf= j.getItem();
1044      if (smallFactorDeg >= liftBound)
1045        result= henselLift (bufEval, result, MOD, diophant, Pi, Mat,
1046                            liftBounds[i -  1], liftBounds[i]);
1047      else if (smallFactorDeg >= degree (buf) + 1)
1048      {
1049        result= henselLift (bufEval, result, MOD, diophant, Pi, Mat,
1050                            liftBounds[i -  1], degree (buf) + 1);
1051
1052        if (Aeval.length() == i + 1)
1053        {
1054          if (!extension)
1055            earlyFactors= earlyFactorDetect
1056                           (buf, result, adaptedLiftBound, earlySuccess,
1057                            degree (buf) + 1, MOD, liftBound);
1058          else
1059            earlyFactors= extEarlyFactorDetect
1060                           (buf, result, adaptedLiftBound, earlySuccess,
1061                            info, evaluation, degree (buf) + 1, MOD, liftBound);
1062        }
1063        else
1064        {
1065          if (!extension)
1066            adaptedLiftBound= liftBoundAdaption
1067                                (buf, result, earlySuccess, degree (buf)
1068                                 + 1,  MOD, liftBound);
1069          else
1070            adaptedLiftBound= extLiftBoundAdaption
1071                                (buf, result, earlySuccess, info, evaluation,
1072                                 degree (buf) + 1, MOD, liftBound);
1073        }
1074
1075        if (!earlySuccess)
1076        {
1077          result.insert (LC (buf, 1));
1078          liftBounds[i]= adaptedLiftBound;
1079          liftBound= adaptedLiftBound;
1080          henselLiftResume (buf, result, degree (buf) + 1, liftBound,
1081                            Pi, diophant, Mat, MOD);
1082        }
1083        else
1084        {
1085          liftBounds[i]= adaptedLiftBound;
1086        }
1087      }
1088      else if (smallFactorDeg < degree (buf) + 1)
1089      {
1090        result= henselLift (bufEval, result, MOD, diophant, Pi, Mat,
1091                            liftBounds[i -  1], smallFactorDeg);
1092
1093        if (Aeval.length() == i + 1)
1094        {
1095          if (!extension)
1096            earlyFactors= earlyFactorDetect
1097                           (buf, result, adaptedLiftBound, earlySuccess,
1098                            smallFactorDeg, MOD, liftBound);
1099          else
1100            earlyFactors= extEarlyFactorDetect
1101                           (buf, result, adaptedLiftBound, earlySuccess,
1102                            info, evaluation, smallFactorDeg, MOD, liftBound);
1103        }
1104        else
1105        {
1106          if (!extension)
1107            adaptedLiftBound= liftBoundAdaption
1108                                (buf, result, earlySuccess,
1109                                 smallFactorDeg, MOD, liftBound);
1110          else
1111            adaptedLiftBound= extLiftBoundAdaption
1112                                (buf, result, earlySuccess, info, evaluation,
1113                                 smallFactorDeg, MOD, liftBound);
1114        }
1115
1116        if (!earlySuccess)
1117        {
1118          result.insert (LC (buf, 1));
1119          henselLiftResume (buf, result, smallFactorDeg,
1120                            degree (buf) + 1, Pi, diophant, Mat, MOD);
1121          if (Aeval.length() == i + 1)
1122          {
1123            if (!extension)
1124              earlyFactors= earlyFactorDetect
1125                             (buf, result, adaptedLiftBound, earlySuccess,
1126                              degree (buf) +  1,  MOD, liftBound);
1127            else
1128              earlyFactors= extEarlyFactorDetect
1129                             (buf, result, adaptedLiftBound, earlySuccess,
1130                              info, evaluation, degree (buf) + 1, MOD,
1131                              liftBound);
1132          }
1133          else
1134          {
1135            if (!extension)
1136              adaptedLiftBound= liftBoundAdaption
1137                                  (buf, result, earlySuccess, degree
1138                                   (buf) +  1,  MOD, liftBound);
1139            else
1140              adaptedLiftBound= extLiftBoundAdaption
1141                                  (buf, result, earlySuccess, info, evaluation,
1142                                   degree (buf) + 1,  MOD, liftBound);
1143          }
1144
1145          if (!earlySuccess)
1146          {
1147            result.insert (LC (buf, 1));
1148            liftBounds[i]= adaptedLiftBound;
1149            liftBound= adaptedLiftBound;
1150            henselLiftResume (buf, result, degree (buf) + 1, liftBound,
1151                              Pi, diophant, Mat, MOD);
1152          }
1153          else
1154            liftBounds[i]= adaptedLiftBound;
1155        }
1156        else
1157          liftBounds[i]= adaptedLiftBound;
1158      }
1159      MOD.append (power (Variable (i + 2), liftBounds[i]));
1160      bufEval.removeFirst();
1161    }
1162    bufFactors= result;
1163  }
1164  else
1165    bufFactors= result;
1166
1167  if (earlySuccess)
1168    A= buf;
1169  return result;
1170}
1171
1172void
1173gcdFreeBasis (CFFList& factors1, CFFList& factors2)
1174{
1175  CanonicalForm g;
1176  int k= factors1.length();
1177  int l= factors2.length();
1178  int n= 1;
1179  int m;
1180  CFFListIterator j;
1181  for (CFFListIterator i= factors1; (n < k && i.hasItem()); i++, n++)
1182  {
1183    m= 1;
1184    for (j= factors2; (m < l && j.hasItem()); j++, m++)
1185    {
1186      g= gcd (i.getItem().factor(), j.getItem().factor());
1187      if (degree (g) > 0)
1188      {
1189        j.getItem()= CFFactor (j.getItem().factor()/g, j.getItem().exp());
1190        i.getItem()= CFFactor (i.getItem().factor()/g, i.getItem().exp());
1191        factors1.append (CFFactor (g, i.getItem().exp()));
1192        factors2.append (CFFactor (g, j.getItem().exp()));
1193      }
1194    }
1195  }
1196}
1197
1198CFList
1199distributeContent (const CFList& L, const CFList* differentSecondVarFactors,
1200                   int length
1201                  )
1202{
1203  CFList l= L;
1204  CanonicalForm content= l.getFirst();
1205
1206  if (content.inCoeffDomain())
1207    return l;
1208
1209  if (l.length() == 1)
1210  {
1211    CFList result;
1212    for (int i= 0; i < length; i++)
1213    {
1214      if (differentSecondVarFactors[i].isEmpty())
1215        continue;
1216      if (result.isEmpty())
1217      {
1218        result= differentSecondVarFactors[i];
1219        for (CFListIterator iter= result; iter.hasItem(); iter++)
1220          content /= iter.getItem();
1221      }
1222      else
1223      {
1224        CFListIterator iter1= result;
1225        for (CFListIterator iter2= differentSecondVarFactors[i]; iter2.hasItem();
1226             iter2++, iter1++)
1227        {
1228          iter1.getItem() *= iter2.getItem();
1229          content /= iter2.getItem();
1230        }
1231      }
1232    }
1233    result.insert (content);
1234    return result;
1235  }
1236
1237  Variable v;
1238  CFListIterator iter1;
1239  CanonicalForm tmp, g;
1240  for (int i= 0; i < length; i++)
1241  {
1242    if (differentSecondVarFactors[i].isEmpty())
1243      continue;
1244    iter1= l;
1245    iter1++;
1246
1247    v= Variable (i + 3);
1248    for (CFListIterator iter2= differentSecondVarFactors[i]; iter2.hasItem();
1249         iter2++, iter1++)
1250    {
1251      if (degree (iter2.getItem(),v) == degree (iter1.getItem(),v))
1252        continue;
1253      tmp= iter1.getItem();
1254      for (int j= tmp.level(); j > 1; j--)
1255      {
1256        if (j == i + 3)
1257          continue;
1258        tmp= tmp (0, j);
1259      }
1260      g= gcd (iter2.getItem(), content);
1261      if (degree (g) > 0)
1262      {
1263        if (!tmp.isZero())
1264          iter2.getItem() /= tmp;
1265        content /= g;
1266        iter1.getItem() *= g;
1267      }
1268    }
1269  }
1270
1271  l.removeFirst();
1272  l.insert (content);
1273  return l;
1274}
1275
1276CFList evaluateAtZero (const CanonicalForm& F)
1277{
1278  CFList result;
1279  CanonicalForm buf= F;
1280  result.insert (buf);
1281  for (int i= F.level(); i > 2; i--)
1282  {
1283    buf= buf (0, i);
1284    result.insert (buf);
1285  }
1286  return result;
1287}
1288
1289CFList evaluateAtEval (const CanonicalForm& F, const CFArray& eval)
1290{
1291  CFList result;
1292  CanonicalForm buf= F;
1293  result.insert (buf);
1294  int k= eval.size();
1295  for (int i= 1; i < k; i++)
1296  {
1297    buf= buf (eval[i], i + 2);
1298    result.insert (buf);
1299  }
1300  return result;
1301}
1302
1303CFList evaluateAtEval (const CanonicalForm& F, const CFList& evaluation, int l)
1304{
1305  CFList result;
1306  CanonicalForm buf= F;
1307  result.insert (buf);
1308  int k= evaluation.length() + l - 1;
1309  CFListIterator j= evaluation;
1310  for (int i= k; j.hasItem() && i > l; i--, j++)
1311  {
1312    if (F.level() < i)
1313      continue;
1314    buf= buf (j.getItem(), i);
1315    result.insert (buf);
1316  }
1317  return result;
1318}
1319
1320int
1321testFactors (const CanonicalForm& G, const CFList& uniFactors,
1322             const Variable& alpha, CanonicalForm& sqrfPartF, CFList& factors,
1323             CFFList*& bufSqrfFactors, CFList& evalSqrfPartF,
1324             const CFArray& evalPoint)
1325{
1326  CanonicalForm tmp;
1327  CFListIterator j;
1328  for (CFListIterator i= uniFactors; i.hasItem(); i++)
1329  {
1330    tmp= i.getItem();
1331    if (i.hasItem())
1332      i++;
1333    else
1334      break;
1335    for (j= i; j.hasItem(); j++)
1336    {
1337      if (tmp == j.getItem())
1338        return 0;
1339    }
1340  }
1341
1342  CanonicalForm F= G;
1343  CFFList sqrfFactorization= squarefreeFactorization (F, alpha);
1344
1345  sqrfPartF= 1;
1346  for (CFFListIterator i= sqrfFactorization; i.hasItem(); i++)
1347    sqrfPartF *= i.getItem().factor();
1348
1349  evalSqrfPartF= evaluateAtEval (sqrfPartF, evalPoint);
1350
1351  CanonicalForm test= evalSqrfPartF.getFirst() (evalPoint[0], 2);
1352
1353  if (degree (test) != degree (sqrfPartF, 1))
1354    return 0;
1355
1356  CFFList sqrfFactors;
1357  CFList tmp2;
1358  int k= 0;
1359  factors= uniFactors;
1360  CFFListIterator iter;
1361  for (CFListIterator i= factors; i.hasItem(); i++, k++)
1362  {
1363    tmp= 1;
1364    sqrfFactors= squarefreeFactorization (i.getItem(), alpha);
1365
1366    for (iter= sqrfFactors; iter.hasItem(); iter++)
1367    {
1368      tmp2.append (iter.getItem().factor());
1369      tmp *= iter.getItem().factor();
1370    }
1371    i.getItem()= tmp/Lc(tmp);
1372    bufSqrfFactors [k]= sqrfFactors;
1373  }
1374
1375  for (int i= 0; i < factors.length() - 1; i++)
1376  {
1377    for (k= i + 1; k < factors.length(); k++)
1378    {
1379      gcdFreeBasis (bufSqrfFactors [i], bufSqrfFactors[k]);
1380    }
1381  }
1382
1383  factors= CFList();
1384  for (int i= 0; i < uniFactors.length(); i++)
1385  {
1386    if (i == 0)
1387    {
1388      for (iter= bufSqrfFactors [i]; iter.hasItem(); iter++)
1389      {
1390        if (iter.getItem().factor().inCoeffDomain())
1391          continue;
1392        iter.getItem()= CFFactor (iter.getItem().factor()/
1393                                  Lc (iter.getItem().factor()),
1394                                  iter.getItem().exp());
1395        factors.append (iter.getItem().factor());
1396      }
1397    }
1398    else
1399    {
1400      for (iter= bufSqrfFactors [i]; iter.hasItem(); iter++)
1401      {
1402        if (iter.getItem().factor().inCoeffDomain())
1403          continue;
1404        iter.getItem()= CFFactor (iter.getItem().factor()/
1405                                  Lc (iter.getItem().factor()),
1406                                  iter.getItem().exp());
1407        if (!find (factors, iter.getItem().factor()))
1408          factors.append (iter.getItem().factor());
1409      }
1410    }
1411  }
1412
1413  test= prod (factors);
1414  tmp= evalSqrfPartF.getFirst() (evalPoint[0],2);
1415  if (test/Lc (test) != tmp/Lc (tmp))
1416    return 0;
1417  else
1418    return 1;
1419}
1420
1421CFList
1422precomputeLeadingCoeff (const CanonicalForm& LCF, const CFList& LCFFactors,
1423                        const Variable& alpha, const CFList& evaluation,
1424                        CFList* & differentSecondVarLCs, int lSecondVarLCs,
1425                        Variable& y
1426                       )
1427{
1428  y= Variable (1);
1429  if (LCF.inCoeffDomain())
1430  {
1431    CFList result;
1432    for (int i= 1; i <= LCFFactors.length() + 1; i++)
1433      result.append (1);
1434    return result;
1435  }
1436
1437  CFMap N;
1438  CanonicalForm F= compress (LCF, N);
1439  if (LCF.isUnivariate())
1440  {
1441    CFList result;
1442    int LCFLevel= LCF.level();
1443    bool found= false;
1444    if (LCFLevel == 2)
1445    {
1446    //bivariate leading coefficients are already the true leading coefficients
1447      result= LCFFactors;
1448      found= true;
1449    }
1450    else
1451    {
1452      CFListIterator j;
1453      for (int i= 0; i < lSecondVarLCs; i++)
1454      {
1455        for (j= differentSecondVarLCs[i]; j.hasItem(); j++)
1456        {
1457          if (j.getItem().level() == LCFLevel)
1458          {
1459            found= true;
1460            break;
1461          }
1462        }
1463        if (found)
1464        {
1465          result= differentSecondVarLCs [i];
1466          break;
1467        }
1468      }
1469      if (!found)
1470        result= LCFFactors;
1471    }
1472    if (found)
1473      result.insert (Lc (LCF));
1474    else
1475      result.append (LCF);
1476    return result;
1477  }
1478
1479  CFList factors= LCFFactors;
1480
1481  CFMap dummy;
1482  for (CFListIterator i= factors; i.hasItem(); i++)
1483    i.getItem()= compress (i.getItem(), dummy);
1484
1485  CanonicalForm sqrfPartF;
1486  CFFList * bufSqrfFactors= new CFFList [factors.length()];
1487  CFList evalSqrfPartF, bufFactors;
1488  CFArray evalPoint= CFArray (evaluation.length() - 1);
1489  CFListIterator iter= evaluation;
1490  for (int i= evaluation.length() - 2; i > -1; i--, iter++)
1491    evalPoint[i]= iter.getItem();
1492
1493  int pass= testFactors (F, factors, alpha, sqrfPartF,
1494                         bufFactors, bufSqrfFactors, evalSqrfPartF, evalPoint);
1495
1496  bool foundDifferent= false;
1497  Variable z, x= y;
1498  int j= 0;
1499  if (!pass)
1500  {
1501    int lev= 0;
1502    // LCF is non-constant here
1503    for (int i= 1; i <= LCF.level(); i++)
1504    {
1505      if(degree (LCF, i) > 0)
1506      {
1507        lev= i - 1;
1508        break;
1509      }
1510    }
1511    CFList bufBufFactors;
1512    CanonicalForm bufF, swap;
1513    CFArray buf;
1514    for (int i= 0; i < lSecondVarLCs; i++)
1515    {
1516      if (!differentSecondVarLCs [i].isEmpty())
1517      {
1518        bool allConstant= true;
1519        for (iter= differentSecondVarLCs[i]; iter.hasItem(); iter++)
1520        {
1521          if (!iter.getItem().inCoeffDomain())
1522          {
1523            allConstant= false;
1524            y= Variable (iter.getItem().level());
1525          }
1526        }
1527        if (allConstant)
1528          continue;
1529
1530        bufFactors= differentSecondVarLCs [i];
1531        for (iter= bufFactors; iter.hasItem(); iter++)
1532          iter.getItem()= swapvar (iter.getItem(), x, y);
1533        bufF= F;
1534        z= Variable (y.level() - lev);
1535        bufF= swapvar (bufF, x, z);
1536        bufBufFactors= bufFactors;
1537        evalPoint= CFArray (evaluation.length() - 1);
1538        buf= CFArray (evaluation.length());
1539        iter= evaluation;
1540        int k= evaluation.length() - 1;
1541        for (; iter.hasItem(); iter++, k--)
1542          buf[k]= iter.getItem();
1543        swap= buf[z.level() - 1];
1544        buf[z.level() - 1]= buf[0];
1545        buf[0]= 0;
1546        int l= 0;
1547        for (k= 0; k < evaluation.length(); k++)
1548        {
1549          if (buf[k].isZero())
1550            continue;
1551          evalPoint[l]= buf[k];
1552          l++;
1553        }
1554        pass= testFactors (bufF, bufBufFactors, alpha, sqrfPartF, bufFactors,
1555                           bufSqrfFactors, evalSqrfPartF, evalPoint);
1556        if (pass)
1557        {
1558          foundDifferent= true;
1559          F= bufF;
1560          CFList l= factors;
1561          for (iter= l; iter.hasItem(); iter++)
1562            iter.getItem()= swapvar (iter.getItem(), x, y);
1563          differentSecondVarLCs [i]= l;
1564          j= i;
1565          break;
1566        }
1567        if (!pass && i == lSecondVarLCs - 1)
1568        {
1569          CFList result;
1570          result.append (LCF);
1571          for (int k= 1; k <= factors.length(); k++)
1572            result.append (LCF);
1573          y= Variable (1);
1574          delete [] bufSqrfFactors;
1575          return result;
1576        }
1577      }
1578    }
1579  }
1580  if (!pass)
1581  {
1582    CFList result;
1583    result.append (LCF);
1584    for (int k= 1; k <= factors.length(); k++)
1585      result.append (LCF);
1586    y= Variable (1);
1587    delete [] bufSqrfFactors;
1588    return result;
1589  }
1590  else
1591    factors= bufFactors;
1592
1593  bufFactors= factors;
1594  CFList evaluation2;
1595  if (y == x)
1596    evaluation2= evaluation;
1597  else
1598  {
1599    CanonicalForm tmp;
1600    evaluation2= evaluation;
1601    int i= evaluation.length() + 1;
1602    for (CFListIterator iter= evaluation2; iter.hasItem(); iter++, i--)
1603    {
1604      if (i == y.level())
1605      {
1606        tmp= iter.getItem();
1607        iter.getItem()= evaluation2.getLast();
1608        evaluation2.removeLast();
1609        evaluation2.append (tmp);
1610        break;
1611      }
1612    }
1613  }
1614
1615  CFList interMedResult;
1616  CanonicalForm oldSqrfPartF= sqrfPartF;
1617  sqrfPartF= shift2Zero (sqrfPartF, evalSqrfPartF, evaluation2, 1);
1618  if (factors.length() > 1)
1619  {
1620    CanonicalForm LC1= LC (oldSqrfPartF, 1);
1621    CFList leadingCoeffs;
1622    for (int i= 0; i < factors.length(); i++)
1623      leadingCoeffs.append (LC1);
1624
1625    CFList LC1eval= evaluateAtEval (LC1, evaluation2, 1);
1626    CFList oldFactors= factors;
1627    for (CFListIterator i= oldFactors; i.hasItem(); i++)
1628      i.getItem() *= LC1eval.getFirst()/Lc (i.getItem());
1629
1630    bool success= false;
1631    CanonicalForm oldSqrfPartFPowLC= oldSqrfPartF*power(LC1,factors.length()-1);
1632    if (size (oldSqrfPartFPowLC)/getNumVars (oldSqrfPartFPowLC) < 500 &&
1633        LucksWangSparseHeuristic (oldSqrfPartFPowLC,
1634                                  oldFactors, 1, leadingCoeffs, factors))
1635    {
1636      interMedResult= recoverFactors (oldSqrfPartF, factors);
1637      if (oldFactors.length() == interMedResult.length())
1638        success= true;
1639    }
1640    if (!success)
1641    {
1642      LC1= LC (evalSqrfPartF.getFirst(), 1);
1643
1644      CFArray leadingCoeffs= CFArray (factors.length());
1645      for (int i= 0; i < factors.length(); i++)
1646        leadingCoeffs[i]= LC1;
1647
1648      for (CFListIterator i= factors; i.hasItem(); i++)
1649      {
1650        i.getItem()= i.getItem() (x + evaluation2.getLast(), x);
1651        i.getItem() *= LC1 (0,2)/Lc (i.getItem());
1652      }
1653      factors.insert (1);
1654
1655      CanonicalForm
1656      newSqrfPartF= evalSqrfPartF.getFirst()*power (LC1, factors.length() - 2);
1657
1658      int liftBound= degree (newSqrfPartF,2) + 1;
1659
1660      CFMatrix M= CFMatrix (liftBound, factors.length() - 1);
1661      CFArray Pi;
1662      CFList diophant;
1663      nonMonicHenselLift12 (newSqrfPartF, factors, liftBound, Pi, diophant, M,
1664                            leadingCoeffs, false);
1665
1666      if (sqrfPartF.level() > 2)
1667      {
1668        int* liftBounds= new int [sqrfPartF.level() - 1];
1669        liftBounds [0]= liftBound;
1670        bool noOneToOne= false;
1671        CFList *leadingCoeffs2= new CFList [sqrfPartF.level()-2];
1672        LC1= LC (evalSqrfPartF.getLast(), 1);
1673        CFList LCs;
1674        for (int i= 0; i < factors.length(); i++)
1675          LCs.append (LC1);
1676        leadingCoeffs2 [sqrfPartF.level() - 3]= LCs;
1677        for (int i= sqrfPartF.level() - 1; i > 2; i--)
1678        {
1679          for (CFListIterator j= LCs; j.hasItem(); j++)
1680            j.getItem()= j.getItem() (0, i + 1);
1681          leadingCoeffs2 [i - 3]= LCs;
1682        }
1683        sqrfPartF *= power (LC1, factors.length()-1);
1684
1685        int liftBoundsLength= sqrfPartF.level() - 1;
1686        for (int i= 1; i < liftBoundsLength; i++)
1687          liftBounds [i]= degree (sqrfPartF, i + 2) + 1;
1688        evalSqrfPartF= evaluateAtZero (sqrfPartF);
1689        evalSqrfPartF.removeFirst();
1690        factors= nonMonicHenselLift (evalSqrfPartF, factors, leadingCoeffs2,
1691                 diophant, Pi, liftBounds, sqrfPartF.level() - 1, noOneToOne);
1692        delete [] leadingCoeffs2;
1693        delete [] liftBounds;
1694      }
1695      for (CFListIterator iter= factors; iter.hasItem(); iter++)
1696        iter.getItem()= reverseShift (iter.getItem(), evaluation2, 1);
1697
1698      interMedResult=
1699      recoverFactors (reverseShift(evalSqrfPartF.getLast(),evaluation2,1),
1700                      factors);
1701    }
1702  }
1703  else
1704  {
1705    factors= CFList (oldSqrfPartF);
1706    interMedResult= recoverFactors (oldSqrfPartF, factors);
1707  }
1708
1709  CFList result;
1710  CFFListIterator k;
1711  for (int i= 0; i < LCFFactors.length(); i++)
1712  {
1713    CanonicalForm tmp= 1;
1714    for (k= bufSqrfFactors[i]; k.hasItem(); k++)
1715    {
1716      int pos= findItem (bufFactors, k.getItem().factor());
1717      if (pos)
1718        tmp *= power (getItem (interMedResult, pos), k.getItem().exp());
1719    }
1720    result.append (tmp);
1721  }
1722
1723  for (CFListIterator i= result; i.hasItem(); i++)
1724  {
1725    F /= i.getItem();
1726    if (foundDifferent)
1727      i.getItem()= swapvar (i.getItem(), x, z);
1728    i.getItem()= N (i.getItem());
1729  }
1730
1731  if (foundDifferent)
1732  {
1733    CFList l= differentSecondVarLCs [j];
1734    for (CFListIterator i= l; i.hasItem(); i++)
1735      i.getItem()= swapvar (i.getItem(), y, z);
1736    differentSecondVarLCs [j]= l;
1737    F= swapvar (F, x, z);
1738  }
1739
1740  result.insert (N (F));
1741
1742  result= distributeContent (result, differentSecondVarLCs, lSecondVarLCs);
1743
1744  if (!result.getFirst().inCoeffDomain())
1745  {
1746    CFListIterator i= result;
1747    CanonicalForm tmp;
1748    if (foundDifferent)
1749      i.getItem()= swapvar (i.getItem(), Variable (2), y);
1750
1751    tmp= i.getItem();
1752
1753    i++;
1754    for (; i.hasItem(); i++)
1755    {
1756      if (foundDifferent)
1757        i.getItem()= swapvar (i.getItem(), Variable (2), y)*tmp;
1758      else
1759        i.getItem() *= tmp;
1760    }
1761  }
1762  else
1763    y= Variable (1);
1764
1765  delete [] bufSqrfFactors;
1766
1767  return result;
1768}
1769
1770void
1771evaluationWRTDifferentSecondVars (CFList*& Aeval, const CFList& evaluation,
1772                                  const CanonicalForm& A)
1773{
1774  CanonicalForm tmp;
1775  CFList tmp2;
1776  CFListIterator iter;
1777  for (int i= A.level(); i > 2; i--)
1778  {
1779    tmp= A;
1780    tmp2= CFList();
1781    iter= evaluation;
1782    bool preserveDegree= true;
1783    for (int j= A.level(); j > 1; j--, iter++)
1784    {
1785      if (j == i)
1786        continue;
1787      else
1788      {
1789        tmp= tmp (iter.getItem(), j);
1790        tmp2.insert (tmp);
1791        if ((degree (tmp, i) != degree (A, i)) ||
1792            (degree (tmp, 1) != degree (A, 1)))
1793        {
1794          preserveDegree= false;
1795          break;
1796        }
1797        if (!content(tmp).inCoeffDomain() || !content(tmp,1).inCoeffDomain())
1798        {
1799          preserveDegree= false;
1800          break;
1801        }
1802      }
1803    }
1804    if (preserveDegree)
1805      Aeval [i - 3]= tmp2;
1806    else
1807      Aeval [i - 3]= CFList();
1808  }
1809}
1810
1811#endif
1812
1813static inline
1814CanonicalForm prodEval (const CFList& l, const CanonicalForm& evalPoint,
1815                        const Variable& v)
1816{
1817  CanonicalForm result= 1;
1818  for (CFListIterator i= l; i.hasItem(); i++)
1819    result *= i.getItem() (evalPoint, v);
1820  return result;
1821}
1822
1823//recombine bivariate factors in case one bivariate factorization yields less
1824// factors than the other
1825CFList
1826recombination (const CFList& factors1, const CFList& factors2, int s, int thres,
1827               const CanonicalForm& evalPoint, const Variable& x)
1828{
1829  CFList T, S;
1830
1831  T= factors1;
1832  CFList result;
1833  CanonicalForm buf;
1834  int * v= new int [T.length()];
1835  for (int i= 0; i < T.length(); i++)
1836    v[i]= 0;
1837  bool nosubset= false;
1838  CFArray TT;
1839  TT= copy (factors1);
1840  while (T.length() >= 2*s && s <= thres)
1841  {
1842    while (nosubset == false) 
1843    {
1844      if (T.length() == s) 
1845      {
1846        delete [] v;
1847        result.append (prod (T));
1848        return result;
1849      }
1850      S= subset (v, s, TT, nosubset);
1851      if (nosubset) break;
1852      buf= prodEval (S, evalPoint, x);
1853      buf /= Lc (buf);
1854      if (find (factors2, buf))
1855      {
1856        T= Difference (T, S);
1857        result.append (prod (S));
1858        TT= copy (T);
1859        indexUpdate (v, s, T.length(), nosubset);
1860        if (nosubset) break;
1861      }
1862    }
1863    s++;
1864    if (T.length() < 2*s || T.length() == s) 
1865    {
1866      delete [] v;
1867      result.append (prod (T));
1868      return result;
1869    }
1870    for (int i= 0; i < T.length(); i++)
1871      v[i]= 0;
1872    nosubset= false;
1873  }
1874
1875  delete [] v;
1876  if (T.length() < 2*s)
1877  {
1878    result.append (prod (T));
1879    return result;
1880  }
1881
1882  return result;
1883}
1884
1885#ifdef HAVE_NTL
1886void
1887factorizationWRTDifferentSecondVars (const CanonicalForm& A, CFList*& Aeval,
1888                                     const ExtensionInfo& info,
1889                                     int& minFactorsLength, bool& irred)
1890{
1891  Variable x= Variable (1);
1892  minFactorsLength= 0;
1893  irred= false;
1894  CFList factors;
1895  Variable v;
1896  for (int j= 0; j < A.level() - 2; j++)
1897  {
1898    if (!Aeval[j].isEmpty())
1899    {
1900      v= Variable (Aeval[j].getFirst().level());
1901      if (CFFactory::gettype() == GaloisFieldDomain)
1902        factors= GFBiSqrfFactorize (Aeval[j].getFirst());
1903      else if (info.getAlpha().level() == 1)
1904        factors= FpBiSqrfFactorize (Aeval[j].getFirst());
1905      else
1906        factors= FqBiSqrfFactorize (Aeval[j].getFirst(), info.getAlpha());
1907
1908      factors.removeFirst();
1909      if (minFactorsLength == 0)
1910        minFactorsLength= factors.length();
1911      else
1912        minFactorsLength= tmin (minFactorsLength, factors.length());
1913
1914      if (factors.length() == 1)
1915      {
1916        irred= true;
1917        return;
1918      }
1919      sortList (factors, x);
1920      Aeval [j]= factors;
1921    }
1922  }
1923}
1924
1925void getLeadingCoeffs (const CanonicalForm& A, CFList*& Aeval,
1926                       const CFList& uniFactors, const CFList& evaluation
1927                      )
1928{
1929  CanonicalForm evalPoint;
1930  int i;
1931  CFListIterator iter, iter2;
1932  Variable v;
1933  CFList l, LCs, buf;
1934  int pos;
1935  for (int j= 0; j < A.level() - 2; j++)
1936  {
1937    if (!Aeval[j].isEmpty())
1938    {
1939      i= A.level();
1940      for (iter= evaluation; iter.hasItem(); iter++, i--)
1941      {
1942        if (i == Aeval[j].getFirst().level())
1943        {
1944          evalPoint= iter.getItem();
1945          break;
1946        }
1947      }
1948
1949      v= Variable (i);
1950      if (Aeval[j].length() > uniFactors.length())
1951        Aeval[j]= recombination (Aeval[j], uniFactors, 1,
1952                                 Aeval[j].length() - uniFactors.length() + 1,
1953                                 evalPoint, v);
1954
1955      l= CFList();
1956      buf= buildUniFactors (Aeval[j], evalPoint, v);
1957      for (iter= uniFactors; iter.hasItem(); iter++)
1958      {
1959        pos= findItem (buf, iter.getItem());
1960        if (pos)
1961          l.append (getItem (Aeval[j], pos));
1962      }
1963      Aeval [j]= l;
1964
1965      LCs= CFList();
1966      for (iter= Aeval[j]; iter.hasItem(); iter++)
1967        LCs.append (LC (iter.getItem(), 1));
1968      normalize (LCs);
1969      Aeval[j]= LCs;
1970    }
1971  }
1972}
1973
1974CFList
1975buildUniFactors (const CFList& biFactors, const CanonicalForm& evalPoint,
1976                 const Variable& y)
1977{
1978  CFList result;
1979  CanonicalForm tmp;
1980  for (CFListIterator i= biFactors; i.hasItem(); i++)
1981  {
1982    tmp= mod (i.getItem(), y - evalPoint);
1983    tmp /= Lc (tmp);
1984    result.append (tmp);
1985  }
1986  return result;
1987}
1988
1989void refineBiFactors (const CanonicalForm& A, CFList& biFactors,
1990                      CFList* const& Aeval, const CFList& evaluation,
1991                      int minFactorsLength)
1992{
1993  CFListIterator iter;
1994  CanonicalForm evalPoint;
1995  int i;
1996  Variable v;
1997  Variable y= Variable (2);
1998  CFList list;
1999  for (int j= 0; j < A.level() - 2; j++)
2000  {
2001    if (Aeval[j].length() == minFactorsLength)
2002    {
2003      i= A.level();
2004
2005      for (iter= evaluation; iter.hasItem(); iter++, i--)
2006      {
2007        if (i == Aeval[j].getFirst().level())
2008        {
2009          evalPoint= iter.getItem();
2010          break;
2011        }
2012      }
2013
2014      v= Variable (i);
2015      list= buildUniFactors (Aeval[j], evalPoint, v);
2016
2017      biFactors= recombination (biFactors, list, 1,
2018                                biFactors.length() - list.length() + 1,
2019                                evaluation.getLast(), y);
2020      return;
2021    }
2022  }
2023}
2024
2025void prepareLeadingCoeffs (CFList*& LCs, int n, const CFList& leadingCoeffs,
2026                           const CFList& biFactors, const CFList& evaluation)
2027{
2028  CFList l= leadingCoeffs;
2029  LCs [n-3]= l;
2030  CFListIterator j;
2031  CFListIterator iter= evaluation;
2032  for (int i= n - 1; i > 2; i--, iter++)
2033  {
2034    for (j= l; j.hasItem(); j++)
2035      j.getItem()= j.getItem() (iter.getItem(), i + 1);
2036    LCs [i - 3]= l;
2037  }
2038  l= LCs [0];
2039  for (CFListIterator i= l; i.hasItem(); i++)
2040    i.getItem()= i.getItem() (iter.getItem(), 3);
2041  CFListIterator ii= biFactors;
2042  CFList normalizeFactor;
2043  for (CFListIterator i= l; i.hasItem(); i++, ii++)
2044    normalizeFactor.append (Lc (LC (ii.getItem(), 1))/Lc (i.getItem()));
2045  for (int i= 0; i < n-2; i++)
2046  {
2047    ii= normalizeFactor;
2048    for (j= LCs [i]; j.hasItem(); j++, ii++)
2049      j.getItem() *= ii.getItem();
2050  }
2051}
2052
2053CFList recoverFactors (const CanonicalForm& F, const CFList& factors)
2054{
2055  CFList result;
2056  CanonicalForm tmp, tmp2;
2057  CanonicalForm G= F;
2058  for (CFListIterator i= factors; i.hasItem(); i++)
2059  {
2060    tmp= i.getItem()/content (i.getItem(), 1);
2061    if (fdivides (tmp, G, tmp2))
2062    {
2063      G= tmp2;
2064      result.append (tmp);
2065    }
2066  }
2067  return result;
2068}
2069
2070CFList recoverFactors (const CanonicalForm& F, const CFList& factors,
2071                       const CFList& evaluation)
2072{
2073  CFList result;
2074  CanonicalForm tmp, tmp2;
2075  CanonicalForm G= F;
2076  for (CFListIterator i= factors; i.hasItem(); i++)
2077  {
2078    tmp= reverseShift (i.getItem(), evaluation);
2079    tmp /= content (tmp, 1);
2080    if (fdivides (tmp, G, tmp2))
2081    {
2082      G= tmp2;
2083      result.append (tmp);
2084    }
2085  }
2086  return result;
2087}
2088
2089CFList
2090extNonMonicFactorRecombination (const CFList& factors, const CanonicalForm& F,
2091                                const ExtensionInfo& info)
2092{
2093  Variable alpha= info.getAlpha();
2094  Variable beta= info.getBeta();
2095  CanonicalForm gamma= info.getGamma();
2096  CanonicalForm delta= info.getDelta();
2097  int k= info.getGFDegree();
2098  CFList source, dest;
2099
2100  int degMipoBeta= 1;
2101  if (!k && beta != Variable(1))
2102    degMipoBeta= degree (getMipo (beta));
2103
2104  CFList T, S;
2105  T= factors;
2106  int s= 1;
2107  CFList result;
2108  CanonicalForm quot, buf= F;
2109
2110  CanonicalForm g;
2111  CanonicalForm buf2;
2112  int * v= new int [T.length()];
2113  for (int i= 0; i < T.length(); i++)
2114    v[i]= 0;
2115  bool noSubset= false;
2116  CFArray TT;
2117  TT= copy (factors);
2118  bool recombination= false;
2119  bool trueFactor= false;
2120  while (T.length() >= 2*s)
2121  {
2122    while (noSubset == false)
2123    {
2124      if (T.length() == s)
2125      {
2126        delete [] v;
2127        if (recombination)
2128        {
2129          g= prod (T);
2130          T.removeFirst();
2131          result.append (g/myContent (g));
2132          g /= Lc (g);
2133          appendTestMapDown (result, g, info, source, dest);
2134          return result;
2135        }
2136        else
2137          return CFList (buf/myContent(buf));
2138      }
2139
2140      S= subset (v, s, TT, noSubset);
2141      if (noSubset) break;
2142
2143      g= prod (S);
2144      g /= myContent (g);
2145      if (fdivides (g, buf, quot))
2146      {
2147        buf2= g;
2148        buf2 /= Lc (buf2);
2149        if (!k && beta.level() == 1)
2150        {
2151          if (degree (buf2, alpha) < degMipoBeta)
2152          {
2153            appendTestMapDown (result, buf2, info, source, dest);
2154            buf= quot;
2155            recombination= true;
2156            trueFactor= true;
2157          }
2158        }
2159        else
2160        {
2161          if (!isInExtension (buf2, gamma, k, delta, source, dest))
2162          {
2163            appendTestMapDown (result, buf2, info, source, dest);
2164            buf= quot;
2165            recombination= true;
2166            trueFactor= true;
2167          }
2168        }
2169        if (trueFactor)
2170        {
2171          T= Difference (T, S);
2172
2173          if (T.length() < 2*s || T.length() == s)
2174          {
2175            delete [] v;
2176            buf /= myContent (buf);
2177            buf /= Lc (buf);
2178            appendTestMapDown (result, buf, info, source, dest);
2179            return result;
2180          }
2181          trueFactor= false;
2182          TT= copy (T);
2183          indexUpdate (v, s, T.length(), noSubset);
2184          if (noSubset) break;
2185        }
2186      }
2187    }
2188    s++;
2189    if (T.length() < 2*s || T.length() == s)
2190    {
2191      delete [] v;
2192      appendTestMapDown (result, buf/myContent(buf), info, source, dest);
2193      return result;
2194    }
2195    for (int i= 0; i < T.length(); i++)
2196      v[i]= 0;
2197    noSubset= false;
2198  }
2199  if (T.length() < 2*s)
2200    appendMapDown (result, F/myContent(F), info, source, dest);
2201
2202  delete [] v;
2203  return result;
2204}
2205
2206CFList
2207extFactorize (const CanonicalForm& F, const ExtensionInfo& info);
2208
2209CFList
2210multiFactorize (const CanonicalForm& F, const ExtensionInfo& info)
2211{
2212
2213  if (F.inCoeffDomain())
2214    return CFList (F);
2215
2216  // compress and find main Variable
2217  CFMap N;
2218  CanonicalForm A= myCompress (F, N);
2219
2220  A /= Lc (A); // make monic
2221
2222  Variable alpha= info.getAlpha();
2223  Variable beta= info.getBeta();
2224  CanonicalForm gamma= info.getGamma();
2225  CanonicalForm delta= info.getDelta();
2226  bool extension= info.isInExtension();
2227  bool GF= (CFFactory::gettype() == GaloisFieldDomain);
2228  //univariate case
2229  if (F.isUnivariate())
2230  {
2231    if (extension == false)
2232      return uniFactorizer (F, alpha, GF);
2233    else
2234    {
2235      CFList source, dest;
2236      A= mapDown (F, info, source, dest);
2237      return uniFactorizer (A, beta, GF);
2238    }
2239  }
2240
2241  //bivariate case
2242  if (A.level() == 2)
2243  {
2244    CFList buf= biFactorize (F, info);
2245    return buf;
2246  }
2247
2248  Variable x= Variable (1);
2249  Variable y= Variable (2);
2250
2251  // remove content
2252  CFList contentAi;
2253  CanonicalForm lcmCont= lcmContent (A, contentAi);
2254  A /= lcmCont;
2255
2256  // trivial after content removal
2257  CFList contentAFactors;
2258  if (A.inCoeffDomain())
2259  {
2260    for (CFListIterator i= contentAi; i.hasItem(); i++)
2261    {
2262      if (i.getItem().inCoeffDomain())
2263        continue;
2264      else
2265      {
2266        lcmCont /= i.getItem();
2267        contentAFactors=
2268        Union (multiFactorize (lcmCont, info),
2269               multiFactorize (i.getItem(), info));
2270        break;
2271      }
2272    }
2273    decompress (contentAFactors, N);
2274    normalize (contentAFactors);
2275    return contentAFactors;
2276  }
2277
2278  // factorize content
2279  contentAFactors= multiFactorize (lcmCont, info);
2280
2281  // univariate after content removal
2282  CFList factors;
2283  if (A.isUnivariate ())
2284  {
2285    factors= uniFactorizer (A, alpha, GF);
2286    append (factors, contentAFactors);
2287    decompress (factors, N);
2288    return factors;
2289  }
2290
2291  // check main variable
2292  int swapLevel= 0;
2293  CanonicalForm derivZ;
2294  CanonicalForm gcdDerivZ;
2295  CanonicalForm bufA= A;
2296  Variable z;
2297  for (int i= 1; i <= A.level(); i++)
2298  {
2299    z= Variable (i);
2300    derivZ= deriv (bufA, z);
2301    if (derivZ.isZero())
2302    {
2303      if (i == 1)
2304        swapLevel= 1;
2305      else
2306        continue;
2307    }
2308    else
2309    {
2310      if (swapLevel == 1)
2311      {
2312        swapLevel= i;
2313        bufA= swapvar (A, x, z);
2314      }
2315      gcdDerivZ= gcd (bufA, derivZ);
2316      if (degree (gcdDerivZ) > 0 && !derivZ.isZero())
2317      {
2318        CanonicalForm g= bufA/gcdDerivZ;
2319        CFList factorsG=
2320        Union (multiFactorize (g, info),
2321               multiFactorize (gcdDerivZ, info));
2322        appendSwapDecompress (factorsG, contentAFactors, N, swapLevel, x);
2323        normalize (factorsG);
2324        return factorsG;
2325      }
2326      else
2327      {
2328        A= bufA;
2329        break;
2330      }
2331    }
2332  }
2333
2334
2335  CFList Aeval, list, evaluation, bufEvaluation, bufAeval;
2336  bool fail= false;
2337  int swapLevel2= 0;
2338  int level;
2339  int factorNums= 3;
2340  CanonicalForm bivarEval;
2341  CFList biFactors, bufBiFactors;
2342  CanonicalForm evalPoly;
2343  int lift, bufLift;
2344  double logarithm= (double) ilog2 (totaldegree (A));
2345  logarithm /= log2exp;
2346  logarithm= ceil (logarithm);
2347  if (factorNums < (int) logarithm)
2348    factorNums= (int) logarithm;
2349  CFList* bufAeval2= new CFList [A.level() - 2];
2350  CFList* Aeval2= new CFList [A.level() - 2];
2351  int counter;
2352  int differentSecondVar= 0;
2353  // several bivariate factorizations
2354  for (int i= 0; i < factorNums; i++)
2355  {
2356    counter= 0;
2357    bufA= A;
2358    bufAeval= CFList();
2359    bufEvaluation= evalPoints (bufA, bufAeval, alpha, list, GF, fail);
2360    evalPoly= 0;
2361
2362    if (fail && (i == 0))
2363    {
2364      if (!swapLevel)
2365        level= 2;
2366      else
2367        level= swapLevel + 1;
2368
2369      CanonicalForm g;
2370      swapLevel2= newMainVariableSearch (A, Aeval, evaluation, alpha, level, g);
2371
2372      if (!swapLevel2) // need to pass to an extension
2373      {
2374        factors= extFactorize (A, info);
2375        appendSwapDecompress (factors, contentAFactors, N, swapLevel, x);
2376        normalize (factors);
2377        delete [] bufAeval2;
2378        delete [] Aeval2;
2379        return factors;
2380      }
2381      else
2382      {
2383        if (swapLevel2 == -1)
2384        {
2385          CFList factorsG=
2386          Union (multiFactorize (g, info),
2387                 multiFactorize (A/g, info));
2388          appendSwapDecompress (factorsG, contentAFactors, N, swapLevel, x);
2389          normalize (factorsG);
2390          delete [] bufAeval2;
2391          delete [] Aeval2;
2392          return factorsG;
2393        }
2394        fail= false;
2395        bufAeval= Aeval;
2396        bufA= A;
2397        bufEvaluation= evaluation;
2398      }
2399    }
2400    else if (fail && (i > 0))
2401      break;
2402
2403    bivarEval= bufEvaluation.getLast();
2404
2405    evaluationWRTDifferentSecondVars (bufAeval2, bufEvaluation, A);
2406
2407    for (int j= 0; j < A.level() - 2; j++)
2408    {
2409      if (!bufAeval2[j].isEmpty())
2410        counter++;
2411    }
2412
2413    bufLift= degree (A, y) + 1 + degree (LC(A, x), y);
2414
2415    TIMING_START (fac_bi_factorizer);
2416    if (!GF && alpha.level() == 1)
2417      bufBiFactors= FpBiSqrfFactorize (bufAeval.getFirst());
2418    else if (GF)
2419      bufBiFactors= GFBiSqrfFactorize (bufAeval.getFirst());
2420    else
2421      bufBiFactors= FqBiSqrfFactorize (bufAeval.getFirst(), alpha);
2422    TIMING_END_AND_PRINT (fac_bi_factorizer,
2423                          "time for bivariate factorization: ");
2424    bufBiFactors.removeFirst();
2425
2426    if (bufBiFactors.length() == 1)
2427    {
2428      if (extension)
2429      {
2430        CFList source, dest;
2431        A= mapDown (A, info, source, dest);
2432      }
2433      factors.append (A);
2434      appendSwapDecompress (factors, contentAFactors, N, swapLevel,
2435                            swapLevel2, x);
2436      normalize (factors);
2437      delete [] bufAeval2;
2438      delete [] Aeval2;
2439      return factors;
2440    }
2441
2442    if (i == 0)
2443    {
2444      Aeval= bufAeval;
2445      evaluation= bufEvaluation;
2446      biFactors= bufBiFactors;
2447      lift= bufLift;
2448      for (int j= 0; j < A.level() - 2; j++)
2449        Aeval2 [j]= bufAeval2 [j];
2450      differentSecondVar= counter;
2451    }
2452    else
2453    {
2454      if (bufBiFactors.length() < biFactors.length() ||
2455          ((bufLift < lift) && (bufBiFactors.length() == biFactors.length())) ||
2456          counter > differentSecondVar)
2457      {
2458        Aeval= bufAeval;
2459        evaluation= bufEvaluation;
2460        biFactors= bufBiFactors;
2461        lift= bufLift;
2462        for (int j= 0; j < A.level() - 2; j++)
2463          Aeval2 [j]= bufAeval2 [j];
2464        differentSecondVar= counter;
2465      }
2466    }
2467    int k= 0;
2468    for (CFListIterator j= bufEvaluation; j.hasItem(); j++, k++)
2469      evalPoly += j.getItem()*power (x, k);
2470    list.append (evalPoly);
2471  }
2472
2473  delete [] bufAeval2;
2474
2475  sortList (biFactors, x);
2476
2477  int minFactorsLength;
2478  bool irred= false;
2479  factorizationWRTDifferentSecondVars (A, Aeval2, info, minFactorsLength, irred);
2480
2481  if (irred)
2482  {
2483    if (extension)
2484    {
2485      CFList source, dest;
2486      A= mapDown (A, info, source, dest);
2487    }
2488    factors.append (A);
2489    appendSwapDecompress (factors, contentAFactors, N, swapLevel,
2490                          swapLevel2, x);
2491    normalize (factors);
2492    delete [] Aeval2;
2493    return factors;
2494  }
2495
2496  if (minFactorsLength == 0)
2497    minFactorsLength= biFactors.length();
2498  else if (biFactors.length() > minFactorsLength)
2499    refineBiFactors (A, biFactors, Aeval2, evaluation, minFactorsLength);
2500  minFactorsLength= tmin (minFactorsLength, biFactors.length());
2501
2502  if (differentSecondVar == A.level() - 2)
2503  {
2504    bool zeroOccured= false;
2505    for (CFListIterator iter= evaluation; iter.hasItem(); iter++)
2506    {
2507      if (iter.getItem().isZero())
2508      {
2509        zeroOccured= true;
2510        break;
2511      }
2512    }
2513    if (!zeroOccured)
2514    {
2515      factors= sparseHeuristic (A, biFactors, Aeval2, evaluation, minFactorsLength);
2516      if (factors.length() == biFactors.length())
2517      {
2518        if (extension)
2519          factors= extNonMonicFactorRecombination (factors, A, info);
2520
2521        appendSwapDecompress (factors, contentAFactors, N, swapLevel,
2522                              swapLevel2, x);
2523        normalize (factors);
2524        delete [] Aeval2;
2525        return factors;
2526      }
2527      else
2528        factors= CFList();
2529      //TODO case where factors.length() > 0
2530    }
2531  }
2532
2533  CFList uniFactors= buildUniFactors (biFactors, evaluation.getLast(), y);
2534
2535  CFList * oldAeval= new CFList [A.level() - 2]; //TODO use bufAeval2 for this
2536  for (int i= 0; i < A.level() - 2; i++)
2537    oldAeval[i]= Aeval2[i];
2538
2539  getLeadingCoeffs (A, Aeval2, uniFactors, evaluation);
2540
2541  CFList biFactorsLCs;
2542  for (CFListIterator i= biFactors; i.hasItem(); i++)
2543    biFactorsLCs.append (LC (i.getItem(), 1));
2544
2545  Variable v;
2546  CFList leadingCoeffs= precomputeLeadingCoeff (LC (A, 1), biFactorsLCs, alpha,
2547                                          evaluation, Aeval2, A.level() - 2, v);
2548
2549  if (v.level() != 1)
2550  {
2551    A= swapvar (A, y, v);
2552    for (int i= 0; i < A.level() - 2; i++)
2553    {
2554      if (oldAeval[i].isEmpty())
2555        continue;
2556      if (oldAeval[i].getFirst().level() == v.level())
2557      {
2558        biFactors= CFList();
2559        for (CFListIterator iter= oldAeval [i]; iter.hasItem(); iter++)
2560          biFactors.append (swapvar (iter.getItem(), v, y));
2561      }
2562    }
2563    int i= A.level();
2564    CanonicalForm evalPoint;
2565    for (CFListIterator iter= evaluation; iter.hasItem(); iter++, i--)
2566    {
2567      if (i == v.level())
2568      {
2569        evalPoint= iter.getItem();
2570        iter.getItem()= evaluation.getLast();
2571        evaluation.removeLast();
2572        evaluation.append (evalPoint);
2573        break;
2574      }
2575    }
2576  }
2577
2578  CFListIterator iter;
2579  CanonicalForm oldA= A;
2580  CFList oldBiFactors= biFactors;
2581  if (!leadingCoeffs.getFirst().inCoeffDomain())
2582  {
2583    CanonicalForm tmp= power (leadingCoeffs.getFirst(), biFactors.length() - 1);
2584    A *= tmp;
2585    tmp= leadingCoeffs.getFirst();
2586    iter= evaluation;
2587    for (int i= A.level(); i > 2; i--, iter++)
2588      tmp= tmp (iter.getItem(), i);
2589    if (!tmp.inCoeffDomain())
2590    {
2591      for (CFListIterator i= biFactors; i.hasItem(); i++)
2592      {
2593        i.getItem() *= tmp/LC (i.getItem(), 1);
2594        i.getItem() /= Lc (i.getItem());
2595      }
2596    }
2597  }
2598
2599  leadingCoeffs.removeFirst();
2600
2601  //prepare leading coefficients
2602  CFList* leadingCoeffs2= new CFList [A.level() - 2];
2603  prepareLeadingCoeffs (leadingCoeffs2, A.level(), leadingCoeffs, biFactors,
2604                        evaluation);
2605
2606  Aeval= evaluateAtEval (A, evaluation, 2);
2607  CanonicalForm hh= Lc (Aeval.getFirst());
2608  for (iter= Aeval; iter.hasItem(); iter++)
2609    iter.getItem() /= hh;
2610
2611  A /= hh;
2612
2613  if (size (A)/getNumVars (A) < 500 && 
2614      LucksWangSparseHeuristic (A, biFactors, 2, leadingCoeffs2 [A.level() - 3],
2615                                factors))
2616  {
2617    int check= factors.length();
2618    factors= recoverFactors (A, factors);
2619
2620    if (check == factors.length())
2621    {
2622      if (extension)
2623        factors= extNonMonicFactorRecombination (factors, A, info);
2624
2625      if (v.level() != 1)
2626      {
2627        for (CFListIterator iter= factors; iter.hasItem(); iter++)
2628          iter.getItem()= swapvar (iter.getItem(), v, y);
2629      }
2630      appendSwapDecompress (factors, contentAFactors, N, swapLevel,
2631                            swapLevel2, x);
2632      normalize (factors);
2633      delete [] Aeval2;
2634      return factors;
2635    }
2636    else
2637      factors= CFList();
2638    //TODO handle this case
2639  }
2640
2641  A= shift2Zero (A, Aeval, evaluation);
2642
2643  for (iter= biFactors; iter.hasItem(); iter++)
2644    iter.getItem()= iter.getItem () (y + evaluation.getLast(), y);
2645
2646  for (int i= 0; i < A.level() - 2; i++)
2647  {
2648    if (i != A.level() - 3)
2649      leadingCoeffs2[i]= CFList();
2650  }
2651  for (iter= leadingCoeffs2[A.level() - 3]; iter.hasItem(); iter++)
2652  {
2653    iter.getItem()= shift2Zero (iter.getItem(), list, evaluation);
2654    for (int i= A.level() - 4; i > -1; i--)
2655    {
2656      if (i + 1 == A.level() - 3)
2657        leadingCoeffs2[i].append (iter.getItem() (0, i + 4));
2658      else
2659        leadingCoeffs2[i].append (leadingCoeffs2[i+1].getLast() (0, i + 4));
2660    }
2661  }
2662
2663  CFArray Pi;
2664  CFList diophant;
2665  int* liftBounds= new int [A.level() - 1];
2666  int liftBoundsLength= A.level() - 1;
2667  for (int i= 0; i < liftBoundsLength; i++)
2668    liftBounds [i]= degree (A, i + 2) + 1;
2669
2670  Aeval.removeFirst();
2671  bool noOneToOne= false;
2672  factors= nonMonicHenselLift (Aeval, biFactors, leadingCoeffs2, diophant,
2673                               Pi, liftBounds, liftBoundsLength, noOneToOne);
2674
2675  if (!noOneToOne)
2676  {
2677    int check= factors.length();
2678    A= oldA;
2679    factors= recoverFactors (A, factors, evaluation);
2680    if (check != factors.length())
2681      noOneToOne= true;
2682
2683    if (extension && !noOneToOne)
2684      factors= extNonMonicFactorRecombination (factors, A, info);
2685  }
2686  if (noOneToOne)
2687  {
2688    A= shift2Zero (oldA, Aeval, evaluation);
2689    biFactors= oldBiFactors;
2690    for (iter= biFactors; iter.hasItem(); iter++)
2691      iter.getItem()= iter.getItem () (y + evaluation.getLast(), y);
2692    CanonicalForm LCA= LC (Aeval.getFirst(), 1);
2693    CanonicalForm yToLift= power (y, lift);
2694    CFListIterator i= biFactors;
2695    lift= degree (i.getItem(), 2) + degree (LC (i.getItem(), 1)) + 1;
2696    i++;
2697
2698    for (; i.hasItem(); i++)
2699      lift= tmax (lift, degree (i.getItem(), 2) + degree (LC (i.getItem(), 1)) + 1);
2700
2701    lift= tmax (degree (Aeval.getFirst() , 2) + 1, lift);
2702
2703    i= biFactors;
2704    yToLift= power (y, lift);
2705    CanonicalForm dummy;
2706    for (; i.hasItem(); i++)
2707    {
2708      LCA= LC (i.getItem(), 1);
2709      extgcd (LCA, yToLift, LCA, dummy);
2710      i.getItem()= mod (i.getItem()*LCA, yToLift);
2711    }
2712
2713    liftBoundsLength= F.level() - 1;
2714    liftBounds= liftingBounds (A, lift);
2715
2716    CFList MOD;
2717    bool earlySuccess;
2718    CFList earlyFactors, liftedFactors;
2719    TIMING_START (fac_hensel_lift);
2720    liftedFactors= henselLiftAndEarly
2721                   (A, MOD, liftBounds, earlySuccess, earlyFactors,
2722                    Aeval, biFactors, evaluation, info);
2723    TIMING_END_AND_PRINT (fac_hensel_lift, "time for hensel lifting: ");
2724
2725    if (!extension)
2726    {
2727      TIMING_START (fac_factor_recombination);
2728      factors= factorRecombination (A, liftedFactors, MOD);
2729      TIMING_END_AND_PRINT (fac_factor_recombination,
2730                            "time for factor recombination: ");
2731    }
2732    else
2733    {
2734      TIMING_START (fac_factor_recombination);
2735      factors= extFactorRecombination (liftedFactors, A, MOD, info, evaluation);
2736      TIMING_END_AND_PRINT (fac_factor_recombination,
2737                            "time for factor recombination: ");
2738    }
2739
2740    if (earlySuccess)
2741      factors= Union (factors, earlyFactors);
2742    if (!extension)
2743    {
2744      for (CFListIterator i= factors; i.hasItem(); i++)
2745      {
2746        int kk= Aeval.getLast().level();
2747        for (CFListIterator j= evaluation; j.hasItem(); j++, kk--)
2748        {
2749          if (i.getItem().level() < kk)
2750            continue;
2751          i.getItem()= i.getItem() (Variable (kk) - j.getItem(), kk);
2752        }
2753      }
2754    }
2755  }
2756
2757  if (v.level() != 1)
2758  {
2759    for (CFListIterator iter= factors; iter.hasItem(); iter++)
2760      iter.getItem()= swapvar (iter.getItem(), v, y);
2761  }
2762
2763  swap (factors, swapLevel, swapLevel2, x);
2764  append (factors, contentAFactors);
2765  decompress (factors, N);
2766  normalize (factors);
2767
2768  delete[] liftBounds;
2769
2770  return factors;
2771}
2772
2773/// multivariate factorization over an extension of the initial field
2774CFList
2775extFactorize (const CanonicalForm& F, const ExtensionInfo& info)
2776{
2777  CanonicalForm A= F;
2778
2779  Variable alpha= info.getAlpha();
2780  Variable beta= info.getBeta();
2781  int k= info.getGFDegree();
2782  char cGFName= info.getGFName();
2783  CanonicalForm delta= info.getDelta();
2784  bool GF= (CFFactory::gettype() == GaloisFieldDomain);
2785  Variable w= Variable (1);
2786
2787  CFList factors;
2788  if (!GF && alpha == w)  // we are in F_p
2789  {
2790    CFList factors;
2791    bool extension= true;
2792    int p= getCharacteristic();
2793    if (p*p < (1<<16)) // pass to GF if possible
2794    {
2795      setCharacteristic (getCharacteristic(), 2, 'Z');
2796      ExtensionInfo info= ExtensionInfo (extension);
2797      A= A.mapinto();
2798      factors= multiFactorize (A, info);
2799
2800      CanonicalForm mipo= gf_mipo;
2801      setCharacteristic (getCharacteristic());
2802      Variable vBuf= rootOf (mipo.mapinto());
2803      for (CFListIterator j= factors; j.hasItem(); j++)
2804        j.getItem()= GF2FalphaRep (j.getItem(), vBuf);
2805    }
2806    else  // not able to pass to GF, pass to F_p(\alpha)
2807    {
2808      CanonicalForm mipo= randomIrredpoly (2, w);
2809      Variable v= rootOf (mipo);
2810      ExtensionInfo info= ExtensionInfo (v);
2811      factors= multiFactorize (A, info);
2812    }
2813    return factors;
2814  }
2815  else if (!GF && (alpha != w)) // we are in F_p(\alpha)
2816  {
2817    if (k == 1) // need factorization over F_p
2818    {
2819      int extDeg= degree (getMipo (alpha));
2820      extDeg++;
2821      CanonicalForm mipo= randomIrredpoly (extDeg + 1, w);
2822      Variable v= rootOf (mipo);
2823      ExtensionInfo info= ExtensionInfo (v);
2824      factors= multiFactorize (A, info);
2825    }
2826    else
2827    {
2828      if (beta == w)
2829      {
2830        Variable v= chooseExtension (alpha, beta, k);
2831        CanonicalForm primElem, imPrimElem;
2832        bool primFail= false;
2833        Variable vBuf;
2834        primElem= primitiveElement (alpha, vBuf, primFail);
2835        ASSERT (!primFail, "failure in integer factorizer");
2836        if (primFail)
2837          ; //ERROR
2838        else
2839          imPrimElem= mapPrimElem (primElem, vBuf, v);
2840
2841        CFList source, dest;
2842        CanonicalForm bufA= mapUp (A, alpha, v, primElem, imPrimElem,
2843                                   source, dest);
2844        ExtensionInfo info= ExtensionInfo (v, alpha, imPrimElem, primElem);
2845        factors= multiFactorize (bufA, info);
2846      }
2847      else
2848      {
2849        Variable v= chooseExtension (alpha, beta, k);
2850        CanonicalForm primElem, imPrimElem;
2851        bool primFail= false;
2852        Variable vBuf;
2853        ASSERT (!primFail, "failure in integer factorizer");
2854        if (primFail)
2855          ; //ERROR
2856        else
2857          imPrimElem= mapPrimElem (delta, beta, v);
2858
2859        CFList source, dest;
2860        CanonicalForm bufA= mapDown (A, info, source, dest);
2861        source= CFList();
2862        dest= CFList();
2863        bufA= mapUp (bufA, beta, v, delta, imPrimElem, source, dest);
2864        ExtensionInfo info= ExtensionInfo (v, beta, imPrimElem, delta);
2865        factors= multiFactorize (bufA, info);
2866      }
2867    }
2868    return factors;
2869  }
2870  else // we are in GF (p^k)
2871  {
2872    int p= getCharacteristic();
2873    int extensionDeg= getGFDegree();
2874    bool extension= true;
2875    if (k == 1) // need factorization over F_p
2876    {
2877      extensionDeg++;
2878      if (pow ((double) p, (double) extensionDeg) < (1<<16))
2879      // pass to GF(p^k+1)
2880      {
2881        CanonicalForm mipo= gf_mipo;
2882        setCharacteristic (p);
2883        Variable vBuf= rootOf (mipo.mapinto());
2884        A= GF2FalphaRep (A, vBuf);
2885        setCharacteristic (p, extensionDeg, 'Z');
2886        ExtensionInfo info= ExtensionInfo (extension);
2887        factors= multiFactorize (A.mapinto(), info);
2888      }
2889      else // not able to pass to another GF, pass to F_p(\alpha)
2890      {
2891        CanonicalForm mipo= gf_mipo;
2892        setCharacteristic (p);
2893        Variable vBuf= rootOf (mipo.mapinto());
2894        A= GF2FalphaRep (A, vBuf);
2895        Variable v= chooseExtension (vBuf, beta, k);
2896        ExtensionInfo info= ExtensionInfo (v, extension);
2897        factors= multiFactorize (A, info);
2898      }
2899    }
2900    else // need factorization over GF (p^k)
2901    {
2902      if (pow ((double) p, (double) 2*extensionDeg) < (1<<16))
2903      // pass to GF(p^2k)
2904      {
2905        setCharacteristic (p, 2*extensionDeg, 'Z');
2906        ExtensionInfo info= ExtensionInfo (k, cGFName, extension);
2907        factors= multiFactorize (GFMapUp (A, extensionDeg), info);
2908        setCharacteristic (p, extensionDeg, cGFName);
2909      }
2910      else // not able to pass to GF (p^2k), pass to F_p (\alpha)
2911      {
2912        CanonicalForm mipo= gf_mipo;
2913        setCharacteristic (p);
2914        Variable v1= rootOf (mipo.mapinto());
2915        A= GF2FalphaRep (A, v1);
2916        Variable v2= chooseExtension (v1, v1, k);
2917        CanonicalForm primElem, imPrimElem;
2918        bool primFail= false;
2919        Variable vBuf;
2920        primElem= primitiveElement (v1, v1, primFail);
2921        if (primFail)
2922          ; //ERROR
2923        else
2924          imPrimElem= mapPrimElem (primElem, v1, v2);
2925        CFList source, dest;
2926        CanonicalForm bufA= mapUp (A, v1, v2, primElem, imPrimElem,
2927                                     source, dest);
2928        ExtensionInfo info= ExtensionInfo (v2, v1, imPrimElem, primElem);
2929        factors= multiFactorize (bufA, info);
2930        setCharacteristic (p, k, cGFName);
2931        for (CFListIterator i= factors; i.hasItem(); i++)
2932          i.getItem()= Falpha2GFRep (i.getItem());
2933      }
2934    }
2935    return factors;
2936  }
2937}
2938
2939#endif
2940/* HAVE_NTL */
2941
Note: See TracBrowser for help on using the repository browser.