source: git/factory/facAbsFact.cc @ 73c222

spielwiese
Last change on this file since 73c222 was 73c222, checked in by Martin Lee <martinlee84@…>, 11 years ago
chg: moved some code from .h to .cc chg: absBiFactorizeMain now also returns other factors
  • Property mode set to 100644
File size: 20.7 KB
Line 
1/*****************************************************************************\
2 * Computer Algebra System SINGULAR
3\*****************************************************************************/
4/** @file facAbsFact.cc
5 *
6 * @author Martin Lee
7 *
8 **/
9/*****************************************************************************/
10
11#ifdef HAVE_CONFIG_H
12#include "config.h"
13#endif /* HAVE_CONFIG_H */
14
15#include "timing.h"
16#include "debug.h"
17
18#include "facAbsFact.h"
19#include "facBivar.h"
20#include "facFqBivar.h"
21#include "cf_reval.h"
22#include "cf_primes.h"
23#include "cf_algorithm.h"
24#ifdef HAVE_FLINT
25#include "FLINTconvert.h"
26#include <flint/fmpz_poly_factor.h>
27#endif
28#ifdef HAVE_NTL
29#include "NTLconvert.h"
30#include <NTL/LLL.h>
31#endif
32
33#ifdef HAVE_NTL
34
35TIMING_DEFINE_PRINT(fac_Qa_factorize)
36TIMING_DEFINE_PRINT(fac_evalpoint)
37
38CFAFList uniAbsFactorize (const CanonicalForm& F)
39{
40  CFFList rationalFactors= factorize (F);
41  CFFListIterator i= rationalFactors;
42  i++;
43  Variable alpha;
44  CFAFList result;
45  CFFList QaFactors;
46  CFFListIterator iter;
47  for (; i.hasItem(); i++)
48  {
49    if (degree (i.getItem().factor()) == 1)
50    {
51      result.append (CFAFactor (i.getItem().factor(), 1, i.getItem().exp()));
52      continue;
53    }
54    alpha= rootOf (i.getItem().factor());
55    QaFactors= factorize (i.getItem().factor(), alpha);
56    iter= QaFactors;
57    if (iter.getItem().factor().inCoeffDomain())
58      iter++;
59    for (;iter.hasItem(); iter++)
60    {
61      if (degree (iter.getItem().factor()) == 1)
62      {
63        result.append (CFAFactor (iter.getItem().factor(), getMipo (alpha),
64                                  i.getItem().exp()));
65        break;
66      }
67    }
68  }
69  result.insert (CFAFactor (rationalFactors.getFirst().factor(), 1, 1));
70  return result;
71}
72
73
74/// absolute factorization of bivariate poly over Q
75///
76/// @return absFactorize returns a list whose entries contain three entities:
77///         an absolute irreducible factor, an irreducible univariate polynomial
78///         that defines the minimal field extension over which the irreducible
79///         factor is defined and the multiplicity of the absolute irreducible
80///         factor
81CFAFList absBiFactorize (const CanonicalForm& G ///<[in] bivariate poly over Q
82                        )
83{
84  //TODO handle homogeneous input
85  ASSERT (getNumVars (G) <= 2, "expected bivariate input");
86  ASSERT (getCharacteristic() == 0, "expected poly over Q");
87
88  CFMap N;
89  CanonicalForm F= compress (G, N);
90  bool isRat= isOn (SW_RATIONAL);
91  if (isRat)
92    F *= bCommonDen (F);
93
94  Off (SW_RATIONAL);
95  F /= icontent (F);
96  if (isRat)
97    On (SW_RATIONAL);
98
99  CanonicalForm contentX= content (F, 1);
100  CanonicalForm contentY= content (F, 2);
101  F /= (contentX*contentY);
102  CFAFList contentXFactors, contentYFactors;
103  contentXFactors= uniAbsFactorize (contentX);
104  contentYFactors= uniAbsFactorize (contentY);
105
106  if (contentXFactors.getFirst().factor().inCoeffDomain())
107    contentXFactors.removeFirst();
108  if (contentYFactors.getFirst().factor().inCoeffDomain())
109    contentYFactors.removeFirst();
110  if (F.inCoeffDomain())
111  {
112    CFAFList result;
113    for (CFAFListIterator i= contentXFactors; i.hasItem(); i++)
114      result.append (CFAFactor (N (i.getItem().factor()), i.getItem().minpoly(),
115                                i.getItem().exp()));
116    for (CFAFListIterator i= contentYFactors; i.hasItem(); i++)
117      result.append (CFAFactor (N (i.getItem().factor()),i.getItem().minpoly(),
118                                i.getItem().exp()));
119    normalize (result);
120    result.insert (CFAFactor (Lc (G), 1, 1));
121    return result;
122  }
123  CFFList rationalFactors= factorize (F);
124
125  CFAFList result, resultBuf;
126
127  CFAFListIterator iter;
128  CFFListIterator i= rationalFactors;
129  i++;
130  for (; i.hasItem(); i++)
131  {
132    resultBuf= absBiFactorizeMain (i.getItem().factor());
133    for (iter= resultBuf; iter.hasItem(); iter++)
134      iter.getItem()= CFAFactor (iter.getItem().factor(),
135                                 iter.getItem().minpoly(), i.getItem().exp());
136    result= Union (result, resultBuf);
137  }
138
139  for (CFAFListIterator i= result; i.hasItem(); i++)
140    i.getItem()= CFAFactor (N (i.getItem().factor()), i.getItem().minpoly(),
141                            i.getItem().exp());
142  for (CFAFListIterator i= contentXFactors; i.hasItem(); i++)
143    result.append (CFAFactor (N(i.getItem().factor()), i.getItem().minpoly(),
144                              i.getItem().exp()));
145  for (CFAFListIterator i= contentYFactors; i.hasItem(); i++)
146    result.append (CFAFactor (N(i.getItem().factor()), i.getItem().minpoly(),
147                              i.getItem().exp()));
148  normalize (result);
149  result.insert (CFAFactor (Lc(G), 1, 1));
150
151  return result;
152}
153
154//TODO optimize choice of p -> choose p as large as possible (better than small p since factorization mod p does not require field extension, also less lifting)
155int choosePoint (const CanonicalForm& F, int tdegF, CFArray& eval, bool rec)
156{
157  REvaluation E1 (1, 1, IntRandom (25));
158  REvaluation E2 (2, 2, IntRandom (25));
159  if (rec)
160  {
161    E1.nextpoint();
162    E2.nextpoint();
163  }
164  CanonicalForm f, f1, f2, Fp;
165  int i, p;
166  eval=CFArray (2);
167  while (1)
168  {
169    f1= E1(F);
170    if (!f1.isZero() && factorize (f1).length() == 2)
171    {
172      Off (SW_RATIONAL);
173      f= E2(f1);
174      f2= E2 (F);
175      if ((!f.isZero()) && (abs(f)>cf_getSmallPrime (cf_getNumSmallPrimes()-1)))
176      {
177        for (i= cf_getNumPrimes()-1; i >= 0; i--)
178        {
179          if (f % CanonicalForm (cf_getPrime (i)) == 0)
180          {
181            p= cf_getPrime(i);
182            Fp= mod (F,p);
183            if (totaldegree (Fp) == tdegF &&
184                degree (mod (f2,p), 1) == degree (F,1) &&
185                degree (mod (f1, p),2) == degree (F,2))
186            {
187              eval[0]= E1[1];
188              eval[1]= E2[2];
189              return p;
190            }
191          }
192        }
193      }
194      else if (!f.isZero())
195      {
196        for (i= cf_getNumSmallPrimes()-1; i >= 0; i--)
197        {
198          if (f % CanonicalForm (cf_getSmallPrime (i)) == 0)
199          {
200            p= cf_getSmallPrime (i);
201            Fp= mod (F,p);
202            if (totaldegree (Fp) == tdegF &&
203                degree (mod (f2, p),1) == degree (F,1) &&
204                degree (mod (f1,p),2) == degree (F,2))
205            {
206              eval[0]= E1[1];
207              eval[1]= E2[2];
208              return p;
209            }
210          }
211        }
212      }
213      E2.nextpoint();
214      On (SW_RATIONAL);
215    }
216    E1.nextpoint();
217  }
218  return 0;
219}
220
221//G is assumed to be bivariate, irreducible over Q, primitive wrt x and y, compressed
222CFAFList absBiFactorizeMain (const CanonicalForm& G, bool full)
223{
224  CanonicalForm F= bCommonDen (G)*G;
225  Off (SW_RATIONAL);
226  F /= icontent (F);
227  On (SW_RATIONAL);
228  CFArray eval;
229  int minTdeg, tdegF= totaldegree (F);
230  CanonicalForm Fp, smallestFactor;
231  int p;
232  CFFList factors;
233  Variable alpha;
234  bool rec= false;
235  Variable x= Variable (1);
236  Variable y= Variable (2);
237  CanonicalForm bufF= F;
238  CFFListIterator iter;
239differentevalpoint:
240  while (1)
241  {
242    TIMING_START (fac_evalpoint);
243    p= choosePoint (F, tdegF, eval, rec);
244    TIMING_END_AND_PRINT (fac_evalpoint, "time to find eval point: ");
245
246    setCharacteristic (p);
247    Fp=F.mapinto();
248    factors= factorize (Fp);
249
250    if (factors.getFirst().factor().inCoeffDomain())
251      factors.removeFirst();
252
253    if (factors.length() == 1 && factors.getFirst().exp() == 1)
254    {
255      if (absIrredTest (Fp))
256      {
257        setCharacteristic(0);
258        return CFAFList (CFAFactor (G, 1, 1));
259      }
260      else
261      {
262        setCharacteristic (0);
263        if (modularIrredTestWithShift (F))
264        {
265          return CFAFList (CFAFactor (G, 1, 1));
266        }
267        rec= true;
268        continue;
269      }
270    }
271    iter= factors;
272    smallestFactor= iter.getItem().factor();
273    while (smallestFactor.isUnivariate() && iter.hasItem())
274    {
275      iter++;
276      if (!iter.hasItem())
277        break;
278      smallestFactor= iter.getItem().factor();
279    }
280
281    minTdeg= totaldegree (smallestFactor);
282    if (iter.hasItem())
283      iter++;
284    for (; iter.hasItem(); iter++)
285    {
286      if (!iter.getItem().factor().isUnivariate() &&
287          (totaldegree (iter.getItem().factor()) < minTdeg))
288      {
289        smallestFactor= iter.getItem().factor();
290        minTdeg= totaldegree (smallestFactor);
291      }
292    }
293    if (tdegF % minTdeg == 0)
294      break;
295    setCharacteristic(0);
296    rec=true;
297  }
298  CanonicalForm Gp= Fp/smallestFactor;
299  Gp= Gp /Lc (Gp);
300
301  CanonicalForm Gpy= Gp (eval[0].mapinto(), 1);
302  CanonicalForm smallestFactorEvaly= smallestFactor (eval[0].mapinto(),1);
303  CanonicalForm Gpx= Gp (eval[1].mapinto(), 2);
304  CanonicalForm smallestFactorEvalx= smallestFactor (eval[1].mapinto(),2);
305
306  bool xValid= !(Gpx.inCoeffDomain() || smallestFactorEvalx.inCoeffDomain() ||
307               !gcd (Gpx, smallestFactorEvalx).inCoeffDomain());
308  bool yValid= !(Gpy.inCoeffDomain() || smallestFactorEvaly.inCoeffDomain() ||
309               !gcd (Gpy, smallestFactorEvaly).inCoeffDomain());
310  if (!xValid && !yValid)
311  {
312    rec= true;
313    setCharacteristic (0);
314    goto differentevalpoint;
315  }
316
317  setCharacteristic (0);
318
319  CanonicalForm mipo;
320
321  int loop, i;
322  if (xValid && yValid)
323  {
324    loop= 3;
325    i=1;
326  }
327  else if (xValid)
328  {
329    loop= 3;
330    i=2;
331  }
332  else
333  {
334    loop= 2;
335    i=1;
336  }
337
338  CFArray mipos= CFArray (loop-i);
339  for (; i < loop; i++)
340  {
341    CanonicalForm Fi= F(eval[i-1],i);
342
343    int s= tdegF/minTdeg + 1;
344    CanonicalForm bound= power (maxNorm (Fi), 2*(s-1));
345    bound *= power (CanonicalForm (s),s-1);
346    bound *= power (CanonicalForm (2), ((s-1)*(s-1))/2); //possible int overflow
347
348    CanonicalForm B = p;
349    long k = 1L;
350    while ( B < bound ) {
351      B *= p;
352      k++;
353    }
354
355    //TODO take floor (log2(k))
356    k= k+1;
357#ifdef HAVE_FLINT
358    fmpz_poly_t FLINTFi;
359    convertFacCF2Fmpz_poly_t (FLINTFi, Fi);
360    setCharacteristic (p);
361    nmod_poly_t FLINTFpi, FLINTGpi;
362    if (i == 2)
363    {
364      convertFacCF2nmod_poly_t (FLINTFpi,
365                                smallestFactorEvalx/lc (smallestFactorEvalx));
366      convertFacCF2nmod_poly_t (FLINTGpi, Gpx/lc (Gpx));
367    }
368    else
369    {
370      convertFacCF2nmod_poly_t (FLINTFpi,
371                                smallestFactorEvaly/lc (smallestFactorEvaly));
372      convertFacCF2nmod_poly_t (FLINTGpi, Gpy/lc (Gpy));
373    }
374    nmod_poly_factor_t nmodFactors;
375    nmod_poly_factor_init (nmodFactors);
376    nmod_poly_factor_insert (nmodFactors, FLINTFpi, 1L);
377    nmod_poly_factor_insert (nmodFactors, FLINTGpi, 1L);
378
379    long * link= new long [2];
380    fmpz_poly_t *v= new fmpz_poly_t[2];
381    fmpz_poly_t *w= new fmpz_poly_t[2];
382    fmpz_poly_init(v[0]);
383    fmpz_poly_init(v[1]);
384    fmpz_poly_init(w[0]);
385    fmpz_poly_init(w[1]);
386
387    fmpz_poly_factor_t liftedFactors;
388    fmpz_poly_factor_init (liftedFactors);
389    _fmpz_poly_hensel_start_lift (liftedFactors, link, v, w, FLINTFi,
390                                  nmodFactors, k);
391
392    nmod_poly_factor_clear (nmodFactors);
393    nmod_poly_clear (FLINTFpi);
394    nmod_poly_clear (FLINTGpi);
395
396    setCharacteristic(0);
397    CanonicalForm liftedSmallestFactor=
398    convertFmpz_poly_t2FacCF ((fmpz_poly_t &)liftedFactors->p[0],x);
399
400    CanonicalForm otherFactor=
401    convertFmpz_poly_t2FacCF ((fmpz_poly_t &)liftedFactors->p[1],x);
402    modpk pk= modpk (p, k);
403#else
404    modpk pk= modpk (p, k);
405    ZZX NTLFi=convertFacCF2NTLZZX (pk (Fi*pk.inverse (lc(Fi))));
406    setCharacteristic (p);
407
408    if (fac_NTL_char != p)
409    {
410      fac_NTL_char= p;
411      zz_p::init (p);
412    }
413    zz_pX NTLFpi, NTLGpi;
414    if (i == 2)
415    {
416      NTLFpi= convertFacCF2NTLzzpX (smallestFactorEvalx/lc (smallestFactorEvalx));
417      NTLGpi= convertFacCF2NTLzzpX (Gpx/lc (Gpx));
418    }
419    else
420    {
421      NTLFpi= convertFacCF2NTLzzpX (smallestFactorEvaly/lc (smallestFactorEvaly));
422      NTLGpi= convertFacCF2NTLzzpX (Gpy/lc (Gpy));
423    }
424    vec_zz_pX modFactors;
425    modFactors.SetLength(2);
426    modFactors[0]= NTLFpi;
427    modFactors[1]= NTLGpi;
428    vec_ZZX liftedFactors;
429    MultiLift (liftedFactors, modFactors, NTLFi, (long) k);
430    setCharacteristic(0);
431    CanonicalForm liftedSmallestFactor=
432                  convertNTLZZX2CF (liftedFactors[0], x);
433
434    CanonicalForm otherFactor=
435                  convertNTLZZX2CF (liftedFactors[1], x);
436#endif
437
438    Off (SW_SYMMETRIC_FF);
439    liftedSmallestFactor= pk (liftedSmallestFactor);
440    if (i == 2)
441      liftedSmallestFactor= liftedSmallestFactor (eval[0],1);
442    else
443      liftedSmallestFactor= liftedSmallestFactor (eval[1],1);
444
445    On (SW_SYMMETRIC_FF);
446    CFMatrix M= CFMatrix (s, s);
447    M(s,s)= power (CanonicalForm (p), k);
448    for (int j= 1; j < s; j++)
449    {
450      M (j,j)= 1;
451      M (j+1,j)= -liftedSmallestFactor;
452    }
453
454    mat_ZZ NTLM= *convertFacCFMatrix2NTLmat_ZZ (M);
455
456    ZZ det;
457
458    transpose (NTLM, NTLM);
459    (void) LLL (det, NTLM, 1L, 1L); //use floating point LLL ?
460    transpose (NTLM, NTLM);
461    M= *convertNTLmat_ZZ2FacCFMatrix (NTLM);
462
463    mipo= 0;
464    for (int j= 1; j <= s; j++)
465      mipo += M (j,1)*power (x,s-j);
466
467    CFFList mipoFactors= factorize (mipo);
468    mipoFactors.removeFirst();
469
470#ifdef HAVE_FLINT
471    fmpz_poly_clear (v[0]);
472    fmpz_poly_clear (v[1]);
473    fmpz_poly_clear (w[0]);
474    fmpz_poly_clear (w[1]);
475    delete [] v;
476    delete [] w;
477    delete [] link;
478    fmpz_poly_factor_clear (liftedFactors);
479#endif
480
481    if (mipoFactors.length() > 1 ||
482        (mipoFactors.length() == 1 && mipoFactors.getFirst().exp() > 1))
483    {
484      if (i+1 >= loop && ((loop-i == 1) || (loop-i==2 && mipos[0].isZero())))
485      {
486        rec=true;
487        goto differentevalpoint;
488      }
489    }
490    else
491      mipos[loop-i-1]= mipo;
492  }
493
494  On (SW_RATIONAL);
495  if (xValid && yValid && !mipos[0].isZero() && !mipos[1].isZero())
496  {
497    if (maxNorm (mipos[0]) < maxNorm (mipos[1]))
498      alpha= rootOf (mipos[0]);
499    else
500      alpha= rootOf (mipos[1]);
501  }
502  else if (xValid && yValid)
503  {
504    if (mipos[0].isZero())
505      alpha= rootOf (mipos[1]);
506    else
507      alpha= rootOf (mipos[0]);
508  }
509  else
510    alpha= rootOf (mipo);
511
512  CanonicalForm F1;
513  CFFList QaF1Factors;
514  int wrongMipo= 0;
515  if (xValid && yValid)
516  {
517    if (degree (F,1) > minTdeg)
518      F1= F (eval[1], 2);
519    else
520      F1= F (eval[0], 1);
521  }
522  else if (xValid)
523    F1= F (eval[1], 2);
524  else
525    F1= F (eval[0], 1);
526
527  bool swap= false;
528  if (F1.level() == 2)
529  {
530    swap= true;
531    F1=swapvar (F1, x, y);
532    F= swapvar (F, x, y);
533  }
534
535  QaF1Factors= factorize (F1, alpha);
536  if (QaF1Factors.getFirst().factor().inCoeffDomain())
537    QaF1Factors.removeFirst();
538  for (iter= QaF1Factors; iter.hasItem(); iter++)
539  {
540    if (degree (iter.getItem().factor()) > minTdeg)
541      wrongMipo++;
542  }
543
544  if (wrongMipo == QaF1Factors.length())
545  {
546    if (xValid && yValid && !mipos[0].isZero() && !mipos[1].isZero())
547    {
548      if (maxNorm (mipos[0]) < maxNorm (mipos[1])) //try the other minpoly
549        alpha= rootOf (mipos[1]);
550      else
551        alpha= rootOf (mipos[0]);
552    }
553    else
554    {
555      rec= true;
556      F= bufF;
557      goto differentevalpoint;
558    }
559
560    wrongMipo= 0;
561    QaF1Factors= factorize (F1, alpha);
562    if (QaF1Factors.getFirst().factor().inCoeffDomain())
563      QaF1Factors.removeFirst();
564    for (iter= QaF1Factors; iter.hasItem(); iter++)
565    {
566      if (degree (iter.getItem().factor()) > minTdeg)
567        wrongMipo++;
568    }
569    if (wrongMipo == QaF1Factors.length())
570    {
571      rec= true;
572      F= bufF;
573      goto differentevalpoint;
574    }
575  }
576
577  CanonicalForm evaluation;
578  if (swap)
579    evaluation= eval[0];
580  else
581    evaluation= eval[1];
582
583  F *= bCommonDen (F);
584  F= F (y + evaluation, y);
585
586  int liftBound= degree (F,y) + 1;
587
588  modpk b= modpk();
589
590  CanonicalForm den= 1;
591
592  mipo= getMipo (alpha);
593
594  CFList uniFactors;
595  for (iter=QaF1Factors; iter.hasItem(); iter++)
596    uniFactors.append (iter.getItem().factor());
597
598  F /= Lc (F1);
599  ZZX NTLmipo= convertFacCF2NTLZZX (mipo);
600  ZZX NTLLcf= convertFacCF2NTLZZX (Lc (F*bCommonDen (F)));
601  ZZ NTLf= resultant (NTLmipo, NTLLcf);
602  ZZ NTLD= discriminant (NTLmipo);
603  den= abs (convertZZ2CF (NTLD*NTLf));
604
605  // make factors elements of Z(a)[x] disable for modularDiophant
606  CanonicalForm multiplier= 1;
607  for (CFListIterator i= uniFactors; i.hasItem(); i++)
608  {
609    multiplier *= bCommonDen (i.getItem());
610    i.getItem()= i.getItem()*bCommonDen(i.getItem());
611  }
612  F *= multiplier;
613  F *= bCommonDen (F);
614
615  Off (SW_RATIONAL);
616  int ii= 0;
617  CanonicalForm discMipo= convertZZ2CF (NTLD);
618  findGoodPrime (bufF*discMipo,ii);
619  findGoodPrime (F1*discMipo,ii);
620  findGoodPrime (F*discMipo,ii);
621
622  int pp=cf_getBigPrime(ii);
623  b = coeffBound( F, pp, mipo );
624  modpk bb= coeffBound (F1, pp, mipo);
625  if (bb.getk() > b.getk() ) b=bb;
626    bb= coeffBound (F, pp, mipo);
627  if (bb.getk() > b.getk() ) b=bb;
628
629  ExtensionInfo dummy= ExtensionInfo (alpha, false);
630  DegreePattern degs= DegreePattern (uniFactors);
631
632  bool earlySuccess= false;
633  CFList earlyFactors;
634  uniFactors= henselLiftAndEarly
635              (F, earlySuccess, earlyFactors, degs, liftBound,
636               uniFactors, dummy, evaluation, b, den);
637
638  DEBOUTLN (cerr, "lifted factors= " << uniFactors);
639
640  CanonicalForm MODl= power (y, liftBound);
641
642  On (SW_RATIONAL);
643  F *= bCommonDen (F);
644  Off (SW_RATIONAL);
645
646  CFList biFactors;
647
648  biFactors= factorRecombination (uniFactors, F, MODl, degs, evaluation, 1,
649                                  uniFactors.length()/2, b, den);
650
651  On (SW_RATIONAL);
652
653  if (earlySuccess)
654    biFactors= Union (earlyFactors, biFactors);
655  else if (!earlySuccess && degs.getLength() == 1)
656    biFactors= earlyFactors;
657
658  bool swap2= false;
659  appendSwapDecompress (biFactors, CFList(), CFList(), swap, swap2, CFMap());
660  if (isOn (SW_RATIONAL))
661    normalize (biFactors);
662
663  CFAFList result;
664  bool found= false;
665
666  for (CFListIterator i= biFactors; i.hasItem(); i++)
667  {
668    if (full)
669      result.append (CFAFactor (i.getItem(), getMipo (alpha), 1));
670
671    if (totaldegree (i.getItem()) == minTdeg)
672    {
673      if (!full)
674        result= CFAFList (CFAFactor (i.getItem(), getMipo (alpha), 1));
675      found= true;
676
677      if (!full)
678        break;
679    }
680  }
681
682  if (!found)
683  {
684    rec= true;
685    F= bufF;
686    goto differentevalpoint;
687  }
688
689  return result;
690}
691
692#endif
693
694#ifdef HAVE_NTL
695/// absolute factorization of bivariate poly over Q
696///
697/// @return absFactorize returns a list whose entries contain three entities:
698///         an absolute irreducible factor, an irreducible univariate polynomial
699///         that defines the minimal field extension over which the irreducible
700///         factor is defined and the multiplicity of the absolute irreducible
701///         factor
702CFAFList absFactorize (const CanonicalForm& G ///<[in] bivariate poly over Q
703                      )
704{
705  //TODO handle homogeneous input
706  ASSERT (getNumVars (G) <= 2, "expected bivariate input");
707  ASSERT (getCharacteristic() == 0, "expected poly over Q");
708
709  CFMap N;
710  CanonicalForm F= compress (G, N);
711  bool isRat= isOn (SW_RATIONAL);
712  if (isRat)
713    F *= bCommonDen (F);
714
715  Off (SW_RATIONAL);
716  F /= icontent (F);
717  if (isRat)
718    On (SW_RATIONAL);
719
720  CanonicalForm contentX= content (F, 1);
721  CanonicalForm contentY= content (F, 2);
722  F /= (contentX*contentY);
723  CFAFList contentXFactors, contentYFactors;
724  contentXFactors= uniAbsFactorize (contentX);
725  contentYFactors= uniAbsFactorize (contentY);
726
727  if (contentXFactors.getFirst().factor().inCoeffDomain())
728    contentXFactors.removeFirst();
729  if (contentYFactors.getFirst().factor().inCoeffDomain())
730    contentYFactors.removeFirst();
731  if (F.inCoeffDomain())
732  {
733    CFAFList result;
734    for (CFAFListIterator i= contentXFactors; i.hasItem(); i++)
735      result.append (CFAFactor (N (i.getItem().factor()), i.getItem().minpoly(),
736                                i.getItem().exp()));
737    for (CFAFListIterator i= contentYFactors; i.hasItem(); i++)
738      result.append (CFAFactor (N (i.getItem().factor()),i.getItem().minpoly(),
739                                i.getItem().exp()));
740    normalize (result);
741    result.insert (CFAFactor (Lc (G), 1, 1));
742    return result;
743  }
744  CFFList rationalFactors= factorize (F);
745
746  CFAFList result, resultBuf;
747
748  CFAFListIterator iter;
749  CFFListIterator i= rationalFactors;
750  i++;
751  for (; i.hasItem(); i++)
752  {
753    resultBuf= absFactorizeMain (i.getItem().factor());
754    for (iter= resultBuf; iter.hasItem(); iter++)
755      iter.getItem()= CFAFactor (iter.getItem().factor(),
756                                 iter.getItem().minpoly(), i.getItem().exp());
757    result= Union (result, resultBuf);
758  }
759
760  for (CFAFListIterator i= result; i.hasItem(); i++)
761    i.getItem()= CFAFactor (N (i.getItem().factor()), i.getItem().minpoly(),
762                            i.getItem().exp());
763  for (CFAFListIterator i= contentXFactors; i.hasItem(); i++)
764    result.append (CFAFactor (N(i.getItem().factor()), i.getItem().minpoly(),
765                              i.getItem().exp()));
766  for (CFAFListIterator i= contentYFactors; i.hasItem(); i++)
767    result.append (CFAFactor (N(i.getItem().factor()), i.getItem().minpoly(),
768                              i.getItem().exp()));
769  normalize (result);
770  result.insert (CFAFactor (Lc(G), 1, 1));
771
772  return result;
773}
774#endif
775
776
Note: See TracBrowser for help on using the repository browser.