source: git/factory/facFactorize.cc @ 86faff

spielwiese
Last change on this file since 86faff was 86faff, checked in by Martin Lee <martinlee84@…>, 11 years ago
chg: more checks in heurstic to distribute leading coefficients
  • Property mode set to 100644
File size: 37.3 KB
RevLine 
[ab547e2]1/*****************************************************************************\
2 * Computer Algebra System SINGULAR
3\*****************************************************************************/
[b4fef8]4/** @file facFactorize.cc
[ab547e2]5 *
6 * multivariate factorization over Q(a)
7 *
8 * @author Martin Lee
9 *
10 **/
11/*****************************************************************************/
12
[e4fe2b]13#include "config.h"
[ab547e2]14
15#include "assert.h"
16#include "debug.h"
17#include "timing.h"
18
[d990001]19#include "cf_algorithm.h"
[ab547e2]20#include "facFqFactorizeUtil.h"
21#include "facFactorize.h"
22#include "facFqFactorize.h"
23#include "cf_random.h"
24#include "facHensel.h"
25#include "cf_gcd_smallp.h"
26#include "cf_map_ext.h"
27#include "algext.h"
28#include "cf_reval.h"
[91788c0]29#include "facSparseHensel.h"
[ab547e2]30
[0851b0]31TIMING_DEFINE_PRINT(fac_bi_factorizer)
[ec664a2]32TIMING_DEFINE_PRINT(fac_hensel_lift)
33TIMING_DEFINE_PRINT(fac_factor_recombination)
[0851b0]34TIMING_DEFINE_PRINT(fac_shift_to_zero)
35TIMING_DEFINE_PRINT(fac_precompute_leadcoeff)
36TIMING_DEFINE_PRINT(fac_evaluation)
37TIMING_DEFINE_PRINT(fac_recover_factors)
38TIMING_DEFINE_PRINT(fac_preprocess_and_content)
39TIMING_DEFINE_PRINT(fac_bifactor_total)
40TIMING_DEFINE_PRINT(fac_luckswang)
41TIMING_DEFINE_PRINT(fac_lcheuristic)
42TIMING_DEFINE_PRINT(fac_cleardenom)
43TIMING_DEFINE_PRINT(fac_content)
44TIMING_DEFINE_PRINT(fac_compress)
[0b618a7]45
[d990001]46#ifdef HAVE_NTL
[ab547e2]47CFList evalPoints (const CanonicalForm& F, CFList& eval, Evaluation& E)
48{
49  CFList result;
50  Variable x= Variable (1);
51
52  bool found= false;
[c8fff0]53  bool allZero= true;
54  bool foundZero= false;
[ab547e2]55  CanonicalForm deriv_x, gcd_deriv;
[38ffb7]56  CFListIterator iter;
[ab547e2]57  do
58  {
59    eval.insert (F);
60    bool bad= false;
61    for (int i= E.max(); i >= E.min(); i--)
62    {
63      eval.insert (eval.getFirst()( E [i], i));
64      result.append (E[i]);
[c8fff0]65      if (!E[i].isZero())
66        allZero= false;
67      else
68        foundZero= true;
69      if (!allZero && foundZero)
70      {
71        result= CFList();
72        eval= CFList();
73        bad= true;
74        foundZero= false;
75        break;
76      }
[ab547e2]77      if (degree (eval.getFirst(), i - 1) != degree (F, i - 1))
78      {
79        result= CFList();
80        eval= CFList();
81        bad= true;
82        break;
83      }
84    }
85
86    if (bad)
87    {
88      E.nextpoint();
89      continue;
90    }
91
92    if (degree (eval.getFirst()) != degree (F, 1))
93    {
94      result= CFList();
95      eval= CFList();
96      E.nextpoint();
97      continue;
98    }
99
100    deriv_x= deriv (eval.getFirst(), x);
101    gcd_deriv= gcd (eval.getFirst(), deriv_x);
102    if (degree (gcd_deriv) > 0)
103    {
104      result= CFList();
105      eval= CFList();
106      E.nextpoint();
107      continue;
108    }
[38ffb7]109    iter= eval;
110    iter++;
111    CanonicalForm contentx= content (iter.getItem(), x);
[ab547e2]112    if (degree (contentx) > 0)
113    {
114      result= CFList();
115      eval= CFList();
116      E.nextpoint();
117      continue;
118    }
119    found= true;
120  }
121  while (!found);
122
123  if (!eval.isEmpty())
124    eval.removeFirst();
125  return result;
126}
127
128void
129factorizationWRTDifferentSecondVars (const CanonicalForm& A, CFList*& Aeval,
130                                     int& minFactorsLength, bool& irred,
131                                     const Variable& w)
132{
133  Variable x= Variable (1);
134  minFactorsLength= 0;
135  irred= false;
[38ffb7]136  Variable v;
137  CFList factors;
[afb6a6]138  CanonicalForm LCA= LC (A,1);
[c5230d1]139  for (int j= 0; j < A.level() - 2; j++)
[ab547e2]140  {
[c5230d1]141    if (!Aeval[j].isEmpty())
[ab547e2]142    {
[c5230d1]143      v= Variable (Aeval[j].getFirst().level());
[ab547e2]144
[c5230d1]145      factors= ratBiSqrfFactorize (Aeval[j].getFirst(), w);
146      if (factors.getFirst().inCoeffDomain())
147        factors.removeFirst();
[ab547e2]148
[c5230d1]149      if (minFactorsLength == 0)
150        minFactorsLength= factors.length();
151      else
152        minFactorsLength= tmin (minFactorsLength, factors.length());
[ab547e2]153
[c5230d1]154      if (factors.length() == 1)
[ab547e2]155      {
[c5230d1]156        irred= true;
157        return;
[ab547e2]158      }
[c5230d1]159      sortList (factors, x);
160      Aeval [j]= factors;
[ab547e2]161    }
[afb6a6]162  }
[ab547e2]163}
164
165CFList
166multiFactorize (const CanonicalForm& F, const Variable& v)
167{
168  if (F.inCoeffDomain())
169    return CFList (F);
170
[0851b0]171  TIMING_START (fac_preprocess_and_content);
[ab547e2]172  // compress and find main Variable
173  CFMap N;
[0851b0]174  TIMING_START (fac_compress)
[ab547e2]175  CanonicalForm A= myCompress (F, N);
[0851b0]176  TIMING_END_AND_PRINT (fac_compress, "time to compress poly over Q: ")
[ab547e2]177
178  //univariate case
179  if (F.isUnivariate())
180  {
[0e6668]181    CFList result;
182    if (v.level() != 1)
183      result= conv (factorize (F, v));
184    else
185      result= conv (factorize (F,true));
[ab547e2]186    if (result.getFirst().inCoeffDomain())
187      result.removeFirst();
188    return result;
189  }
190
191  //bivariate case
192  if (A.level() == 2)
193  {
194    CFList buf= biFactorize (F, v);
195    return buf;
196  }
197
198  Variable x= Variable (1);
199  Variable y= Variable (2);
200
201  // remove content
[0851b0]202  TIMING_START (fac_content);
[ab547e2]203  CFList contentAi;
204  CanonicalForm lcmCont= lcmContent (A, contentAi);
205  A /= lcmCont;
[0851b0]206  TIMING_END_AND_PRINT (fac_content, "time to extract content over Q: ");
[ab547e2]207
208  // trivial after content removal
209  CFList contentAFactors;
210  if (A.inCoeffDomain())
211  {
212    for (CFListIterator i= contentAi; i.hasItem(); i++)
213    {
214      if (i.getItem().inCoeffDomain())
215        continue;
216      else
217      {
218        lcmCont /= i.getItem();
219        contentAFactors=
220        Union (multiFactorize (lcmCont, v),
221               multiFactorize (i.getItem(), v));
222        break;
223      }
224    }
225    decompress (contentAFactors, N);
226    if (isOn (SW_RATIONAL))
227      normalize (contentAFactors);
228    return contentAFactors;
229  }
230
231  // factorize content
[0851b0]232  TIMING_START (fac_content);
[ab547e2]233  contentAFactors= multiFactorize (lcmCont, v);
[0851b0]234  TIMING_END_AND_PRINT (fac_content, "time to factor content over Q: ");
[ab547e2]235
236  // univariate after content removal
237  CFList factors;
238  if (A.isUnivariate ())
239  {
[0e6668]240    if (v.level() != 1)
241      factors= conv (factorize (A, v));
242    else
243      factors= conv (factorize (A,true));
[ab547e2]244    append (factors, contentAFactors);
245    decompress (factors, N);
246    return factors;
247  }
248
[f55f39]249  A *= bCommonDen (A);
[ab547e2]250  CFList Aeval, list, evaluation, bufEvaluation, bufAeval;
[bf4c6e]251  int factorNums= 1;
[ab547e2]252  CFList biFactors, bufBiFactors;
253  CanonicalForm evalPoly;
[5dad7c7]254  int lift, bufLift, lengthAeval2= A.level()-2;
255  CFList* bufAeval2= new CFList [lengthAeval2];
256  CFList* Aeval2= new CFList [lengthAeval2];
[ab547e2]257  int counter;
258  int differentSecondVar= 0;
259  CanonicalForm bufA;
[0851b0]260  TIMING_END_AND_PRINT (fac_preprocess_and_content,
261                       "time to preprocess poly and extract content over Q: ");
262
[ab547e2]263  // several bivariate factorizations
[0851b0]264  TIMING_START (fac_bifactor_total);
[ab547e2]265  REvaluation E (2, A.level(), IntRandom (25));
266  for (int i= 0; i < factorNums; i++)
267  {
268    counter= 0;
269    bufA= A;
270    bufAeval= CFList();
[0851b0]271    TIMING_START (fac_evaluation);
[ab547e2]272    bufEvaluation= evalPoints (bufA, bufAeval, E);
[0851b0]273    TIMING_END_AND_PRINT (fac_evaluation,
274                          "time to find evaluation point over Q: ");
[ab547e2]275    E.nextpoint();
276    evalPoly= 0;
277
[0851b0]278    TIMING_START (fac_evaluation);
[ab547e2]279    evaluationWRTDifferentSecondVars (bufAeval2, bufEvaluation, A);
[0851b0]280    TIMING_END_AND_PRINT (fac_evaluation,
281                          "time to eval wrt diff second vars over Q: ");
[ab547e2]282
[c5230d1]283    for (int j= 0; j < lengthAeval2; j++)
[ab547e2]284    {
285      if (!bufAeval2[j].isEmpty())
286        counter++;
287    }
288
289    bufLift= degree (A, y) + 1 + degree (LC(A, x), y);
290
[0851b0]291    TIMING_START (fac_bi_factorizer);
[ab547e2]292    bufBiFactors= ratBiSqrfFactorize (bufAeval.getFirst(), v);
[0851b0]293    TIMING_END_AND_PRINT (fac_bi_factorizer,
[ab547e2]294                          "time for bivariate factorization: ");
295    bufBiFactors.removeFirst();
296
297    if (bufBiFactors.length() == 1)
298    {
299      factors.append (A);
300      appendSwapDecompress (factors, contentAFactors, N, 0, 0, x);
301      if (isOn (SW_RATIONAL))
302        normalize (factors);
303      delete [] bufAeval2;
304      delete [] Aeval2;
305      return factors;
306    }
307
308    if (i == 0)
309    {
310      Aeval= bufAeval;
311      evaluation= bufEvaluation;
312      biFactors= bufBiFactors;
313      lift= bufLift;
[5dad7c7]314      for (int j= 0; j < lengthAeval2; j++)
[ab547e2]315        Aeval2 [j]= bufAeval2 [j];
316      differentSecondVar= counter;
317    }
318    else
319    {
320      if (bufBiFactors.length() < biFactors.length() ||
321          ((bufLift < lift) && (bufBiFactors.length() == biFactors.length())) ||
322          counter > differentSecondVar)
323      {
324        Aeval= bufAeval;
325        evaluation= bufEvaluation;
326        biFactors= bufBiFactors;
327        lift= bufLift;
[5dad7c7]328        for (int j= 0; j < lengthAeval2; j++)
[ab547e2]329          Aeval2 [j]= bufAeval2 [j];
330        differentSecondVar= counter;
331      }
332    }
333    int k= 0;
334    for (CFListIterator j= bufEvaluation; j.hasItem(); j++, k++)
335      evalPoly += j.getItem()*power (x, k);
336    list.append (evalPoly);
337  }
338
339  delete [] bufAeval2;
340
341  sortList (biFactors, x);
342
343  int minFactorsLength;
344  bool irred= false;
[0851b0]345  TIMING_START (fac_bi_factorizer);
[ab547e2]346  factorizationWRTDifferentSecondVars (A, Aeval2, minFactorsLength, irred, v);
[0851b0]347  TIMING_END_AND_PRINT (fac_bi_factorizer,
348             "time for bivariate factorization wrt diff second vars over Q: ");
[ab547e2]349
[0851b0]350  TIMING_END_AND_PRINT (fac_bifactor_total,
351                        "total time for eval and bivar factors over Q: ");
[ab547e2]352  if (irred)
353  {
354    factors.append (A);
355    appendSwapDecompress (factors, contentAFactors, N, 0, 0, x);
356    if (isOn (SW_RATIONAL))
357      normalize (factors);
358    delete [] Aeval2;
359    return factors;
360  }
361
362  if (minFactorsLength == 0)
363    minFactorsLength= biFactors.length();
364  else if (biFactors.length() > minFactorsLength)
365    refineBiFactors (A, biFactors, Aeval2, evaluation, minFactorsLength);
366
[c5230d1]367  if (differentSecondVar == lengthAeval2)
[6036d90]368  {
369    bool zeroOccured= false;
370    for (CFListIterator iter= evaluation; iter.hasItem(); iter++)
371    {
372      if (iter.getItem().isZero())
373      {
374        zeroOccured= true;
375        break;
376      }
377    }
378    if (!zeroOccured)
379    {
[8256fd]380      factors= sparseHeuristic (A, biFactors, Aeval2, evaluation,
381                                minFactorsLength);
[6036d90]382      if (factors.length() == biFactors.length())
383      {
384        appendSwapDecompress (factors, contentAFactors, N, 0, 0, x);
385        normalize (factors);
386        delete [] Aeval2;
387        return factors;
388      }
389      else
390        factors= CFList();
391      //TODO case where factors.length() > 0
392    }
393  }
394
[ab547e2]395  CFList uniFactors= buildUniFactors (biFactors, evaluation.getLast(), y);
396
[5dad7c7]397  sortByUniFactors (Aeval2, lengthAeval2, uniFactors, evaluation);
[fce807]398
[5dad7c7]399  CFList * oldAeval= new CFList [lengthAeval2];
400  for (int i= 0; i < lengthAeval2; i++)
[ab547e2]401    oldAeval[i]= Aeval2[i];
402
[a0adc3]403  getLeadingCoeffs (A, Aeval2);
[ab547e2]404
405  CFList biFactorsLCs;
406  for (CFListIterator i= biFactors; i.hasItem(); i++)
407    biFactorsLCs.append (LC (i.getItem(), 1));
408
409  Variable w;
[0851b0]410  TIMING_START (fac_precompute_leadcoeff);
[afbebe]411  CFList leadingCoeffs= precomputeLeadingCoeff (LC (A, 1), biFactorsLCs, x,
[5dad7c7]412                                          evaluation, Aeval2, lengthAeval2, w);
[ab547e2]413
414  if (w.level() != 1)
415  {
416    A= swapvar (A, y, w);
417    int i= A.level();
418    CanonicalForm evalPoint;
419    for (CFListIterator iter= evaluation; iter.hasItem(); iter++, i--)
420    {
421      if (i == w.level())
422      {
423        evalPoint= iter.getItem();
424        iter.getItem()= evaluation.getLast();
425        evaluation.removeLast();
426        evaluation.append (evalPoint);
427        break;
428      }
429    }
[5dad7c7]430    for (i= 0; i < lengthAeval2; i++)
[8256fd]431    {
432      if (oldAeval[i].isEmpty())
433        continue;
434      if (oldAeval[i].getFirst().level() == w.level())
435      {
436        CFArray tmp= copy (oldAeval[i]);
[753fb61]437        oldAeval[i]= biFactors; 
438        for (CFListIterator iter= oldAeval[i]; iter.hasItem(); iter++)
439          iter.getItem()= swapvar (iter.getItem(), w, y);
[8256fd]440        for (int ii= 0; ii < tmp.size(); ii++)
441          tmp[ii]= swapvar (tmp[ii], w, y);
442        CFArray tmp2= CFArray (tmp.size());
443        CanonicalForm buf;
444        for (int ii= 0; ii < tmp.size(); ii++)
445        {
446          buf= tmp[ii] (evaluation.getLast(),y);
447          buf /= Lc (buf);
448          tmp2[findItem (uniFactors, buf)-1]=tmp[ii];
449        }
450        biFactors= CFList();
451        for (int j= 0; j < tmp2.size(); j++)
452          biFactors.append (tmp2[j]);
453      }
454    }
[ab547e2]455  }
456
[eefc3a]457  CFListIterator iter;
[ab547e2]458  CanonicalForm oldA= A;
459  CFList oldBiFactors= biFactors;
460  if (!leadingCoeffs.getFirst().inCoeffDomain())
461  {
462    CanonicalForm tmp= power (leadingCoeffs.getFirst(), biFactors.length() - 1);
463    A *= tmp;
464    tmp= leadingCoeffs.getFirst();
[eefc3a]465    iter= evaluation;
466    for (int i= A.level(); i > 2; i--, iter++)
467      tmp= tmp (iter.getItem(), i);
[ab547e2]468    if (!tmp.inCoeffDomain())
469    {
470      for (CFListIterator i= biFactors; i.hasItem(); i++)
471      {
472        i.getItem() *= tmp/LC (i.getItem(), 1);
473        i.getItem() /= Lc (i.getItem());
474      }
475    }
476  }
477
[5295f9]478  CanonicalForm LCmultiplier= leadingCoeffs.getFirst();
[fcd296]479  bool LCmultiplierIsConst= LCmultiplier.inCoeffDomain();
[ab547e2]480  leadingCoeffs.removeFirst();
481
482  //prepare leading coefficients
[5dad7c7]483  CFList* leadingCoeffs2= new CFList [lengthAeval2];
[eefc3a]484  prepareLeadingCoeffs (leadingCoeffs2, A.level(), leadingCoeffs, biFactors,
485                        evaluation);
486
487
488  Aeval= evaluateAtEval (A, evaluation, 2);
489
[f7e9c6]490  CanonicalForm hh= 1/Lc (Aeval.getFirst());
[eefc3a]491
[464b18]492  for (iter= Aeval; iter.hasItem(); iter++)
[f7e9c6]493    iter.getItem() *= hh;
[eefc3a]494
[f7e9c6]495  A *= hh;
[eefc3a]496
[6af78bb]497  CFListIterator iter2;
[5dad7c7]498  CFList bufLeadingCoeffs2= leadingCoeffs2[lengthAeval2-1];
[fcd296]499  bufBiFactors= biFactors;
500  bufA= A;
[c1ea88]501  CanonicalForm bufLCmultiplier= LCmultiplier;
[5dad7c7]502  CanonicalForm testVars;
503  if (!LCmultiplierIsConst)
504  {
505    testVars= Variable (2);
506    for (int i= 0; i < lengthAeval2; i++)
507    {
508      if (!oldAeval[i].isEmpty())
509        testVars *= oldAeval[i].getFirst().mvar();
510    }
511  }
[0851b0]512  TIMING_END_AND_PRINT(fac_precompute_leadcoeff,
513                       "time to precompute LC over Q: ");
514
515  TIMING_START (fac_luckswang);
[5079887]516  CFList bufFactors= CFList();
[fcd296]517  bool LCheuristic= false;
[5dad7c7]518  if (LucksWangSparseHeuristic (A, biFactors, 2, leadingCoeffs2[lengthAeval2-1],
[64b824]519                                factors))
[4fe8a3]520  {
[5079887]521    int check= biFactors.length();
522    int * index= new int [factors.length()];
523    CFList oldFactors= factors;
524    factors= recoverFactors (A, factors, index);
[4fe8a3]525
526    if (check == factors.length())
527    {
[8256fd]528      if (w.level() != 1)
529      {
[5079887]530        for (iter= factors; iter.hasItem(); iter++)
[8256fd]531          iter.getItem()= swapvar (iter.getItem(), w, y);
532      }
533
[4fe8a3]534      appendSwapDecompress (factors, contentAFactors, N, 0, 0, x);
535      normalize (factors);
[5079887]536      delete [] index;
[4fe8a3]537      delete [] Aeval2;
[0851b0]538      TIMING_END_AND_PRINT (fac_luckswang,
539                            "time for successful LucksWang over Q: ");
[4fe8a3]540      return factors;
541    }
[5079887]542    else if (factors.length() > 0)
543    {
544      int oneCount= 0;
545      CFList l;
546      for (int i= 0; i < check; i++)
547      {
548        if (index[i] == 1)
549        {
550          iter=biFactors;
551          for (int j=1; j <= i-oneCount; j++)
552            iter++;
553          iter.remove (1);
[5dad7c7]554          for (int j= 0; j < lengthAeval2; j++)
[5079887]555          {
556            l= leadingCoeffs2[j];
557            iter= l;
558            for (int k=1; k <= i-oneCount; k++)
559              iter++;
560            iter.remove (1);
561            leadingCoeffs2[j]=l;
562          }
563          oneCount++;
564        }
565      }
566      bufFactors= factors;
567      factors= CFList();
568    }
[5dad7c7]569    else if (!LCmultiplierIsConst && factors.length() == 0)
[5295f9]570    {
[fcd296]571      LCheuristic= true;
[5295f9]572      factors= oldFactors;
573      CanonicalForm cont;
574      CFList contents, LCs;
575      int index=1;
576      bool foundTrueMultiplier= false;
577      for (iter= factors; iter.hasItem(); iter++, index++)
578      {
579        cont= content (iter.getItem(), 1);
580        cont= gcd (cont , LCmultiplier);
581        contents.append (cont);
582        if (cont.inCoeffDomain()) // trivial content->LCmultiplier needs to go there
583        {
584          foundTrueMultiplier= true;
585          int index2= 1;
[5dad7c7]586          for (iter2= leadingCoeffs2[lengthAeval2-1]; iter2.hasItem(); iter2++,
[5295f9]587                                                                    index2++)
588          {
589            if (index2 == index)
590              continue;
591            iter2.getItem() /= LCmultiplier;
592          }
[397c12]593          A= oldA;
[5dad7c7]594          leadingCoeffs= leadingCoeffs2[lengthAeval2-1];
595          for (int i= lengthAeval2-1; i > -1; i--)
[5295f9]596            leadingCoeffs2[i]= CFList();
597          prepareLeadingCoeffs (leadingCoeffs2, A.level(), leadingCoeffs,
598                                biFactors, evaluation );
599          Aeval= evaluateAtEval (A, evaluation, 2);
600
[f7e9c6]601          hh= 1/Lc (Aeval.getFirst());
[5295f9]602
603          for (iter2= Aeval; iter2.hasItem(); iter2++)
[f7e9c6]604            iter2.getItem() *= hh;
[5295f9]605
[f7e9c6]606          A *= hh;
[5295f9]607          break;
608        }
609        else
610          LCs.append (LC (iter.getItem()/cont, 1));
611      }
612      if (!foundTrueMultiplier)
613      {
614        index= 1;
615        iter2= factors;
616        bool foundMultiplier= false;
617        for (iter= contents; iter.hasItem(); iter++, iter2++, index++)
618        {
619          if (fdivides (iter.getItem(), LCmultiplier))
620          {
621            if ((LCmultiplier/iter.getItem()).inCoeffDomain() &&
622                !isOnlyLeadingCoeff(iter2.getItem())) //content divides LCmultiplier completely and factor consists of more terms than just the leading coeff
623            {
[86faff]624              Variable xx= Variable (2);
625              CanonicalForm vars;
626              vars= power (xx, degree (LC (getItem(oldBiFactors, index),1),
627                                        xx));
628              for (int i= 0; i < lengthAeval2; i++)
[5295f9]629              {
[86faff]630                if (oldAeval[i].isEmpty())
631                  continue;
632                xx= oldAeval[i].getFirst().mvar();
633                vars *= power (xx, degree (LC (getItem(oldAeval[i], index),1),
634                                           xx));
635              }
636              if (vars.level() <= 2)
637              {
638                int index2= 1;
639                for (CFListIterator iter3= leadingCoeffs2[lengthAeval2-1];
640                     iter3.hasItem(); iter3++, index2++)
[5295f9]641                {
[86faff]642                  if (index2 == index)
643                  {
644                    iter3.getItem() /= LCmultiplier;
645                    break;
646                  }
[5295f9]647                }
[86faff]648                A /= LCmultiplier;
649                foundMultiplier= true;
650                iter.getItem()= 1;
[5295f9]651              }
652            }
653          }
654        }
655        // coming from above: divide out more LCmultiplier if possible
656        if (foundMultiplier)
657        {
658          foundMultiplier= false;
659          index=1;
660          iter2= factors;
661          for (iter= contents; iter.hasItem(); iter++, iter2++, index++)
662          {
[4cf4dc]663            if (!iter.getItem().isOne() &&
[5295f9]664                fdivides (iter.getItem(), LCmultiplier))
665            {
666              if (!isOnlyLeadingCoeff (iter2.getItem())) // factor is more than just leading coeff
667              {
668                int index2= 1;
[5dad7c7]669                for (iter2= leadingCoeffs2[lengthAeval2-1]; iter2.hasItem();
[5295f9]670                     iter2++, index2++)
671                {
672                  if (index2 == index)
673                  {
674                    iter2.getItem() /= iter.getItem();
675                    foundMultiplier= true;
676                    break;
677                  }
678                }
679                A /= iter.getItem();
680                LCmultiplier /= iter.getItem();
681                iter.getItem()= 1;
682              }
[5dad7c7]683              else if (fdivides (getVars (LCmultiplier), testVars))//factor consists of just leading coeff
[5295f9]684              {
[4cf4dc]685                Variable xx= Variable (2);
686                CanonicalForm vars;
687                vars= power (xx, degree (LC (getItem(oldBiFactors, index),1),
688                                          xx));
[5dad7c7]689                for (int i= 0; i < lengthAeval2; i++)
[5295f9]690                {
691                  if (oldAeval[i].isEmpty())
692                    continue;
693                  xx= oldAeval[i].getFirst().mvar();
[4cf4dc]694                  vars *= power (xx, degree (LC (getItem(oldAeval[i], index),1),
695                                             xx));
[5295f9]696                }
[5dad7c7]697                if (myGetVars(content(getItem(leadingCoeffs2[lengthAeval2-1],index),1))
[4cf4dc]698                    /myGetVars (LCmultiplier) == vars)
[5295f9]699                {
700                  int index2= 1;
[5dad7c7]701                  for (iter2= leadingCoeffs2[lengthAeval2-1]; iter2.hasItem();
[5295f9]702                       iter2++, index2++)
703                  {
704                    if (index2 == index)
705                    {
[4cf4dc]706                      iter2.getItem() /= LCmultiplier;
[5295f9]707                      foundMultiplier= true;
708                      break;
709                    }
710                  }
[4cf4dc]711                  A /= LCmultiplier;
[5295f9]712                  iter.getItem()= 1;
713                }
714              }
715            }
716          }
717        }
718        else
719        {
720          CanonicalForm pLCs= prod (LCs);
721          if (fdivides (pLCs, LC (oldA,1)) && (LC(oldA,1)/pLCs).inCoeffDomain()) // check if the product of the lead coeffs of the primitive factors equals the lead coeff of the old A
722          {
723            A= oldA;
[5dad7c7]724            iter2= leadingCoeffs2[lengthAeval2-1];
[5295f9]725            for (iter= contents; iter.hasItem(); iter++, iter2++)
726              iter2.getItem() /= iter.getItem();
[f3ff22]727            foundMultiplier= true;
728          }
[5dad7c7]729          if (!foundMultiplier && fdivides (getVars (LCmultiplier), testVars))
[f3ff22]730          {
731            Variable xx;
732            CFList vars1;
733            CFFList sqrfMultiplier= sqrFree (LCmultiplier);
734            if (sqrfMultiplier.getFirst().factor().inCoeffDomain())
735              sqrfMultiplier.removeFirst();
736            sqrfMultiplier= sortCFFListByNumOfVars (sqrfMultiplier);
737            xx= Variable (2);
738            for (iter= oldBiFactors; iter.hasItem(); iter++)
739              vars1.append (power (xx, degree (LC (iter.getItem(),1), xx)));
[5dad7c7]740            for (int i= 0; i < lengthAeval2; i++)
[f3ff22]741            {
742              if (oldAeval[i].isEmpty())
743                continue;
744              xx= oldAeval[i].getFirst().mvar();
745              iter2= vars1;
746              for (iter= oldAeval[i]; iter.hasItem(); iter++, iter2++)
747                iter2.getItem() *= power(xx,degree (LC (iter.getItem(),1), xx));
748            }
749            CanonicalForm tmp;
750            iter2= vars1;
[5dad7c7]751            for (iter= leadingCoeffs2[lengthAeval2-1]; iter.hasItem(); iter++,
[f3ff22]752                                                                    iter2++)
753            {
754              tmp= iter.getItem()/LCmultiplier;
755              for (int i=1; i <= tmp.level(); i++)
756              {
[bed38b]757                if (degree(tmp,i) > 0 &&
758                    (degree(iter2.getItem(),i) > degree (tmp,i)))
[f3ff22]759                  iter2.getItem() /= power (Variable (i), degree (tmp,i));
760              }
761            }
762            int multi;
763            for (CFFListIterator ii= sqrfMultiplier; ii.hasItem(); ii++)
764            {
765              multi= 0;
766              for (iter= vars1; iter.hasItem(); iter++)
767              {
768                tmp= iter.getItem();
769                while (fdivides (myGetVars (ii.getItem().factor()), tmp))
770                {
771                  multi++;
772                  tmp /= myGetVars (ii.getItem().factor());
773                }
774              }
775              if (multi == ii.getItem().exp())
776              {
777                index= 1;
778                for (iter= vars1; iter.hasItem(); iter++, index++)
779                {
780                  while (fdivides (myGetVars(ii.getItem().factor()),
781                                   iter.getItem()
782                                  )
783                        )
784                  {
785                    int index2= 1;
[5dad7c7]786                    for (iter2= leadingCoeffs2[lengthAeval2-1]; iter2.hasItem();
[f3ff22]787                         iter2++, index2++)
788                    {
789                      if (index2 == index)
790                        continue;
791                      else
792                      {
793                        tmp= ii.getItem().factor();
794                        iter2.getItem() /= tmp;
795                        CFListIterator iter3= evaluation;
796                        for (int jj= A.level(); jj > 2; jj--, iter3++)
797                          tmp= tmp (iter3.getItem(), jj);
798                        if (!tmp.inCoeffDomain())
799                        {
800                          int index3= 1;
801                          for (iter3= biFactors; iter3.hasItem(); iter3++,
802                                                                  index3++)
803                          {
804                            if (index3 == index2)
805                            {
806                              iter3.getItem() /= tmp;
807                              iter3.getItem() /= Lc (iter3.getItem());
808                              break;
809                            }
810                          }
811                        }
812                        A /= ii.getItem().factor();
813                      }
814                    }
815                    iter.getItem() /= getVars (ii.getItem().factor());
816                  }
817                }
818              }
819              else
820              {
821                index= 1;
822                for (iter= vars1; iter.hasItem(); iter++, index++)
823                {
824                  if (!fdivides (myGetVars (ii.getItem().factor()),
825                                 iter.getItem()
826                                )
827                     )
828                  {
829                    int index2= 1;
[5dad7c7]830                    for (iter2= leadingCoeffs2[lengthAeval2-1]; iter2.hasItem();
[f3ff22]831                         iter2++, index2++)
832                    {
833                      if (index2 == index)
834                      {
835                        tmp= power (ii.getItem().factor(), ii.getItem().exp());
836                        iter2.getItem() /= tmp;
837                        A /= tmp;
838                        CFListIterator iter3= evaluation;
839                        for (int jj= A.level(); jj > 2; jj--, iter3++)
840                          tmp= tmp (iter3.getItem(), jj);
841                        if (!tmp.inCoeffDomain())
842                        {
843                          int index3= 1;
844                          for (iter3= biFactors; iter3.hasItem(); iter3++,
845                                                                  index3++)
846                          {
847                            if (index3 == index2)
848                            {
849                              iter3.getItem() /= tmp;
850                              iter3.getItem() /= Lc (iter3.getItem());
851                              break;
852                            }
853                          }
854                        }
855                      }
856                    }
857                  }
858                }
859              }
860            }
[5295f9]861          }
862        }
863
864        // patch everything together again
[5dad7c7]865        leadingCoeffs= leadingCoeffs2[lengthAeval2-1];
866        for (int i= lengthAeval2-1; i > -1; i--)
[5295f9]867          leadingCoeffs2[i]= CFList();
868        prepareLeadingCoeffs (leadingCoeffs2,A.level(),leadingCoeffs, biFactors,
869                              evaluation);
870        Aeval= evaluateAtEval (A, evaluation, 2);
871
[f7e9c6]872        hh= 1/Lc (Aeval.getFirst());
[5295f9]873
874        for (CFListIterator i= Aeval; i.hasItem(); i++)
[f7e9c6]875          i.getItem() *= hh;
[5295f9]876
[f7e9c6]877        A *= hh;
[5295f9]878      }
879      factors= CFList();
[5dad7c7]880      if (!fdivides (LC (oldA,1),prod (leadingCoeffs2[lengthAeval2-1])))
[c1ea88]881      {
882        LCheuristic= false;
883        A= bufA;
884        biFactors= bufBiFactors;
[5dad7c7]885        leadingCoeffs2[lengthAeval2-1]= bufLeadingCoeffs2;
[c1ea88]886        LCmultiplier= bufLCmultiplier;
887      }
[5295f9]888    }
[4fe8a3]889    else
890      factors= CFList();
[5079887]891    delete [] index;
[f3ff22]892  }
[0851b0]893  TIMING_END_AND_PRINT (fac_luckswang, "time for LucksWang over Q: ");
894
895  TIMING_START (fac_lcheuristic);
[5dad7c7]896  if (!LCheuristic && !LCmultiplierIsConst && bufFactors.isEmpty()
897      && fdivides (getVars (LCmultiplier), testVars))
[f3ff22]898  {
[41e77d]899    LCheuristic= true;
[f3ff22]900    int index;
901    Variable xx;
902    CFList vars1;
903    CFFList sqrfMultiplier= sqrFree (LCmultiplier);
904    if (sqrfMultiplier.getFirst().factor().inCoeffDomain())
905      sqrfMultiplier.removeFirst();
906    sqrfMultiplier= sortCFFListByNumOfVars (sqrfMultiplier);
907    xx= Variable (2);
908    for (iter= oldBiFactors; iter.hasItem(); iter++)
909      vars1.append (power (xx, degree (LC (iter.getItem(),1), xx)));
[5dad7c7]910    for (int i= 0; i < lengthAeval2; i++)
[f3ff22]911    {
912      if (oldAeval[i].isEmpty())
913        continue;
914      xx= oldAeval[i].getFirst().mvar();
915      iter2= vars1;
916      for (iter= oldAeval[i]; iter.hasItem(); iter++, iter2++)
917        iter2.getItem() *= power (xx, degree (LC (iter.getItem(),1), xx));
918    }
919    CanonicalForm tmp;
920    iter2= vars1;
[5dad7c7]921    for (iter= leadingCoeffs2[lengthAeval2-1]; iter.hasItem(); iter++, iter2++)
[f3ff22]922    {
923      tmp= iter.getItem()/LCmultiplier;
924      for (int i=1; i <= tmp.level(); i++)
925      {
[bed38b]926        if (degree(tmp,i) > 0 && (degree(iter2.getItem(),i) > degree (tmp,i)))
[f3ff22]927          iter2.getItem() /= power (Variable (i), degree (tmp,i));
928      }
929    }
930    int multi;
931    for (CFFListIterator ii= sqrfMultiplier; ii.hasItem(); ii++)
932    {
933      multi= 0;
934      for (iter= vars1; iter.hasItem(); iter++)
935      {
936        tmp= iter.getItem();
937        while (fdivides (myGetVars (ii.getItem().factor()), tmp))
938        {
939          multi++;
940          tmp /= myGetVars (ii.getItem().factor());
941        }
942      }
943      if (multi == ii.getItem().exp())
944      {
945        index= 1;
946        for (iter= vars1; iter.hasItem(); iter++, index++)
947        {
948          while (fdivides (myGetVars (ii.getItem().factor()), iter.getItem()))
949          {
950            int index2= 1;
[5dad7c7]951            for (iter2= leadingCoeffs2[lengthAeval2-1];iter2.hasItem();iter2++,
[f3ff22]952                                                                      index2++)
953            {
954              if (index2 == index)
955                continue;
956              else
957              {
958                tmp= ii.getItem().factor();
959                iter2.getItem() /= tmp;
960                CFListIterator iter3= evaluation;
961                for (int jj= A.level(); jj > 2; jj--, iter3++)
962                  tmp= tmp (iter3.getItem(), jj);
963                if (!tmp.inCoeffDomain())
964                {
965                  int index3= 1;
966                  for (iter3= biFactors; iter3.hasItem(); iter3++, index3++)
967                  {
968                    if (index3 == index2)
969                    {
970                      iter3.getItem() /= tmp;
971                      iter3.getItem() /= Lc (iter3.getItem());
972                      break;
973                    }
974                  }
975                }
976                A /= ii.getItem().factor();
977              }
978            }
979            iter.getItem() /= getVars (ii.getItem().factor());
980          }
981        }
982      }
983      else
984      {
985        index= 1;
986        for (iter= vars1; iter.hasItem(); iter++, index++)
987        {
988          if (!fdivides (myGetVars (ii.getItem().factor()), iter.getItem()))
989          {
990            int index2= 1;
[5dad7c7]991            for (iter2= leadingCoeffs2[lengthAeval2-1];iter2.hasItem();iter2++,
[f3ff22]992                                                                      index2++)
993            {
994              if (index2 == index)
995              {
996                tmp= power (ii.getItem().factor(), ii.getItem().exp());
997                iter2.getItem() /= tmp;
998                A /= tmp;
999                CFListIterator iter3= evaluation;
1000                for (int jj= A.level(); jj > 2; jj--, iter3++)
1001                  tmp= tmp (iter3.getItem(), jj);
1002                if (!tmp.inCoeffDomain())
1003                {
1004                  int index3= 1;
1005                  for (iter3= biFactors; iter3.hasItem(); iter3++, index3++)
1006                  {
1007                    if (index3 == index2)
1008                    {
1009                      iter3.getItem() /= tmp;
1010                      iter3.getItem() /= Lc (iter3.getItem());
1011                      break;
1012                    }
1013                  }
1014                }
1015              }
1016            }
1017          }
1018        }
1019      }
1020    }
1021
[5dad7c7]1022    leadingCoeffs= leadingCoeffs2[lengthAeval2-1];
1023    for (int i= lengthAeval2-1; i > -1; i--)
[f3ff22]1024      leadingCoeffs2[i]= CFList();
1025    prepareLeadingCoeffs (leadingCoeffs2,A.level(),leadingCoeffs, biFactors,
1026                          evaluation);
1027    Aeval= evaluateAtEval (A, evaluation, 2);
1028
[f7e9c6]1029    hh= 1/Lc (Aeval.getFirst());
[f3ff22]1030
1031    for (CFListIterator i= Aeval; i.hasItem(); i++)
[f7e9c6]1032      i.getItem() *= hh;
[f3ff22]1033
[f7e9c6]1034    A *= hh;
[83e5de]1035
1036    if (!fdivides (LC (oldA,1),prod (leadingCoeffs2[lengthAeval2-1])))
1037    {
1038      LCheuristic= false;
1039      A= bufA;
1040      biFactors= bufBiFactors;
1041      leadingCoeffs2[lengthAeval2-1]= bufLeadingCoeffs2;
1042      LCmultiplier= bufLCmultiplier;
1043    }
[64b824]1044  }
[0851b0]1045  TIMING_END_AND_PRINT (fac_lcheuristic, "time for Lc heuristic over Q: ");
[4fe8a3]1046
[fcd296]1047tryAgainWithoutHeu:
[eefc3a]1048  //shifting to zero
[0851b0]1049  TIMING_START (fac_shift_to_zero);
[e5fdcc4]1050  CanonicalForm denA= bCommonDen (A);
1051  A *= denA;
[eefc3a]1052  A= shift2Zero (A, Aeval, evaluation);
[e5fdcc4]1053  A /= denA;
[eefc3a]1054
1055  for (iter= biFactors; iter.hasItem(); iter++)
1056    iter.getItem()= iter.getItem () (y + evaluation.getLast(), y);
1057
[5dad7c7]1058  for (int i= 0; i < lengthAeval2-1; i++)
[464b18]1059    leadingCoeffs2[i]= CFList();
[5dad7c7]1060  for (iter= leadingCoeffs2[lengthAeval2-1]; iter.hasItem(); iter++)
[eefc3a]1061  {
1062    iter.getItem()= shift2Zero (iter.getItem(), list, evaluation);
1063    for (int i= A.level() - 4; i > -1; i--)
1064    {
[5dad7c7]1065      if (i + 1 == lengthAeval2-1)
[eefc3a]1066        leadingCoeffs2[i].append (iter.getItem() (0, i + 4));
1067      else
1068        leadingCoeffs2[i].append (leadingCoeffs2[i+1].getLast() (0, i + 4));
1069    }
1070  }
[0851b0]1071  TIMING_END_AND_PRINT (fac_shift_to_zero,
1072                        "time to shift evaluation point to zero: ");
[ab547e2]1073
1074  CFArray Pi;
1075  CFList diophant;
1076  int* liftBounds= new int [A.level() - 1];
1077  int liftBoundsLength= A.level() - 1;
1078  for (int i= 0; i < liftBoundsLength; i++)
1079    liftBounds [i]= degree (A, i + 2) + 1;
1080
1081  Aeval.removeFirst();
1082  bool noOneToOne= false;
1083
[0851b0]1084  TIMING_START (fac_cleardenom);
[f55f39]1085  CFList commonDenominators;
[6af78bb]1086  for (iter=biFactors; iter.hasItem(); iter++)
[f55f39]1087    commonDenominators.append (bCommonDen (iter.getItem()));
1088  CanonicalForm tmp1, tmp2, tmp3=1;
[5dad7c7]1089  for (int i= 0; i < lengthAeval2; i++)
[f55f39]1090  {
1091    iter2= commonDenominators;
[6af78bb]1092    for (iter= leadingCoeffs2[i]; iter.hasItem(); iter++, iter2++)
[f55f39]1093    {
[6af78bb]1094      tmp1= bCommonDen (iter.getItem());
[f55f39]1095      Off (SW_RATIONAL);
1096      iter2.getItem()= lcm (iter2.getItem(), tmp1);
1097      On (SW_RATIONAL);
1098    }
1099  }
1100  tmp1= prod (commonDenominators);
[6af78bb]1101  for (iter= Aeval; iter.hasItem(); iter++)
[f55f39]1102  {
[e5fdcc4]1103    tmp2= bCommonDen (iter.getItem()/denA);
[f55f39]1104    Off (SW_RATIONAL);
1105    tmp3= lcm (tmp2,tmp3);
1106    On (SW_RATIONAL);
1107  }
1108  CanonicalForm multiplier;
1109  multiplier= tmp3/tmp1;
1110  iter2= commonDenominators;
[6af78bb]1111  for (iter=biFactors; iter.hasItem(); iter++, iter2++)
1112    iter.getItem() *= iter2.getItem()*multiplier;
[f55f39]1113
[6af78bb]1114  for (iter= Aeval; iter.hasItem(); iter++)
[e5fdcc4]1115    iter.getItem() *= tmp3*power (multiplier, biFactors.length() - 1)/denA;
[f55f39]1116
[5dad7c7]1117  for (int i= 0; i < lengthAeval2; i++)
[f55f39]1118  {
1119    iter2= commonDenominators;
[6af78bb]1120    for (iter= leadingCoeffs2[i]; iter.hasItem(); iter++, iter2++)
1121      iter.getItem() *= iter2.getItem()*multiplier;
[f55f39]1122  }
[0851b0]1123  TIMING_END_AND_PRINT (fac_cleardenom, "time to clear denominators: ");
[f55f39]1124
[0851b0]1125  TIMING_START (fac_hensel_lift);
[ab547e2]1126  factors= nonMonicHenselLift (Aeval, biFactors, leadingCoeffs2, diophant,
1127                               Pi, liftBounds, liftBoundsLength, noOneToOne);
[0851b0]1128  TIMING_END_AND_PRINT (fac_hensel_lift,
1129                        "time for non monic hensel lifting over Q: ");
[ab547e2]1130
1131  if (!noOneToOne)
1132  {
1133    int check= factors.length();
[eefc3a]1134    A= oldA;
[0851b0]1135    TIMING_START (fac_recover_factors);
[eefc3a]1136    factors= recoverFactors (A, factors, evaluation);
[0851b0]1137    TIMING_END_AND_PRINT (fac_recover_factors,
1138                          "time to recover factors over Q: ");
[ab547e2]1139    if (check != factors.length())
1140      noOneToOne= true;
[5079887]1141    else
1142      factors= Union (factors, bufFactors);
[ab547e2]1143  }
1144  if (noOneToOne)
1145  {
[fcd296]1146    if (!LCmultiplierIsConst && LCheuristic)
1147    {
1148      A= bufA;
1149      biFactors= bufBiFactors;
[5dad7c7]1150      leadingCoeffs2[lengthAeval2-1]= bufLeadingCoeffs2;
[fcd296]1151      delete [] liftBounds;
1152      LCheuristic= false;
1153      goto tryAgainWithoutHeu;
1154      //something probably went wrong in the heuristic
1155    }
1156
[eefc3a]1157    A= shift2Zero (oldA, Aeval, evaluation);
[ab547e2]1158    biFactors= oldBiFactors;
[eefc3a]1159    for (iter= biFactors; iter.hasItem(); iter++)
1160      iter.getItem()= iter.getItem () (y + evaluation.getLast(), y);
[ab547e2]1161    CanonicalForm LCA= LC (Aeval.getFirst(), 1);
1162    CanonicalForm yToLift= power (y, lift);
1163    CFListIterator i= biFactors;
1164    lift= degree (i.getItem(), 2) + degree (LC (i.getItem(), 1)) + 1;
1165    i++;
1166
1167    for (; i.hasItem(); i++)
1168      lift= tmax (lift, degree (i.getItem(), 2)+degree (LC (i.getItem(), 1))+1);
1169
1170    lift= tmax (degree (Aeval.getFirst() , 2) + 1, lift);
1171
1172    i= biFactors;
1173    yToLift= power (y, lift);
1174    CanonicalForm dummy;
1175    for (; i.hasItem(); i++)
1176    {
1177      LCA= LC (i.getItem(), 1);
1178      extgcd (LCA, yToLift, LCA, dummy);
1179      i.getItem()= mod (i.getItem()*LCA, yToLift);
1180    }
1181
1182    liftBoundsLength= F.level() - 1;
1183    liftBounds= liftingBounds (A, lift);
1184
1185    CFList MOD;
1186    bool earlySuccess;
1187    CFList earlyFactors;
1188    ExtensionInfo info= ExtensionInfo (false);
[0b618a7]1189    CFList liftedFactors;
[ec664a2]1190    TIMING_START (fac_hensel_lift);
[0b618a7]1191    liftedFactors= henselLiftAndEarly
1192                   (A, MOD, liftBounds, earlySuccess, earlyFactors,
1193                    Aeval, biFactors, evaluation, info);
[ec664a2]1194    TIMING_END_AND_PRINT (fac_hensel_lift, "time for hensel lifting: ");
[0b618a7]1195
[ec664a2]1196    TIMING_START (fac_factor_recombination);
[ab547e2]1197    factors= factorRecombination (A, liftedFactors, MOD);
[ec664a2]1198    TIMING_END_AND_PRINT (fac_factor_recombination,
[ab547e2]1199                          "time for factor recombination: ");
1200
1201    if (earlySuccess)
1202      factors= Union (factors, earlyFactors);
1203
[c79a9d]1204    for (CFListIterator i= factors; i.hasItem(); i++)
[ab547e2]1205    {
[c79a9d]1206      int kk= Aeval.getLast().level();
1207      for (CFListIterator j= evaluation; j.hasItem(); j++, kk--)
1208      {
1209        if (i.getItem().level() < kk)
1210          continue;
1211       i.getItem()= i.getItem() (Variable (kk) - j.getItem(), kk);
1212      }
[ab547e2]1213    }
1214  }
1215
1216  if (w.level() != 1)
1217  {
1218    for (CFListIterator iter= factors; iter.hasItem(); iter++)
1219      iter.getItem()= swapvar (iter.getItem(), w, y);
1220  }
1221
1222  swap (factors, 0, 0, x);
1223  append (factors, contentAFactors);
1224  decompress (factors, N);
1225  if (isOn (SW_RATIONAL))
1226    normalize (factors);
1227
[229530]1228  delete [] leadingCoeffs2;
1229  delete [] oldAeval;
1230  delete [] Aeval2;
[ab547e2]1231  delete[] liftBounds;
1232
1233  return factors;
1234}
[d990001]1235
1236#endif
Note: See TracBrowser for help on using the repository browser.