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