Context Navigation

source: git/factory/facAbsBiFact.cc @ 4946e3

spielwiese

Last change on this file since 4946e3 was 4946e3, checked in by Hans Schoenemann <hannes@…>, 3 years ago
factory: LL via FLINT
Property mode set to `100644`
File size: 21.5 KB

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

Note: See TracBrowser for help on using the repository browser.

Download in other formats: