My Project
Loading...
Searching...
No Matches
cfModGcd.cc
Go to the documentation of this file.
1// -*- c++ -*-
2//*****************************************************************************
3/** @file cfModGcd.cc
4 *
5 * This file implements the GCD of two polynomials over \f$ F_{p} \f$ ,
6 * \f$ F_{p}(\alpha ) \f$, GF or Z based on Alg. 7.1. and 7.2. as described in
7 * "Algorithms for Computer Algebra" by Geddes, Czapor, Labahn via modular
8 * computations. And sparse modular variants as described in Zippel
9 * "Effective Polynomial Computation", deKleine, Monagan, Wittkopf
10 * "Algorithms for the non-monic case of the sparse modular GCD algorithm" and
11 * Javadi "A new solution to the normalization problem"
12 *
13 * @author Martin Lee
14 * @date 22.10.2009
15 *
16 * @par Copyright:
17 * (c) by The SINGULAR Team, see LICENSE file
18 *
19**/
20//*****************************************************************************
21
22
23#include "config.h"
24
25#include <math.h>
26
27#include "cf_assert.h"
28#include "debug.h"
29#include "timing.h"
30
31#include "canonicalform.h"
32#include "cfGcdUtil.h"
33#include "cf_map.h"
34#include "cf_util.h"
35#include "cf_irred.h"
37#include "cf_random.h"
38#include "cf_reval.h"
39#include "facHensel.h"
40#include "cf_iter.h"
41#include "cfNewtonPolygon.h"
42#include "cf_algorithm.h"
43#include "cf_primes.h"
44
45// inline helper functions:
46#include "cf_map_ext.h"
47
48#ifdef HAVE_NTL
49#include "NTLconvert.h"
50#endif
51
52#ifdef HAVE_FLINT
53#include "FLINTconvert.h"
54#endif
55
56#include "cfModGcd.h"
57
65
66bool
70{
72 if (LCCand*abs (LC (coF)) == abs (LC (F)))
73 {
74 if (LCCand*abs (LC (coG)) == abs (LC (G)))
75 {
76 if (abs (cand)*abs (coF) == abs (F))
77 {
78 if (abs (cand)*abs (coG) == abs (G))
79 return true;
80 }
81 return false;
82 }
83 return false;
84 }
85 return false;
86}
87
88#if defined(HAVE_NTL)|| defined(HAVE_FLINT)
89
90/// compressing two polynomials F and G, M is used for compressing,
91/// N to reverse the compression
92int myCompress (const CanonicalForm& F, const CanonicalForm& G, CFMap & M,
93 CFMap & N, bool topLevel)
94{
95 int n= tmax (F.level(), G.level());
96 int * degsf= NEW_ARRAY(int,n + 1);
97 int * degsg= NEW_ARRAY(int,n + 1);
98
99 for (int i = n; i >= 0; i--)
100 degsf[i]= degsg[i]= 0;
101
102 degsf= degrees (F, degsf);
103 degsg= degrees (G, degsg);
104
105 int both_non_zero= 0;
106 int f_zero= 0;
107 int g_zero= 0;
108 int both_zero= 0;
109
110 if (topLevel)
111 {
112 for (int i= 1; i <= n; i++)
113 {
114 if (degsf[i] != 0 && degsg[i] != 0)
115 {
117 continue;
118 }
119 if (degsf[i] == 0 && degsg[i] != 0 && i <= G.level())
120 {
121 f_zero++;
122 continue;
123 }
124 if (degsg[i] == 0 && degsf[i] && i <= F.level())
125 {
126 g_zero++;
127 continue;
128 }
129 }
130
131 if (both_non_zero == 0)
132 {
135 return 0;
136 }
137
138 // map Variables which do not occur in both polynomials to higher levels
139 int k= 1;
140 int l= 1;
141 for (int i= 1; i <= n; i++)
142 {
143 if (degsf[i] != 0 && degsg[i] == 0 && i <= F.level())
144 {
145 if (k + both_non_zero != i)
146 {
147 M.newpair (Variable (i), Variable (k + both_non_zero));
148 N.newpair (Variable (k + both_non_zero), Variable (i));
149 }
150 k++;
151 }
152 if (degsf[i] == 0 && degsg[i] != 0 && i <= G.level())
153 {
154 if (l + g_zero + both_non_zero != i)
155 {
156 M.newpair (Variable (i), Variable (l + g_zero + both_non_zero));
157 N.newpair (Variable (l + g_zero + both_non_zero), Variable (i));
158 }
159 l++;
160 }
161 }
162
163 // sort Variables x_{i} in increasing order of
164 // min(deg_{x_{i}}(f),deg_{x_{i}}(g))
165 int m= tmax (F.level(), G.level());
166 int min_max_deg;
168 l= 0;
169 int i= 1;
170 while (k > 0)
171 {
172 if (degsf [i] != 0 && degsg [i] != 0)
174 else
175 min_max_deg= 0;
176 while (min_max_deg == 0)
177 {
178 i++;
179 if (degsf [i] != 0 && degsg [i] != 0)
181 else
182 min_max_deg= 0;
183 }
184 for (int j= i + 1; j <= m; j++)
185 {
186 if (degsf[j] != 0 && degsg [j] != 0 &&
187 tmax (degsf[j],degsg[j]) <= min_max_deg)
188 {
190 l= j;
191 }
192 }
193 if (l != 0)
194 {
195 if (l != k)
196 {
197 M.newpair (Variable (l), Variable(k));
198 N.newpair (Variable (k), Variable(l));
199 degsf[l]= 0;
200 degsg[l]= 0;
201 l= 0;
202 }
203 else
204 {
205 degsf[l]= 0;
206 degsg[l]= 0;
207 l= 0;
208 }
209 }
210 else if (l == 0)
211 {
212 if (i != k)
213 {
214 M.newpair (Variable (i), Variable (k));
215 N.newpair (Variable (k), Variable (i));
216 degsf[i]= 0;
217 degsg[i]= 0;
218 }
219 else
220 {
221 degsf[i]= 0;
222 degsg[i]= 0;
223 }
224 i++;
225 }
226 k--;
227 }
228 }
229 else
230 {
231 //arrange Variables such that no gaps occur
232 for (int i= 1; i <= n; i++)
233 {
234 if (degsf[i] == 0 && degsg[i] == 0)
235 {
236 both_zero++;
237 continue;
238 }
239 else
240 {
241 if (both_zero != 0)
242 {
243 M.newpair (Variable (i), Variable (i - both_zero));
244 N.newpair (Variable (i - both_zero), Variable (i));
245 }
246 }
247 }
248 }
249
252
253 return 1;
254}
255
256static inline CanonicalForm
257uni_content (const CanonicalForm & F);
258
261{
262 if (F.inCoeffDomain())
263 return F.genOne();
264 if (F.level() == x.level() && F.isUnivariate())
265 return F;
266 if (F.level() != x.level() && F.isUnivariate())
267 return F.genOne();
268
269 if (x.level() != 1)
270 {
271 CanonicalForm f= swapvar (F, x, Variable (1));
273 return swapvar (result, x, Variable (1));
274 }
275 else
276 return uni_content (F);
277}
278
279/// compute the content of F, where F is considered as an element of
280/// \f$ R[x_{1}][x_{2},\ldots ,x_{n}] \f$
281static inline CanonicalForm
283{
284 if (F.inBaseDomain())
285 return F.genOne();
286 if (F.level() == 1 && F.isUnivariate())
287 return F;
288 if (F.level() != 1 && F.isUnivariate())
289 return F.genOne();
290 if (degree (F,1) == 0) return F.genOne();
291
292 int l= F.level();
293 if (l == 2)
294 return content(F);
295 else
296 {
297 CanonicalForm pol, c = 0;
298 CFIterator i = F;
299 for (; i.hasTerms(); i++)
300 {
301 pol= i.coeff();
303 c= gcd (c, pol);
304 if (c.isOne())
305 return c;
306 }
307 return c;
308 }
309}
310
314 CanonicalForm& ppF, CanonicalForm& ppG, const int d)
315{
317 contentF= 1;
318 contentG= 1;
319 ppF= F;
320 ppG= G;
322 for (int i= 1; i <= d; i++)
323 {
329 ppF /= uniContentF;
330 ppG /= uniContentG;
331 result *= gcdcFcG;
332 }
333 return result;
334}
335
336/// compute the leading coefficient of F, where F is considered as an element
337/// of \f$ R[x_{1}][x_{2},\ldots ,x_{n}] \f$, order on
338/// \f$ Mon (x_{2},\ldots ,x_{n}) \f$ is dp.
339static inline
341{
342 if (F.level() > 1)
343 {
344 Variable x= Variable (2);
345 int deg= totaldegree (F, x, F.mvar());
346 for (CFIterator i= F; i.hasTerms(); i++)
347 {
348 if (i.exp() + totaldegree (i.coeff(), x, i.coeff().mvar()) == deg)
349 return uni_lcoeff (i.coeff());
350 }
351 }
352 return F;
353}
354
355/// Newton interpolation - Incremental algorithm.
356/// Given a list of values alpha_i and a list of polynomials u_i, 1 <= i <= n,
357/// computes the interpolation polynomial assuming that
358/// the polynomials in u are the results of evaluating the variabe x
359/// of the unknown polynomial at the alpha values.
360/// This incremental version receives only the values of alpha_n and u_n and
361/// the previous interpolation polynomial for points 1 <= i <= n-1, and computes
362/// the polynomial interpolating in all the points.
363/// newtonPoly must be equal to (x - alpha_1) * ... * (x - alpha_{n-1})
364static inline CanonicalForm
375
376/// compute a random element a of \f$ F_{p}(\alpha ) \f$ ,
377/// s.t. F(a) \f$ \neq 0 \f$ , F is a univariate polynomial, returns
378/// fail if there are no field elements left which have not been used before
379static inline CanonicalForm
380randomElement (const CanonicalForm & F, const Variable & alpha, CFList & list,
381 bool & fail)
382{
383 fail= false;
384 Variable x= F.mvar();
388 mipo= getMipo (alpha);
389 int p= getCharacteristic ();
390 int d= degree (mipo);
391 double bound= pow ((double) p, (double) d);
392 do
393 {
394 if (list.length() == bound)
395 {
396 fail= true;
397 break;
398 }
399 if (list.length() < p)
400 {
401 random= genFF.generate();
402 while (find (list, random))
403 random= genFF.generate();
404 }
405 else
406 {
407 random= genAlgExt.generate();
408 while (find (list, random))
409 random= genAlgExt.generate();
410 }
411 if (F (random, x) == 0)
412 {
413 list.append (random);
414 continue;
415 }
416 } while (find (list, random));
417 return random;
418}
419
420static inline
422{
423 int i, m;
424 // extension of F_p needed
425 if (alpha.level() == 1)
426 {
427 i= 1;
428 m= 2;
429 } //extension of F_p(alpha)
430 if (alpha.level() != 1)
431 {
432 i= 4;
433 m= degree (getMipo (alpha));
434 }
435 #ifdef HAVE_FLINT
441 #else
443 {
445 zz_p::init (getCharacteristic());
446 }
450 #endif
451 return rootOf (newMipo);
452}
453
454#if defined(HAVE_NTL) || defined(HAVE_FLINT)
456modGCDFq (const CanonicalForm& F, const CanonicalForm& G,
458 Variable & alpha, CFList& l, bool& topLevel);
459#endif
460
461#if defined(HAVE_NTL) || defined(HAVE_FLINT)
464 Variable & alpha, CFList& l, bool& topLevel)
465{
468 topLevel);
469 return result;
470}
471#endif
472
473/// GCD of F and G over \f$ F_{p}(\alpha ) \f$ ,
474/// l and topLevel are only used internally, output is monic
475/// based on Alg. 7.2. as described in "Algorithms for
476/// Computer Algebra" by Geddes, Czapor, Labahn
477#if defined(HAVE_NTL) || defined(HAVE_FLINT)
481 Variable & alpha, CFList& l, bool& topLevel)
482{
483 CanonicalForm A= F;
485 if (F.isZero() && degree(G) > 0)
486 {
487 coF= 0;
488 coG= Lc (G);
489 return G/Lc(G);
490 }
491 else if (G.isZero() && degree (F) > 0)
492 {
493 coF= Lc (F);
494 coG= 0;
495 return F/Lc(F);
496 }
497 else if (F.isZero() && G.isZero())
498 {
499 coF= coG= 0;
500 return F.genOne();
501 }
502 if (F.inBaseDomain() || G.inBaseDomain())
503 {
504 coF= F;
505 coG= G;
506 return F.genOne();
507 }
508 if (F.isUnivariate() && fdivides(F, G, coG))
509 {
510 coF= Lc (F);
511 return F/Lc(F);
512 }
513 if (G.isUnivariate() && fdivides(G, F, coF))
514 {
515 coG= Lc (G);
516 return G/Lc(G);
517 }
518 if (F == G)
519 {
520 coF= coG= Lc (F);
521 return F/Lc(F);
522 }
523
524 CFMap M,N;
525 int best_level= myCompress (A, B, M, N, topLevel);
526
527 if (best_level == 0)
528 {
529 coF= F;
530 coG= G;
531 return B.genOne();
532 }
533
534 A= M(A);
535 B= M(B);
536
537 Variable x= Variable(1);
538
539 //univariate case
540 if (A.isUnivariate() && B.isUnivariate())
541 {
543 coF= N (A/result);
544 coG= N (B/result);
545 return N (result);
546 }
547
548 CanonicalForm cA, cB; // content of A and B
549 CanonicalForm ppA, ppB; // primitive part of A and B
551
552 cA = uni_content (A);
553 cB = uni_content (B);
554 gcdcAcB= gcd (cA, cB);
555 ppA= A/cA;
556 ppB= B/cB;
557
558 CanonicalForm lcA, lcB; // leading coefficients of A and B
560
561 lcA= uni_lcoeff (ppA);
562 lcB= uni_lcoeff (ppB);
563
564 gcdlcAlcB= gcd (lcA, lcB);
565
566 int d= totaldegree (ppA, Variable(2), Variable (ppA.level()));
567
568 if (d == 0)
569 {
570 coF= N (ppA*(cA/gcdcAcB));
571 coG= N (ppB*(cB/gcdcAcB));
572 return N(gcdcAcB);
573 }
574
575 int d0= totaldegree (ppB, Variable(2), Variable (ppB.level()));
576 if (d0 < d)
577 d= d0;
578 if (d == 0)
579 {
580 coF= N (ppA*(cA/gcdcAcB));
581 coG= N (ppB*(cB/gcdcAcB));
582 return N(gcdcAcB);
583 }
584
587 coG_m, ppCoF, ppCoG;
588
589 newtonPoly= 1;
590 m= gcdlcAlcB;
591 G_m= 0;
592 coF= 0;
593 coG= 0;
594 H= 0;
595 bool fail= false;
596 topLevel= false;
597 bool inextension= false;
601 CFList source, dest;
602 int bound1= degree (ppA, 1);
603 int bound2= degree (ppB, 1);
604 do
605 {
607 if (fail)
608 {
609 source= CFList();
610 dest= CFList();
611
614 bool prim_fail= false;
617 if (V_buf4 == alpha)
619
620 if (V_buf3 != V_buf4)
621 {
627 source, dest);
631 source, dest);
634 for (CFListIterator i= l; i.hasItem(); i++)
635 i.getItem()= mapDown (i.getItem(), prim_elem, im_prim_elem, V_buf4,
636 source, dest);
637 }
638
639 ASSERT (!prim_fail, "failure in integer factorizer");
640 if (prim_fail)
641 ; //ERROR
642 else
644
645 if (V_buf4 == alpha)
647 else
649 im_prim_elem, source, dest);
650 DEBOUTLN (cerr, "getMipo (V_buf4)= " << getMipo (V_buf4));
651 DEBOUTLN (cerr, "getMipo (V_buf2)= " << getMipo (V_buf2));
652 inextension= true;
653 for (CFListIterator i= l; i.hasItem(); i++)
654 i.getItem()= mapUp (i.getItem(), V_buf4, V_buf, prim_elem,
655 im_prim_elem, source, dest);
661 source, dest);
665 source, dest);
668
669 fail= false;
671 DEBOUTLN (cerr, "fail= " << fail);
672 CFList list;
677 list, topLevel);
679 "time for recursive call: ");
680 DEBOUTLN (cerr, "G_random_element= " << G_random_element);
681 V_buf4= V_buf;
682 }
683 else
684 {
685 CFList list;
690 list, topLevel);
692 "time for recursive call: ");
693 DEBOUTLN (cerr, "G_random_element= " << G_random_element);
694 }
695
696 if (!G_random_element.inCoeffDomain())
698 Variable (G_random_element.level()));
699 else
700 d0= 0;
701
702 if (d0 == 0)
703 {
704 if (inextension)
705 {
706 CFList u, v;
709 prune1 (alpha);
710 }
711 coF= N (ppA*(cA/gcdcAcB));
712 coG= N (ppB*(cB/gcdcAcB));
713 return N(gcdcAcB);
714 }
715 if (d0 > d)
716 {
717 if (!find (l, random_element))
719 continue;
720 }
721
729
730 if (!G_random_element.inCoeffDomain())
732 Variable (G_random_element.level()));
733 else
734 d0= 0;
735
736 if (d0 < d)
737 {
738 m= gcdlcAlcB;
739 newtonPoly= 1;
740 G_m= 0;
741 d= d0;
742 coF_m= 0;
743 coG_m= 0;
744 }
745
751 "time for newton interpolation: ");
752
753 //termination test
754 if ((uni_lcoeff (H) == gcdlcAlcB) || (G_m == H))
755 {
757 if (gcdlcAlcB.isOne())
758 cH= 1;
759 else
760 cH= uni_content (H);
761 ppH= H/cH;
762 ppH /= Lc (ppH);
766 ppCoF= coF/ccoF;
767 ppCoG= coG/ccoG;
768 if (inextension)
769 {
770 if (((degree (ppCoF,1)+degree (ppH,1) == bound1) &&
771 (degree (ppCoG,1)+degree (ppH,1) == bound2) &&
773 (fdivides (ppH, ppA, ppCoF) && fdivides (ppH, ppB, ppCoG)))
774 {
775 CFList u, v;
776 DEBOUTLN (cerr, "ppH before mapDown= " << ppH);
780 DEBOUTLN (cerr, "ppH after mapDown= " << ppH);
781 coF= N ((cA/gcdcAcB)*ppCoF);
782 coG= N ((cB/gcdcAcB)*ppCoG);
784 "time for successful termination test Fq: ");
785 prune1 (alpha);
786 return N(gcdcAcB*ppH);
787 }
788 }
789 else if (((degree (ppCoF,1)+degree (ppH,1) == bound1) &&
790 (degree (ppCoG,1)+degree (ppH,1) == bound2) &&
792 (fdivides (ppH, ppA, ppCoF) && fdivides (ppH, ppB, ppCoG)))
793 {
794 coF= N ((cA/gcdcAcB)*ppCoF);
795 coG= N ((cB/gcdcAcB)*ppCoG);
797 "time for successful termination test Fq: ");
798 return N(gcdcAcB*ppH);
799 }
801 "time for unsuccessful termination test Fq: ");
802 }
803
804 G_m= H;
805 coF_m= coF;
806 coG_m= coG;
808 m= m*(x - random_element);
809 if (!find (l, random_element))
811 } while (1);
812}
813#endif
814
815/// compute a random element a of GF, s.t. F(a) \f$ \neq 0 \f$ , F is a
816/// univariate polynomial, returns fail if there are no field elements left
817/// which have not been used before
818static inline
820GFRandomElement (const CanonicalForm& F, CFList& list, bool& fail)
821{
822 fail= false;
823 Variable x= F.mvar();
826 int p= getCharacteristic();
827 int d= getGFDegree();
828 int bound= ipower (p, d);
829 do
830 {
831 if (list.length() == bound)
832 {
833 fail= true;
834 break;
835 }
836 if (list.length() < 1)
837 random= 0;
838 else
839 {
840 random= genGF.generate();
841 while (find (list, random))
842 random= genGF.generate();
843 }
844 if (F (random, x) == 0)
845 {
846 list.append (random);
847 continue;
848 }
849 } while (find (list, random));
850 return random;
851}
852
854modGCDGF (const CanonicalForm& F, const CanonicalForm& G,
856 CFList& l, bool& topLevel);
857
860 bool& topLevel)
861{
864 return result;
865}
866
867/// GCD of F and G over GF, based on Alg. 7.2. as described in "Algorithms for
868/// Computer Algebra" by Geddes, Czapor, Labahn
869/// Usually this algorithm will be faster than modGCDFq since GF has
870/// faster field arithmetics, however it might fail if the input is large since
871/// the size of the base field is bounded by 2^16, output is monic
875 CFList& l, bool& topLevel)
876{
877 CanonicalForm A= F;
879 if (F.isZero() && degree(G) > 0)
880 {
881 coF= 0;
882 coG= Lc (G);
883 return G/Lc(G);
884 }
885 else if (G.isZero() && degree (F) > 0)
886 {
887 coF= Lc (F);
888 coG= 0;
889 return F/Lc(F);
890 }
891 else if (F.isZero() && G.isZero())
892 {
893 coF= coG= 0;
894 return F.genOne();
895 }
896 if (F.inBaseDomain() || G.inBaseDomain())
897 {
898 coF= F;
899 coG= G;
900 return F.genOne();
901 }
902 if (F.isUnivariate() && fdivides(F, G, coG))
903 {
904 coF= Lc (F);
905 return F/Lc(F);
906 }
907 if (G.isUnivariate() && fdivides(G, F, coF))
908 {
909 coG= Lc (G);
910 return G/Lc(G);
911 }
912 if (F == G)
913 {
914 coF= coG= Lc (F);
915 return F/Lc(F);
916 }
917
918 CFMap M,N;
919 int best_level= myCompress (A, B, M, N, topLevel);
920
921 if (best_level == 0)
922 {
923 coF= F;
924 coG= G;
925 return B.genOne();
926 }
927
928 A= M(A);
929 B= M(B);
930
931 Variable x= Variable(1);
932
933 //univariate case
934 if (A.isUnivariate() && B.isUnivariate())
935 {
937 coF= N (A/result);
938 coG= N (B/result);
939 return N (result);
940 }
941
942 CanonicalForm cA, cB; // content of A and B
943 CanonicalForm ppA, ppB; // primitive part of A and B
945
946 cA = uni_content (A);
947 cB = uni_content (B);
948 gcdcAcB= gcd (cA, cB);
949 ppA= A/cA;
950 ppB= B/cB;
951
952 CanonicalForm lcA, lcB; // leading coefficients of A and B
954
955 lcA= uni_lcoeff (ppA);
956 lcB= uni_lcoeff (ppB);
957
958 gcdlcAlcB= gcd (lcA, lcB);
959
960 int d= totaldegree (ppA, Variable(2), Variable (ppA.level()));
961 if (d == 0)
962 {
963 coF= N (ppA*(cA/gcdcAcB));
964 coG= N (ppB*(cB/gcdcAcB));
965 return N(gcdcAcB);
966 }
967 int d0= totaldegree (ppB, Variable(2), Variable (ppB.level()));
968 if (d0 < d)
969 d= d0;
970 if (d == 0)
971 {
972 coF= N (ppA*(cA/gcdcAcB));
973 coG= N (ppB*(cB/gcdcAcB));
974 return N(gcdcAcB);
975 }
976
979 coG_m, ppCoF, ppCoG;
980
981 newtonPoly= 1;
982 m= gcdlcAlcB;
983 G_m= 0;
984 coF= 0;
985 coG= 0;
986 H= 0;
987 bool fail= false;
988 topLevel= false;
989 bool inextension= false;
990 int p=-1;
991 int k= getGFDegree();
992 int kk;
993 int expon;
994 char gf_name_buf= gf_name;
995 int bound1= degree (ppA, 1);
996 int bound2= degree (ppB, 1);
997 do
998 {
1000 if (fail)
1001 {
1003 expon= 2;
1004 kk= getGFDegree();
1005 if (ipower (p, kk*expon) < (1 << 16))
1006 setCharacteristic (p, kk*(int)expon, 'b');
1007 else
1008 {
1009 expon= (int) floor((log ((double)((1<<16) - 1)))/(log((double)p)*kk));
1010 ASSERT (expon >= 2, "not enough points in modGCDGF");
1011 setCharacteristic (p, (int)(kk*expon), 'b');
1012 }
1013 inextension= true;
1014 fail= false;
1015 for (CFListIterator i= l; i.hasItem(); i++)
1016 i.getItem()= GFMapUp (i.getItem(), kk);
1017 m= GFMapUp (m, kk);
1018 G_m= GFMapUp (G_m, kk);
1020 coF_m= GFMapUp (coF_m, kk);
1021 coG_m= GFMapUp (coG_m, kk);
1022 ppA= GFMapUp (ppA, kk);
1023 ppB= GFMapUp (ppB, kk);
1025 lcA= GFMapUp (lcA, kk);
1026 lcB= GFMapUp (lcB, kk);
1028 DEBOUTLN (cerr, "fail= " << fail);
1029 CFList list;
1033 list, topLevel);
1035 "time for recursive call: ");
1036 DEBOUTLN (cerr, "G_random_element= " << G_random_element);
1037 }
1038 else
1039 {
1040 CFList list;
1044 list, topLevel);
1046 "time for recursive call: ");
1047 DEBOUTLN (cerr, "G_random_element= " << G_random_element);
1048 }
1049
1050 if (!G_random_element.inCoeffDomain())
1052 Variable (G_random_element.level()));
1053 else
1054 d0= 0;
1055
1056 if (d0 == 0)
1057 {
1058 if (inextension)
1059 {
1060 ppA= GFMapDown (ppA, k);
1061 ppB= GFMapDown (ppB, k);
1063 }
1064 coF= N (ppA*(cA/gcdcAcB));
1065 coG= N (ppB*(cB/gcdcAcB));
1066 return N(gcdcAcB);
1067 }
1068 if (d0 > d)
1069 {
1070 if (!find (l, random_element))
1072 continue;
1073 }
1074
1078
1083
1084 if (!G_random_element.inCoeffDomain())
1086 Variable (G_random_element.level()));
1087 else
1088 d0= 0;
1089
1090 if (d0 < d)
1091 {
1092 m= gcdlcAlcB;
1093 newtonPoly= 1;
1094 G_m= 0;
1095 d= d0;
1096 coF_m= 0;
1097 coG_m= 0;
1098 }
1099
1105 "time for newton interpolation: ");
1106
1107 //termination test
1108 if ((uni_lcoeff (H) == gcdlcAlcB) || (G_m == H))
1109 {
1111 if (gcdlcAlcB.isOne())
1112 cH= 1;
1113 else
1114 cH= uni_content (H);
1115 ppH= H/cH;
1116 ppH /= Lc (ppH);
1120 ppCoF= coF/ccoF;
1121 ppCoG= coG/ccoG;
1122 if (inextension)
1123 {
1124 if (((degree (ppCoF,1)+degree (ppH,1) == bound1) &&
1125 (degree (ppCoG,1)+degree (ppH,1) == bound2) &&
1127 (fdivides (ppH, ppA, ppCoF) && fdivides (ppH, ppB, ppCoG)))
1128 {
1129 DEBOUTLN (cerr, "ppH before mapDown= " << ppH);
1130 ppH= GFMapDown (ppH, k);
1131 ppCoF= GFMapDown (ppCoF, k);
1132 ppCoG= GFMapDown (ppCoG, k);
1133 DEBOUTLN (cerr, "ppH after mapDown= " << ppH);
1134 coF= N ((cA/gcdcAcB)*ppCoF);
1135 coG= N ((cB/gcdcAcB)*ppCoG);
1138 "time for successful termination GF: ");
1139 return N(gcdcAcB*ppH);
1140 }
1141 }
1142 else
1143 {
1144 if (((degree (ppCoF,1)+degree (ppH,1) == bound1) &&
1145 (degree (ppCoG,1)+degree (ppH,1) == bound2) &&
1147 (fdivides (ppH, ppA, ppCoF) && fdivides (ppH, ppB, ppCoG)))
1148 {
1149 coF= N ((cA/gcdcAcB)*ppCoF);
1150 coG= N ((cB/gcdcAcB)*ppCoG);
1152 "time for successful termination GF: ");
1153 return N(gcdcAcB*ppH);
1154 }
1155 }
1157 "time for unsuccessful termination GF: ");
1158 }
1159
1160 G_m= H;
1161 coF_m= coF;
1162 coG_m= coG;
1164 m= m*(x - random_element);
1165 if (!find (l, random_element))
1167 } while (1);
1168}
1169
1170static inline
1173{
1174 fail= false;
1175 Variable x= F.mvar();
1178 int p= getCharacteristic();
1179 int bound= p;
1180 do
1181 {
1182 if (list.length() == bound)
1183 {
1184 fail= true;
1185 break;
1186 }
1187 if (list.length() < 1)
1188 random= 0;
1189 else
1190 {
1191 random= genFF.generate();
1192 while (find (list, random))
1193 random= genFF.generate();
1194 }
1195 if (F (random, x) == 0)
1196 {
1197 list.append (random);
1198 continue;
1199 }
1200 } while (find (list, random));
1201 return random;
1202}
1203
1204#if defined(HAVE_NTL) || defined(HAVE_FLINT)
1206modGCDFp (const CanonicalForm& F, const CanonicalForm& G,
1208 bool& topLevel, CFList& l);
1209#endif
1210
1211#if defined(HAVE_NTL) || defined(HAVE_FLINT)
1214 bool& topLevel, CFList& l)
1215{
1218 return result;
1219}
1220#endif
1221
1222#if defined(HAVE_NTL) || defined(HAVE_FLINT)
1226 bool& topLevel, CFList& l)
1227{
1228 CanonicalForm A= F;
1229 CanonicalForm B= G;
1230 if (F.isZero() && degree(G) > 0)
1231 {
1232 coF= 0;
1233 coG= Lc (G);
1234 return G/Lc(G);
1235 }
1236 else if (G.isZero() && degree (F) > 0)
1237 {
1238 coF= Lc (F);
1239 coG= 0;
1240 return F/Lc(F);
1241 }
1242 else if (F.isZero() && G.isZero())
1243 {
1244 coF= coG= 0;
1245 return F.genOne();
1246 }
1247 if (F.inBaseDomain() || G.inBaseDomain())
1248 {
1249 coF= F;
1250 coG= G;
1251 return F.genOne();
1252 }
1253 if (F.isUnivariate() && fdivides(F, G, coG))
1254 {
1255 coF= Lc (F);
1256 return F/Lc(F);
1257 }
1258 if (G.isUnivariate() && fdivides(G, F, coF))
1259 {
1260 coG= Lc (G);
1261 return G/Lc(G);
1262 }
1263 if (F == G)
1264 {
1265 coF= coG= Lc (F);
1266 return F/Lc(F);
1267 }
1268 CFMap M,N;
1269 int best_level= myCompress (A, B, M, N, topLevel);
1270
1271 if (best_level == 0)
1272 {
1273 coF= F;
1274 coG= G;
1275 return B.genOne();
1276 }
1277
1278 A= M(A);
1279 B= M(B);
1280
1281 Variable x= Variable (1);
1282
1283 //univariate case
1284 if (A.isUnivariate() && B.isUnivariate())
1285 {
1287 coF= N (A/result);
1288 coG= N (B/result);
1289 return N (result);
1290 }
1291
1292 CanonicalForm cA, cB; // content of A and B
1293 CanonicalForm ppA, ppB; // primitive part of A and B
1295
1296 cA = uni_content (A);
1297 cB = uni_content (B);
1298 gcdcAcB= gcd (cA, cB);
1299 ppA= A/cA;
1300 ppB= B/cB;
1301
1302 CanonicalForm lcA, lcB; // leading coefficients of A and B
1304 lcA= uni_lcoeff (ppA);
1305 lcB= uni_lcoeff (ppB);
1306
1307 gcdlcAlcB= gcd (lcA, lcB);
1308
1309 int d= totaldegree (ppA, Variable (2), Variable (ppA.level()));
1310 int d0;
1311
1312 if (d == 0)
1313 {
1314 coF= N (ppA*(cA/gcdcAcB));
1315 coG= N (ppB*(cB/gcdcAcB));
1316 return N(gcdcAcB);
1317 }
1318
1319 d0= totaldegree (ppB, Variable (2), Variable (ppB.level()));
1320
1321 if (d0 < d)
1322 d= d0;
1323
1324 if (d == 0)
1325 {
1326 coF= N (ppA*(cA/gcdcAcB));
1327 coG= N (ppB*(cB/gcdcAcB));
1328 return N(gcdcAcB);
1329 }
1330
1334
1335 newtonPoly= 1;
1336 m= gcdlcAlcB;
1337 H= 0;
1338 coF= 0;
1339 coG= 0;
1340 G_m= 0;
1342 bool fail= false;
1343 bool inextension= false;
1344 topLevel= false;
1345 CFList source, dest;
1346 int bound1= degree (ppA, 1);
1347 int bound2= degree (ppB, 1);
1348 do
1349 {
1350 if (inextension)
1352 else
1354
1355 if (!fail && !inextension)
1356 {
1357 CFList list;
1362 list);
1364 "time for recursive call: ");
1365 DEBOUTLN (cerr, "G_random_element= " << G_random_element);
1366 }
1367 else if (!fail && inextension)
1368 {
1369 CFList list;
1374 list, topLevel);
1376 "time for recursive call: ");
1377 DEBOUTLN (cerr, "G_random_element= " << G_random_element);
1378 }
1379 else if (fail && !inextension)
1380 {
1381 source= CFList();
1382 dest= CFList();
1383 CFList list;
1385 int deg= 2;
1386 bool initialized= false;
1387 do
1388 {
1389 mipo= randomIrredpoly (deg, x);
1390 if (initialized)
1391 setMipo (alpha, mipo);
1392 else
1393 alpha= rootOf (mipo);
1394 inextension= true;
1395 initialized= true;
1396 fail= false;
1398 deg++;
1399 } while (fail);
1400 list= CFList();
1401 V_buf= alpha;
1402 cleanUp= alpha;
1407 list, topLevel);
1409 "time for recursive call: ");
1410 DEBOUTLN (cerr, "G_random_element= " << G_random_element);
1411 }
1412 else if (fail && inextension)
1413 {
1414 source= CFList();
1415 dest= CFList();
1416
1419 bool prim_fail= false;
1423
1424 if (V_buf3 != alpha)
1425 {
1431 source, dest);
1435 dest);
1438 for (CFListIterator i= l; i.hasItem(); i++)
1439 i.getItem()= mapDown (i.getItem(), prim_elem, im_prim_elem, alpha,
1440 source, dest);
1441 }
1442
1443 ASSERT (!prim_fail, "failure in integer factorizer");
1444 if (prim_fail)
1445 ; //ERROR
1446 else
1448
1449 DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (alpha));
1450 DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (V_buf2));
1451
1452 for (CFListIterator i= l; i.hasItem(); i++)
1453 i.getItem()= mapUp (i.getItem(), alpha, V_buf, prim_elem,
1454 im_prim_elem, source, dest);
1460 source, dest);
1464 source, dest);
1467 fail= false;
1469 DEBOUTLN (cerr, "fail= " << fail);
1470 CFList list;
1475 list, topLevel);
1477 "time for recursive call: ");
1478 DEBOUTLN (cerr, "G_random_element= " << G_random_element);
1479 alpha= V_buf;
1480 }
1481
1482 if (!G_random_element.inCoeffDomain())
1484 Variable (G_random_element.level()));
1485 else
1486 d0= 0;
1487
1488 if (d0 == 0)
1489 {
1490 if (inextension)
1491 prune (cleanUp);
1492 coF= N (ppA*(cA/gcdcAcB));
1493 coG= N (ppB*(cB/gcdcAcB));
1494 return N(gcdcAcB);
1495 }
1496
1497 if (d0 > d)
1498 {
1499 if (!find (l, random_element))
1501 continue;
1502 }
1503
1506
1511
1512 if (!G_random_element.inCoeffDomain())
1514 Variable (G_random_element.level()));
1515 else
1516 d0= 0;
1517
1518 if (d0 < d)
1519 {
1520 m= gcdlcAlcB;
1521 newtonPoly= 1;
1522 G_m= 0;
1523 d= d0;
1524 coF_m= 0;
1525 coG_m= 0;
1526 }
1527
1533 "time for newton_interpolation: ");
1534
1535 //termination test
1536 if ((uni_lcoeff (H) == gcdlcAlcB) || (G_m == H))
1537 {
1539 if (gcdlcAlcB.isOne())
1540 cH= 1;
1541 else
1542 cH= uni_content (H);
1543 ppH= H/cH;
1544 ppH /= Lc (ppH);
1548 ppCoF= coF/ccoF;
1549 ppCoG= coG/ccoG;
1550 DEBOUTLN (cerr, "ppH= " << ppH);
1551 if (((degree (ppCoF,1)+degree (ppH,1) == bound1) &&
1552 (degree (ppCoG,1)+degree (ppH,1) == bound2) &&
1554 (fdivides (ppH, ppA, ppCoF) && fdivides (ppH, ppB, ppCoG)))
1555 {
1556 if (inextension)
1557 prune (cleanUp);
1558 coF= N ((cA/gcdcAcB)*ppCoF);
1559 coG= N ((cB/gcdcAcB)*ppCoG);
1561 "time for successful termination Fp: ");
1562 return N(gcdcAcB*ppH);
1563 }
1565 "time for unsuccessful termination Fp: ");
1566 }
1567
1568 G_m= H;
1569 coF_m= coF;
1570 coG_m= coG;
1572 m= m*(x - random_element);
1573 if (!find (l, random_element))
1575 } while (1);
1576}
1577#endif
1578
1579CFArray
1581{
1582 int r= M.size();
1583 ASSERT (A.size() == r, "vector does not have right size");
1584
1585 if (r == 1)
1586 {
1587 CFArray result= CFArray (1);
1588 result [0]= A [0] / M [0];
1589 return result;
1590 }
1591 // check solvability
1592 bool notDistinct= false;
1593 for (int i= 0; i < r - 1; i++)
1594 {
1595 for (int j= i + 1; j < r; j++)
1596 {
1597 if (M [i] == M [j])
1598 {
1599 notDistinct= true;
1600 break;
1601 }
1602 }
1603 }
1604 if (notDistinct)
1605 return CFArray();
1606
1608 Variable x= Variable (1);
1609 for (int i= 0; i < r; i++)
1610 master *= x - M [i];
1611 CFList Pj;
1613 for (int i= 0; i < r; i++)
1614 {
1615 tmp= master/(x - M [i]);
1616 tmp /= tmp (M [i], 1);
1617 Pj.append (tmp);
1618 }
1619 CFArray result= CFArray (r);
1620
1622 for (int i= 1; i <= r; i++, j++)
1623 {
1624 tmp= 0;
1625 for (int l= 0; l < A.size(); l++)
1626 tmp += A[l]*j.getItem()[l];
1627 result[i - 1]= tmp;
1628 }
1629 return result;
1630}
1631
1632CFArray
1634{
1635 int r= M.size();
1636 ASSERT (A.size() == r, "vector does not have right size");
1637 if (r == 1)
1638 {
1639 CFArray result= CFArray (1);
1640 result [0]= A[0] / M [0];
1641 return result;
1642 }
1643 // check solvability
1644 bool notDistinct= false;
1645 for (int i= 0; i < r - 1; i++)
1646 {
1647 for (int j= i + 1; j < r; j++)
1648 {
1649 if (M [i] == M [j])
1650 {
1651 notDistinct= true;
1652 break;
1653 }
1654 }
1655 }
1656 if (notDistinct)
1657 return CFArray();
1658
1660 Variable x= Variable (1);
1661 for (int i= 0; i < r; i++)
1662 master *= x - M [i];
1663 master *= x;
1664 CFList Pj;
1666 for (int i= 0; i < r; i++)
1667 {
1668 tmp= master/(x - M [i]);
1669 tmp /= tmp (M [i], 1);
1670 Pj.append (tmp);
1671 }
1672
1673 CFArray result= CFArray (r);
1674
1676 for (int i= 1; i <= r; i++, j++)
1677 {
1678 tmp= 0;
1679
1680 for (int l= 1; l <= A.size(); l++)
1681 tmp += A[l - 1]*j.getItem()[l];
1682 result[i - 1]= tmp;
1683 }
1684 return result;
1685}
1686
1687/// M in row echolon form, rk rank of M
1688CFArray
1689readOffSolution (const CFMatrix& M, const long rk)
1690{
1693 for (int i= rk; i >= 1; i--)
1694 {
1695 tmp3= 0;
1696 tmp1= M (i, M.columns());
1697 for (int j= M.columns() - 1; j >= 1; j--)
1698 {
1699 tmp2= M (i, j);
1700 if (j == i)
1701 break;
1702 else
1703 tmp3 += tmp2*result[j - 1];
1704 }
1705 result[i - 1]= (tmp1 - tmp3)/tmp2;
1706 }
1707 return result;
1708}
1709
1710CFArray
1712{
1713 CFArray result= CFArray (M.rows());
1715 int k;
1716 for (int i= M.rows(); i >= 1; i--)
1717 {
1718 tmp3= 0;
1719 tmp1= L[i - 1];
1720 k= 0;
1721 for (int j= M.columns(); j >= 1; j--, k++)
1722 {
1723 tmp2= M (i, j);
1724 if (j == i)
1725 break;
1726 else
1727 {
1728 if (k > partialSol.size() - 1)
1729 tmp3 += tmp2*result[j - 1];
1730 else
1731 tmp3 += tmp2*partialSol[partialSol.size() - k - 1];
1732 }
1733 }
1734 result [i - 1]= (tmp1 - tmp3)/tmp2;
1735 }
1736 return result;
1737}
1738
1739long
1741{
1742 ASSERT (L.size() <= M.rows(), "dimension exceeded");
1743 CFMatrix *N;
1744 N= new CFMatrix (M.rows(), M.columns() + 1);
1745
1746 for (int i= 1; i <= M.rows(); i++)
1747 for (int j= 1; j <= M.columns(); j++)
1748 (*N) (i, j)= M (i, j);
1749
1750 int j= 1;
1751 for (int i= 0; i < L.size(); i++, j++)
1752 (*N) (j, M.columns() + 1)= L[i];
1753#ifdef HAVE_FLINT
1756 long rk= nmod_mat_rref (FLINTN);
1757
1758 delete N;
1761#else
1762 int p= getCharacteristic ();
1763 if (fac_NTL_char != p)
1764 {
1765 fac_NTL_char= p;
1766 zz_p::init (p);
1767 }
1769 delete N;
1770 long rk= gauss (*NTLN);
1771
1773 delete NTLN;
1774#endif
1775
1776 L= CFArray (M.rows());
1777 for (int i= 0; i < M.rows(); i++)
1778 L[i]= (*N) (i + 1, M.columns() + 1);
1779 M= (*N) (1, M.rows(), 1, M.columns());
1780 delete N;
1781 return rk;
1782}
1783
1784long
1786{
1787 ASSERT (L.size() <= M.rows(), "dimension exceeded");
1788 CFMatrix *N;
1789 N= new CFMatrix (M.rows(), M.columns() + 1);
1790
1791 for (int i= 1; i <= M.rows(); i++)
1792 for (int j= 1; j <= M.columns(); j++)
1793 (*N) (i, j)= M (i, j);
1794
1795 int j= 1;
1796 for (int i= 0; i < L.size(); i++, j++)
1797 (*N) (j, M.columns() + 1)= L[i];
1798 #ifdef HAVE_FLINT
1799 // convert mipo
1805 // convert matrix
1808 // rank
1809 long rk= fq_nmod_mat_rref (FLINTN,ctx);
1810 // clean up
1813 #elif defined(HAVE_NTL)
1814 int p= getCharacteristic ();
1815 if (fac_NTL_char != p)
1816 {
1817 fac_NTL_char= p;
1818 zz_p::init (p);
1819 }
1821 zz_pE::init (NTLMipo);
1823 long rk= gauss (*NTLN);
1825 delete NTLN;
1826 #else
1827 factoryError("NTL/FLINT missing: gaussianElimFq");
1828 #endif
1829 delete N;
1830
1831 M= (*N) (1, M.rows(), 1, M.columns());
1832 L= CFArray (M.rows());
1833 for (int i= 0; i < M.rows(); i++)
1834 L[i]= (*N) (i + 1, M.columns() + 1);
1835
1836 delete N;
1837 return rk;
1838}
1839
1840CFArray
1842{
1843 ASSERT (L.size() <= M.rows(), "dimension exceeded");
1844 CFMatrix *N;
1845 N= new CFMatrix (M.rows(), M.columns() + 1);
1846
1847 for (int i= 1; i <= M.rows(); i++)
1848 for (int j= 1; j <= M.columns(); j++)
1849 (*N) (i, j)= M (i, j);
1850
1851 int j= 1;
1852 for (int i= 0; i < L.size(); i++, j++)
1853 (*N) (j, M.columns() + 1)= L[i];
1854
1855#ifdef HAVE_FLINT
1858 long rk= nmod_mat_rref (FLINTN);
1859#else
1860 int p= getCharacteristic ();
1861 if (fac_NTL_char != p)
1862 {
1863 fac_NTL_char= p;
1864 zz_p::init (p);
1865 }
1867 long rk= gauss (*NTLN);
1868#endif
1869 delete N;
1870 if (rk != M.columns())
1871 {
1872#ifdef HAVE_FLINT
1874#else
1875 delete NTLN;
1876#endif
1877 return CFArray();
1878 }
1879#ifdef HAVE_FLINT
1882#else
1884 delete NTLN;
1885#endif
1887
1888 delete N;
1889 return A;
1890}
1891
1892CFArray
1893solveSystemFq (const CFMatrix& M, const CFArray& L, const Variable& alpha)
1894{
1895 ASSERT (L.size() <= M.rows(), "dimension exceeded");
1896 CFMatrix *N;
1897 N= new CFMatrix (M.rows(), M.columns() + 1);
1898
1899 for (int i= 1; i <= M.rows(); i++)
1900 for (int j= 1; j <= M.columns(); j++)
1901 (*N) (i, j)= M (i, j);
1902 int j= 1;
1903 for (int i= 0; i < L.size(); i++, j++)
1904 (*N) (j, M.columns() + 1)= L[i];
1905 #ifdef HAVE_FLINT
1906 // convert mipo
1912 // convert matrix
1915 // rank
1916 long rk= fq_nmod_mat_rref (FLINTN,ctx);
1917 #elif defined(HAVE_NTL)
1918 int p= getCharacteristic ();
1919 if (fac_NTL_char != p)
1920 {
1921 fac_NTL_char= p;
1922 zz_p::init (p);
1923 }
1925 zz_pE::init (NTLMipo);
1927 long rk= gauss (*NTLN);
1928 #else
1929 factoryError("NTL/FLINT missing: solveSystemFq");
1930 #endif
1931
1932 delete N;
1933 if (rk != M.columns())
1934 {
1935 #if defined(HAVE_NTL) && !defined(HAVE_FLINT)
1936 delete NTLN;
1937 #endif
1938 return CFArray();
1939 }
1940 #ifdef HAVE_FLINT
1941 // convert and clean up
1945 #elif defined(HAVE_NTL)
1947 delete NTLN;
1948 #endif
1949
1951
1952 delete N;
1953 return A;
1954}
1955#endif
1956
1957CFArray
1959{
1960 if (F.inCoeffDomain())
1961 {
1962 CFArray result= CFArray (1);
1963 result [0]= 1;
1964 return result;
1965 }
1966 if (F.isUnivariate())
1967 {
1968 CFArray result= CFArray (size(F));
1969 int j= 0;
1970 for (CFIterator i= F; i.hasTerms(); i++, j++)
1971 result[j]= power (F.mvar(), i.exp());
1972 return result;
1973 }
1974 int numMon= size (F);
1976 int j= 0;
1978 Variable x= F.mvar();
1980 for (CFIterator i= F; i.hasTerms(); i++)
1981 {
1982 powX= power (x, i.exp());
1983 recResult= getMonoms (i.coeff());
1984 for (int k= 0; k < recResult.size(); k++)
1985 result[j+k]= powX*recResult[k];
1986 j += recResult.size();
1987 }
1988 return result;
1989}
1990
1991#if defined(HAVE_NTL) || defined(HAVE_FLINT)
1992CFArray
1994{
1995 if (F.inCoeffDomain())
1996 {
1997 CFArray result= CFArray (1);
1998 result [0]= F;
1999 return result;
2000 }
2001 if (F.isUnivariate())
2002 {
2003 ASSERT (evalPoints.length() == 1,
2004 "expected an eval point with only one component");
2005 CFArray result= CFArray (size(F));
2006 int j= 0;
2008 for (CFIterator i= F; i.hasTerms(); i++, j++)
2009 result[j]= power (evalPoint, i.exp());
2010 return result;
2011 }
2012 int numMon= size (F);
2014 int j= 0;
2017 buf.removeLast();
2020 for (CFIterator i= F; i.hasTerms(); i++)
2021 {
2022 powEvalPoint= power (evalPoint, i.exp());
2023 recResult= evaluateMonom (i.coeff(), buf);
2024 for (int k= 0; k < recResult.size(); k++)
2026 j += recResult.size();
2027 }
2028 return result;
2029}
2030
2031CFArray
2033{
2034 CFArray result= A.size();
2036 int k;
2037 for (int i= 0; i < A.size(); i++)
2038 {
2039 tmp= A[i];
2040 k= 1;
2041 for (CFListIterator j= evalPoints; j.hasItem(); j++, k++)
2042 tmp= tmp (j.getItem(), k);
2043 result[i]= tmp;
2044 }
2045 return result;
2046}
2047
2048CFList
2051 const CanonicalForm& LCF, const bool& GF,
2052 const Variable& alpha, bool& fail, CFList& list
2053 )
2054{
2055 int k= tmax (F.level(), G.level()) - 1;
2056 Variable x= Variable (1);
2057 CFList result;
2060 int p= getCharacteristic ();
2061 double bound;
2062 if (alpha != Variable (1))
2063 {
2064 bound= pow ((double) p, (double) degree (getMipo(alpha)));
2065 bound= pow (bound, (double) k);
2066 }
2067 else if (GF)
2068 {
2069 bound= pow ((double) p, (double) getGFDegree());
2070 bound= pow ((double) bound, (double) k);
2071 }
2072 else
2073 bound= pow ((double) p, (double) k);
2074
2076 int j;
2077 bool zeroOneOccured= false;
2078 bool allEqual= false;
2080 do
2081 {
2082 random= 0;
2083 // possible overflow if list.length() does not fit into a int
2084 if (list.length() >= bound)
2085 {
2086 fail= true;
2087 break;
2088 }
2089 for (int i= 0; i < k; i++)
2090 {
2091 if (GF)
2092 {
2093 result.append (genGF.generate());
2094 random += result.getLast()*power (x, i);
2095 }
2096 else if (alpha.level() != 1)
2097 {
2099 result.append (genAlgExt.generate());
2100 random += result.getLast()*power (x, i);
2101 }
2102 else
2103 {
2104 result.append (genFF.generate());
2105 random += result.getLast()*power (x, i);
2106 }
2107 if (result.getLast().isOne() || result.getLast().isZero())
2108 zeroOneOccured= true;
2109 }
2110 if (find (list, random))
2111 {
2112 zeroOneOccured= false;
2113 allEqual= false;
2114 result= CFList();
2115 continue;
2116 }
2117 if (zeroOneOccured)
2118 {
2119 list.append (random);
2120 zeroOneOccured= false;
2121 allEqual= false;
2122 result= CFList();
2123 continue;
2124 }
2125 // no zero at this point
2126 if (k > 1)
2127 {
2128 allEqual= true;
2130 buf= iter.coeff();
2131 iter++;
2132 for (; iter.hasTerms(); iter++)
2133 if (buf != iter.coeff())
2134 allEqual= false;
2135 }
2136 if (allEqual)
2137 {
2138 list.append (random);
2139 allEqual= false;
2140 zeroOneOccured= false;
2141 result= CFList();
2142 continue;
2143 }
2144
2145 Feval= F;
2146 Geval= G;
2148 j= 1;
2149 for (CFListIterator i= result; i.hasItem(); i++, j++)
2150 {
2151 Feval= Feval (i.getItem(), j);
2152 Geval= Geval (i.getItem(), j);
2153 LCeval= LCeval (i.getItem(), j);
2154 }
2155
2156 if (LCeval.isZero())
2157 {
2158 if (!find (list, random))
2159 list.append (random);
2160 zeroOneOccured= false;
2161 allEqual= false;
2162 result= CFList();
2163 continue;
2164 }
2165
2166 if (list.length() >= bound)
2167 {
2168 fail= true;
2169 break;
2170 }
2171 } while (find (list, random));
2172
2173 return result;
2174}
2175
2176/// multiply two lists componentwise
2177void mult (CFList& L1, const CFList& L2)
2178{
2179 ASSERT (L1.length() == L2.length(), "lists of the same size expected");
2180
2182 for (CFListIterator i= L1; i.hasItem(); i++, j++)
2183 i.getItem() *= j.getItem();
2184}
2185
2187 CanonicalForm& Beval, const CFList& L)
2188{
2189 Aeval= A;
2190 Beval= B;
2191 int j= 1;
2192 for (CFListIterator i= L; i.hasItem(); i++, j++)
2193 {
2194 Aeval= Aeval (i.getItem(), j);
2195 Beval= Beval (i.getItem(), j);
2196 }
2197}
2198
2201 const CanonicalForm& skeleton, const Variable& alpha,
2203 )
2204{
2205 CanonicalForm A= F;
2206 CanonicalForm B= G;
2207 if (F.isZero() && degree(G) > 0) return G/Lc(G);
2208 else if (G.isZero() && degree (F) > 0) return F/Lc(F);
2209 else if (F.isZero() && G.isZero()) return F.genOne();
2210 if (F.inBaseDomain() || G.inBaseDomain()) return F.genOne();
2211 if (F.isUnivariate() && fdivides(F, G)) return F/Lc(F);
2212 if (G.isUnivariate() && fdivides(G, F)) return G/Lc(G);
2213 if (F == G) return F/Lc(F);
2214
2215 ASSERT (degree (A, 1) == 0, "expected degree (F, 1) == 0");
2216 ASSERT (degree (B, 1) == 0, "expected degree (G, 1) == 0");
2217
2218 CFMap M,N;
2219 int best_level= myCompress (A, B, M, N, false);
2220
2221 if (best_level == 0)
2222 return B.genOne();
2223
2224 A= M(A);
2225 B= M(B);
2226
2227 Variable x= Variable (1);
2228
2229 //univariate case
2230 if (A.isUnivariate() && B.isUnivariate())
2231 return N (gcd (A, B));
2232
2234 CanonicalForm cA, cB; // content of A and B
2235 CanonicalForm ppA, ppB; // primitive part of A and B
2237 cA = uni_content (A);
2238 cB = uni_content (B);
2239 gcdcAcB= gcd (cA, cB);
2240 ppA= A/cA;
2241 ppB= B/cB;
2242
2243 CanonicalForm lcA, lcB; // leading coefficients of A and B
2245 lcA= uni_lcoeff (ppA);
2246 lcB= uni_lcoeff (ppB);
2247
2248 if (fdivides (lcA, lcB))
2249 {
2250 if (fdivides (A, B))
2251 return F/Lc(F);
2252 }
2253 if (fdivides (lcB, lcA))
2254 {
2255 if (fdivides (B, A))
2256 return G/Lc(G);
2257 }
2258
2259 gcdlcAlcB= gcd (lcA, lcB);
2260 int skelSize= size (skel, skel.mvar());
2261
2262 int j= 0;
2263 int biggestSize= 0;
2264
2265 for (CFIterator i= skel; i.hasTerms(); i++, j++)
2266 biggestSize= tmax (biggestSize, size (i.coeff()));
2267
2269
2271 bool evalFail= false;
2272 CFList list;
2273 bool GF= false;
2274 CanonicalForm LCA= LC (A);
2279 CFList source, dest;
2282 for (int i= 0; i < biggestSize; i++)
2283 {
2284 if (i == 0)
2286 list);
2287 else
2288 {
2290 eval (A, B, Aeval, Beval, evalPoints);
2291 }
2292
2293 if (evalFail)
2294 {
2295 if (V_buf.level() != 1)
2296 {
2297 do
2298 {
2301 source= CFList();
2302 dest= CFList();
2303
2304 bool prim_fail= false;
2307 if (V_buf4 == alpha && alpha.level() != 1)
2309
2310 ASSERT (!prim_fail, "failure in integer factorizer");
2311 if (prim_fail)
2312 ; //ERROR
2313 else
2315
2316 DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (V_buf));
2317 DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (V_buf2));
2318
2319 if (V_buf4 == alpha && alpha.level() != 1)
2321 else if (alpha.level() != 1)
2323 prim_elem, im_prim_elem, source, dest);
2324
2325 for (CFListIterator j= list; j.hasItem(); j++)
2326 j.getItem()= mapUp (j.getItem(), V_buf4, V_buf, prim_elem,
2327 im_prim_elem, source, dest);
2328 for (int k= 0; k < i; k++)
2329 {
2330 for (CFListIterator j= pEvalPoints[k]; j.hasItem(); j++)
2331 j.getItem()= mapUp (j.getItem(), V_buf4, V_buf, prim_elem,
2332 im_prim_elem, source, dest);
2334 source, dest);
2335 }
2336
2337 if (alpha.level() != 1)
2338 {
2341 }
2342 V_buf4= V_buf;
2343 evalFail= false;
2345 evalFail, list);
2346 } while (evalFail);
2347 }
2348 else
2349 {
2351 int deg= 2;
2352 bool initialized= false;
2353 do
2354 {
2355 mipo= randomIrredpoly (deg, x);
2356 if (initialized)
2357 setMipo (V_buf, mipo);
2358 else
2359 V_buf= rootOf (mipo);
2360 evalFail= false;
2361 initialized= true;
2363 evalFail, list);
2364 deg++;
2365 } while (evalFail);
2366 V_buf4= V_buf;
2367 }
2368 }
2369
2370 g= gcd (Aeval, Beval);
2371 g /= Lc (g);
2372
2373 if (degree (g) != degree (skel) || (skelSize != size (g)))
2374 {
2375 delete[] pEvalPoints;
2376 fail= true;
2377 if (alpha.level() != 1 && V_buf != alpha)
2378 prune1 (alpha);
2379 return 0;
2380 }
2381 CFIterator l= skel;
2382 for (CFIterator k= g; k.hasTerms(); k++, l++)
2383 {
2384 if (k.exp() != l.exp())
2385 {
2386 delete[] pEvalPoints;
2387 fail= true;
2388 if (alpha.level() != 1 && V_buf != alpha)
2389 prune1 (alpha);
2390 return 0;
2391 }
2392 }
2394 gcds[i]= g;
2395
2396 tmp= 0;
2397 int j= 0;
2398 for (CFListIterator k= evalPoints; k.hasItem(); k++, j++)
2399 tmp += k.getItem()*power (x, j);
2400 list.append (tmp);
2401 }
2402
2403 if (Monoms.size() == 0)
2405
2407 j= 0;
2408 for (CFIterator i= skel; i.hasTerms(); i++, j++)
2409 coeffMonoms[j]= getMonoms (i.coeff());
2410
2411 CFArray* pL= new CFArray [skelSize];
2412 CFArray* pM= new CFArray [skelSize];
2413 for (int i= 0; i < biggestSize; i++)
2414 {
2415 CFIterator l= gcds [i];
2417 for (int k= 0; k < skelSize; k++, l++)
2418 {
2419 if (i == 0)
2420 pL[k]= CFArray (biggestSize);
2421 pL[k] [i]= l.coeff();
2422
2423 if (i == 0)
2425 }
2426 }
2427
2430 int ind= 0;
2432 CFMatrix Mat;
2433 for (int k= 0; k < skelSize; k++)
2434 {
2435 if (biggestSize != coeffMonoms[k].size())
2436 {
2438 for (int i= 0; i < coeffMonoms[k].size(); i++)
2439 bufArray [i]= pL[k] [i];
2441 }
2442 else
2444
2445 if (solution.size() == 0)
2446 {
2447 delete[] pEvalPoints;
2448 delete[] pM;
2449 delete[] pL;
2450 delete[] coeffMonoms;
2451 fail= true;
2452 if (alpha.level() != 1 && V_buf != alpha)
2453 prune1 (alpha);
2454 return 0;
2455 }
2456 for (int l= 0; l < solution.size(); l++)
2457 result += solution[l]*Monoms [ind + l];
2458 ind += solution.size();
2459 }
2460
2461 delete[] pEvalPoints;
2462 delete[] pM;
2463 delete[] pL;
2464 delete[] coeffMonoms;
2465
2466 if (alpha.level() != 1 && V_buf != alpha)
2467 {
2468 CFList u, v;
2470 prune1 (alpha);
2471 }
2472
2473 result= N(result);
2474 if (fdivides (result, F) && fdivides (result, G))
2475 return result;
2476 else
2477 {
2478 fail= true;
2479 return 0;
2480 }
2481}
2482
2485 const CanonicalForm& skeleton, const Variable& alpha,
2487 )
2488{
2489 CanonicalForm A= F;
2490 CanonicalForm B= G;
2491 if (F.isZero() && degree(G) > 0) return G/Lc(G);
2492 else if (G.isZero() && degree (F) > 0) return F/Lc(F);
2493 else if (F.isZero() && G.isZero()) return F.genOne();
2494 if (F.inBaseDomain() || G.inBaseDomain()) return F.genOne();
2495 if (F.isUnivariate() && fdivides(F, G)) return F/Lc(F);
2496 if (G.isUnivariate() && fdivides(G, F)) return G/Lc(G);
2497 if (F == G) return F/Lc(F);
2498
2499 ASSERT (degree (A, 1) == 0, "expected degree (F, 1) == 0");
2500 ASSERT (degree (B, 1) == 0, "expected degree (G, 1) == 0");
2501
2502 CFMap M,N;
2503 int best_level= myCompress (A, B, M, N, false);
2504
2505 if (best_level == 0)
2506 return B.genOne();
2507
2508 A= M(A);
2509 B= M(B);
2510
2511 Variable x= Variable (1);
2512
2513 //univariate case
2514 if (A.isUnivariate() && B.isUnivariate())
2515 return N (gcd (A, B));
2516
2518
2519 CanonicalForm cA, cB; // content of A and B
2520 CanonicalForm ppA, ppB; // primitive part of A and B
2522 cA = uni_content (A);
2523 cB = uni_content (B);
2524 gcdcAcB= gcd (cA, cB);
2525 ppA= A/cA;
2526 ppB= B/cB;
2527
2528 CanonicalForm lcA, lcB; // leading coefficients of A and B
2530 lcA= uni_lcoeff (ppA);
2531 lcB= uni_lcoeff (ppB);
2532
2533 if (fdivides (lcA, lcB))
2534 {
2535 if (fdivides (A, B))
2536 return F/Lc(F);
2537 }
2538 if (fdivides (lcB, lcA))
2539 {
2540 if (fdivides (B, A))
2541 return G/Lc(G);
2542 }
2543
2544 gcdlcAlcB= gcd (lcA, lcB);
2545 int skelSize= size (skel, skel.mvar());
2546
2547 int j= 0;
2548 int biggestSize= 0;
2549 int bufSize;
2550 int numberUni= 0;
2551 for (CFIterator i= skel; i.hasTerms(); i++, j++)
2552 {
2553 bufSize= size (i.coeff());
2555 numberUni += bufSize;
2556 }
2557 numberUni--;
2558 numberUni= (int) ceil ( (double) numberUni / (skelSize - 1));
2560
2561 numberUni= biggestSize + size (LC(skel)) - 1;
2563
2565
2568 bool evalFail= false;
2569 CFList list;
2570 bool GF= false;
2571 CanonicalForm LCA= LC (A);
2576 CFList source, dest;
2579 for (int i= 0; i < biggestSize; i++)
2580 {
2581 if (i == 0)
2582 {
2583 if (getCharacteristic() > 3)
2585 evalFail, list);
2586 else
2587 evalFail= true;
2588
2589 if (evalFail)
2590 {
2591 if (V_buf.level() != 1)
2592 {
2593 do
2594 {
2597 source= CFList();
2598 dest= CFList();
2599
2600 bool prim_fail= false;
2603 if (V_buf4 == alpha && alpha.level() != 1)
2605
2606 ASSERT (!prim_fail, "failure in integer factorizer");
2607 if (prim_fail)
2608 ; //ERROR
2609 else
2611
2612 DEBOUTLN (cerr, "getMipo (V_buf)= " << getMipo (V_buf));
2613 DEBOUTLN (cerr, "getMipo (V_buf2)= " << getMipo (V_buf2));
2614
2615 if (V_buf4 == alpha && alpha.level() != 1)
2617 else if (alpha.level() != 1)
2619 prim_elem, im_prim_elem, source, dest);
2620
2621 for (CFListIterator i= list; i.hasItem(); i++)
2622 i.getItem()= mapUp (i.getItem(), V_buf4, V_buf, prim_elem,
2623 im_prim_elem, source, dest);
2624 if (alpha.level() != 1)
2625 {
2628 }
2629 evalFail= false;
2630 V_buf4= V_buf;
2632 evalFail, list);
2633 } while (evalFail);
2634 }
2635 else
2636 {
2638 int deg= 2;
2639 bool initialized= false;
2640 do
2641 {
2642 mipo= randomIrredpoly (deg, x);
2643 if (initialized)
2644 setMipo (V_buf, mipo);
2645 else
2646 V_buf= rootOf (mipo);
2647 evalFail= false;
2648 initialized= true;
2650 evalFail, list);
2651 deg++;
2652 } while (evalFail);
2653 V_buf4= V_buf;
2654 }
2655 }
2656 }
2657 else
2658 {
2660 eval (A, B, Aeval, Beval, evalPoints);
2661 }
2662
2663 g= gcd (Aeval, Beval);
2664 g /= Lc (g);
2665
2666 if (degree (g) != degree (skel) || (skelSize != size (g)))
2667 {
2668 delete[] pEvalPoints;
2669 fail= true;
2670 if (alpha.level() != 1 && V_buf != alpha)
2671 prune1 (alpha);
2672 return 0;
2673 }
2674 CFIterator l= skel;
2675 for (CFIterator k= g; k.hasTerms(); k++, l++)
2676 {
2677 if (k.exp() != l.exp())
2678 {
2679 delete[] pEvalPoints;
2680 fail= true;
2681 if (alpha.level() != 1 && V_buf != alpha)
2682 prune1 (alpha);
2683 return 0;
2684 }
2685 }
2687 gcds[i]= g;
2688
2689 tmp= 0;
2690 int j= 0;
2691 for (CFListIterator k= evalPoints; k.hasItem(); k++, j++)
2692 tmp += k.getItem()*power (x, j);
2693 list.append (tmp);
2694 }
2695
2696 if (Monoms.size() == 0)
2698
2700
2701 j= 0;
2702 for (CFIterator i= skel; i.hasTerms(); i++, j++)
2703 coeffMonoms[j]= getMonoms (i.coeff());
2704
2706 if (skelSize > 1)
2708 else
2710 int minimalColumns=-1;
2711
2712 CFArray* pM= new CFArray [skelSize];
2713 CFMatrix Mat;
2714 // find the Matrix with minimal number of columns
2715 for (int i= 0; i < skelSize; i++)
2716 {
2717 pM[i]= CFArray (coeffMonoms[i].size());
2718 if (i == 1)
2719 minimalColumns= coeffMonoms[i].size();
2720 if (i > 1)
2721 {
2723 {
2724 minimalColumns= coeffMonoms[i].size();
2726 }
2727 }
2728 }
2729 CFMatrix* pMat= new CFMatrix [2];
2732 CFArray* pL= new CFArray [skelSize];
2733 for (int i= 0; i < biggestSize; i++)
2734 {
2735 CFIterator l= gcds [i];
2737 for (int k= 0; k < skelSize; k++, l++)
2738 {
2739 if (i == 0)
2740 pL[k]= CFArray (biggestSize);
2741 pL[k] [i]= l.coeff();
2742
2743 if (i == 0 && k != 0 && k != minimalColumnsIndex)
2745 else if (i == 0 && (k == 0 || k == minimalColumnsIndex))
2747 if (i > 0 && (k == 0 || k == minimalColumnsIndex))
2749
2750 if (k == 0)
2751 {
2752 if (pMat[k].rows() >= i + 1)
2753 {
2754 for (int ind= 1; ind < coeffEval.size() + 1; ind++)
2755 pMat[k] (i + 1, ind)= coeffEval[ind - 1];
2756 }
2757 }
2758 if (k == minimalColumnsIndex)
2759 {
2760 if (pMat[1].rows() >= i + 1)
2761 {
2762 for (int ind= 1; ind < coeffEval.size() + 1; ind++)
2763 pMat[1] (i + 1, ind)= coeffEval[ind - 1];
2764 }
2765 }
2766 }
2767 }
2768
2771 int ind= 1;
2772 int matRows, matColumns;
2773 matRows= pMat[1].rows();
2774 matColumns= pMat[0].columns() - 1;
2775 matColumns += pMat[1].columns();
2776
2778 for (int i= 1; i <= matRows; i++)
2779 for (int j= 1; j <= pMat[1].columns(); j++)
2780 Mat (i, j)= pMat[1] (i, j);
2781
2782 for (int j= pMat[1].columns() + 1; j <= matColumns;
2783 j++, ind++)
2784 {
2785 for (int i= 1; i <= matRows; i++)
2786 Mat (i, j)= (-pMat [0] (i, ind + 1))*pL[minimalColumnsIndex] [i - 1];
2787 }
2788
2789 CFArray firstColumn= CFArray (pMat[0].rows());
2790 for (int i= 0; i < pMat[0].rows(); i++)
2791 firstColumn [i]= pMat[0] (i + 1, 1);
2792
2794
2795 for (int i= 0; i < biggestSize; i++)
2797
2798 CFMatrix bufMat= pMat[1];
2799 pMat[1]= Mat;
2800
2801 if (V_buf.level() != 1)
2804 else
2807
2808 if (solution.size() == 0)
2809 { //Javadi's method failed, try deKleine, Monagan, Wittkopf instead
2810 CFMatrix bufMat0= pMat[0];
2811 delete [] pMat;
2812 pMat= new CFMatrix [skelSize];
2817 {
2820 bufGcds= gcds;
2822 for (int i= 0; i < biggestSize; i++)
2823 {
2825 gcds[i]= bufGcds[i];
2826 }
2827 for (int i= 0; i < biggestSize2 - biggestSize; i++)
2828 {
2830 eval (A, B, Aeval, Beval, evalPoints);
2831 g= gcd (Aeval, Beval);
2832 g /= Lc (g);
2833 if (degree (g) != degree (skel) || (skelSize != size (g)))
2834 {
2835 delete[] pEvalPoints;
2836 delete[] pMat;
2837 delete[] pL;
2838 delete[] coeffMonoms;
2839 delete[] pM;
2840 if (bufpEvalPoints != NULL)
2841 delete [] bufpEvalPoints;
2842 fail= true;
2843 if (alpha.level() != 1 && V_buf != alpha)
2844 prune1 (alpha);
2845 return 0;
2846 }
2847 CFIterator l= skel;
2848 for (CFIterator k= g; k.hasTerms(); k++, l++)
2849 {
2850 if (k.exp() != l.exp())
2851 {
2852 delete[] pEvalPoints;
2853 delete[] pMat;
2854 delete[] pL;
2855 delete[] coeffMonoms;
2856 delete[] pM;
2857 if (bufpEvalPoints != NULL)
2858 delete [] bufpEvalPoints;
2859 fail= true;
2860 if (alpha.level() != 1 && V_buf != alpha)
2861 prune1 (alpha);
2862 return 0;
2863 }
2864 }
2866 gcds[i + biggestSize]= g;
2867 }
2868 }
2869 for (int i= 0; i < biggestSize; i++)
2870 {
2871 CFIterator l= gcds [i];
2873 for (int k= 1; k < skelSize; k++, l++)
2874 {
2875 if (i == 0)
2877 if (k == minimalColumnsIndex)
2878 continue;
2880 if (pMat[k].rows() >= i + 1)
2881 {
2882 for (int ind= 1; ind < coeffEval.size() + 1; ind++)
2883 pMat[k] (i + 1, ind)= coeffEval[ind - 1];
2884 }
2885 }
2886 }
2887 Mat= bufMat0;
2889 for (int j= 1; j <= Mat.rows(); j++)
2890 for (int k= 1; k <= Mat.columns(); k++)
2891 pMat [0] (j,k)= Mat (j,k);
2892 Mat= bufMat;
2893 for (int j= 1; j <= Mat.rows(); j++)
2894 for (int k= 1; k <= Mat.columns(); k++)
2895 pMat [minimalColumnsIndex] (j,k)= Mat (j,k);
2896 // write old matrix entries into new matrices
2897 for (int i= 0; i < skelSize; i++)
2898 {
2899 bufArray= pL[i];
2901 for (int j= 0; j < bufArray.size(); j++)
2902 pL[i] [j]= bufArray [j];
2903 }
2904 //write old vector entries into new and add new entries to old matrices
2905 for (int i= 0; i < biggestSize2 - biggestSize; i++)
2906 {
2909 for (int k= 0; k < skelSize; k++, l++)
2910 {
2911 pL[k] [i + biggestSize]= l.coeff();
2913 if (pMat[k].rows() >= i + biggestSize + 1)
2914 {
2915 for (int ind= 1; ind < coeffEval.size() + 1; ind++)
2916 pMat[k] (i + biggestSize + 1, ind)= coeffEval[ind - 1];
2917 }
2918 }
2919 }
2920 // begin new
2921 for (int i= 0; i < skelSize; i++)
2922 {
2923 if (pL[i].size() > 1)
2924 {
2925 for (int j= 2; j <= pMat[i].rows(); j++)
2926 pMat[i] (j, coeffMonoms[i].size() + j - 1)=
2927 -pL[i] [j - 1];
2928 }
2929 }
2930
2932 matRows= 0;
2933 for (int i= 0; i < skelSize; i++)
2934 {
2935 if (V_buf.level() == 1)
2936 (void) gaussianElimFp (pMat[i], pL[i]);
2937 else
2938 (void) gaussianElimFq (pMat[i], pL[i], V_buf);
2939
2940 if (pMat[i] (coeffMonoms[i].size(), coeffMonoms[i].size()) == 0)
2941 {
2942 delete[] pEvalPoints;
2943 delete[] pMat;
2944 delete[] pL;
2945 delete[] coeffMonoms;
2946 delete[] pM;
2947 if (bufpEvalPoints != NULL)
2948 delete [] bufpEvalPoints;
2949 fail= true;
2950 if (alpha.level() != 1 && V_buf != alpha)
2951 prune1 (alpha);
2952 return 0;
2953 }
2954 matRows += pMat[i].rows() - coeffMonoms[i].size();
2955 }
2956
2959 ind= 0;
2962 for (int i= 0; i < skelSize; i++)
2963 {
2964 if (coeffMonoms[i].size() + 1 >= pMat[i].rows() || coeffMonoms[i].size() + 1 >= pMat[i].columns())
2965 {
2966 delete[] pEvalPoints;
2967 delete[] pMat;
2968 delete[] pL;
2969 delete[] coeffMonoms;
2970 delete[] pM;
2971 if (bufpEvalPoints != NULL)
2972 delete [] bufpEvalPoints;
2973 fail= true;
2974 if (alpha.level() != 1 && V_buf != alpha)
2975 prune1 (alpha);
2976 return 0;
2977 }
2978 bufMat= pMat[i] (coeffMonoms[i].size() + 1, pMat[i].rows(),
2979 coeffMonoms[i].size() + 1, pMat[i].columns());
2980
2981 for (int j= 1; j <= bufMat.rows(); j++)
2982 for (int k= 1; k <= bufMat.columns(); k++)
2983 Mat (j + ind, k)= bufMat(j, k);
2984 bufArray2= coeffMonoms[i].size();
2985 for (int j= 1; j <= pMat[i].rows(); j++)
2986 {
2987 if (j > coeffMonoms[i].size())
2988 bufArray [j-coeffMonoms[i].size() + ind - 1]= pL[i] [j - 1];
2989 else
2990 bufArray2 [j - 1]= pL[i] [j - 1];
2991 }
2992 pL[i]= bufArray2;
2993 ind += bufMat.rows();
2994 pMat[i]= pMat[i] (1, coeffMonoms[i].size(), 1, pMat[i].columns());
2995 }
2996
2997 if (V_buf.level() != 1)
2999 else
3001
3002 if (solution.size() == 0)
3003 {
3004 delete[] pEvalPoints;
3005 delete[] pMat;
3006 delete[] pL;
3007 delete[] coeffMonoms;
3008 delete[] pM;
3009 if (bufpEvalPoints != NULL)
3010 delete [] bufpEvalPoints;
3011 fail= true;
3012 if (alpha.level() != 1 && V_buf != alpha)
3013 prune1 (alpha);
3014 return 0;
3015 }
3016
3017 ind= 0;
3018 result= 0;
3020 for (int i= 0; i < skelSize; i++)
3021 {
3023 for (int i= 0; i < bufSolution.size(); i++, ind++)
3025 }
3026 if (alpha.level() != 1 && V_buf != alpha)
3027 {
3028 CFList u, v;
3030 prune1 (alpha);
3031 }
3032 result= N(result);
3033 delete[] pEvalPoints;
3034 delete[] pMat;
3035 delete[] pL;
3036 delete[] coeffMonoms;
3037 delete[] pM;
3038
3039 if (bufpEvalPoints != NULL)
3040 delete [] bufpEvalPoints;
3041 if (fdivides (result, F) && fdivides (result, G))
3042 return result;
3043 else
3044 {
3045 fail= true;
3046 return 0;
3047 }
3048 } // end of deKleine, Monagan & Wittkopf
3049
3050 result += Monoms[0];
3051 int ind2= 0, ind3= 2;
3052 ind= 0;
3053 for (int l= 0; l < minimalColumnsIndex; l++)
3054 ind += coeffMonoms[l].size();
3055 for (int l= solution.size() - pMat[0].columns() + 1; l < solution.size();
3056 l++, ind2++, ind3++)
3057 {
3058 result += solution[l]*Monoms [1 + ind2];
3059 for (int i= 0; i < pMat[0].rows(); i++)
3060 firstColumn[i] += solution[l]*pMat[0] (i + 1, ind3);
3061 }
3062 for (int l= 0; l < solution.size() + 1 - pMat[0].columns(); l++)
3063 result += solution[l]*Monoms [ind + l];
3064 ind= coeffMonoms[0].size();
3065 for (int k= 1; k < skelSize; k++)
3066 {
3067 if (k == minimalColumnsIndex)
3068 {
3069 ind += coeffMonoms[k].size();
3070 continue;
3071 }
3072 if (k != minimalColumnsIndex)
3073 {
3074 for (int i= 0; i < biggestSize; i++)
3075 pL[k] [i] *= firstColumn [i];
3076 }
3077
3079 {
3081 for (int i= 0; i < bufArray.size(); i++)
3082 bufArray [i]= pL[k] [i];
3084 }
3085 else
3087
3088 if (solution.size() == 0)
3089 {
3090 delete[] pEvalPoints;
3091 delete[] pMat;
3092 delete[] pL;
3093 delete[] coeffMonoms;
3094 delete[] pM;
3095 fail= true;
3096 if (alpha.level() != 1 && V_buf != alpha)
3097 prune1 (alpha);
3098 return 0;
3099 }
3100 if (k != minimalColumnsIndex)
3101 {
3102 for (int l= 0; l < solution.size(); l++)
3103 result += solution[l]*Monoms [ind + l];
3104 ind += solution.size();
3105 }
3106 }
3107
3108 delete[] pEvalPoints;
3109 delete[] pMat;
3110 delete[] pL;
3111 delete[] pM;
3112 delete[] coeffMonoms;
3113
3114 if (alpha.level() != 1 && V_buf != alpha)
3115 {
3116 CFList u, v;
3118 prune1 (alpha);
3119 }
3120 result= N(result);
3121
3122 if (fdivides (result, F) && fdivides (result, G))
3123 return result;
3124 else
3125 {
3126 fail= true;
3127 return 0;
3128 }
3129}
3130
3132 const Variable & alpha, CFList& l, bool& topLevel)
3133{
3134 CanonicalForm A= F;
3135 CanonicalForm B= G;
3136 if (F.isZero() && degree(G) > 0) return G/Lc(G);
3137 else if (G.isZero() && degree (F) > 0) return F/Lc(F);
3138 else if (F.isZero() && G.isZero()) return F.genOne();
3139 if (F.inBaseDomain() || G.inBaseDomain()) return F.genOne();
3140 if (F.isUnivariate() && fdivides(F, G)) return F/Lc(F);
3141 if (G.isUnivariate() && fdivides(G, F)) return G/Lc(G);
3142 if (F == G) return F/Lc(F);
3143
3144 CFMap M,N;
3145 int best_level= myCompress (A, B, M, N, topLevel);
3146
3147 if (best_level == 0) return B.genOne();
3148
3149 A= M(A);
3150 B= M(B);
3151
3152 Variable x= Variable (1);
3153
3154 //univariate case
3155 if (A.isUnivariate() && B.isUnivariate())
3156 return N (gcd (A, B));
3157
3158 CanonicalForm cA, cB; // content of A and B
3159 CanonicalForm ppA, ppB; // primitive part of A and B
3161
3162 cA = uni_content (A);
3163 cB = uni_content (B);
3164 gcdcAcB= gcd (cA, cB);
3165 ppA= A/cA;
3166 ppB= B/cB;
3167
3168 CanonicalForm lcA, lcB; // leading coefficients of A and B
3170 lcA= uni_lcoeff (ppA);
3171 lcB= uni_lcoeff (ppB);
3172
3173 if (fdivides (lcA, lcB))
3174 {
3175 if (fdivides (A, B))
3176 return F/Lc(F);
3177 }
3178 if (fdivides (lcB, lcA))
3179 {
3180 if (fdivides (B, A))
3181 return G/Lc(G);
3182 }
3183
3184 gcdlcAlcB= gcd (lcA, lcB);
3185
3186 int d= totaldegree (ppA, Variable (2), Variable (ppA.level()));
3187 int d0;
3188
3189 if (d == 0)
3190 return N(gcdcAcB);
3191 d0= totaldegree (ppB, Variable (2), Variable (ppB.level()));
3192
3193 if (d0 < d)
3194 d= d0;
3195
3196 if (d == 0)
3197 return N(gcdcAcB);
3198
3201 m= gcdlcAlcB;
3202 G_m= 0;
3203 H= 0;
3204 bool fail= false;
3205 topLevel= false;
3206 bool inextension= false;
3210 CFList source, dest;
3211 do // first do
3212 {
3214 if (random_element == 0 && !fail)
3215 {
3216 if (!find (l, random_element))
3218 continue;
3219 }
3220 if (fail)
3221 {
3222 source= CFList();
3223 dest= CFList();
3224
3227 bool prim_fail= false;
3230 if (V_buf4 == alpha)
3232
3233 if (V_buf3 != V_buf4)
3234 {
3238 source, dest);
3242 dest);
3243 for (CFListIterator i= l; i.hasItem(); i++)
3244 i.getItem()= mapDown (i.getItem(), prim_elem, im_prim_elem, V_buf4,
3245 source, dest);
3246 }
3247
3248 ASSERT (!prim_fail, "failure in integer factorizer");
3249 if (prim_fail)
3250 ; //ERROR
3251 else
3253
3254 if (V_buf4 == alpha)
3256 else
3258 im_prim_elem, source, dest);
3259
3260 DEBOUTLN (cerr, "getMipo (V_buf4)= " << getMipo (V_buf4));
3261 DEBOUTLN (cerr, "getMipo (V_buf2)= " << getMipo (V_buf2));
3262 inextension= true;
3263 for (CFListIterator i= l; i.hasItem(); i++)
3264 i.getItem()= mapUp (i.getItem(), V_buf4, V_buf, prim_elem,
3265 im_prim_elem, source, dest);
3269 source, dest);
3273 source, dest);
3274
3275 fail= false;
3277 DEBOUTLN (cerr, "fail= " << fail);
3278 CFList list;
3282 list, topLevel);
3284 "time for recursive call: ");
3285 DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3286 V_buf4= V_buf;
3287 }
3288 else
3289 {
3290 CFList list;
3294 list, topLevel);
3296 "time for recursive call: ");
3297 DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3298 }
3299
3300 if (!G_random_element.inCoeffDomain())
3302 Variable (G_random_element.level()));
3303 else
3304 d0= 0;
3305
3306 if (d0 == 0)
3307 {
3308 if (inextension)
3309 prune1 (alpha);
3310 return N(gcdcAcB);
3311 }
3312 if (d0 > d)
3313 {
3314 if (!find (l, random_element))
3316 continue;
3317 }
3318
3322
3324 if (!G_random_element.inCoeffDomain())
3326 Variable (G_random_element.level()));
3327 else
3328 d0= 0;
3329
3330 if (d0 < d)
3331 {
3332 m= gcdlcAlcB;
3333 newtonPoly= 1;
3334 G_m= 0;
3335 d= d0;
3336 }
3337
3339 if (uni_lcoeff (H) == gcdlcAlcB)
3340 {
3341 cH= uni_content (H);
3342 ppH= H/cH;
3343 if (inextension)
3344 {
3345 CFList u, v;
3346 //maybe it's better to test if ppH is an element of F(\alpha) before
3347 //mapping down
3348 if (fdivides (ppH, ppA) && fdivides (ppH, ppB))
3349 {
3350 DEBOUTLN (cerr, "ppH before mapDown= " << ppH);
3352 ppH /= Lc(ppH);
3353 DEBOUTLN (cerr, "ppH after mapDown= " << ppH);
3354 prune1 (alpha);
3355 return N(gcdcAcB*ppH);
3356 }
3357 }
3358 else if (fdivides (ppH, ppA) && fdivides (ppH, ppB))
3359 return N(gcdcAcB*ppH);
3360 }
3361 G_m= H;
3363 m= m*(x - random_element);
3364 if (!find (l, random_element))
3366
3367 if (getCharacteristic () > 3 && size (skeleton) < 100)
3368 {
3371 do //second do
3372 {
3374 if (random_element == 0 && !fail)
3375 {
3376 if (!find (l, random_element))
3378 continue;
3379 }
3380 if (fail)
3381 {
3382 source= CFList();
3383 dest= CFList();
3384
3387 bool prim_fail= false;
3390 if (V_buf4 == alpha)
3392
3393 if (V_buf3 != V_buf4)
3394 {
3398 source, dest);
3402 source, dest);
3403 for (CFListIterator i= l; i.hasItem(); i++)
3404 i.getItem()= mapDown (i.getItem(), prim_elem, im_prim_elem, V_buf4,
3405 source, dest);
3406 }
3407
3408 ASSERT (!prim_fail, "failure in integer factorizer");
3409 if (prim_fail)
3410 ; //ERROR
3411 else
3413
3414 if (V_buf4 == alpha)
3416 else
3418 prim_elem, im_prim_elem, source, dest);
3419
3420 DEBOUTLN (cerr, "getMipo (V_buf4)= " << getMipo (V_buf4));
3421 DEBOUTLN (cerr, "getMipo (V_buf2)= " << getMipo (V_buf2));
3422 inextension= true;
3423 for (CFListIterator i= l; i.hasItem(); i++)
3424 i.getItem()= mapUp (i.getItem(), V_buf4, V_buf, prim_elem,
3425 im_prim_elem, source, dest);
3429 source, dest);
3432
3434 source, dest);
3435
3436 fail= false;
3438 DEBOUTLN (cerr, "fail= " << fail);
3439 CFList list;
3441
3442 V_buf4= V_buf;
3443
3444 //sparseInterpolation
3445 bool sparseFail= false;
3446 if (LC (skeleton).inCoeffDomain())
3450 else
3454 Monoms);
3455 if (sparseFail)
3456 break;
3457
3459 "time for recursive call: ");
3460 DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3461 }
3462 else
3463 {
3464 CFList list;
3466 bool sparseFail= false;
3467 if (LC (skeleton).inCoeffDomain())
3471 else
3475 Monoms);
3476 if (sparseFail)
3477 break;
3478
3480 "time for recursive call: ");
3481 DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3482 }
3483
3484 if (!G_random_element.inCoeffDomain())
3486 Variable (G_random_element.level()));
3487 else
3488 d0= 0;
3489
3490 if (d0 == 0)
3491 {
3492 if (inextension)
3493 prune1 (alpha);
3494 return N(gcdcAcB);
3495 }
3496 if (d0 > d)
3497 {
3498 if (!find (l, random_element))
3500 continue;
3501 }
3502
3506
3507 if (!G_random_element.inCoeffDomain())
3509 Variable (G_random_element.level()));
3510 else
3511 d0= 0;
3512
3513 if (d0 < d)
3514 {
3515 m= gcdlcAlcB;
3516 newtonPoly= 1;
3517 G_m= 0;
3518 d= d0;
3519 }
3520
3524 "time for newton interpolation: ");
3525
3526 //termination test
3527 if (uni_lcoeff (H) == gcdlcAlcB)
3528 {
3529 cH= uni_content (H);
3530 ppH= H/cH;
3531 if (inextension)
3532 {
3533 CFList u, v;
3534 //maybe it's better to test if ppH is an element of F(\alpha) before
3535 //mapping down
3536 if (fdivides (ppH, ppA) && fdivides (ppH, ppB))
3537 {
3538 DEBOUTLN (cerr, "ppH before mapDown= " << ppH);
3540 ppH /= Lc(ppH);
3541 DEBOUTLN (cerr, "ppH after mapDown= " << ppH);
3542 prune1 (alpha);
3543 return N(gcdcAcB*ppH);
3544 }
3545 }
3546 else if (fdivides (ppH, ppA) && fdivides (ppH, ppB))
3547 {
3548 return N(gcdcAcB*ppH);
3549 }
3550 }
3551
3552 G_m= H;
3554 m= m*(x - random_element);
3555 if (!find (l, random_element))
3557
3558 } while (1);
3559 }
3560 } while (1);
3561}
3562
3564 bool& topLevel, CFList& l)
3565{
3566 CanonicalForm A= F;
3567 CanonicalForm B= G;
3568 if (F.isZero() && degree(G) > 0) return G/Lc(G);
3569 else if (G.isZero() && degree (F) > 0) return F/Lc(F);
3570 else if (F.isZero() && G.isZero()) return F.genOne();
3571 if (F.inBaseDomain() || G.inBaseDomain()) return F.genOne();
3572 if (F.isUnivariate() && fdivides(F, G)) return F/Lc(F);
3573 if (G.isUnivariate() && fdivides(G, F)) return G/Lc(G);
3574 if (F == G) return F/Lc(F);
3575
3576 CFMap M,N;
3577 int best_level= myCompress (A, B, M, N, topLevel);
3578
3579 if (best_level == 0) return B.genOne();
3580
3581 A= M(A);
3582 B= M(B);
3583
3584 Variable x= Variable (1);
3585
3586 //univariate case
3587 if (A.isUnivariate() && B.isUnivariate())
3588 return N (gcd (A, B));
3589
3590 CanonicalForm cA, cB; // content of A and B
3591 CanonicalForm ppA, ppB; // primitive part of A and B
3593
3594 cA = uni_content (A);
3595 cB = uni_content (B);
3596 gcdcAcB= gcd (cA, cB);
3597 ppA= A/cA;
3598 ppB= B/cB;
3599
3600 CanonicalForm lcA, lcB; // leading coefficients of A and B
3602 lcA= uni_lcoeff (ppA);
3603 lcB= uni_lcoeff (ppB);
3604
3605 if (fdivides (lcA, lcB))
3606 {
3607 if (fdivides (A, B))
3608 return F/Lc(F);
3609 }
3610 if (fdivides (lcB, lcA))
3611 {
3612 if (fdivides (B, A))
3613 return G/Lc(G);
3614 }
3615
3616 gcdlcAlcB= gcd (lcA, lcB);
3617
3618 int d= totaldegree (ppA, Variable (2), Variable (ppA.level()));
3619 int d0;
3620
3621 if (d == 0)
3622 return N(gcdcAcB);
3623
3624 d0= totaldegree (ppB, Variable (2), Variable (ppB.level()));
3625
3626 if (d0 < d)
3627 d= d0;
3628
3629 if (d == 0)
3630 return N(gcdcAcB);
3631
3634 m= gcdlcAlcB;
3635 G_m= 0;
3636 H= 0;
3637 bool fail= false;
3638 topLevel= false;
3639 bool inextension= false;
3642 CFList source, dest;
3643 do //first do
3644 {
3645 if (inextension)
3647 else
3649 if (random_element == 0 && !fail)
3650 {
3651 if (!find (l, random_element))
3653 continue;
3654 }
3655
3656 if (!fail && !inextension)
3657 {
3658 CFList list;
3662 list);
3664 "time for recursive call: ");
3665 DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3666 }
3667 else if (!fail && inextension)
3668 {
3669 CFList list;
3673 list, topLevel);
3675 "time for recursive call: ");
3676 DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3677 }
3678 else if (fail && !inextension)
3679 {
3680 source= CFList();
3681 dest= CFList();
3682 CFList list;
3684 int deg= 2;
3685 bool initialized= false;
3686 do
3687 {
3688 mipo= randomIrredpoly (deg, x);
3689 if (initialized)
3690 setMipo (alpha, mipo);
3691 else
3692 alpha= rootOf (mipo);
3693 inextension= true;
3694 fail= false;
3695 initialized= true;
3697 deg++;
3698 } while (fail);
3699 cleanUp= alpha;
3700 V_buf= alpha;
3701 list= CFList();
3705 list, topLevel);
3707 "time for recursive call: ");
3708 DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3709 }
3710 else if (fail && inextension)
3711 {
3712 source= CFList();
3713 dest= CFList();
3714
3717 bool prim_fail= false;
3721
3722 if (V_buf3 != alpha)
3723 {
3727 dest);
3731 dest);
3732 for (CFListIterator i= l; i.hasItem(); i++)
3733 i.getItem()= mapDown (i.getItem(), prim_elem, im_prim_elem, alpha,
3734 source, dest);
3735 }
3736
3737 ASSERT (!prim_fail, "failure in integer factorizer");
3738 if (prim_fail)
3739 ; //ERROR
3740 else
3742
3743 DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (alpha));
3744 DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (V_buf2));
3745
3746 for (CFListIterator i= l; i.hasItem(); i++)
3747 i.getItem()= mapUp (i.getItem(), alpha, V_buf, prim_elem,
3748 im_prim_elem, source, dest);
3752 source, dest);
3756 source, dest);
3757 fail= false;
3759 DEBOUTLN (cerr, "fail= " << fail);
3760 CFList list;
3764 list, topLevel);
3766 "time for recursive call: ");
3767 DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3768 alpha= V_buf;
3769 }
3770
3771 if (!G_random_element.inCoeffDomain())
3773 Variable (G_random_element.level()));
3774 else
3775 d0= 0;
3776
3777 if (d0 == 0)
3778 {
3779 if (inextension)
3780 prune (cleanUp);
3781 return N(gcdcAcB);
3782 }
3783 if (d0 > d)
3784 {
3785 if (!find (l, random_element))
3787 continue;
3788 }
3789
3793
3795
3796 if (!G_random_element.inCoeffDomain())
3798 Variable (G_random_element.level()));
3799 else
3800 d0= 0;
3801
3802 if (d0 < d)
3803 {
3804 m= gcdlcAlcB;
3805 newtonPoly= 1;
3806 G_m= 0;
3807 d= d0;
3808 }
3809
3811
3812 if (uni_lcoeff (H) == gcdlcAlcB)
3813 {
3814 cH= uni_content (H);
3815 ppH= H/cH;
3816 ppH /= Lc (ppH);
3817 DEBOUTLN (cerr, "ppH= " << ppH);
3818
3819 if (fdivides (ppH, ppA) && fdivides (ppH, ppB))
3820 {
3821 if (inextension)
3822 prune (cleanUp);
3823 return N(gcdcAcB*ppH);
3824 }
3825 }
3826 G_m= H;
3828 m= m*(x - random_element);
3829 if (!find (l, random_element))
3831
3832 if ((getCharacteristic() > 3 && size (skeleton) < 200))
3833 {
3836
3837 do //second do
3838 {
3839 if (inextension)
3841 else
3843 if (random_element == 0 && !fail)
3844 {
3845 if (!find (l, random_element))
3847 continue;
3848 }
3849
3850 bool sparseFail= false;
3851 if (!fail && !inextension)
3852 {
3853 CFList list;
3855 if (LC (skeleton).inCoeffDomain())
3859 Monoms);
3860 else
3866 "time for recursive call: ");
3867 DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3868 }
3869 else if (!fail && inextension)
3870 {
3871 CFList list;
3873 if (LC (skeleton).inCoeffDomain())
3877 Monoms);
3878 else
3882 Monoms);
3884 "time for recursive call: ");
3885 DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3886 }
3887 else if (fail && !inextension)
3888 {
3889 source= CFList();
3890 dest= CFList();
3891 CFList list;
3893 int deg= 2;
3894 bool initialized= false;
3895 do
3896 {
3897 mipo= randomIrredpoly (deg, x);
3898 if (initialized)
3899 setMipo (alpha, mipo);
3900 else
3901 alpha= rootOf (mipo);
3902 inextension= true;
3903 fail= false;
3904 initialized= true;
3906 deg++;
3907 } while (fail);
3908 cleanUp= alpha;
3909 V_buf= alpha;
3910 list= CFList();
3912 if (LC (skeleton).inCoeffDomain())
3916 Monoms);
3917 else
3921 Monoms);
3923 "time for recursive call: ");
3924 DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3925 }
3926 else if (fail && inextension)
3927 {
3928 source= CFList();
3929 dest= CFList();
3930
3933 bool prim_fail= false;
3937
3938 if (V_buf3 != alpha)
3939 {
3943 source, dest);
3947 source, dest);
3948 for (CFListIterator i= l; i.hasItem(); i++)
3949 i.getItem()= mapDown (i.getItem(), prim_elem, im_prim_elem, alpha,
3950 source, dest);
3951 }
3952
3953 ASSERT (!prim_fail, "failure in integer factorizer");
3954 if (prim_fail)
3955 ; //ERROR
3956 else
3958
3959 DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (alpha));
3960 DEBOUTLN (cerr, "getMipo (alpha)= " << getMipo (V_buf2));
3961
3962 for (CFListIterator i= l; i.hasItem(); i++)
3963 i.getItem()= mapUp (i.getItem(), alpha, V_buf, prim_elem,
3964 im_prim_elem, source, dest);
3968 source, dest);
3972 source, dest);
3973 fail= false;
3975 DEBOUTLN (cerr, "fail= " << fail);
3976 CFList list;
3978 if (LC (skeleton).inCoeffDomain())
3982 Monoms);
3983 else
3989 "time for recursive call: ");
3990 DEBOUTLN (cerr, "G_random_element= " << G_random_element);
3991 alpha= V_buf;
3992 }
3993
3994 if (sparseFail)
3995 break;
3996
3997 if (!G_random_element.inCoeffDomain())
3999 Variable (G_random_element.level()));
4000 else
4001 d0= 0;
4002
4003 if (d0 == 0)
4004 {
4005 if (inextension)
4006 prune (cleanUp);
4007 return N(gcdcAcB);
4008 }
4009 if (d0 > d)
4010 {
4011 if (!find (l, random_element))
4013 continue;
4014 }
4015
4019
4020 if (!G_random_element.inCoeffDomain())
4022 Variable (G_random_element.level()));
4023 else
4024 d0= 0;
4025
4026 if (d0 < d)
4027 {
4028 m= gcdlcAlcB;
4029 newtonPoly= 1;
4030 G_m= 0;
4031 d= d0;
4032 }
4033
4037 "time for newton interpolation: ");
4038
4039 //termination test
4040 if (uni_lcoeff (H) == gcdlcAlcB)
4041 {
4042 cH= uni_content (H);
4043 ppH= H/cH;
4044 ppH /= Lc (ppH);
4045 DEBOUTLN (cerr, "ppH= " << ppH);
4046 if (fdivides (ppH, ppA) && fdivides (ppH, ppB))
4047 {
4048 if (inextension)
4049 prune (cleanUp);
4050 return N(gcdcAcB*ppH);
4051 }
4052 }
4053
4054 G_m= H;
4056 m= m*(x - random_element);
4057 if (!find (l, random_element))
4059
4060 } while (1); //end of second do
4061 }
4062 else
4063 {
4064 if (inextension)
4065 prune (cleanUp);
4066 return N(gcdcAcB*modGCDFp (ppA, ppB));
4067 }
4068 } while (1); //end of first do
4069}
4070
4073
4074/// modular gcd algorithm, see Keith, Czapor, Geddes "Algorithms for Computer
4075/// Algebra", Algorithm 7.1
4077{
4078 CanonicalForm f, g, cl, q(0), Dp, newD, D, newq, newqh;
4079 int p, i, dp_deg, d_deg=-1;
4080
4085 cf= icontent (f);
4086 f /= cf;
4087 //cd = bCommonDen( f ); f *=cd;
4088 //f /=vcontent(f,Variable(1));
4089
4092 cg= icontent (g);
4093 g /= cg;
4094 //cd = bCommonDen( g ); g *=cd;
4095 //g /=vcontent(g,Variable(1));
4096
4099 lcf= f.lc();
4100 lcg= g.lc();
4101 cl = gcd (f.lc(),g.lc());
4106 for (i= tmax (f.level(), g.level()); i > 0; i--)
4107 {
4108 if (degree (f, i) <= 0 || degree (g, i) <= 0)
4109 continue;
4110 else
4111 {
4112 minCommonDeg= tmin (degree (g, i), degree (f, i));
4113 break;
4114 }
4115 }
4116 if (i == 0)
4117 return gcdcfcg;
4118 for (; i > 0; i--)
4119 {
4120 if (degree (f, i) <= 0 || degree (g, i) <= 0)
4121 continue;
4122 else
4124 }
4125 b= 2*tmin (maxNorm (f), maxNorm (g))*abs (cl)*
4127 bool equal= false;
4128 i = cf_getNumBigPrimes() - 1;
4129
4132 //Off (SW_RATIONAL);
4133 while ( true )
4134 {
4135 p = cf_getBigPrime( i );
4136 i--;
4137 while ( i >= 0 && mod( cl*(lc(f)/cl)*(lc(g)/cl), p ) == 0 )
4138 {
4139 p = cf_getBigPrime( i );
4140 i--;
4141 }
4142 //printf("try p=%d\n",p);
4144 fp= mapinto (f);
4145 gp= mapinto (g);
4147#if defined(HAVE_NTL) || defined(HAVE_FLINT)
4148 if (size (fp)/maxNumVars > 500 && size (gp)/maxNumVars > 500)
4149 Dp = modGCDFp (fp, gp, cofp, cogp);
4150 else
4151 {
4152 Dp= gcd_poly (fp, gp);
4153 cofp= fp/Dp;
4154 cogp= gp/Dp;
4155 }
4156#else
4157 Dp= gcd_poly (fp, gp);
4158 cofp= fp/Dp;
4159 cogp= gp/Dp;
4160#endif
4162 "time for gcd mod p in modular gcd: ");
4163 Dp /=Dp.lc();
4164 Dp *= mapinto (cl);
4165 cofp /= lc (cofp);
4166 cofp *= mapinto (lcf);
4167 cogp /= lc (cogp);
4168 cogp *= mapinto (lcg);
4169 setCharacteristic( 0 );
4171 if ( dp_deg == 0 )
4172 {
4173 //printf(" -> 1\n");
4174 return CanonicalForm(gcdcfcg);
4175 }
4176 if ( q.isZero() )
4177 {
4178 D = mapinto( Dp );
4179 cof= mapinto (cofp);
4180 cog= mapinto (cogp);
4181 d_deg=dp_deg;
4182 q = p;
4183 Dn= balance_p (D, p);
4184 cofn= balance_p (cof, p);
4185 cogn= balance_p (cog, p);
4186 }
4187 else
4188 {
4189 if ( dp_deg == d_deg )
4190 {
4191 chineseRemainder( D, q, mapinto( Dp ), p, newD, newq );
4194 cof= newCof;
4195 cog= newCog;
4196 newqh= newq/2;
4197 Dn= balance_p (newD, newq, newqh);
4200 if (test != Dn) //balance_p (newD, newq))
4201 test= balance_p (newD, newq);
4202 else
4203 equal= true;
4204 q = newq;
4205 D = newD;
4206 }
4207 else if ( dp_deg < d_deg )
4208 {
4209 //printf(" were all bad, try more\n");
4210 // all previous p's are bad primes
4211 q = p;
4212 D = mapinto( Dp );
4213 cof= mapinto (cof);
4214 cog= mapinto (cog);
4215 d_deg=dp_deg;
4216 test= 0;
4217 equal= false;
4218 Dn= balance_p (D, p);
4219 cofn= balance_p (cof, p);
4220 cogn= balance_p (cog, p);
4221 }
4222 else
4223 {
4224 //printf(" was bad, try more\n");
4225 }
4226 //else dp_deg > d_deg: bad prime
4227 }
4228 if ( i >= 0 )
4229 {
4230 cDn= icontent (Dn);
4231 Dn /= cDn;
4232 cofn /= cl/cDn;
4233 //cofn /= icontent (cofn);
4234 cogn /= cl/cDn;
4235 //cogn /= icontent (cogn);
4237 if ((terminationTest (f,g, cofn, cogn, Dn)) ||
4238 ((equal || q > b) && fdivides (Dn, f) && fdivides (Dn, g)))
4239 {
4241 "time for successful termination in modular gcd: ");
4242 //printf(" -> success\n");
4243 return Dn*gcdcfcg;
4244 }
4246 "time for unsuccessful termination in modular gcd: ");
4247 equal= false;
4248 //else: try more primes
4249 }
4250 else
4251 { // try other method
4252 //printf("try other gcd\n");
4254 D=gcd_poly( f, g );
4256 return D*gcdcfcg;
4257 }
4258 }
4259}
4260#endif
CFMatrix * convertNmod_mat_t2FacCFMatrix(const nmod_mat_t m)
conversion of a FLINT matrix over Z/p to a factory matrix
void convertFacCFMatrix2nmod_mat_t(nmod_mat_t M, const CFMatrix &m)
conversion of a factory matrix over Z/p to a nmod_mat_t
void convertFacCFMatrix2Fq_nmod_mat_t(fq_nmod_mat_t M, const fq_nmod_ctx_t fq_con, const CFMatrix &m)
conversion of a factory matrix over F_q to a fq_nmod_mat_t
CanonicalForm convertnmod_poly_t2FacCF(const nmod_poly_t poly, const Variable &x)
conversion of a FLINT poly over Z/p to CanonicalForm
CFMatrix * convertFq_nmod_mat_t2FacCFMatrix(const fq_nmod_mat_t m, const fq_nmod_ctx_t &fq_con, const Variable &alpha)
conversion of a FLINT matrix over F_q to a factory matrix
This file defines functions for conversion to FLINT (www.flintlib.org) and back.
Rational pow(const Rational &a, int e)
Definition GMPrat.cc:411
Rational abs(const Rational &a)
Definition GMPrat.cc:436
CFMatrix * convertNTLmat_zz_p2FacCFMatrix(const mat_zz_p &m)
CFMatrix * convertNTLmat_zz_pE2FacCFMatrix(const mat_zz_pE &m, const Variable &alpha)
CanonicalForm convertNTLzzpX2CF(const zz_pX &poly, const Variable &x)
mat_zz_pE * convertFacCFMatrix2NTLmat_zz_pE(const CFMatrix &m)
zz_pX convertFacCF2NTLzzpX(const CanonicalForm &f)
mat_zz_p * convertFacCFMatrix2NTLmat_zz_p(const CFMatrix &m)
VAR long fac_NTL_char
Definition NTLconvert.cc:46
Conversion to and from NTL.
void On(int sw)
switches
void Off(int sw)
switches
CanonicalForm power(const CanonicalForm &f, int n)
exponentiation
Header for factory's main class CanonicalForm.
CanonicalForm FACTORY_PUBLIC content(const CanonicalForm &)
CanonicalForm content ( const CanonicalForm & f )
Definition cf_gcd.cc:603
CanonicalForm mapinto(const CanonicalForm &f)
int size(const CanonicalForm &f, const Variable &v)
int size ( const CanonicalForm & f, const Variable & v )
Definition cf_ops.cc:600
int getNumVars(const CanonicalForm &f)
int getNumVars ( const CanonicalForm & f )
Definition cf_ops.cc:314
CanonicalForm lc(const CanonicalForm &f)
CanonicalForm FACTORY_PUBLIC icontent(const CanonicalForm &f)
CanonicalForm icontent ( const CanonicalForm & f )
Definition cf_gcd.cc:74
int * degrees(const CanonicalForm &f, int *degs=0)
int * degrees ( const CanonicalForm & f, int * degs )
Definition cf_ops.cc:493
CF_NO_INLINE FACTORY_PUBLIC CanonicalForm mod(const CanonicalForm &, const CanonicalForm &)
int degree(const CanonicalForm &f)
int getGFDegree()
Definition cf_char.cc:75
Array< CanonicalForm > CFArray
void FACTORY_PUBLIC setCharacteristic(int c)
Definition cf_char.cc:28
Matrix< CanonicalForm > CFMatrix
CanonicalForm FACTORY_PUBLIC swapvar(const CanonicalForm &, const Variable &, const Variable &)
swapvar() - swap variables x1 and x2 in f.
Definition cf_ops.cc:168
int totaldegree(const CanonicalForm &f)
int totaldegree ( const CanonicalForm & f )
Definition cf_ops.cc:523
CanonicalForm Lc(const CanonicalForm &f)
CanonicalForm FACTORY_PUBLIC gcd_poly(const CanonicalForm &f, const CanonicalForm &g)
CanonicalForm gcd_poly ( const CanonicalForm & f, const CanonicalForm & g )
Definition cf_gcd.cc:492
List< CanonicalForm > CFList
CanonicalForm LC(const CanonicalForm &f)
int FACTORY_PUBLIC getCharacteristic()
Definition cf_char.cc:70
const CanonicalForm CFMap CFMap & N
Definition cfEzgcd.cc:56
int * degsf
Definition cfEzgcd.cc:59
int f_zero
Definition cfEzgcd.cc:69
int l
Definition cfEzgcd.cc:100
int m
Definition cfEzgcd.cc:128
const CanonicalForm CFMap CFMap int & both_non_zero
Definition cfEzgcd.cc:57
int g_zero
Definition cfEzgcd.cc:70
int k
Definition cfEzgcd.cc:99
int * degsg
Definition cfEzgcd.cc:60
const CanonicalForm CFMap CFMap bool topLevel
int both_zero
coprimality check and change of representation mod n
CanonicalForm sparseGCDFq(const CanonicalForm &F, const CanonicalForm &G, const Variable &alpha, CFList &l, bool &topLevel)
Definition cfModGcd.cc:3131
CFArray readOffSolution(const CFMatrix &M, const long rk)
M in row echolon form, rk rank of M.
Definition cfModGcd.cc:1689
CanonicalForm modGCDFq(const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &coF, CanonicalForm &coG, Variable &alpha, CFList &l, bool &topLevel)
GCD of F and G over , l and topLevel are only used internally, output is monic based on Alg....
Definition cfModGcd.cc:479
static CanonicalForm GFRandomElement(const CanonicalForm &F, CFList &list, bool &fail)
compute a random element a of GF, s.t. F(a) , F is a univariate polynomial, returns fail if there ar...
Definition cfModGcd.cc:820
CFArray evaluateMonom(const CanonicalForm &F, const CFList &evalPoints)
Definition cfModGcd.cc:1993
const CanonicalForm const CanonicalForm const CanonicalForm & coG
Definition cfModGcd.cc:68
CFArray getMonoms(const CanonicalForm &F)
extract monomials of F, parts in algebraic variable are considered coefficients
Definition cfModGcd.cc:1958
void mult(CFList &L1, const CFList &L2)
multiply two lists componentwise
Definition cfModGcd.cc:2177
const CanonicalForm const CanonicalForm const CanonicalForm const CanonicalForm & cand
Definition cfModGcd.cc:70
CanonicalForm cofp
Definition cfModGcd.cc:4130
long gaussianElimFq(CFMatrix &M, CFArray &L, const Variable &alpha)
Definition cfModGcd.cc:1785
int myCompress(const CanonicalForm &F, const CanonicalForm &G, CFMap &M, CFMap &N, bool topLevel)
compressing two polynomials F and G, M is used for compressing, N to reverse the compression
Definition cfModGcd.cc:92
Variable x
Definition cfModGcd.cc:4083
CanonicalForm lcg
Definition cfModGcd.cc:4098
static Variable chooseExtension(const Variable &alpha)
Definition cfModGcd.cc:421
int dp_deg
Definition cfModGcd.cc:4079
CanonicalForm newCog
Definition cfModGcd.cc:4130
CFArray evaluate(const CFArray &A, const CFList &evalPoints)
Definition cfModGcd.cc:2032
CanonicalForm monicSparseInterpol(const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &skeleton, const Variable &alpha, bool &fail, CFArray *&coeffMonoms, CFArray &Monoms)
Definition cfModGcd.cc:2200
long gaussianElimFp(CFMatrix &M, CFArray &L)
Definition cfModGcd.cc:1740
CanonicalForm cogn
Definition cfModGcd.cc:4130
int p
Definition cfModGcd.cc:4079
static CanonicalForm uni_content(const CanonicalForm &F)
compute the content of F, where F is considered as an element of
Definition cfModGcd.cc:282
CFArray solveGeneralVandermonde(const CFArray &M, const CFArray &A)
Definition cfModGcd.cc:1633
cl
Definition cfModGcd.cc:4101
f
Definition cfModGcd.cc:4082
CanonicalForm extractContents(const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &contentF, CanonicalForm &contentG, CanonicalForm &ppF, CanonicalForm &ppG, const int d)
Definition cfModGcd.cc:312
CanonicalForm lcf
Definition cfModGcd.cc:4098
static CanonicalForm uni_lcoeff(const CanonicalForm &F)
compute the leading coefficient of F, where F is considered as an element of , order on is dp.
Definition cfModGcd.cc:340
g
Definition cfModGcd.cc:4091
int maxNumVars
Definition cfModGcd.cc:4131
CanonicalForm fp
Definition cfModGcd.cc:4103
CFArray solveVandermonde(const CFArray &M, const CFArray &A)
Definition cfModGcd.cc:1580
int i
Definition cfModGcd.cc:4079
CanonicalForm cogp
Definition cfModGcd.cc:4130
const CanonicalForm const CanonicalForm & coF
Definition cfModGcd.cc:68
CanonicalForm cofn
Definition cfModGcd.cc:4130
CanonicalForm cof
Definition cfModGcd.cc:4130
const CanonicalForm & GG
Definition cfModGcd.cc:4077
CanonicalForm cd(bCommonDen(FF))
Definition cfModGcd.cc:4090
CanonicalForm cog
Definition cfModGcd.cc:4130
CanonicalForm nonMonicSparseInterpol(const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &skeleton, const Variable &alpha, bool &fail, CFArray *&coeffMonoms, CFArray &Monoms)
Definition cfModGcd.cc:2484
int minCommonDeg
Definition cfModGcd.cc:4105
CanonicalForm modGCDFp(const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &coF, CanonicalForm &coG, bool &topLevel, CFList &l)
Definition cfModGcd.cc:1224
CFArray solveSystemFp(const CFMatrix &M, const CFArray &L)
Definition cfModGcd.cc:1841
const CanonicalForm & G
Definition cfModGcd.cc:67
CanonicalForm sparseGCDFp(const CanonicalForm &F, const CanonicalForm &G, bool &topLevel, CFList &l)
Definition cfModGcd.cc:3563
static CanonicalForm randomElement(const CanonicalForm &F, const Variable &alpha, CFList &list, bool &fail)
compute a random element a of , s.t. F(a) , F is a univariate polynomial, returns fail if there are...
Definition cfModGcd.cc:380
CanonicalForm gcdcfcg
Definition cfModGcd.cc:4102
CanonicalForm cf
Definition cfModGcd.cc:4084
CanonicalForm Dn
Definition cfModGcd.cc:4097
CanonicalForm newCof
Definition cfModGcd.cc:4130
static CanonicalForm FpRandomElement(const CanonicalForm &F, CFList &list, bool &fail)
Definition cfModGcd.cc:1172
CanonicalForm gp
Definition cfModGcd.cc:4103
CanonicalForm modGCDGF(const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &coF, CanonicalForm &coG, CFList &l, bool &topLevel)
GCD of F and G over GF, based on Alg. 7.2. as described in "Algorithms for Computer Algebra" by Gedde...
Definition cfModGcd.cc:873
static CanonicalForm newtonInterp(const CanonicalForm &alpha, const CanonicalForm &u, const CanonicalForm &newtonPoly, const CanonicalForm &oldInterPoly, const Variable &x)
Newton interpolation - Incremental algorithm. Given a list of values alpha_i and a list of polynomial...
Definition cfModGcd.cc:365
bool equal
Definition cfModGcd.cc:4127
CFArray solveSystemFq(const CFMatrix &M, const CFArray &L, const Variable &alpha)
Definition cfModGcd.cc:1893
CanonicalForm test
Definition cfModGcd.cc:4097
CanonicalForm b
Definition cfModGcd.cc:4104
CanonicalForm cg
Definition cfModGcd.cc:4084
CanonicalForm cDn
Definition cfModGcd.cc:4130
int d_deg
Definition cfModGcd.cc:4079
CFList evaluationPoints(const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &Feval, CanonicalForm &Geval, const CanonicalForm &LCF, const bool &GF, const Variable &alpha, bool &fail, CFList &list)
Definition cfModGcd.cc:2049
modular and sparse modular GCD algorithms over finite fields and Z.
bool terminationTest(const CanonicalForm &F, const CanonicalForm &G, const CanonicalForm &coF, const CanonicalForm &coG, const CanonicalForm &cand)
CanonicalForm modGCDZ(const CanonicalForm &FF, const CanonicalForm &GG)
modular GCD over Z
static void evalPoint(const CanonicalForm &F, const CanonicalForm &G, CanonicalForm &FEval, CanonicalForm &GEval, CFGenerator &evalPoint)
This file provides functions to compute the Newton polygon of a bivariate polynomial.
CanonicalForm bCommonDen(const CanonicalForm &f)
CanonicalForm bCommonDen ( const CanonicalForm & f )
CanonicalForm maxNorm(const CanonicalForm &f)
CanonicalForm maxNorm ( const CanonicalForm & f )
bool fdivides(const CanonicalForm &f, const CanonicalForm &g)
bool fdivides ( const CanonicalForm & f, const CanonicalForm & g )
declarations of higher level algorithms.
void FACTORY_PUBLIC chineseRemainder(const CanonicalForm &x1, const CanonicalForm &q1, const CanonicalForm &x2, const CanonicalForm &q2, CanonicalForm &xnew, CanonicalForm &qnew)
void chineseRemainder ( const CanonicalForm & x1, const CanonicalForm & q1, const CanonicalForm & x2,...
Definition cf_chinese.cc:57
assertions for Factory
#define ASSERT(expression, message)
Definition cf_assert.h:99
static const int SW_USE_CHINREM_GCD
set to 1 to use modular gcd over Z
Definition cf_defs.h:41
#define DELETE_ARRAY(P)
Definition cf_defs.h:65
#define NEW_ARRAY(T, N)
Definition cf_defs.h:64
static CanonicalForm balance_p(const CanonicalForm &f, const CanonicalForm &q, const CanonicalForm &qh)
Definition cf_gcd.cc:149
CanonicalForm randomIrredpoly(int i, const Variable &x)
computes a random monic irreducible univariate polynomial in x over Fp of degree i via NTL/FLINT
Definition cf_irred.cc:26
generate random irreducible univariate polynomials
Iterators for CanonicalForm's.
static CanonicalForm bound(const CFMatrix &M)
Definition cf_linsys.cc:460
map polynomials
CanonicalForm mapPrimElem(const CanonicalForm &primElem, const Variable &alpha, const Variable &beta)
compute the image of a primitive element of in . We assume .
CanonicalForm GFMapDown(const CanonicalForm &F, int k)
maps a polynomial over to a polynomial over , d needs to be a multiple of k
CanonicalForm primitiveElement(const Variable &alpha, Variable &beta, bool &fail)
determine a primitive element of , is a primitive element of a field which is isomorphic to
static CanonicalForm mapDown(const CanonicalForm &F, const Variable &alpha, const CanonicalForm &G, CFList &source, CFList &dest)
the CanonicalForm G is the output of map_up, returns F considered as an element over ,...
static CanonicalForm mapUp(const Variable &alpha, const Variable &beta)
and is a primitive element, returns the image of
Definition cf_map_ext.cc:71
CanonicalForm GFMapUp(const CanonicalForm &F, int k)
maps a polynomial over to a polynomial over , d needs to be a multiple of k
This file implements functions to map between extensions of finite fields.
int cf_getBigPrime(int i)
Definition cf_primes.cc:39
int cf_getNumBigPrimes()
Definition cf_primes.cc:45
access to prime tables
GLOBAL_VAR flint_rand_t FLINTrandom
Definition cf_random.cc:25
generate random integers, random elements of finite fields
generate random evaluation points
VAR void(* factoryError)(const char *s)
Definition cf_util.cc:80
int ipower(int b, int m)
int ipower ( int b, int m )
Definition cf_util.cc:27
generate random elements in F_p(alpha)
Definition cf_random.h:70
int size() const
class to iterate through CanonicalForm's
Definition cf_iter.h:44
class CFMap
Definition cf_map.h:85
factory's main class
CanonicalForm genOne() const
CF_NO_INLINE bool isZero() const
Variable mvar() const
mvar() returns the main variable of CO or Variable() if CO is in a base domain.
CF_NO_INLINE bool isOne() const
bool inCoeffDomain() const
int level() const
level() returns the level of CO.
bool inBaseDomain() const
CanonicalForm lc() const
CanonicalForm CanonicalForm::lc (), Lc (), LC (), LC ( v ) const.
bool isUnivariate() const
generate random elements in F_p
Definition cf_random.h:44
generate random elements in GF
Definition cf_random.h:32
int length() const
void append(const T &)
T getLast() const
void removeLast()
factory's class for variables
Definition factory.h:127
int level() const
Definition factory.h:143
functions to print debug output
#define DEBOUTLN(stream, objects)
Definition debug.h:49
CFFListIterator iter
Variable alpha
return result
CanonicalForm Feval
Definition facAbsFact.cc:60
CanonicalForm H
Definition facAbsFact.cc:60
CanonicalForm mipo
Definition facAlgExt.cc:57
b *CanonicalForm B
Definition facBivar.cc:52
const Variable & v
< [in] a sqrfree bivariate poly
Definition facBivar.h:39
CanonicalForm LCF
CFList & eval
CFList *& Aeval
<[in] poly
CFList tmp1
Definition facFqBivar.cc:75
CFList tmp2
Definition facFqBivar.cc:75
CFList evalPoints(const CanonicalForm &F, CFList &eval, const Variable &alpha, CFList &list, const bool &GF, bool &fail)
evaluation point search for multivariate factorization, looks for a (F.level() - 1)-tuple such that t...
int j
Definition facHensel.cc:110
fq_nmod_ctx_clear(fq_con)
nmod_poly_init(FLINTmipo, getCharacteristic())
fq_nmod_ctx_init_modulus(fq_con, FLINTmipo, "Z")
convertFacCF2nmod_poly_t(FLINTmipo, M)
nmod_poly_clear(FLINTmipo)
This file defines functions for Hensel lifting.
Variable FACTORY_PUBLIC rootOf(const CanonicalForm &, char name='@')
returns a symbolic root of polynomial with name name Use it to define algebraic variables
Definition variable.cc:162
void prune1(const Variable &alpha)
Definition variable.cc:291
CanonicalForm getMipo(const Variable &alpha, const Variable &x)
Definition variable.cc:207
void FACTORY_PUBLIC prune(Variable &alpha)
Definition variable.cc:261
void setMipo(const Variable &alpha, const CanonicalForm &mipo)
Definition variable.cc:219
#define const
Definition fegetopt.c:39
some useful template functions.
template CanonicalForm tmax(const CanonicalForm &, const CanonicalForm &)
template CanonicalForm tmin(const CanonicalForm &, const CanonicalForm &)
template bool find(const List< CanonicalForm > &, const CanonicalForm &)
#define D(A)
Definition gentable.cc:131
VAR char gf_name
Definition gfops.cc:52
gmp_float log(const gmp_float &a)
#define NULL
Definition omList.c:12
int status int void * buf
Definition si_signals.h:59
#define A
Definition sirandom.c:24
#define M
Definition sirandom.c:25
#define TIMING_DEFINE_PRINT(t)
Definition timing.h:95
#define TIMING_START(t)
Definition timing.h:92
#define TIMING_END_AND_PRINT(t, msg)
Definition timing.h:94
int gcd(int a, int b)